Gamma-Poisson Mixing and the Negative-Binomial Law

From a Random Rate to a Count Distribution

A Poisson model assumes that its rate is fixed. In some settings, however, different observational units have different underlying rates, even after accounting for measured covariates. A useful model represents that heterogeneity by making the Poisson rate random. When the random rate has a gamma distribution, integrating it out produces a negative-binomial count distribution.

The derivation matters for more than recognizing a familiar name. It specifies exactly which negative-binomial parameterization is being used and shows how rate variation changes the mean and variance. In particular, a random rate adds variability beyond the Poisson variance.

Definition: Let $\Lambda\sim\operatorname{Gamma}(\alpha,\beta)$ in the shape-rate parameterization, where $\alpha>0$ and $\beta>0$. Its density is $$ f_\Lambda(\lambda)=\frac{\beta^\alpha}{\Gamma(\alpha)}\lambda^{\alpha-1}e^{-\beta\lambda},\qquad \lambda>0, $$ where $\Gamma(a)=\int_0^\infty t^{a-1}e^{-t}\,dt$ for $a>0$. Conditional on $\Lambda=\lambda$, suppose $X\mid\Lambda=\lambda\sim\operatorname{Poisson}(\lambda)$, so $P(X=k\mid\Lambda=\lambda)=e^{-\lambda}\lambda^k/k!$ for $k=0,1,2,\ldots$. This construction is called a gamma-Poisson mixture.

Deriving the Negative-Binomial Mass Function

For each nonnegative integer $k$, the unconditional probability is obtained by integrating the conditional Poisson probability over all possible rates. The gamma density contributes a factor $e^{-\beta\lambda}$, while the Poisson probability contributes $e^{-\lambda}$; these combine to give the rate $\beta+1$ inside the integral.

Theorem: If $\Lambda\sim\operatorname{Gamma}(\alpha,\beta)$ and $X\mid\Lambda\sim\operatorname{Poisson}(\Lambda)$, then, for $k=0,1,2,\ldots$, $$ P(X=k)=\frac{\Gamma(\alpha+k)}{\Gamma(\alpha)k!} \left(\frac{\beta}{\beta+1}\right)^\alpha \left(\frac{1}{\beta+1}\right)^k. $$ Equivalently, $X$ has a negative-binomial distribution with shape $r=\alpha$ and success probability $p=\beta/(\beta+1)$, using the convention that the random variable counts failures before $r$ successes.

To verify the formula, start from the mixture integral:

\begin{aligned} P(X=k) &=\int_0^\infty P(X=k\mid\Lambda=\lambda)f_\Lambda(\lambda)\,d\lambda\\ &=\int_0^\infty \frac{e^{-\lambda}\lambda^k}{k!} \frac{\beta^\alpha}{\Gamma(\alpha)} \lambda^{\alpha-1}e^{-\beta\lambda}\,d\lambda\\ &=\frac{\beta^\alpha}{\Gamma(\alpha)k!} \int_0^\infty \lambda^{\alpha+k-1}e^{-(\beta+1)\lambda}\,d\lambda. \end{aligned}

The integral is finite because $\alpha+k>0$ and $\beta+1>0$. Substituting $t=(\beta+1)\lambda$ gives

\int_0^\infty \lambda^{\alpha+k-1}e^{-(\beta+1)\lambda}\,d\lambda =\frac{\Gamma(\alpha+k)}{(\beta+1)^{\alpha+k}}.

Thus the mass function above follows. It is a probability distribution: summing the masses over $k\geq 0$ gives the binomial series $(1-q)^{-\alpha}$, with $q=1/(\beta+1)<1$, multiplied by $(\beta/(\beta+1))^\alpha=(1-q)^\alpha$, so the sum is one.

Under the stated negative-binomial convention, the mass function with parameters $r>0$ and $0<p<1$ is

P(X=k)=\frac{\Gamma(r+k)}{\Gamma(r)k!}p^r(1-p)^k,\qquad k=0,1,2,\ldots.

Matching terms gives $r=\alpha$, $p=\beta/(\beta+1)$, and $1-p=1/(\beta+1)$. The gamma rate is therefore $\beta=p/(1-p)$. If a gamma distribution is instead specified by a scale, that scale is $1/\beta$, not $\beta$.

Worked Example: A Gamma Rate of 2

Suppose $\Lambda\sim\operatorname{Gamma}(3,2)$ in the shape-rate convention and $X\mid\Lambda\sim\operatorname{Poisson}(\Lambda)$. Find the marginal distribution and calculate $P(X=2)$, the mean, and the variance.

