On 2016 Jan 30, Wolfgang Huber commented:

The method is superseded by our more recent work: Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Love MI, Huber W, Anders S. Genome Biology 2014 15:550. Love MI, 2014

On 2013 Nov 11, James Hadfield commented:

This paper describes the DESeq method, for differential RNA-Seq, ChIP-Seq and other analyses. DEseq was developed to work on low replicate numbers and indeed many people cannot generate high numbers of replicates. But I would challenge the community to consider that the costs of NGS have dropped very significantly since these methods were conceived and that increasing replicate numbers to higher levels is now inexcusable in many scenarios.

Both of the papers referred to in the comments so far reference multiple RNA-seq, and/or other, datasets that were used to test the methods from which their conclusions are drawn. Wolfgang Huber mentions the constraints of samples-size in his comments and also has a section on working without replicates in the Anders/Huber paper above, in it they discuss the impact that within and between group sample variability have on the results.

Some very real difficulties in appraising which approach (DESeq2 or SamSeq) is best include the limited amount of time the community has been testing the different approaches, that the approaches themselves are still very much in development, and that very different datasets are used in each study.

This last issue is made more of a problem since the experimental methods section in many NGS papers is generally not clear enough. It would help to have clear guidelines on the number samples used and their relationship e.g. biological or technical replicates, and if technical at which stage is replication being made; the number and type of reads generated at a per sample and per group level would also be useful. Getting this information can be painful as evidenced by digging through the DESeq2 and SamSeq references:

DESeq2

Wilczynski: very difficult to determine from the data provided or online. Engstrom: mRNATag-seq, 5 samples (3 & 2 replicates per group), no indication of reads per sample.

Nagalakshmi: mRNA-seq, 4 samples (2 technical & 2 biological replicates per group), 7M reads per sample (possibly).

Kasowski: ChIP-seq, 10 biological samples, 2 groups, 660M reads 33M reads per sample.

SamSeq

Hoen: mRNA-seq, 4 replicates each of around 2.5M reads.

Marioni: mRNA-seq, 2 groups (liver & kidney), 7 technical replicates (at the lane level), 85M reads per group, 12M reads per replicate.

Witten: miRNA RNA-seq, 29 biological samples per group (Tumour vs normal), average 0.75M reads per sample.

Perhaps a simple format could be agreed on by the community as a table to be added in to each publication as a supplemental?

On 2013 Jul 02, Wolfgang Huber commented:

Robert Tibshirani raises an important point: how does the method react to outliers, e.g. to cases where a gene hovers around the same value for most samples but has extremely different values for a few isolated samples. In some cases, such genes might be of interest, but in most biological applications, users seem to be more interested in genes that show a consistent shift in expression for all or most samples in a condition.

We have recently adapted DESeq2 to be able to flag (and ignore) such outliers, see Appendix F of the packages vignette, which is available from http://www.bioconductor.org.

Regarding SamSeq (or related methods), these are very promising in large-sample situations. Section 3.3 of the SamSeq paper Li J, 2013 is titled "Data with small sample size" and considers a comparison of 5 versus 5 samples. However, sample sizes in some applications of RNA-Seq are smaller, controlled experiments with 2 versus 2 samples are not uncommmon. They can be reasonably analysed with parametric methods such as DESeq/DESeq2, and requiring a 250%-fold experimental effort in order to be able to apply a non-parametric method could be uneconomic.

There is a pressing need for benchmarks; A few years ago, Rafael Irizarry's AffyComp Irizarry RA, 2006 was an excellent example of how to do this, at the time, on Affymetrix arrays applied to designed experiments. It would be interesting to see similar benchmarks for RNA-Seq data, for various experimental or study designs.

On 2013 Jun 22, Robert Tibshirani commented:

DESeq represents a promising approach for detecting differential abundance in RNA-seq and other count-based biological data. However it is a parametric method, and my student Jun Li found that as a result, it can be very non-robust to outlying points. The same is true of EdgeR, another well-known method, that like DESeq, uses the negative binomial distribution. And such outlying points can be common in this kind of data--- there are often some extremely large counts.

