DESeq2 compare all levels using wald and lrt test

Question

[aiswaryabioinfo]

I am currently working in plant pathogen interaction and is trying to make a gene co-expression network using RNA-seq data. I used HTseq-count to get the read counts.

My dataset consist of Resis_treated (1st day, 3rd day and 7th day) samples taken under 3 time point with 2 replicated for each time point. I have a factor with 3 levels

condition
Day1
Day1
Day3
Day3
Day7
Day7

I am trying to get the DEGs between Day3_Vs_Day1, Day7_VsDay1 and Day7_Vs_Day3. The code I used is as follows

library(DESeq2)

directory<-'D:/test_deseq/tempora'

sampleFiles<-grep('Day',list.files(directory),value=TRUE)

sampleCondition<-c('Day1','Day1','Day3','Day3', 'Day7','Day7')

sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles, condition=sampleCondition)

ddsHTSeq<-DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory=directory, design=~condition)

colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c('Day1','Day3','Day7'))

colData(ddsHTSeq)

ddsHTSeq <- ddsHTSeq[ rowSums(counts(ddsHTSeq)) > 1, ]

ddsHTSeq

dds <- DESeq(ddsHTSeq)

resultsNames(dds)

[1] "Intercept"              "condition_Day3_vs_Day1" "condition_Day7_vs_Day1"

But, how do I get all pairwise combinations of all levels? ie; condition_Day7_vs_Day3 along with the other two combinations

Will LRT test give a good all_vs_all level comparison, If so how I should design the DEseq2 condition for it.

Answer 1

[Devon Ryan]

If you want all of the pairwise comparisons then request all of them. Namely, use lfcShrink() with the contrast= option 3 times.

An LRT is mostly useful for asking the general "is there a change according to day?" question.

[aiswaryabioinfo]

Hi Devon,

Will I get the condition_Day7_vs_Day3 combinations through the following command

day73 <- results(dds, contrast=c("condition", "day7", "day3"))
day73 <- lfcShrink(dds, contrast=c("condition", "day3", "day3"), res=day73)` #(Use lfcShrink() with the contrast= option 3 times.)

I want to see the change according to days across all the three changes. For that what full and reduced design I should use

[Devon Ryan]

day73 <- lfcShrink(dds, contrast=c("condition", "day7", "day3"))
day71 <- lfcShrink(dds, contrast=c("condition", "day7", "day1"))
day31 <- lfcShrink(dds, contrast=c("condition", "day3", "day1"))

[aiswaryabioinfo]

What is the difference in using LRT in this context (All_vs_all comparisons)

[Devon Ryan]

An LRT will tell you whether there's a change in day, but not whether there's a difference between any particular days. It compares the model where day is included against one without it (i.e., it's what a one-way ANOVA is doing).

[aiswaryabioinfo]

Thank you sir for your help

[aiswaryabioinfo]

After running the code, I tried to contrast day7_Vs_Day3 data. But it showed an error as follows

day73 <- lfcShrink(dds, contrast=c("condition", "day7", "day3")) Error in cleanContrast(object, contrast, expanded = isExpanded, listValues = listValues, : day7 and day3 should be levels of condition such that condition_day7_vs_Day1 and condition_day3_vs_Day1 are contained in 'resultsNames(object)'

The Sample condition I have provided is

sampleCondition<-c('Day1','Day1','Day3','Day3', 'Day7','Day7')

sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles, condition=sampleCondition)

ddsHTSeq<-DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory=directory, design=~condition)

colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c('Day1','Day3','Day7'))

Kindly help me to resolve the error

[Devon Ryan]

Capitalize the D in Day.

[aiswaryabioinfo]

Hi Devon, I got the error as typo

Is it necessary to use lfcshrink() in all_vs_all design comparison. Will it affect in identification of significant genes?

What affect lfcshrink () will have in the expression values of genes?

I tried the all_vs_all design comparison in DESeq2 (Version 2.11.40.2) provided in online Galaxy toolbox and later tried in R (Version DESeq2_1.22.2). After applying padj and log2fc cutoff to both the results, I noticed more number of significant DEGs in the result online Galaxy toolbox. Is the difference in number of genes due to the difference in DESeq2 version.

[Devon Ryan]

You should always use lfcShrink in all cases, the results are more reliable. The version of DESeq2 used in the Galaxy wrapper are likely quite old, which accounts for differences. Please see the DESeq2 documentation for the other questions.