Sample from a negative binomial in terms of mean and precision parameters

neg_binomial_2_rng.py

import numpy as np

def neg_binomial_2_rng(mu, phi, size=None):
    """
    RNG for the negative binomial parametrised in terms
    of the mean and the precision.
    mean = mu
    variance = mu + (mu**2)/phi
    """
    tau = np.random.gamma(phi, mu/phi, size=size)
    return np.random.poisson(tau)