Dear all,

I have performed a Differential Expression Analysis using

dds <- DESeqDataSetFromMatrix(countData = cts, colData = colData, design = ~ Type+ batch)

dds<- DESeq(dds)

I understand that even though I included "batch" in the formula, but I still see a batch effect.

I read this:

"Why after VST are there still batches in the PCA plot? The transformations implemented in DESeq2, vst and rlog , compute a variance stabilizing transformation which is roughly similar to putting the data on the log2 scale, while also dealing with the sampling variability of low counts. It uses the design formula to calculate the within-group variability (if blind=FALSE ) or the acrossall-samples variability (if blind=TRUE ). It does not use the design to remove variation in the data. It therefore does not remove variation that can be associated with batch or other covariates (nor does DESeq2 have a way to specify which covariates are nuisance and which are of interest). It is possible to visualize the transformed data with batch variation removed, using the removeBatchEffect function from limma. This simply removes any shifts in the log2-scale expression data that can be explained by batch. The paradigm for this operation for designs with balanced batches would be: For unbalanced batches (e.g. the condition groups are not distributed balanced across batches), the design argument should be used, see ?removeBatchEffect in the limma package for details."

So I tried to run removeBatchEffect after the above formula, but still it seems to be having some issues.

Would that be the right approach to run removeBatchEffect after dds<- DESeq(dds)?

Thank you, Bine

"I still see a batch effect"

Can you explain what you mean by this exactly?

Re: removeBatchEffect(), this is just for visualization of transformed data in RNA-seq, e.g. PCA plots. It is not possible to use in combination with DESeq(), where we recommend to account for batch effects in the design as you have done.