Best statistical method for differential expression analysis of my RNA-seq data?

Question

[amyfm]

Hi,

I have to do differential expression analsysis of RNA-seq data. I have 2 different conditions, and 5 animals in each condition, and I want to see the DEGs between both conditions.

Around 5 years ago, this analysis with my same data was done using the limma and puma libraries in R and conducting a two-way ANOVA2, obtaining more than 100 DEGS.

Now I have to repeat the analysis using the new genome reference. However, after my analysis with DESeq2 (and previously alignment using STAR and counting reads per feature using htseq-count), I only obtain 5 DEGs.

What's going on? Which method do you think is best?

[Devon Ryan]

Try it with limma again but this time using voom() to prepare the data (this is described in the limma vignette). Limma should only be used in conjunction with voom when it comes to RNAseq data. You probably won't get identical results to DESeq2, but they'll likely be similar.

Answer 1

[Johanna Schott]

According to this article by Law et al. in Genome Biology (voom: precision weights unlock linear model analysis tools for RNA-seq read counts, http://www.genomebiology.com/2014/15/2/R29), DESeq is expected to be extremely conservative compared to limma in conjunction with voom. The major difference is that DESeq is count based and uses the negative binomial distribution to model sampling error of read counts plus the biological variability between replicates, while limma was originally developed for microarrays and is therefore not based on counting statistics.

Have you ever performed validation experiments on your old data? If many of the initial hits could be confirmed, I would trust the first analysis and assume that DESeq was too strict.