Hi Dear Scientist, Thank you so much for the platform and as usual, I may have your few minutes to my question? I used both DESeq2 and edgeR to analyze my RNAseq data. However, I found a higher number of significant genes in my DESeq2 analysis compared to edgeR. The difference is like 1000 which is high in my opinion. Here I posted the code I used below and please show me my mistake case I missed something. My variable(sample) is continuous data.
## edgeR
x <- DGEList(counts = muscle, group = Sample)
design <- model.matrix(~Sample)
fit <- estimateDisp(fat, design = design, robust = TRUE)
QL <- glmQLFit(fit, design = design)
table(p.adjust(QL$table$PValue, method="BH")<0.05) #### 5928 genes
### DESeq2
dds <- DESeqDataSetFromMatrix(countData = countData,
colData = colData,
design = ~ Sample)
dds <- DESeq(dds, fitType = "mean")
resultsNames(dds )
Sample <- results(dds)
sum(Sample$padj < 0.05, na.rm = TRUE) #### 6042 genes
Thank you so much! Best, Amare
I'm going to echo "swbarnes2". This has been asked and answered, even recently, on the support site. The methods are different, and it's a lot less interesting or surprising when you note that a method might call a gene DE because adj p = 0.04, and another method might call a gene not DE because adj p = 0.06.
We recommend on the site (and have on many previous threads), pick a tool and use it, but it's not a good idea to alternate through various methods on the dataset you are going to use these methods for analysis. You can certainly do this on another dataset in order to choose which method to use, or just look at the dozens of papers comparing methods systematically on simulated and real datasets.
Thank you so much for your replay. All of the DE genes (adjp < 0.05) I found from edgeR analysis are actually present in the DESeq2 analysis result(adjp<0.05). Thanks!