DESEq2: Two factor design without interaction
2
0
Entering edit mode
hwildha • 0
@hwildha-7562
Last seen 7.3 years ago
Germany

Hello,

i am using DESeq2 to analyse an experiment with two conditions and three groups, as described in "Example 3" in the ?results help page.

I started the analysis by testing for an interaction effect with the following code:

## Example 3: two conditions, three groups, with interaction terms
dds <- makeExampleDESeqDataSet(n=100,m=18)
dds$group <- factor(rep(rep(c("X","Y","Z"),each=3),2)) design(dds) <- ~ group + condition + group:condition dds <- DESeq(dds) resultsNames(dds) ## Interaction effect with LRT test group.x.cond <- nbinomLRT(dds,reduced=~group + condition) group.x.cond.res <- results(group.x.cond) tab.group.x.cond <- table(group.x.cond.res$padj < 0.05)
tab.group.x.cond

In much the same way as here in the example data, i did not found any gene showing a significant interaction effect in my real data set.
Consequently, i set-up a new full model without interaction term:

## New model without interaction term:
design(dds) <- ~ group + condition
dds <- DESeq(dds)
resultsNames(dds)

I have now three questions:
1. Is there a way to test for a condition effect in a specific group (say Z) based on the new model, as it is described in the ?results help for a model with interaction term?
2. In my own data set, if i set-up a full model with interaction term, and test for treatment effects in specific groups, the number of genes showing a significant treatment effect is much lower compared to the number of genes showing a treatment main effect? I am grateful for any comments on possible reasons for this observation.
3. Is there a way to test for a main effect for group based based on contrasts, i.e. not with a LRT test, but with a WALD test?

deseq2 • 5.6k views
0
Entering edit mode
@mikelove
Last seen 1 day ago
United States

1) The "condition effect for a specific group" can be found either i) by including an interaction term, or ii) by making a new factor variable which consists of unique combinations of group and condition, e.g.:

dds$var = factor(paste0(dds$group, dds$condition)) dds = DESeq(dds) results(dds, contrast=c("var", "ZB", "ZA")) (i) makes it easier to test for differences in the condition effect across groups and to pull out the main effect, (ii) makes it easier to test the condition effect in each group. 2) I wonder if you are testing only the interaction term(s). Contrasts involving only the interaction term(s) are tests for whether the condition effect is different in a specific group or across groups. In order to test the condition effect in a specific group, you add the interaction term to the main effect. Remember: there can be a large condition effect, but it is not different across groups. In this case, the main effect is large, but the interaction terms are small. If the contrast only involves the interaction terms you are not asking, what is the effect, but are the effects different across groups. 3) Yes, with the formula with interactions this is easy. Just use the default DESeq() pipeline, (don't put in test="LRT" or the reduced formula), and then results(dds, contrast=c("condition", "B", "A")) or whatever your condition levels are. ADD COMMENT 0 Entering edit mode Dear Michael Love, I got 2 questions about 2-way factorial design on DESeq2. I understand what the interaction term answers. If it is significant, it means the condition effect is different across genotypes. But, I don't understand which pair-wise comparison this " results(dds, name="genotypeII.conditionB") " for? For example, the results as below, it means the condition effect is different across genotypes for these genes, right? But, which comparison these 106 genes are significantly differently expressed between at FDR 0.05? The second question is how should I code the contrast (fix a level for one factor) for each of 6 pair-wise comparisons: conditionAgenotypeI VS conditionAgenotypeII, conditionBgenotypeI VS conditionBgenotypeII, conditionAgenotypeI VS conditionBgenotypeI, conditionAgenotypeII VS conditionBgenotypeII, conditionAgenotypeI VS conditionBgenotypeII, and conditionAgenotypeII VS conditionBgenotypeI? For these 2 questions, I got stuck with applying post-hoc test (https://www.sheffield.ac.uk/polopoly_fs/1.536444!/file/MASH_2way_ANOVA_in_R.pdf) to DESeq2 in a word. Thanks a lot! subset(results(dds, name="genotypeII.conditionB"), padj<0.05) log2 fold change (MLE): genotypeII.conditionB Wald test p-value: genotypeII.conditionB DataFrame with 106 rows and 6 columns baseMean log2FoldChange lfcSE stat pvalue padj <numeric> <numeric> <numeric> <numeric> <numeric> <numeric> gene1 21759.8 -19.2656 3.94362 -4.88525 1.03298e-06 0.00186987 gene2 16509.7 -19.1612 3.93043 -4.87510 1.08751e-06 0.00186987 gene3 13845.8 -18.5270 3.81758 -4.85309 1.21554e-06 0.00195494   ## Example 2: two conditions, two genotypes, with an interaction term dds <- makeExampleDESeqDataSet(n=100,m=12) dds$genotype <- factor(rep(rep(c("I","II"),each=3),2))

design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)

# the condition effect for genotype I (the main effect)
results(dds, contrast=c("condition","B","A"))

# the condition effect for genotype II
# this is, by definition, the main effect *plus* the interaction term
# (the extra condition effect in genotype II compared to genotype I).
results(dds, list( c("condition_B_vs_A","genotypeII.conditionB") ))

# the interaction term, answering: is the condition effect *different* across genotypes?
results(dds, name="genotypeII.conditionB")


0
Entering edit mode

I don't understand which pair-wise comparison...

For statistical consultation on how to set up and interpret designs in linear model, I recommend to work with a local statistician. I only have sufficient time to answer software related support posts.

0
Entering edit mode

My confusion was fixed by reading "Contrast" and "Interaction" parts here https://www.bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#contrasts. Thanks!

0
Entering edit mode
hwildha • 0
@hwildha-7562
Last seen 7.3 years ago
Germany

Dear Michael,

Add1) Do you think it is valid to include an interaction term, if the result of the LRT is telling me that there is no evidence for an interaction term for any gene?

Add 3) In my third question:

"Is there a way to test for a main effect for group based based on contrasts, i.e. not with a LRT test, but with a WALD test?"

i was actually not asking for the main effect for factor "condition", but for factor "group", which has three levels.

Best, Henning

0
Entering edit mode

1) yes you can include interaction terms even if they are not significant in the LRT. It could be that individual terms are significant, although the LRT (which tests all terms together) is not.

3) Yes, you can just substitute group for condition in the results() call, and the two groups you want to compare.