Negative-Binomial Models For Overdispersed Counts in Multiple Imputation

Why Use a Negative-Binomial Model?

A Poisson regression model is a natural starting point for a nonnegative integer count, such as the number of hospital visits during a fixed period. Its key restriction is that, conditional on the predictors, the count variance equals its mean. Real count data can be more variable than this restriction allows: individuals with the same modeled predictors may still differ substantially in their underlying count rates. If a Poisson model is used for imputation in that situation, its predictive draws can be too concentrated around the conditional mean.

A negative-binomial model retains a regression model for the conditional mean while allowing the conditional variance to exceed it. In MICE, this distinction matters because imputed values should reflect the variation in the conditional distribution, not only its average. The negative-binomial model is useful when extra-Poisson variation remains after the predictors and exposure have been represented appropriately.

Definition: For a count \(Y_i\) and predictors \(x_i\), a Poisson regression specifies \(Y_i\mid x_i\sim\operatorname{Poisson}(\mu_i)\), with \(\log(\mu_i)=x_i^\mathsf{T}\beta\) (or a linear predictor that also includes an offset). It implies \(\operatorname{Var}(Y_i\mid x_i)=\mu_i\). A negative-binomial regression specifies the same type of conditional mean model but allows \(\operatorname{Var}(Y_i\mid x_i)=\mu_i+\alpha\mu_i^2\), where \(\alpha>0\) is a dispersion parameter. The Poisson model is recovered as \(\alpha\) approaches zero.

Parameterization and Conditional Model

This tutorial uses the mean-dispersion parameterization. Write \(r=1/\alpha\). Given the predictors, a negative-binomial count has probability mass function

P(Y_i=y\mid x_i)= \frac{\Gamma(y+r)}{\Gamma(r)y!} \left(\frac{r}{r+\mu_i}\right)^r \left(\frac{\mu_i}{r+\mu_i}\right)^y, \qquad y=0,1,2,\ldots

Here \(\Gamma(\cdot)\) is the gamma function, \(\mu_i>0\) is the conditional mean, and \(r>0\) controls dispersion. Equivalently, the two probability factors are \(p_i=r/(r+\mu_i)\) and \(1-p_i=\mu_i/(r+\mu_i)\). This parameterization gives \(E(Y_i\mid x_i)=\mu_i\) and \(\operatorname{Var}(Y_i\mid x_i)=\mu_i+\mu_i^2/r=\mu_i+\alpha\mu_i^2\). A larger \(\alpha\) means greater variance at a fixed mean; it does not change the modeled mean itself.

For regression, a common specification is \(\log(\mu_i)=x_i^\mathsf{T}\beta\). If individuals have different observation times or exposure amounts \(t_i\), include \(\log(t_i)\) as an offset: \(\log(\mu_i)=\log(t_i)+x_i^\mathsf{T}\beta\). This models a rate through the predictors while allowing the expected count to scale with exposure.

Model assumptions: Conditional on the included predictors and any offset, the count is a nonnegative integer with a negative-binomial distribution; its mean follows the specified log-linear form; and the conditional variance follows \(\mu_i+\alpha\mu_i^2\) with a common dispersion parameter unless the model explicitly allows otherwise. Observations are conditionally independent under the fitted model. For an imputation model, predictors should also represent the relevant relationships in the analysis model and information associated with missingness; the imputation procedure must account for uncertainty in model parameters when generating imputations.

Why the Variance Has This Form

A useful justification comes from viewing the negative binomial as a Poisson distribution with an individual-specific rate. Suppose \(\Lambda_i\), the rate for observation \(i\), varies across otherwise similar individuals, with a gamma distribution having shape \(r\) and rate \(r/\mu_i\). Then \(E(\Lambda_i)=\mu_i\) and \(\operatorname{Var}(\Lambda_i)=\mu_i^2/r\). Conditional on \(\Lambda_i\), let \(Y_i\mid\Lambda_i\sim\operatorname{Poisson}(\Lambda_i)\).

Applying the laws of total expectation and total variance gives the result. The first term in the variance is the Poisson variation within a given rate; the second is the variation between individual rates.

E(Y_i\mid x_i)=E\{E(Y_i\mid\Lambda_i,x_i)\}=E(\Lambda_i)=\mu_i

