Hello guys,

I am doing some RNA-seq analysis using DESeq2 and following the vignette from here, where they use time points 0,15,30,60 ... 180 min. and two groups (strains)

http://www.bioconductor.org/help/workflows/rnaseqGene/

In the time series chapter it says:

"Genes with small p values from this test are those which, at one or more time points after time 0 showed a strain-specific effect."

As an example the following is shown:

## log2 fold change (MLE): strainmut.minute180 
## LRT p-value: '~ strain + minute + strain:minute' vs '~ strain + minute' 
## DataFrame with 4 rows and 7 columns
##               baseMean log2FoldChange     lfcSE      stat       pvalue         padj      symbol
##              <numeric>      <numeric> <numeric> <numeric>    <numeric>    <numeric> <character>
## SPBC2F12.09c  174.6712    -2.65763737 0.7498270  99.23199 7.671942e-20 5.184698e-16       atf21

Question 1)

As I am understanding it- this gene has a strain specific expression profile over time. This effect occures over the other time points (from 0 to 180).

If I had picked the strainmut.minute60  then we would have seen a list of DE genes from 0 to 60 minutes, right?

Just for clarification- here are the options of which one can choose:

resultsNames(ddsTC)

##  [1] "Intercept"           "strain_mut_vs_wt"    "minute_15_vs_0"      "minute_30_vs_0"     
##  [5] "minute_60_vs_0"      "minute_120_vs_0"     "minute_180_vs_0"     "strainmut.minute15" 
##  [9] "strainmut.minute30"  "strainmut.minute60"  "strainmut.minute120" "strainmut.minute180"

Question 2)

Imagine a gene that behaves in the same way in all timepoints (between the two groups) but time point 0 to 15.

Will this be detected by deseq?

Cheers,

john

"If I had picked the strainmut.minute60  then we would have seen a list of DE genes from 0 to 60 minutes, right?"

This is a model with interaction terms, and this is an interaction term. .You can tell it is an interaction term because it contains the names of both variables strain and minute. So this term is a test for if the mut vs WT fold change is different at minute 60 than at minute 0. 

To generate the tables of log fold change of 60 minutes vs 0 minutes for the WT strain would be:

results(dds, name="minute_60_vs_0")

To generate the tables of log fold change of 60 minutes vs 0 minutes for the mut strain would be the sum of the WT term above and the interaction term which is an additional effect beyond the effect for the reference level (WT):

results(dds, contrast=list(c("minute_60_vs_0","strainmut.minute60"))

"Imagine a gene that behaves in the same way in all timepoints (between the two groups) but time point 0 to 15. Will this be detected by deseq?"

This is the description of an interaction term in the model. You can find this by testing the interaction term strainmut.minute15:

results(dds, name="strainmut.minute15")

Hello Michael,

You said "So this term is a test for if the mut vs WT fold change is different at minute 60 than at minute 0. "

Alright. A gene that has a low pvalue in strainmut.minute60 AND 30 AND 120 AND 180 is then to be considered DE in all time points (compared to 0) between WT and MUT.

----------------

You said:

"To generate the tables of log fold change of 60 minutes vs 0 minutes for the WT strain would be:

results(dds, name="minute_60_vs_0")"

So this shows me genes that are DE between two time points in one group (WT). Not considering the MUT at all.

---------------

You said: "reference level (WT):"

Why is WT our reference? Where did we specify that WT is the reference and not MUT? Why this way?

-------------

By calling 

results(dds, name="strainmut.minute180")

Deseq finds this gene because it takes 0minutes and 180minutes into account. That's alright.

Is there a way to take all time points into account? In order to find a gene that is DE in all time points between two groups?