This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Correcting for batch effect in RNA-seq data

I used DESeq2 to process RNA-seq data from different sources. And I found harsh batch effect when plotted PCA (different shapes of the figures represent 3 different batches, for example, ctr and PH.7d from different batches cluster apart):

enter image description here

I tried to remove it using limma package as described here:

colData
      sample   condition batch
1         100       PH.7d     1
..........
7          75         ctr     1
8  SRR5035380 hblast.10.5     2
..........
25 SRR5035397 hblast.18.5     2
26 SRR8437299         ctr     3
..........
37 SRR8437324       PH.7d     3

vsd<-vst(dds)
assay(vsd)<-limma::removeBatchEffect(assay(vsd),vsd$data1)
data2<-plotPCA(vsd, intgroup=c('condition','batch'),returnData=T)
data2<-as.data.frame(data2)
percentVar<- round(100*attr(data2,'percentVar'))
plot2<-qplot(PC1,PC2,color=condition,shape=batch,data=data2)

However, there is no changes when I plot the results:

enter image description here

What am I doing wrong?

Also, I tried to remove batch effect using design in DESeq:

ddsB=DESeqDataSetFromMatrix(countData = countData,colData = colData, design = ~batch+condition)

I'm getting this error:

Error in checkFullRank(modelMatrix) : 
the model matrix is not full rank, so the model cannot be fit as specified.
One or more variables or interaction terms in the design formula are linear
combinations of the others and must be removed.

Can somebody help me to solve it?

Thanks in advance!

rna-seq batch-effect
• 5,312 views • link
updated by Ram Ram • written by Rimma • Rimma

Are you sure that vsd$data1 corresponds to the vector encoding the batch variable? Seems to me it should be vsd$batch.

• 0 views • link
written by Friederike Friederike

It looks like batch 2 doesn't contain any of the groups in batch 1 and 3, therefore it is not possible to correct for that batch. Are you sure there is at least one overlapping group in batch 2, that is also found in batch 1 and 3?

• 0 views • link
written by Benn Benn

1 answer

RNA-seq is strongly confounded by the kit and library preparation method from what I've seen. The confounding effect kight dominate the biological variability. The confounding effect probably dominates any kind of biological differences, see here for example a PCA that I made from five independent data sources, processed identically from the in silico side.

Edit: Check if correct use of batch removal attempts as Benn says below can limit the confounding effect.

enter image description here

• 0 views • link
written by ATpoint ATpoint

But you can correct for batch when there are overlapping groups. However, I suspect that OP's batch 2 doesn't have any overlapping group...

• 0 views • link
written by Benn Benn

True, but to what extend. Do you have experience on how well this works. I mean "mild" batch effects like different culture conditions in the lab, samples taken on different days or different sequencing protocols might be correctable, but can you really "regress" out the effect of different kits and laboratories?

• 0 views • link
written by ATpoint ATpoint

I have only experience with removeBatchEffect() from edgeR/limma, they work fine, especially for visualization. Clearly the limma::removeBatchEffect code from OP did not work properly. Like Friederike is already suspecting.

• 0 views • link
written by Benn Benn

Log in to answer this question.

More posts like this