\begin{aligned} \operatorname{Var}(Y_i\mid x_i) &=E\{\operatorname{Var}(Y_i\mid\Lambda_i,x_i)\} +\operatorname{Var}\{E(Y_i\mid\Lambda_i,x_i)\}\\ &=E(\Lambda_i)+\operatorname{Var}(\Lambda_i)\\ &=\mu_i+\frac{\mu_i^2}{r} =\mu_i+\alpha\mu_i^2. \end{aligned}

Thus, extra-Poisson variation can be represented as heterogeneity in latent rates. The gamma mixture is a justification for the distribution and its variance; it does not establish that every source of observed overdispersion is actually caused by gamma-distributed individual rates.

Worked Example: Screening for Extra-Poisson Variation

Suppose a hypothetical group of individuals with broadly similar recorded predictors has a sample mean of 2.4 visits and a sample variance of 8.1 visits squared. As a preliminary check, compare the observed variance with the Poisson variance expected at that mean.

Under the Poisson restriction, the expected variance is the mean, or 2.4. The variance-to-mean ratio is

\frac{s^2}{\bar{x}}=\frac{8.1}{2.4}=3.375.

The observed sample variance is much larger than the Poisson benchmark. A rough method-of-moments estimate of the negative-binomial dispersion, using \(\alpha\approx(s^2-\bar{x})/\bar{x}^2\), is

\widehat{\alpha}_{\mathrm{mom}} =\frac{8.1-2.4}{(2.4)^2} =\frac{5.7}{5.76} \approx 0.9896.

This is evidence to investigate overdispersion, not a final model diagnosis: the calculation is unadjusted for individual predictors and does not account for sampling uncertainty. A regression diagnostic should assess variation remaining conditional on the predictors and exposure. In particular, mixing individuals with different expected counts can itself make the marginal sample variance exceed the marginal mean even when the conditional Poisson model is appropriate.

Using the Model for Count Imputation

For an incomplete count variable in a MICE iteration, fit the conditional count model to cases where that count is observed, using the current available values of its predictors. If the negative-binomial model is supported, estimate or draw the regression coefficients and dispersion, calculate each missing case's conditional mean, and draw a nonnegative integer from its negative-binomial predictive distribution. The count draw, rather than the fitted mean, is inserted as the imputed value.

For valid multiple imputation, the process should reflect uncertainty in the regression coefficients and the dispersion as well as residual variation in the count distribution. The exact parameter-drawing procedure depends on the implementation; use an approach that propagates this uncertainty rather than treating fitted parameters as known. The conditional model is refitted as the chained-equations procedure updates the other incomplete variables.

Worked Example: Calculating a Predictive Distribution

Consider a hypothetical count of events with exposure \(t_i\) and a binary predictor \(x_i\). Suppose a fitted model is \(\log(\mu_i)=\log(t_i)+\beta_0+\beta_1x_i\), with \(\beta_0=\log(1.5)\), \(\beta_1=\log(0.8)\), and \(\alpha=0.5\). For a missing count with \(t_i=2\) and \(x_i=1\), calculate its mean and variance under the fitted parameters.

The conditional mean is

\widehat{\mu}_i =\exp\{\log(2)+\log(1.5)+\log(0.8)\} =2(1.5)(0.8)=2.4\text{ events}.

The dispersion parameter corresponds to \(r=1/\alpha=2\). Therefore the conditional variance is

\widehat{\operatorname{Var}}(Y_i\mid x_i,t_i) =2.4+0.5(2.4)^2 =2.4+2.88 =5.28\text{ events}^2.

For illustration, the probability of imputing zero events at these fixed parameter values is

P(Y_i=0\mid x_i,t_i) =\left(\frac{r}{r+\mu_i}\right)^r =\left(\frac{2}{2+2.4}\right)^2 \approx 0.2066.

This calculation describes one conditional predictive distribution; it is not itself a complete MICE imputation because it holds the fitted parameters fixed. An actual imputation procedure should also reflect parameter uncertainty before drawing the missing count.

Diagnostics, Alternatives, and Limitations

Overdispersion should be assessed after accounting for the mean structure, predictors, and exposure. A Pearson dispersion summary, for example, compares the sum of squared Pearson residuals with residual degrees of freedom. A value substantially above one can flag extra variation, but it is not a universal cutoff or proof that a negative-binomial model is correct. Residual patterns against fitted means or predictors can reveal a misspecified mean structure that a dispersion parameter alone will not repair.

