Question: DESeq2 - Method-of-moments - negative binomial
0
gravatar for Fischer-philipp
3 months ago by
Fischer-philipp20 wrote:

Within the function momentsDisEstimate() of DESeq2 a rough method-of-moments estimate of the mean counts is derived by this equation (bv - xim*bm)/bm^2. (I am aware that this is just an initial set to maximize Cox-Reid adjusted likelihood of the gene-wise dispersion estimat.)

So I was curious and tried to calculate the MoME myself - but I always came to something looking like the inverse of it.

Unfortunately I am not aware how to write mathematical formulas on the support bioconductor page. Thats why I formulated a question on:

https://stats.stackexchange.com/questions/399685/deseq2-method-of-moments-negative-binomial

momentsDispEstimate <- function(object) {
xim <- if (!is.null(normalizationFactors(object))) {
mean(1/colMeans(normalizationFactors(object)))}
else { mean(1/sizeFactors(object)) }
bv <- mcols(object)$baseVar
bm <- mcols(object)$baseMean
(bv - xim*bm)/bm^2 }

Thank you for your time and answer.

ADD COMMENTlink modified 3 months ago by Michael Love24k • written 3 months ago by Fischer-philipp20
Answer: DESeq2 - Method-of-moments - negative binomial
1
gravatar for Michael Love
3 months ago by
Michael Love24k
United States
Michael Love24k wrote:

There are two definitions of “dispersion” parameter for the negative binomial, one is the inverse of the other.

ADD COMMENTlink written 3 months ago by Michael Love24k

Thank you for your quick answers.

so it comes from the relationship sigma^2 = \mu + \alpha \mu^2 ?

It still remains a question to me why it is: (bv - xim*bm)/bm^2 }

and not

(bv - xim*bm)/(xim*bm)^2 }

or

(bv - bm)/bm^2

ADD REPLYlink modified 3 months ago • written 3 months ago by Fischer-philipp20

I don't remember the details, but I believe I profiled both of these options and went with the code that is there. Remember, this value is only used as part of a process to pick an initial point for the line search, so it's ok to be very rough. This particular estimate is used as an upper bound for the much better estimate of dispersion from roughDispEstimate(). I think the latter was sometimes giving implausibly high values and so this estimate we are talking about is to stabilize the initialization process.

ADD REPLYlink written 3 months ago by Michael Love24k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 191 users visited in the last hour