Solved – Negative binomial distribution with Python scipy.stats

count-datanegative-binomial-distributionpoisson distributionpythonscipy

I have a dataset of counts on which I tried to fit a Poisson distribution, but my variance is larger than the average so I decided to use a negative binomial distribution.

I use these formulas enter image description here

to estimate r and p based on the mean and variance of my dataset. However, the nbinom.pmf function requires n and p as parameters. How can I estimate n based on r? The plot is not right if I use r as n.

Best Answer

def convert_params(mu, alpha):
    """ 
    Convert mean/dispersion parameterization of a negative binomial to the ones scipy supports

    Parameters
    ----------
    mu : float 
       Mean of NB distribution.
    alpha : float
       Overdispersion parameter used for variance calculation.

    See https://en.wikipedia.org/wiki/Negative_binomial_distribution#Alternative_formulations
    """
    var = mu + alpha * mu ** 2
    p = (var - mu) / var
    r = mu ** 2 / (var - mu)
    return r, p