Plotting PCA plot in DESeq2 on SVA corrected counts

Question

[Mozart]

Hello everyone, I am trying to understand how to show a PCA plot on SVA-corrected counts by using DESeq2 (I know that expert guys tend to prefer some other tools but I need this for a presentation). Many posts out there (like this, this and this) but none of them seem to address exactly what I need. Specifically, looking at the following post and applying an rlog on it, PCA plot doesn't show up any difference from the uncorrected version. Any useful hints?

Thanks

Answer 1

[ATpoint]

I prefer to perform PCA outside of DESeq2 (but inspired by the code): Assume we start from a count matrix counts that should be on log scale. That could be your matrix of batch-corrected counts. I did not read the entire question because it was too extensive so tell me if I got you wrong.

# Find most variable genes, say top 500:
rv <- rowVars(as.matrix(counts))
select <- order(rv, decreasing = TRUE)[seq_len(min(500, length(rv)))]

# perform PCA on these genes:
pca <- prcomp(t(counts[select, ]))

# calculate the explained variation per principal component:
percentVar <- pca$sdev^2/sum(pca$sdev^2)

# feed into ggplot:
d <- data.frame(PC1 = pca$x[, 1], PC2 = pca$x[, 2], stringsAsFactors = FALSE)

# Most simplistic plot:
ggplot(data=d, aes(x = PC1, y = PC2)) + geom_point()

You can add legends etc. to indicate which data point is which sample etc.

[Mozart]

Thanks ATpoint, that basically solved my problem. It really sounds so reassuring how important this forum could be in the learning process. It's good to have a good string of code to use for all the pipelines, to be honest.

Thanks a million.

[ATpoint]

Sure think. Most of what I know today I learned from online platforms such as this one here, blogs as well as from the many Bioconductor manuals and packages.

[Mozart]

Thanks for the useful hints ATpoint

Not sure why applying this further normalisation step before plotting results? Any ideas/recommendations? Difficult to catch up with her during this period of time, unfortunately.

[ATpoint]

I personally cannot comment since I am not an active user of this package, therefore no experience with it.