Here $\alpha=3$ and $\beta=2$. The negative-binomial parameters are $r=3$ and $p=2/(2+1)=2/3$, so $1-p=1/3$. For $k=2$,

\begin{aligned} P(X=2) &=\frac{\Gamma(3+2)}{\Gamma(3)2!} \left(\frac{2}{3}\right)^3 \left(\frac{1}{3}\right)^2\\ &=\frac{24}{2\cdot 2}\cdot\frac{8}{27}\cdot\frac{1}{9} =6\cdot\frac{8}{243} =\frac{48}{243} =\frac{16}{81}\approx 0.1975. \end{aligned}

The mean and variance, derived below, are

E[X]=\frac{\alpha}{\beta}=\frac{3}{2}, \qquad \operatorname{Var}(X)=\frac{\alpha}{\beta}+\frac{\alpha}{\beta^2} =\frac{3}{2}+\frac{3}{4}=\frac{9}{4}.

The variance exceeds the mean by $3/4$, reflecting the variability of the random rate.

Probability-Generating Function and Moments

The probability-generating function gives another compact description of the mixture and makes its moment structure visible. For $0\leq z\leq 1$, conditional expectation gives $E[z^X\mid\Lambda]=e^{\Lambda(z-1)}$. Integrating this expression against the gamma density yields

G_X(z)=E[z^X] =\left(\frac{\beta}{\beta+1-z}\right)^\alpha =\left(\frac{p}{1-(1-p)z}\right)^r.

The integral is finite on this range because $\beta+1-z\geq\beta>0$. The expression is also differentiable in an open interval containing $z=1$, since $\beta>0$. Thus its derivatives at $1$ give ordinary factorial moments.

Theorem: For the gamma-Poisson mixture above, $$ E[X]=\frac{\alpha}{\beta},\qquad \operatorname{Var}(X)=\frac{\alpha}{\beta}+\frac{\alpha}{\beta^2}. $$ More generally, for every integer $m\geq 1$, the $m$th falling-factorial moment is $$ E[(X)_m]=\frac{\Gamma(\alpha+m)}{\Gamma(\alpha)\beta^m}, $$ where $(X)_m=X(X-1)\cdots(X-m+1)$.

Here is a direct proof of the mean and variance. The gamma density gives $E[\Lambda]=\alpha/\beta$ and $\operatorname{Var}(\Lambda)=\alpha/\beta^2$. These follow by integrating powers of $\lambda$ against the density and using $\Gamma(a+1)=a\Gamma(a)$. Conditional on $\Lambda$, a Poisson variable has mean and variance $\Lambda$. The laws of total expectation and total variance therefore give

E[X]=E[E[X\mid\Lambda]]=E[\Lambda]=\frac{\alpha}{\beta},

and

\begin{aligned} \operatorname{Var}(X) &=E[\operatorname{Var}(X\mid\Lambda)] +\operatorname{Var}(E[X\mid\Lambda])\\ &=E[\Lambda]+\operatorname{Var}(\Lambda) =\frac{\alpha}{\beta}+\frac{\alpha}{\beta^2}. \end{aligned}

For the factorial moments, differentiate $G_X(z)=(\beta/(\beta+1-z))^\alpha$ $m$ times. For each integer $m\geq 1$,

G_X^{(m)}(z) =\frac{\Gamma(\alpha+m)}{\Gamma(\alpha)} \frac{\beta^\alpha}{(\beta+1-z)^{\alpha+m}}.

The power-series definition of the probability-generating function gives $G_X^{(m)}(1)=E[(X)_m]$; the differentiation is valid because the generating function is finite in a neighborhood of $1$. Substitution of $z=1$ gives the stated formula. In particular, $E[X(X-1)]=\alpha(\alpha+1)/\beta^2$; using $X^2=X(X-1)+X$ recovers the variance formula. $\square$

Worked Example: Translating Negative-Binomial Parameters

A negative-binomial count $X$ has shape $r=4$ and success probability $p=0.6$, with $X$ counting failures before $r$ successes. Find the corresponding gamma mixing distribution and calculate the mean and variance.

Matching $r=\alpha$ and $p=\beta/(\beta+1)$ gives $\alpha=4$ and

\beta=\frac{p}{1-p}=\frac{0.6}{0.4}=\frac{3}{2}.

Thus this negative-binomial law is the marginal distribution obtained by mixing a Poisson count over $\Lambda\sim\operatorname{Gamma}(4,3/2)$. Its mean and variance are

