Skip to content
Snippets Groups Projects
Commit 295fe9b8 authored by Hugo  VARET's avatar Hugo VARET
Browse files

fixedBugs

parent e982f643
Branches
No related tags found
No related merge requests found
...@@ -5,7 +5,7 @@ ...@@ -5,7 +5,7 @@
#' @name countsBoxplots #' @name countsBoxplots
#' @param results results produced by DESeq2 or limma with ChIPuana #' @param results results produced by DESeq2 or limma with ChIPuana
#' @param conditions factor vector of the condition from which each sample belongs #' @param conditions factor vector of the condition from which each sample belongs
#' @param method either "DEseq2" or "Limma" #' @param method either "DESeq2" or "Limma"
#' @param col colors of the boxplots (one per biological condition) #' @param col colors of the boxplots (one per biological condition)
#' @param outfile TRUE to export the figure in a png file #' @param outfile TRUE to export the figure in a png file
#' @return A boxplots of the raw and normalized counts #' @return A boxplots of the raw and normalized counts
...@@ -38,14 +38,14 @@ countsBoxplots <- function(results, conditions, method, ...@@ -38,14 +38,14 @@ countsBoxplots <- function(results, conditions, method,
theme(axis.text.x=element_text(angle=90, hjust=1, vjust=0.5)) theme(axis.text.x=element_text(angle=90, hjust=1, vjust=0.5))
# norm counts # norm counts
if(method=="DEseq2") { if(method=="DESeq2") {
norm.count <- log2(norm.counts+1) norm.count <- log2(norm.counts+1)
} }
norm.count <- removeNull(norm.counts) norm.count <- removeNull(norm.counts)
d <- stack(as.data.frame(norm.counts)) d <- stack(as.data.frame(norm.counts))
d$conditions <- rep(conditions, each=nrow(norm.counts)) d$conditions <- rep(conditions, each=nrow(norm.counts))
p2 <- ggplot(d) + p2 <- ggplot(d) +
geom_boxplot(aes(x=.data$ind, y=.data$values+1, fill=.data$conditions), show.legend=TRUE) + geom_boxplot(aes(x=gsub("norm.", "", .data$ind), y=.data$values+1, fill=.data$conditions), show.legend=TRUE) +
labs(fill="") + labs(fill="") +
scale_y_continuous(trans = log10_trans(), scale_y_continuous(trans = log10_trans(),
breaks = trans_breaks("log10", function(x) 10^x), breaks = trans_breaks("log10", function(x) 10^x),
......
...@@ -49,6 +49,6 @@ run.DESeq2 <- function(counts, conditions, batch=NULL, pAdjustMethod="BH", spike ...@@ -49,6 +49,6 @@ run.DESeq2 <- function(counts, conditions, batch=NULL, pAdjustMethod="BH", spike
results[[paste0(levelTest,"_vs_",levelRef)]] <- resAnDif results[[paste0(levelTest,"_vs_",levelRef)]] <- resAnDif
cat(paste("Comparison", levelTest, "vs", levelRef, "done\n")) cat(paste("Comparison", levelTest, "vs", levelRef, "done\n"))
} }
return(list(dds=dds, results=results, sf=sizeFactors(dds))) return(list(dds=dds, results=results, sf=sizeFactors(dds), counts.trans=assay(vst(dds))))
} }
...@@ -29,25 +29,27 @@ counts <- as.matrix(data[,-c(1:6)]) ...@@ -29,25 +29,27 @@ counts <- as.matrix(data[,-c(1:6)])
rownames(counts) <- data[,1] rownames(counts) <- data[,1]
``` ```
```{r echo=FALSE, results="asis", results="hide"} ```{r echo=FALSE, results="asis", results="hide", message=FALSE}
############### LIMMA ######################################## ############### LIMMA ########################################
if (method=="Limma") { if (method=="Limma") {
library(limma) library(limma)
if (is.null(spikes)){ if (is.null(spikes)){
# Diferential analysis # Diferential analysis
resAnDif <- run.Limma(counts = counts, conditions = Conditions, resAnDif <- run.Limma(counts = counts, conditions = Conditions,
normalize.method = normalisation, pAdjustMethod = pAdjustMethod, outfile = FALSE) normalize.method = normalisation,
batch=batch, pAdjustMethod = pAdjustMethod, outfile = FALSE)
} else { } else {
# Diferential analysis # Diferential analysis
resAnDif <- run.Limma(counts = counts, conditions = Conditions, resAnDif <- run.Limma(counts = counts, conditions = Conditions,
normalize.method = "none", pAdjustMethod = pAdjustMethod, spikes = spikes, outfile = FALSE) normalize.method = "none", pAdjustMethod = pAdjustMethod,
batch = batch, spikes = spikes, outfile = FALSE)
# Normalisation factor # Normalisation factor
sf <- as.data.frame(spikes/max(spikes)) sf <- as.data.frame(spikes/max(spikes))
colnames(sf) <- "Size factor" colnames(sf) <- "Size factor"
} }
counts.trans <- resAnDif$voom$E
############### DESeq2 ######################################## ############### DESeq2 ########################################
} else if(method=="DEseq2") { } else if(method=="DESeq2") {
library(DESeq2) library(DESeq2)
if (is.null(spikes)){ if (is.null(spikes)){
normalisation = "geometrical" normalisation = "geometrical"
...@@ -65,6 +67,7 @@ if (method=="Limma") { ...@@ -65,6 +67,7 @@ if (method=="Limma") {
sf <- as.data.frame(spikes/max(spikes)) sf <- as.data.frame(spikes/max(spikes))
colnames(sf) <- "Size factor" colnames(sf) <- "Size factor"
} }
counts.trans <- resAnDif$counts.trans
} }
# Export result table # Export result table
...@@ -135,14 +138,14 @@ Figure 3 sample clustering based on normalized data. An euclidean distance is co ...@@ -135,14 +138,14 @@ Figure 3 sample clustering based on normalized data. An euclidean distance is co
```{r clusterplot, echo=FALSE,fig.align="center",fig.cap="Figure 3: Sample clustering based on normalized data.", out.width="600px", warning=FALSE, message=FALSE} ```{r clusterplot, echo=FALSE,fig.align="center",fig.cap="Figure 3: Sample clustering based on normalized data.", out.width="600px", warning=FALSE, message=FALSE}
# Cluster plot # Cluster plot
clusterPlot(counts.trans = resAnDif$voom$E,conditions = Conditions, outfile = FALSE) clusterPlot(counts.trans = counts.trans,conditions = Conditions, outfile = FALSE)
``` ```
Another way of visualizing the experiment variability is to look at the first principal components of the PCA, as shown on the figure 4. On this figure, the first principal component (PC1) is expected to separate samples from the different biological conditions, meaning that the biological variability is the main source of variance in the data. Another way of visualizing the experiment variability is to look at the first principal components of the PCA, as shown on the figure 4. On this figure, the first principal component (PC1) is expected to separate samples from the different biological conditions, meaning that the biological variability is the main source of variance in the data.
```{r PCA, echo=FALSE,fig.align="center",fig.cap="Figure 4: First three components of a Principal Component Analysis, with percentages of variance associated with each axis.", out.width="1200px"} ```{r PCA, echo=FALSE,fig.align="center",fig.cap="Figure 4: First three components of a Principal Component Analysis, with percentages of variance associated with each axis.", out.width="1200px"}
# PCA plot # PCA plot
PCAPlot(counts.trans = resAnDif$voom$E, conditions = Conditions, outfile = FALSE) PCAPlot(counts.trans = counts.trans, conditions = Conditions, outfile = FALSE)
``` ```
......
file <- "K4Me3_MOUSE_CuU_Peak.mtx" file <- "K4Me3_MOUSE_CuU_Peak.mtx"
Conditions <- c("C","C","U","U") Conditions <- c("C","C","U","U")
Replicates <- c("Rep1", "Rep2", "Rep1", "Rep2") Replicates <- c("Rep1", "Rep2", "Rep1", "Rep2")
method <- "Limma" #Limma method <- "DESeq2" # Limma or DESeq2
normalisation <- "scale" # "scale" normalisation <- "scale" # "scale"
spikes <- NULL spikes <- NULL
pAdjustMethod <- "BH" pAdjustMethod <- "BH"
......
...@@ -18,7 +18,7 @@ countsBoxplots( ...@@ -18,7 +18,7 @@ countsBoxplots(
\item{conditions}{factor vector of the condition from which each sample belongs} \item{conditions}{factor vector of the condition from which each sample belongs}
\item{method}{either "DEseq2" or "Limma"} \item{method}{either "DESeq2" or "Limma"}
\item{col}{colors of the boxplots (one per biological condition)} \item{col}{colors of the boxplots (one per biological condition)}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment