Negative Binomial Response For Event Counts

Why Use a Negative Binomial Response?

An event count can be nonnegative and have no known fixed upper limit: for example, the number of service calls received by a unit during a specified period. A Poisson response is a natural starting point for such data, but it makes a strong conditional assumption: its variance equals its mean. Counts often vary more than that assumption allows. The negative binomial response keeps the count support while adding a dispersion parameter that permits greater variance.

In a generalized linear mixed model (GLMM), the response distribution is specified conditional on covariates and random effects. The negative binomial dispersion therefore describes extra response variability at that conditional level. It is not interchangeable with random-effect variation, which also induces dependence among observations sharing a random effect.

Definition: A negative binomial response $Y$ with conditional mean $\mu>0$ and dispersion $\alpha>0$ has support $y=0,1,2,\ldots$ and probability mass function $$ P(Y=y)= \frac{\Gamma(y+1/\alpha)}{\Gamma(1/\alpha)y!} \left(\frac{1}{1+\alpha\mu}\right)^{1/\alpha} \left(\frac{\alpha\mu}{1+\alpha\mu}\right)^y . $$ Its mean is $\mu$, and its variance is $\mu+\alpha\mu^2$. Here $\Gamma(\cdot)$ is the gamma function, which extends the factorial because $\Gamma(m)=(m-1)!$ for positive integers $m$.

This parameterization uses $\alpha$ as a dispersion parameter: larger $\alpha$ means greater variance at a fixed mean. Some software instead reports a size parameter $r=1/\alpha$. When using software output, check its parameterization rather than assuming that a reported dispersion number is $\alpha$.

How the Extra Dispersion Arises

A useful construction is a Poisson distribution whose rate varies across otherwise comparable observations. Let a latent rate $\Lambda$ follow a gamma distribution with shape $1/\alpha$ and scale $\alpha\mu$, and, conditional on that rate, let $Y$ be Poisson:

\Lambda\sim \operatorname{Gamma}(1/\alpha,\ \text{scale}=\alpha\mu), \qquad Y\mid\Lambda\sim \operatorname{Poisson}(\Lambda).

This mixture has a direct interpretation: observations can have different underlying event rates, even when they share the same modeled mean. Integrating over those rates produces a negative binomial count distribution. The mixture is one way to understand the family; it does not mean that every negative binomial model must literally include a separately estimated rate for each observation.

Derivation: A gamma variable with shape $k$ and scale $\theta$ has mean $k\theta$ and variance $k\theta^2$. Thus, for the stated gamma distribution, $E(\Lambda)=\mu$ and $\operatorname{Var}(\Lambda)=\alpha\mu^2$. Using the conditional Poisson mean and variance, followed by the laws of total expectation and variance, gives $$ E(Y)=E\{E(Y\mid\Lambda)\}=E(\Lambda)=\mu, $$ $$ \operatorname{Var}(Y) =E\{\operatorname{Var}(Y\mid\Lambda)\} +\operatorname{Var}\{E(Y\mid\Lambda)\} =E(\Lambda)+\operatorname{Var}(\Lambda) =\mu+\alpha\mu^2. $$ For completeness, integrating the conditional Poisson probability against the gamma density gives $$ \begin{aligned} P(Y=y) &=\int_0^\infty \frac{e^{-\lambda}\lambda^y}{y!} \frac{\lambda^{1/\alpha-1}e^{-\lambda/(\alpha\mu)}} {\Gamma(1/\alpha)(\alpha\mu)^{1/\alpha}}\,d\lambda\\ &=\frac{\Gamma(y+1/\alpha)}{\Gamma(1/\alpha)y!} \left(\frac{1}{1+\alpha\mu}\right)^{1/\alpha} \left(\frac{\alpha\mu}{1+\alpha\mu}\right)^y, \end{aligned} $$ where the second equality follows by integrating a gamma-kernel with rate $1+1/(\alpha\mu)$. This is the probability mass function in the definition.

The Poisson distribution is recovered as a limiting case when $\alpha$ approaches zero: the gamma rate variation vanishes because $\operatorname{Var}(\Lambda)=\alpha\mu^2$ approaches zero, leaving a fixed rate $\mu$. For positive $\alpha$, the variance exceeds the mean by $\alpha\mu^2$. The extra variance grows with both the dispersion and the squared mean.

