I decided to split DE genes over time with DESeq2 into several parts, as I figured out that it is too long. I hope it was a good idea.
I'm doing a time course analysis with deseq2, I have 8 different points and multiple replica.
when I am looking at
resultsNames(dds) in my analysis I see
>resultsNames(dds) #  "Intercept" "replica_2_vs_1" "replica_3_vs_1" "time_16h_vs_0h" "time_1d_vs_0h" #  "time_24h_vs_0h" "time_30h_vs_0h" "time_48h_vs_0h" "time_72h_vs_0h" "time_90h_vs_0h"
but when I do
results(dds), I get only the last time point
> res log2 fold change (MLE): time 90h vs 0h LRT p-value: '~ replica + time' vs '~ replica' DataFrame with 17558 rows and 6 columns ...
Does it means that with this coefficient I get all the genes which are changing over all TPs?
How can I interpret the other coefficients? Are the changes over time up to a specific TP (e.g. 16h_vs_0h) means all genes changing from 0h to 16h?
Is it possible to extract with this kind of analysis also all pair-wise comparisons (like 1d vs. 16h, or 24h vs. 48h), or do i need to run a "normal" wald test for the pair-wise comparisons?
the script I run is below:
#reading in the phenotypic data from a file phenotype <- read.delim2("phenotypeData.txt", colClasses = c(rep("factor", 3))) # create the DESeq object from a matrix. dds<-DESeqDataSetFromMatrix(countData=countTable, colData=phenotype, design= ~ replica + time ) dds = DESeq(dds, test="LRT", reduced=~replica) res<-results(dds)