[Gordon Smyth]

DESeq is basically edgeR with some improvements, so if you want common sense, that seems to be the winner. Since DESeq and edgeR use the same distribution while DEGseq uses a different one, they naturally get more similar results, and that's not a sensible way to conclude that they're better. However, both of the negative-binomial methods' authors provide good evidence that DEGseq's Poisson assumption is invalid.

Here is the DESeq paper: http://genomebiology.com/2010/11/10/R106

My common sense suggests that one needs to estimate the variability for each transcript, in such a way that more variable transcripts get higher dispersion estimates, and only edgeR does that. I'm senior author for the edgeR papers, so I'm biased, but other authors agree, for example Hardcastle et al (BMC Bioinformatics, 2010). Hardcastle et al show edgeR and bayseq to have better performance than other possibilities.

[A Oshlack]

It's interesting that your number of up and down genes is so different with the two approaches. I wonder if it is actually a "normalization" difference. If I were you I would look at some MA plots or smear plots coloured for each of the methods and see if you can see anything the might push you to choose one over the other. Nothing better than visualising your data in my opinion

[tianyub836]

It's interesting that your number of up and down genes is so different with the two approaches. I wonder if it is actually a "normalization" difference. If I were you I would look at some MA plots or smear plots coloured for each of the methods and see if you can see anything the might push you to choose one over the other. Nothing better than visualising your data in my opinion

well, thanks, A Oshlack.

i did have MA plots and smear plots generated by those two packages, but i do not ubderstand what you have inferred.

would you please give me more details.

thank you

[A Oshlack]

As I'm sure you know, an MA plot is just a plot of the raw data. In the normalization assumptions, at least for TMM in edgeR, we assume that the "cloud" of points in an MA plot should be centred at M=0. You have used a two fold cut-off as well as the statistical test so I wonder if the difference between methods that you are using is actually more dependent on the normalization that you have chosen. Are the DE gene in the MA plot equal distance up and down from where you think that the normalization line should sit. Often your eyes are very good at assessing this if you generate the right plots. Hope that helps.

[joxcargator73]

I believe that replicates are very important to have good quality results. RNA-seq is becoming cheaper and cheaper but still quite expensive for small labs. In this case I also believe that RNAseq without replicates could be used as screening and then confirm by replicating qRT-PCR and based you conclusion on these results.

[rfilbert]

I think both EdgeR and DESeq are both pretty bad. Not sure how this community came to believe that "borrowing variance" from genes with similar average expression makes sense. It does not make biological sense and it is mathematically provably wrong. Therefore, both EdgeR, DESeq, Limma produce more false positives and more false negatives than a more robust and reliable statistical test. If it made good sense, respected statistical leaders like SAS and Partek would have adopted this methodology, and neither company has.

[chadn737]

I think both EdgeR and DESeq are both pretty bad. Not sure how this community came to believe that "borrowing variance" from genes with similar average expression makes sense. It does not make biological sense and it is mathematically provably wrong. Therefore, both EdgeR, DESeq, Limma produce more false positives and more false negatives than a more robust and reliable statistical test. If it made good sense, respected statistical leaders like SAS and Partek would have adopted this methodology, and neither company has.

Then why does Partek normalize its data using RPKM, which has been shown to be problematic time and again?

https://springerlink3.metapress.com/content/18284k158v887k57/resource-secured/?target=fulltext.html

http://bib.oxfordjournals.org/conten...ib.bbs046.long

http://www.ncbi.nlm.nih.gov/pmc/arti...rtype=abstract.

etc

I don't buy this argument of "SAS and Partek don't do it, therefore its wrong".

[kopi-o]

So ... what does SAS use for differential gene expression analysis, and who uses SAS for this purpose? Genuinely curious ...

[rfilbert]

I'm not sure what SAS's JMP Genomics team is up to lately, but I was just pointing out that they never jumped on the "borrowing variance" bandwagon for microarrays. The idea of borrowing variance is not unique to NGS data.

[chadn737]

Also, I think your claims about how DESeq works are somewhat inaccurate.

The default for DESeq is to calculate both the variance for each gene and the variance by borrowing from other genes of similar expression and then take the maximum of the two. Since one is taking the maximum of the two values, this would reduce the number of false positives, but increase the number of false negatives. Its easy enough to alter DESeq so that it uses one or the other so that if you really feel that passionate about not borrowing variance, than you don't have to.

[kopi-o]

The same goes for edgeR - you don't need to use the moderated dispersions if you don't want to.

I usually prefer SAMSeq (a non-parametric method) for DE analysis when there are enough replicates, but the nice thing with DESeq and edgeR is that you can consider complex designs in your analysis.

[rfilbert]

I would definitely recommend Partek Flow for the following reasons:
1. They fit 5 different distribution assumptions and use the best fit for each gene or transcript. This gives more statistical power and more biological meaning.
2. There are no limit to the number of factors (handles any type of experiment design)
3. It has a really easy to use point & click web-based GUI.
4. They have excellent technical support.

[chadn737]

And how is it normalized?

[rfilbert]

There are multiple options for normalization, but I believe the default option is to simply normalize to the total number of reads for each sample. I don't think any normalization based on the length of the transcript (like RPKM) matters as for this analysis you are comparing the same transcript in different groups of samples.

[chadn737]

1. They fit 5 different distribution assumptions and use the best fit for each gene or transcript. This gives more statistical power and more biological meaning.

I believe that this approach makes it more likely to call a gene as differentially expressed. After all, it sounds like you are simply trying everything until you find the one fit that will give you a positive result. As far as making biological sense.......I'm not at all convinced.