Worked Example: Computing Negative Binomial Probabilities

Suppose a hypothetical unit has a modeled mean of $\mu=2.5$ events per observation period and dispersion $\alpha=0.4$. Calculate its variance, the probability of no events, and the probability of exactly three events.

The variance is $$ \operatorname{Var}(Y)=2.5+0.4(2.5)^2=2.5+2.5=5.0 $$ events squared. To calculate probabilities, set $r=1/\alpha=2.5$ and $p=1/(1+\alpha\mu)=1/(1+1)=0.5$. The mass function can then be written as $$ P(Y=y)=\frac{\Gamma(y+r)}{\Gamma(r)y!}p^r(1-p)^y. $$ For zero events, $$ P(Y=0)=0.5^{2.5}=0.1768. $$ For three events, $$ \begin{aligned} P(Y=3) &=\frac{\Gamma(5.5)}{\Gamma(2.5)3!}(0.5)^{2.5}(0.5)^3\\ &=\frac{(2.5)(3.5)(4.5)}{6}(0.5)^{5.5}\\ &=6.5625(0.02210)\approx 0.1450. \end{aligned} $$ Thus the model assigns about a 17.68% probability to zero events and a 14.50% probability to exactly three. These are model-based probabilities for the stated hypothetical mean and dispersion.

Using the Family in a GLMM

A common GLMM specifies a log link for the conditional mean. For observation $i$ in cluster $j$, one possible structure is $$ \log(\mu_{ij})=\mathbf{x}_{ij}^{\mathsf T}\boldsymbol{\beta}+b_j, \qquad Y_{ij}\mid b_j,\mathbf{x}_{ij}\sim \operatorname{NB}(\mu_{ij},\alpha). $$ Here $\mathbf{x}_{ij}$ is the covariate vector, $\boldsymbol{\beta}$ is the coefficient vector, and $b_j$ is the cluster random effect. Conditional on the covariates and random effects, the response has mean $\mu_{ij}$ and variance $\mu_{ij}+\alpha\mu_{ij}^2$. The negative binomial family specifies the response distribution; the link and linear predictor specify how its mean depends on predictors.

The usual conditional assumptions include nonnegative integer responses, a correctly specified conditional mean structure, and the stated negative binomial conditional distribution. If the model assumes conditional independence, responses are independent given the modeled covariates and random effects. Responses from the same cluster need not be independent after the random effects are integrated out.

Worked Example: A Conditional Mean and Variance in a GLMM

Consider a hypothetical count model with $\log(\mu_{ij})=\beta_0+\beta_1x_{ij}+b_j$. Let $\beta_0=\log(1.5)$, $\beta_1=\log(2)$, $x_{ij}=1$, and $b_j=\log(1.2)$, with $\alpha=0.2$. Find the conditional mean and variance for this observation.

Substitution in the linear predictor gives $$ \log(\mu_{ij})=\log(1.5)+\log(2)+\log(1.2)=\log(3.6), \qquad \mu_{ij}=3.6. $$ The conditional variance is $$ \operatorname{Var}(Y_{ij}\mid b_j,\mathbf{x}_{ij}) =3.6+0.2(3.6)^2 =3.6+2.592 =6.192. $$ These values describe the response conditional on the specified random effect and covariate. Holding the random effect fixed, increasing $x_{ij}$ by one multiplies the conditional mean by $\exp(\beta_1)=2$. That is a conditional mean ratio, not a claim that the predictor causes a change in event counts.

A random effect and a negative binomial dispersion parameter address different features. The random effect represents modeled heterogeneity shared by observations in a cluster and can create within-cluster association. The negative binomial dispersion permits extra variation in individual counts conditional on that random effect. A model can include both, but fitting both does not guarantee that either structure is appropriate.

Worked Example: Comparing Poisson and Negative Binomial Variability

Suppose two hypothetical models assign the same mean, $\mu=4$ events per observation period. One is Poisson; the other is negative binomial with $\alpha=0.25$. Compare their variances and probabilities of zero events.

