I have RNA-seq data from 27 samples to analyse. I have 2 differents treatments (Compound1 and compound 2) + control (DMSO) and I have 3 times. I have 3 replicates for each combination.
I am trying to use DESeq2 to find the differentially expressed genes. First, may I double check with you that the formula looks OK: ~treatment*time ?? (I have fount ~cell + dex on http://www.bioconductor.org/help/workflows/rnaseqGene/ and I am confused. I am wondering why you would use this formula. Is there something special I should now about DESeq2 formulas?)
So, using a pooled analysis with all the samples at once and the design ~treatment*time:
- I tried to plot the dispersion versus the mean normalized counts (using plotDispEsts(dds)). I am not sure what is the red line in this plot because I thought the parameters were calculated per gene using the replicates... But for the low mean normalized counts (where the dispersion is bigger), this red line is well above the black and shrinked blue points. Is this the sign of a problem?
-I also looked at the p-value plots for different comparisons and some of them are not flat at all. The pvalues seem to slowly increase for the low p-values to the high pvalues with even some spike at 0.4 or 0.8. Is this classical for DESeq2 or is this the sign of a problem?
I thought that there weird observations were coming from a high variability between times so I tried to do a similar analysis but per time (formula: ~treatment ). My problem is that doing this I obtain very different lists of differentially expressed genes...
I am not sure which one is the best analysis really... any suggestion please?