2
$\begingroup$

I ran SVA to remove batch effects for my bulk RNAseq experiments, but now I need to somehow correct my data matrix in order to run pca, mds. I am using DESeq2 for the analysis. Here is the code that I got now:

dds <- estimateSizeFactors(dds)
dat = counts(dds, normalized = TRUE)
idx = rowMeans(dat) > 1
dat = dat[idx,]
mod = model.matrix(~Group, colData(dds))
mod0 = model.matrix(~1, colData(dds))

nsv = num.sv(dat, mod)
svseq = svaseq(dat, mod, mod0, n.sv = nsv)

dds_sva = dds
dds_sva$SV1 = svseq$sv[,1]
design(dds_sva) = ~ SV1 + Group
dds_sva = DESeq(dds_sva)

So, I am wondering now, how to use svseq object to plot, say, the corrected version of MDS plot.

$\endgroup$
2
$\begingroup$

If you want to plot the "corrected" expression, you will need to remove the variation introduced by these surrogate variables. Removing the expression affected can introduce some bias too and it is usually not recommended (despite comBat doing so). You should apply linear algebra, you can look at here is an example how to do it:

library("corpcor")
s <- corpcor::fast.svd(t(scale(t(dat), center = TRUE, scale = TRUE)))
pcSds <- s$d
pcSds[1] <- 0
svdexp <- s$u %*% diag(pcSds) %*% t(s$v)
colnames(svdexp) <- colnames(dat)

But if your surrogate variables do have a big enough influence on your MDS plot to be observable (IMHO) there would be too much batch effect to work with the data.

| improve this answer | |
$\endgroup$
  • $\begingroup$ But batch effect often appears as grouping of points, and that is how often I can see it. So, I am a bit confused regarding it. I would expect this grouping to disappear and get a different clustering of points. Are you talking only about MDS or about pca too? $\endgroup$ – Nikita Vlasenko Oct 30 '18 at 14:15
  • 1
    $\begingroup$ Well, without knowing the data is difficult to know if there is really a batch effect or other variables that could be adjusted for. I'm talking about both PCA and MDS, usually there isn't such a big difference between them (unless you use some other metric for distance in the MDS that is not the eucledian distance). $\endgroup$ – llrs Oct 30 '18 at 14:22

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.