Jun and I developed a resampling-based non-parametric method called "SamSeq" http://www.ncbi.nlm.nih.gov/pubmed/22127579 which in our experiments, showed considerably better robustness, with little loss of power when there were no outliers. This method is implemented in the R package "samr" http://cran.r-project.org/web/packages/samr/index.html and in the Excel add-on package SAM http://www-stat.stanford.edu/~tibs/SAM/. Our comparisons were made with the original DESeq method: the new outlier detection features added in DESeq2 might help matters. We haven't tried this comparison and it would be important to do this. A word of caution: I have written only a few papers on this topic, and have not kept up with the latest developments.

On 2013 Jun 18, Wolfgang Huber commented:

This paper lays out the main ideas behind the DESeq method, which can be used to detect differential abundance in count-based biological data such as obtained from RNA-Seq, ChIP-Seq, HiC, mass spectrometry. In the meanwhile, the method has evolved substantiallly, new features include:

• treatment of general NB-GLMs (generalised linear models of the Negative Binomial family) including paired designs and interactions
• better dispersion estimation (especially in the important small sample size setting)
• outlier detection
• specialised tools to detect differential relative exon usage, described in Anders S, 2012 and provided by the DEXSeq package on Bioconductor

For more details on these improvements, see the vignette of the DESeq2 package http://www.bioconductor.org/packages/release/bioc/html/DESeq2.html (which is in the process of replacing the previous implementation, DESeq) and also the Supplement of Anders S, 2012.

In addition, much practical information is available in online fora and mailing list, in particular https://stat.ethz.ch/mailman/listinfo/bioconductor and http://seqanswers.com/forums/forumdisplay.php?f=18. Generic search engines seem to be a good way to search through these.

This approach was developed particularly for experiments with limited numbers of replicates. For larger data sets (e.g. observational studies involving hundreds of samples), the granularity of the counts becomes less important and estimating the variability becomes easier, so that also other, more generic methods become sufficiently performant. In particular, data transformations, such as DESeq2's varianceStabilizingTransformation and rlogTransformation, together with general-purpose dimension reduction, clustering or classification algorithms, or with normal ANOVA in these cases often provide attractive options.

On 2013 Jun 18, Wolfgang Huber commented:

This paper lays out the main ideas behind the DESeq method, which can be used to detect differential abundance in count-based biological data such as obtained from RNA-Seq, ChIP-Seq, HiC, mass spectrometry. In the meanwhile, the method has evolved substantiallly, new features include:

• treatment of general NB-GLMs (generalised linear models of the Negative Binomial family) including paired designs and interactions
• better dispersion estimation (especially in the important small sample size setting)
• outlier detection
• specialised tools to detect differential relative exon usage, described in Anders S, 2012 and provided by the DEXSeq package on Bioconductor

For more details on these improvements, see the vignette of the DESeq2 package http://www.bioconductor.org/packages/release/bioc/html/DESeq2.html (which is in the process of replacing the previous implementation, DESeq) and also the Supplement of Anders S, 2012.

In addition, much practical information is available in online fora and mailing list, in particular https://stat.ethz.ch/mailman/listinfo/bioconductor and http://seqanswers.com/forums/forumdisplay.php?f=18. Generic search engines seem to be a good way to search through these.

This approach was developed particularly for experiments with limited numbers of replicates. For larger data sets (e.g. observational studies involving hundreds of samples), the granularity of the counts becomes less important and estimating the variability becomes easier, so that also other, more generic methods become sufficiently performant. In particular, data transformations, such as DESeq2's varianceStabilizingTransformation and rlogTransformation, together with general-purpose dimension reduction, clustering or classification algorithms, or with normal ANOVA in these cases often provide attractive options.

On 2013 Jun 22, Robert Tibshirani commented:

DESeq represents a promising approach for detecting differential abundance in RNA-seq and other count-based biological data. However it is a parametric method, and my student Jun Li found that as a result, it can be very non-robust to outlying points. The same is true of EdgeR, another well-known method, that like DESeq, uses the negative binomial distribution. And such outlying points can be common in this kind of data--- there are often some extremely large counts.