E[X]=\frac{4}{3/2}=\frac{8}{3}, \qquad \operatorname{Var}(X)=\frac{4}{3/2}+\frac{4}{(3/2)^2} =\frac{8}{3}+\frac{16}{9} =\frac{40}{9}.

The variance-to-mean ratio is $(40/9)/(8/3)=5/3$, greater than one.

Interpreting the Extra Variance

A fixed-rate Poisson variable has variance equal to its mean. In the mixture, the conditional variance still equals the conditional mean, but the rate itself varies across units. The variance decomposition separates these contributions:

\underbrace{\operatorname{Var}(X)}_{\text{total count variation}} = \underbrace{E[\Lambda]}_{\text{Poisson variation at a given rate}} + \underbrace{\operatorname{Var}(\Lambda)}_{\text{variation across rates}}.

For $\alpha,\beta>0$, the second term $\alpha/\beta^2$ is positive. Consequently, the marginal count is overdispersed relative to a Poisson distribution with the same mean: its variance is strictly larger than its mean. This is why a gamma-Poisson mixture is a useful model when count variation exceeds what a single fixed-rate Poisson model allows.

Worked Example: Zero Counts and Overdispersion

Let $\Lambda\sim\operatorname{Gamma}(2,1/4)$, where $1/4$ is the rate, and let $X\mid\Lambda\sim\operatorname{Poisson}(\Lambda)$. Calculate $P(X=0)$, $P(X=1)$, and compare the variance with the mean.

The negative-binomial parameters are $r=2$, $p=(1/4)/(1+1/4)=1/5$, and $1-p=4/5$. Therefore

P(X=0)=\left(\frac{1}{5}\right)^2=\frac{1}{25}=0.04,

and

\begin{aligned} P(X=1) &=\frac{\Gamma(3)}{\Gamma(2)1!} \left(\frac{1}{5}\right)^2\left(\frac{4}{5}\right)\\ &=2\cdot\frac{1}{25}\cdot\frac{4}{5} =\frac{8}{125}=0.064. \end{aligned}

The mean and variance are

E[X]=\frac{2}{1/4}=8,\qquad \operatorname{Var}(X)=8+\frac{2}{(1/4)^2}=8+32=40.

The variance is five times the mean. This difference is not an arithmetic anomaly: the random rate has variance $2/(1/4)^2=32$, and that between-rate variation adds to the conditional Poisson variation of $8$.

Common Parameterization Pitfalls

Negative-binomial notation is not universal, so a parameter pair without a stated convention can be ambiguous. In this tutorial, $r$ is the shape and $p$ is the success probability in a count of failures before $r$ successes. Some sources instead count successes before a fixed number of failures, or use a mean parameter. Always check the mass function rather than relying on the distribution name alone.

Confusing gamma rate and scale. With shape-rate parameters, $E[\Lambda]=\alpha/\beta$. If the gamma is given by shape and scale $s$, then $s=1/\beta$ and $E[\Lambda]=\alpha s$.
Reversing the negative-binomial probabilities. The gamma rate maps to $p=\beta/(\beta+1)$, while the factor raised to $k$ is $1-p=1/(\beta+1)$.
Assuming the marginal count remains Poisson. Conditional on the rate it is Poisson; after averaging over a non-degenerate gamma rate, its variance exceeds its mean.
Using the wrong count convention. Under the failure-count convention, the mass at zero is $p^r$, which provides a quick check on parameter matching.

Key takeaway: Mixing a Poisson rate $\Lambda\sim\operatorname{Gamma}(\alpha,\beta)$ in the shape-rate convention gives a negative-binomial count with $r=\alpha$ and $p=\beta/(\beta+1)$. Its mean is $\alpha/\beta$, and its variance adds the positive rate-heterogeneity term $\alpha/\beta^2$ to the Poisson contribution.

Check Your Understanding

Use the failure-count negative-binomial convention throughout.

If $\Lambda\sim\operatorname{Gamma}(5,3)$ in the shape-rate convention, what are the negative-binomial shape and success probability?
For a gamma-Poisson mixture with $\alpha=2$ and $\beta=1$, calculate the mean and variance, and explain why they differ.
A negative-binomial count has $r=3$ and $p=1/4$. What gamma shape and rate produce it as a Poisson mixture?
Why does a non-degenerate gamma mixing distribution make the marginal count variance exceed its mean?
For parameters $r$ and $p$, what is the probability of zero under the failure-count convention, and how can it check a proposed parameter match?