  RNA-Seq data characteristically exhibits large variances, which need to be
appropriately accounted for in the model. We first explore the effects of this
variability on the maximum likelihood estimator (MLE) of the overdispersion
parameter of the negative binomial distribution, and propose instead the use an
estimator obtained via maximization of the marginal likelihood in a conjugate
Bayesian framework. We show, via simulation studies, that the marginal MLE can
better control this variation and produce a more stable and reliable estimator.
We then formulate a conjugate Bayesian hierarchical model, in which the
estimate of overdispersion is a marginalized estimate and use this estimator to
propose a Bayesian test to detect differentially expressed genes with RNA-Seq
data. We use numerical studies to show that our much simpler approach is
competitive with other negative binomial based procedures, and we use a real
data set to illustrate the implementation and flexibility of the procedure.
