Hello:

Could you please answer some questions about Figure 1 in Love et al (2014), “Moderated estimation of fold change and dispersion”. I am interested in how and why the shrinkage is overridden for the features that exhibit extraordinarily high dispersion.

1) The paper says that such features “do not obey our modeling assumptions”. At the same time the, a similar reasoning is not applied to the features with extraordinarily low dispersion. Why?

2) Suppose the dispersion parameter indeed follows the lognormal prior assumption and no overriding is required. Using the “2 residual standard deviations” rule will result in about 5% of features having overestimated p-values. Don't you think it's a problem, especially in a plausible situation when only a small proportion of features is differentially expressed?

Regards,

Nik

hi Nik,

These are very good questions. For simplicity, I will be speaking of within-group estimates: the estimates of dispersion tend toward zero occur essentially whenever the variance is less than the mean, even with a positive true dispersion for the distribution. This happens often when the counts are low. For example:

varmean <- function(x) var(x)/mean(x)
​sum(replicate(100, varmean(rnbinom(n=3,mu=10,size=10))) < 1)

So the low values are often not a violation of the modeling assumptions, but just occur due to random sampling. 

Raising low dispersion values is a conservative choice, as is not lowering very high values. The real statistical challenge with RNA-seq differential expression is typically small replicate size (n=3-5), and so doing conservative things regarding shrinkage (which is strongest when sample sizes are small and which disappears as sample size grows) makes sense.

Regarding (2), a few points: yes, if the log dispersion values do follow a Normal then we are not shrinking for ~5%. But again we think this is the best choice. The MLE values are another valid estimator for dispersion, and here we are just making a conservative choice when sample size is very small and shrinkage would be strong.

Regarding differential expression (DE), note that dispersion and DE status are not linked. These genes with high dispersion are not the DE genes, they are the genes with extremely high within-group variance.

Of course, there could be DE genes in this set of genes with high dispersion (as there can be DE genes at all levels of dispersion), but these are the hardest genes to include in an FDR-bounded set, because the noise level is so high. The point of a statistical method is to maximize sensitivity while controlling the FDR, which inevitably means not calling such genes DE, so as to control the FDR. Here we are not guaranteeing that we won't call these genes DE, but we are not moderating dispersion and so the DE evidence needs to be strong as the gene-level (without the information sharing across genes).