Question: design versus interaction term similarities and/or differences?
gravatar for tarun2
5 months ago by
tarun20 wrote:

To the developers,

I've been working on differential expression analysis of drought-tolerance in rice. I have 2 genotypes (tolerant and susceptible) and 2 conditions (drought and well-watered) with 4 replications each, essentially a 2x2 factorial experiment with 4 replications.

One of my main objectives is to identify drought-responsive genes. I set-up the codes as follows:

colData <- data.frame(genotype=rep(c("IL","Swarna"),each=8),
rownames(colData) <- colnames(tx.all$counts)
dds <- DESeqDataSetFromTximport(tx.all, colData, formula(~genotype+condition+genotype:condition))
colData(dds)$condition<-relevel(colData(dds)$condition, ref = "Control")
dds$group<-factor(paste0(dds$genotype, dds$condition))
design(dds) <- ~group
dds<-DESeq(dds, betaPrior = TRUE, parallel = TRUE)
#[1] "Intercept"          "groupILControl"     "groupILDrought"     "groupSwarnaControl"
[5] "groupSwarnaDrought

Is adding the design term similar to adding an interaction term? We want to make use of the interaction term because we still don't know the mechanisms underlying drought in this experiment hence, we use several contrast arguments from resultsNames

We also want to perform ANOVA so we can identify genes to be used for WGCNA analysis. Does using the codes below is acceptable enough to address this?

dds.ANOVA<-DESeq(dds, test="LRT", reduced = ~1, parallel = TRUE)
res.anova <- results(dds.ANOVA)
sum(res.anova$padj < 0.05, na.rm = T )

Please advise.



R version 3.3.3 (2017-03-06)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1

[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252   
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] hexbin_1.27.1              vsn_3.42.3                 readr_1.1.1               
 [4] tximport_1.2.0             genefilter_1.56.0          pheatmap_1.0.8            
 [7] RColorBrewer_1.1-2         ggplot2_2.2.1              gplots_3.0.1              
[10] DESeq2_1.14.1              SummarizedExperiment_1.4.0 Biobase_2.34.0            
[13] GenomicRanges_1.26.4       GenomeInfoDb_1.10.3        IRanges_2.8.2             
[16] S4Vectors_0.12.2           BiocGenerics_0.20.0       
ADD COMMENTlink modified 5 months ago by Michael Love19k • written 5 months ago by tarun20
gravatar for Michael Love
5 months ago by
Michael Love19k
United States
Michael Love19k wrote:

"One of my main objectives is to identify drought-responsive genes."

Can you be more specific than this? What hypotheses do you want to test?

ADD COMMENTlink written 5 months ago by Michael Love19k

Hi Michael,

Since we have two genotypes, one is very drought tolerant (IL) and the other is very drought susceptible (Swarna). The drought-tolerant genotype is a cross between an inbred line (tolerant) and the other genotype in the study (susceptible).

There's a major QTL that is introgressed into the tolerant genotype and this QTL makes the tolerant genotype obviously tolerant to drought. 

Now we want to identify genes that are uniquely upregulated and downregulated in both genotypes under drought. Specifically, we want to identify unique genes in the QTL region of the very drought tolerant genotype that are upregulated and downregulated.  We then want to use these genes for functional validation through CRISPR-Cpf1.

One of the major hypothesis that we want to test is that there are differentially expressed genes (DEGs) between the two genotypes in the QTL region under drought and we want to identify them.Specifically, we want to identify DEGs under drought in the tolerant genotype. Since we don't know the mechanisms of drought tolerance at the reproductive-stage, we set-up contrast to identify DEGs between several groups.

We also want to establish the difference between the tolerant and the susceptible genotype under drought on the whole genome level.That is, what are the DEGs especially under drought and the ontologies/pathways that are associated with this group of DEGs.



ADD REPLYlink modified 5 months ago • written 5 months ago by tarun20

It sounds like you could accomplish all of this with pairwise comparisons using results() and contrast with the group variable as you have it above.  This will give you Wald test statistics and adjusted p values, which can be used by downstream software packages. 

ADD REPLYlink written 5 months ago by Michael Love19k

Thank you as always Michael.

I would also like to ask if the following codes will address testing for overall difference across the groups like oneway-ANOVA?

dds.ANOVA<-DESeq(dds, test="LRT", reduced = ~1, parallel = TRUE)
res.anova <- results(dds.ANOVA)

We want to select genes for WGCNA analysis and filtering genes by a differential expression is not recommended. It is recommended that probesets or genes may be filtered by mean expression or variance.



ADD REPLYlink written 5 months ago by tarun20

That will give you something like an ANOVA. But it is also a test of differential expression.

ADD REPLYlink written 5 months ago by Michael Love19k

Thanks, Michael.

So is there any other way to do ANOVA after transformation? Or the codes above would be the one that DESEq2 users can use?

Or is there any other way that DESEq2 can employ filtering by mean expression or variance after transformation?

Please advise.



ADD REPLYlink written 5 months ago by tarun20

Sorry these questions don’t really make sense and I can’t answer them meaningfully.


ADD REPLYlink written 5 months ago by Michael Love19k

Sorry for the confusion and unclear questions Michael.

Probably, what I'm trying to ask is would the above codes be sufficient enough to run a one-way ANOVA?


ADD REPLYlink written 5 months ago by tarun20
Please log in to add an answer.


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