Jun and I developed a resampling-based non-parametric method called "SamSeq" http://www.ncbi.nlm.nih.gov/pubmed/22127579 which in our experiments, showed considerably better robustness, with little loss of power when there were no outliers. This method is implemented in the R package "samr" http://cran.r-project.org/web/packages/samr/index.html and in the Excel add-on package SAM http://www-stat.stanford.edu/~tibs/SAM/. Our comparisons were made with the original DESeq method: the new outlier detection features added in DESeq2 might help matters. We haven't tried this comparison and it would be important to do this. A word of caution: I have written only a few papers on this topic, and have not kept up with the latest developments.

On 2013 Jul 02, Wolfgang Huber commented:

Robert Tibshirani raises an important point: how does the method react to outliers, e.g. to cases where a gene hovers around the same value for most samples but has extremely different values for a few isolated samples. In some cases, such genes might be of interest, but in most biological applications, users seem to be more interested in genes that show a consistent shift in expression for all or most samples in a condition.

We have recently adapted DESeq2 to be able to flag (and ignore) such outliers, see Appendix F of the packages vignette, which is available from http://www.bioconductor.org.

Regarding SamSeq (or related methods), these are very promising in large-sample situations. Section 3.3 of the SamSeq paper Li J, 2013 is titled "Data with small sample size" and considers a comparison of 5 versus 5 samples. However, sample sizes in some applications of RNA-Seq are smaller, controlled experiments with 2 versus 2 samples are not uncommmon. They can be reasonably analysed with parametric methods such as DESeq/DESeq2, and requiring a 250%-fold experimental effort in order to be able to apply a non-parametric method could be uneconomic.

There is a pressing need for benchmarks; A few years ago, Rafael Irizarry's AffyComp Irizarry RA, 2006 was an excellent example of how to do this, at the time, on Affymetrix arrays applied to designed experiments. It would be interesting to see similar benchmarks for RNA-Seq data, for various experimental or study designs.

On 2013 Nov 11, James Hadfield commented:

This paper describes the DESeq method, for differential RNA-Seq, ChIP-Seq and other analyses. DEseq was developed to work on low replicate numbers and indeed many people cannot generate high numbers of replicates. But I would challenge the community to consider that the costs of NGS have dropped very significantly since these methods were conceived and that increasing replicate numbers to higher levels is now inexcusable in many scenarios.

Both of the papers referred to in the comments so far reference multiple RNA-seq, and/or other, datasets that were used to test the methods from which their conclusions are drawn. Wolfgang Huber mentions the constraints of samples-size in his comments and also has a section on working without replicates in the Anders/Huber paper above, in it they discuss the impact that within and between group sample variability have on the results.

Some very real difficulties in appraising which approach (DESeq2 or SamSeq) is best include the limited amount of time the community has been testing the different approaches, that the approaches themselves are still very much in development, and that very different datasets are used in each study.

This last issue is made more of a problem since the experimental methods section in many NGS papers is generally not clear enough. It would help to have clear guidelines on the number samples used and their relationship e.g. biological or technical replicates, and if technical at which stage is replication being made; the number and type of reads generated at a per sample and per group level would also be useful. Getting this information can be painful as evidenced by digging through the DESeq2 and SamSeq references:

DESeq2

Wilczynski: very difficult to determine from the data provided or online. Engstrom: mRNATag-seq, 5 samples (3 & 2 replicates per group), no indication of reads per sample.

Nagalakshmi: mRNA-seq, 4 samples (2 technical & 2 biological replicates per group), 7M reads per sample (possibly).

Kasowski: ChIP-seq, 10 biological samples, 2 groups, 660M reads 33M reads per sample.

SamSeq

Hoen: mRNA-seq, 4 replicates each of around 2.5M reads.

Marioni: mRNA-seq, 2 groups (liver & kidney), 7 technical replicates (at the lane level), 85M reads per group, 12M reads per replicate.

Witten: miRNA RNA-seq, 29 biological samples per group (Tumour vs normal), average 0.75M reads per sample.

Perhaps a simple format could be agreed on by the community as a table to be added in to each publication as a supplemental?

On 2016 Jan 30, Wolfgang Huber commented:

The method is superseded by our more recent work: Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Love MI, Huber W, Anders S. Genome Biology 2014 15:550. Love MI, 2014