The Poisson variance equals its mean, so it is $4$ events squared. For the negative binomial, $$ \operatorname{Var}(Y)=4+0.25(4)^2=4+4=8 $$ events squared. For the zero probability, the Poisson model gives $$ P(Y=0)=e^{-4}\approx 0.01832. $$ For the negative binomial, $r=1/0.25=4$ and $p=1/(1+0.25(4))=0.5$, so $$ P(Y=0)=p^r=(0.5)^4=0.0625. $$ At this mean and dispersion, the negative binomial assigns a higher probability to zero counts and has twice the variance. This comparison illustrates one possible consequence of the extra dispersion; it does not imply that every negative binomial model must assign more probability to zero than a particular Poisson model at every count.

Checking Fit and Choosing Among Approaches

The negative binomial is useful when a Poisson model’s conditional variance restriction is too narrow, but a variance larger than the mean in raw data is not, by itself, proof that the negative binomial is the right model. In clustered data, raw variation can also reflect random effects, differing covariates, omitted predictors, or unequal observation exposure. A model should be assessed after accounting for its mean structure and dependence.

Check the response and sampling structure. The family is for nonnegative integer counts without a fixed upper bound. A known maximum number of opportunities can make a binomial response more suitable; a fixed total allocated among categories calls for a multinomial response.
Examine conditional fit. Compare observed counts with predictions, including the frequency of zeros and the upper tail. Residual patterns across fitted means or covariates can indicate a wrong mean structure or remaining dispersion.
Check clustering and exposure. If observations share a cluster, assess whether the random-effects structure captures their dependence. If observation windows differ, account for exposure in the mean rather than interpreting all counts as if they covered equal time.
Look for remaining structure. A negative binomial distribution does not automatically fix excess zeros, dependence not represented by the model, or covariate misspecification. Persistent discrepancies in predictive checks are evidence to reconsider the model.

The main alternative for unbounded event counts is the Poisson family. It is simpler and appropriate when its conditional mean–variance restriction is credible; the negative binomial relaxes that restriction by adding dispersion. A binomial model is preferable when each count has a known finite number of eligible opportunities. Random effects can model cluster-level variation and dependence, but they do not make the conditional negative binomial dispersion redundant by definition. Choosing among these approaches requires attention to how the observations were generated, not merely which model yields a larger variance.

Limitations and a Common Mistake

The negative binomial adds flexibility, but its variance still has a specified relationship to the mean, and its probability distribution may not match every count process. It does not correct a wrong link, omitted predictors, an unsuitable random-effects structure, or selection that excludes some counts. If observations are truncated or their recording windows differ, those features require explicit treatment. Model diagnostics and knowledge of the observation process remain necessary.

Common mistake: “The counts are clustered, so I should use a negative binomial; its dispersion parameter handles the clustering.” The correction is that negative binomial dispersion describes extra conditional response variability, while a random effect can represent shared cluster-level variation and induce within-cluster association. Depending on the design, a model may need a random effect, negative binomial dispersion, or both; one is not a substitute for automatically specifying the other.

Key takeaway: For a negative binomial response with mean $\mu$ and dispersion $\alpha>0$, the variance is $\mu+\alpha\mu^2$. In a GLMM this is a conditional response variance, distinct from variation and dependence introduced by random effects.

Check Your Understanding

Answer each question, then compare with the solutions below.

For $\mu=3$ and $\alpha=0.2$, calculate the negative binomial variance and state whether it exceeds the mean.
For those same parameter values, calculate the probability of zero events.
In a GLMM, what does the negative binomial dispersion describe, and how does that differ from a cluster random effect?
Why is a raw sample variance larger than a raw sample mean not enough to establish that a negative binomial model is appropriate?

Solutions. (1) $\operatorname{Var}(Y)=3+0.2(3)^2=3+1.8=4.8$ events squared, which exceeds the mean of 3 events. (2) Here $r=1/0.2=5$ and $p=1/(1+0.2(3))=1/1.6=0.625$. Therefore $P(Y=0)=p^r=(0.625)^5\approx0.09537$, or about 9.537%. (3) Dispersion describes extra response variability conditional on covariates and random effects; a cluster random effect represents shared cluster variation and can induce association among responses in that cluster. (4) Raw variation may also reflect differing means, clustering, omitted predictors, or exposure differences. The conditional model and observation structure must be considered before choosing a response family.