Check whether observed and simulated or predicted counts have similar spread and whether the model reproduces important features such as the frequency of zeros and the upper tail. If the model systematically misses these features, inspect omitted predictors, exposure measurement, influential observations, dependence between observations, and the assumed variance form. A high zero frequency may have more than one explanation; overdispersion alone does not establish a separate zero-generating process.

Approach	What it allows	Implication for imputation
Poisson regression	Conditional variance equals the conditional mean.	Appropriate when the conditional Poisson distribution is adequate; can understate predictive spread when extra-Poisson variation remains.
Negative-binomial regression	Conditional variance is \(\mu_i+\alpha\mu_i^2\).	Provides a full count distribution for predictive draws while allowing extra-Poisson variation under its specified variance structure.
Quasi-Poisson regression	Allows a variance proportional to the mean without specifying a complete count probability distribution.	Can adjust variance-based inference, but by itself does not define a full predictive count distribution for imputations.

A negative-binomial model is not a cure for a misspecified conditional mean, strong unmodeled dependence, or a poor predictor set. Its variance form also imposes a particular relationship between mean and variance; it may not capture every count distribution. Finally, a large estimated \(\alpha\) should prompt careful checking rather than automatic acceptance: unusual counts, omitted structure, or a small observed sample can affect the estimate.

Common pitfall: “The sample variance exceeds the sample mean, so the Poisson model is wrong and a negative-binomial model must be used.” This compares marginal summaries and ignores predictors. The correction is to assess dispersion conditional on the modeled predictors and exposure, then check whether the negative-binomial mean and variance structure adequately represents the data. Overdispersion is evidence to investigate, not a model-selection rule by itself.

Specify the conditional mean.
Choose predictors and, when needed, an exposure offset that represent how expected counts vary across cases.

Check for remaining overdispersion.
Use conditional diagnostics and inspect residual patterns rather than relying only on a raw variance-to-mean ratio.

Fit and assess the negative-binomial model.
Estimate the mean-model coefficients and dispersion, and check whether the implied distribution is plausible for the observed counts.

Draw imputations with uncertainty.
Propagate uncertainty in model parameters and generate missing values from the conditional count distribution.

Worked Example: Comparing Poisson and Negative-Binomial Spread

Suppose a hypothetical missing count has fitted conditional mean \(\mu_i=2\), and the estimated negative-binomial dispersion is \(\alpha=0.75\). Compare the two models' conditional variances and probabilities of zero.

The Poisson model gives variance 2. The negative-binomial model gives

\operatorname{Var}(Y_i\mid x_i)=2+0.75(2)^2=5.

For the negative-binomial model, \(r=1/0.75=4/3\), so the probability of zero is

P(Y_i=0\mid x_i) =\left(\frac{r}{r+\mu_i}\right)^r =\left(\frac{4/3}{4/3+2}\right)^{4/3} =(0.4)^{4/3} \approx 0.294.

For comparison, the Poisson probability is \(P(Y_i=0\mid x_i)=e^{-2}\approx0.1353\). At the same conditional mean, the negative-binomial distribution has greater variance and, in this example, a higher zero probability. This illustrates why a distributional choice changes the values that may be imputed; it does not mean that the negative-binomial model is automatically preferable for every count.

Check Your Understanding

Answer each question, then check the brief solutions below.

A Poisson model has fitted conditional mean 3. What conditional variance does it imply? What variance does a negative-binomial model with \(\alpha=0.4\) imply?
For a count with mean 4 and variance 12, calculate the method-of-moments dispersion \(\alpha\approx(s^2-\bar{x})/\bar{x}^2\), treating the supplied values as moment estimates.
Why is a raw sample variance greater than a raw sample mean not, by itself, sufficient to establish conditional overdispersion?
What distributional information does a negative-binomial model provide for an incomplete count that a quasi-Poisson variance specification alone does not provide?

Solutions. (1) The Poisson variance is 3. The negative-binomial variance is \(3+0.4(3^2)=3+3.6=6.6\). (2) \(\alpha\approx(12-4)/4^2=8/16=0.5\). (3) The raw summaries do not condition on predictors or exposure; differing expected counts across cases can raise the marginal variance even if the conditional Poisson model is adequate. (4) The negative-binomial model specifies a full count probability distribution, enabling predictive count draws; a quasi-Poisson variance alone does not specify that distribution.