Question: DESeq2 compare all levels
0
gravatar for rmf
13 months ago by
rmf690
rmf690 wrote:

For RNA-Seq DGE analysis using DESeq2, if I have a factor with let's say 5 levels:

condition
A
A
B
B
C
C
D
D
E
E

The standard workflow

d <- DESeqDataSetFromMatrix(countData=counts,colData=meta,design=as.formula(~condition))
d <- DESeq2::estimateSizeFactors(d,type="ratio")
d <- DESeq2::estimateDispersions(d)
d1 <- nbinomWaldTest(d)
resultsNames(d1)

gives the following contrast combinations:

condition_A_vs_B
condition_A_vs_C
condition_A_vs_D
condition_A_vs_E

and then:

results(d1,name="condition_A_vs_B")

But, how do I get all pairwise combinations of all levels? ie; The first combinations as well as these below?

B - C
B - D
B - E
C - D
C - E
D - E
ADD COMMENTlink modified 13 months ago by Kevin Blighe46k • written 13 months ago by rmf690
4
gravatar for Kevin Blighe
13 months ago by
Kevin Blighe46k
Kevin Blighe46k wrote:

To get all pairwise comparisons, it should be as easy as:

BC <- results(dds, contrast=c("condition", "B", "C"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
BC <- lfcShrink(dds, contrast=c("condition", "B", "C"), res=BC)

BD <- results(dds, contrast=c("condition", "B", "D"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
BD <- lfcShrink(dds, contrast=c("condition", "B", "D"), res=BD)

BE <- results(dds, contrast=c("condition", "B", "E"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
BE <- lfcShrink(dds, contrast=c("condition", "B", "E"), res=BE)

CD <- results(dds, contrast=c("condition", "C", "D"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
CD <- lfcShrink(dds, contrast=c("condition", "C", "D"), res=CD)

CE <- results(dds, contrast=c("condition", "C", "E"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
CE <- lfcShrink(dds, contrast=c("condition", "C", "E"), res=CE)

DE <- results(dds, contrast=c("condition", "D", "E"), independentFiltering=TRUE, alpha=0.05, pAdjustMethod="BH", parallel=TRUE)
DE <- lfcShrink(dds, contrast=c("condition", "D", "E"), res=DE)

Note that I use lfcShrink, which you should also use if you had previously used the function DESeq() with betaPrior=FALSE

You can also do a LRT (likelihood ratio test if you want to compare all levels simultaneously).

What you can do after this is then concatenate all of the statistically significant genes into a unique list and cluster your samples using these, which should bring out good separation in a sample dendrogram (using the regularised log or variance-stabilised counts).

Kevin

ADD COMMENTlink written 13 months ago by Kevin Blighe46k

What you can do after this is then concatenate all of the statistically significant genes into a unique list and cluster your samples >using these, which should bring out good separation in a sample dendrogram (using the regularised log or variance-stabilised counts)

I am very much interested in doing that, but I have no idea how I would go about doing that step-by-step. Could you give me any pointers as to which tools to use for the concatenating and clustering? I would appreciate this very much.

ADD REPLYlink written 5 months ago by paul-georg10
1

Hey paul, you can do something like this. If you have conducted 3 different differential expression comparisons, you can subset the results objects (res1, res2, res3) to include the statistically significant genes. You then subset your original data matrix for these genes.

log2cutoff <- 2
qvaluecutoff <- 0.05

sigGenes <- unique(c(
  rownames(subset(res1, padj<=qvaluecutoff & abs(log2FoldChange)>=log2cutoff)),
  rownames(subset(res2, padj<=qvaluecutoff & abs(log2FoldChange)>=log2cutoff)),
  rownames(subset(res3, padj<=qvaluecutoff & abs(log2FoldChange)>=log2cutoff)),
))

heat <- assay(rld)[sigGenes,]

Makes sense?

ADD REPLYlink modified 5 months ago • written 5 months ago by Kevin Blighe46k
1

Oh, yes! Absolutely. I did not think it could be that easy. Thank you very much!

ADD REPLYlink written 5 months ago by paul-georg10
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1302 users visited in the last hour