Integrated Count Effects and Prediction Workshop

From Regression Effects to Count Predictions

A Poisson regression analysis often answers several related questions in sequence. How does an exposure effect differ between two groups? What is the expected count for a particular covariate profile? And what range of values might a future count take? These questions concern different quantities: a rate-ratio contrast, a conditional mean, and a random future observation. Keeping those quantities distinct is essential for interpreting the calculations correctly.

We work through one model and one set of fitted values throughout. Let \(Y_i\) be a count observed over exposure \(t_i>0\), and let \(A_i,B_i\) be binary covariates taking values zero or one. Assume independent observations with a Poisson distribution and a log link:

Y_i\mid A_i,B_i,t_i \sim \operatorname{Poisson}(\mu_i),\qquad \log(\mu_i)=\log(t_i)+\beta_0+\beta_A A_i+\beta_B B_i+\beta_{AB}A_iB_i.

The offset \(\log(t_i)\) makes the model describe a rate per exposure unit; when \(t_i=1\), the mean count and rate have the same numerical value. The interaction coefficient allows the effect of \(A\) to depend on \(B\). For the worked calculations, take the fitted coefficient vector, in the order \((\hat\beta_0,\hat\beta_A,\hat\beta_B,\hat\beta_{AB})\), to be

\hat{\boldsymbol\beta} = \left(\log 2,\ \log 1.5,\ \log 1.2,\ \log(83/45)\right), \qquad \widehat{\operatorname{Var}}(\hat{\boldsymbol\beta}) = \begin{pmatrix} 0.0400&0&0&0\\ 0&0.0225&0&0\\ 0&0&0.0100&0\\ 0&0&0&0.0625 \end{pmatrix}.

This covariance matrix is an illustrative fitted covariance estimate. Its diagonal entries give standard errors \(0.20,0.15,0.10,0.25\), respectively. Wald intervals below use the usual large-sample normal approximation for the fitted coefficients and the supplied covariance estimate; they are not exact finite-sample intervals.

Interaction-Specific Rate Ratio

Definition: In this model, the rate ratio comparing \(A=1\) with \(A=0\) at a fixed value \(B=b\) is \(\exp(\beta_A+\beta_{AB}b)\). The interaction-specific ratio of these rate ratios, comparing \(B=1\) with \(B=0\), is \(\exp(\beta_{AB})\).

To see why, the rates for \(A=1\) and \(A=0\) at a common exposure and fixed \(B=b\) are proportional to \(\exp(\beta_0+\beta_A+\beta_Bb+\beta_{AB}b)\) and \(\exp(\beta_0+\beta_Bb)\). Their ratio is \(\exp(\beta_A+\beta_{AB}b)\). Taking the ratio of the \(A\)-rate ratios at \(b=1\) and \(b=0\) cancels \(\exp(\beta_A)\), leaving \(\exp(\beta_{AB})\).

Worked Example: Estimating the Interaction Ratio

At \(B=0\), the estimated rate ratio for \(A=1\) versus \(A=0\) is

\exp(\hat\beta_A)=\exp(\log 1.5)=1.5.

At \(B=1\), it is

\exp(\hat\beta_A+\hat\beta_{AB}) =1.5\left(\frac{83}{45}\right) =\frac{83}{30} \approx 2.767.

Thus the interaction-specific ratio of rate ratios is \((83/30)/1.5=83/45\approx1.844\). Equivalently, \(\exp(\hat\beta_{AB})=83/45\). In this fitted model, the rate ratio associated with changing \(A\) from zero to one is about 1.84 times as large when \(B=1\) as when \(B=0\). This is a multiplicative comparison of rate ratios, not a difference of rates.

For a scalar coefficient contrast \(c^\mathsf{T}\hat{\boldsymbol\beta}\), its estimated variance is \(c^\mathsf{T}\widehat{\operatorname{Var}}(\hat{\boldsymbol\beta})c\). Here \(c\) selects \(\hat\beta_{AB}\), so the estimated standard error is \(\sqrt{0.0625}=0.25\). A 95% Wald interval on the log ratio scale is \(\hat\beta_{AB}\pm1.96(0.25)\); exponentiating its endpoints gives an interval on the rate-ratio scale.

Derived Wald interval: If a fitted coefficient contrast \(\hat\psi\) has estimated standard error \(s\), the usual two-sided 95% Wald interval for \(\exp(\psi)\) is \(\big(\exp(\hat\psi-1.96s),\exp(\hat\psi+1.96s)\big)\). This follows by forming the normal-approximation interval on the \(\psi\) scale and applying the increasing function \(\exp(\cdot)\) to both endpoints.

Worked Example: A Confidence Interval for the Interaction

The estimated log interaction is \(\hat\beta_{AB}=\log(83/45)\approx0.6122\), and its standard error is \(0.25\). The lower and upper endpoints on the log scale are

0.6122-1.96(0.25)=0.1222,\qquad 0.6122+1.96(0.25)=1.1022.

Exponentiating,

\left(\exp(0.1222),\exp(1.1022)\right) \approx (1.13,3.01).

The estimated interaction ratio is about 1.84, and the 95% Wald interval is approximately 1.13 to 3.01. Under the large-sample approximation, this interval excludes one, the value corresponding to no multiplicative modification of the \(A\) rate ratio across the two levels of \(B\). The interpretation remains conditional on the stated model and its covariates.

Fitted Mean at a Covariate Profile

For a profile \((A,B,t)\), substitute its values into the fitted linear predictor and exponentiate to obtain the fitted conditional mean. This is a point estimate of \(E[Y\mid A,B,t]\), not a realized count and not, by itself, a statement about the range of future observations.

Worked Example: The Mean for Profile \(A=1,B=1,t=1\)

For the profile \(A=1,B=1,t=1\), the fitted log mean is

\hat\eta =\log(1)+\log 2+\log 1.5+\log 1.2+\log(83/45) =\log\left(2\times1.5\times1.2\times\frac{83}{45}\right).

The product inside the logarithm is \(3.6(83/45)=6.64\), so \(\hat\eta=\log(6.64)\) and the fitted mean is \(\hat\mu=\exp(\hat\eta)=6.64\) counts per exposure unit. The arithmetic verifies the result directly from the fitted coefficients.

The coefficient covariance estimate also gives an approximate standard error for the fitted log mean at this profile. The linear predictor uses all four coefficients with coefficient one, so its estimated variance is

0.0400+0.0225+0.0100+0.0625=0.1350, \qquad \operatorname{SE}(\hat\eta)=\sqrt{0.1350}\approx0.3674.

A 95% Wald interval for the mean is obtained by exponentiating the interval for \(\eta\):

\left( \exp\{\log(6.64)-1.96(0.3674)\}, \exp\{\log(6.64)+1.96(0.3674)\} \right) \approx (3.23,13.64).

This is an interval for the model’s conditional mean at that profile. It reflects estimated coefficient uncertainty under the Wald approximation; it is not an interval that is intended to contain a particular future count with 95% probability.

Predicting a Future Count

Suppose a new, conditionally independent count is observed at the same profile and its model parameters are treated as fixed at the fitted values. Its plug-in predictive distribution is \(\operatorname{Poisson}(6.64)\). A prediction interval must account for the discreteness of this distribution, rather than merely using a symmetric interval around its mean.

Definition: For a discrete random variable with cumulative distribution function \(F(k)=P(Y\leq k)\), an equal-tail interval with nominal tail probability \(\alpha/2\) uses integer endpoints \(L,U\) chosen so that \(F(L-1)\leq\alpha/2\), \(F(L)\geq\alpha/2\), \(F(U-1)<1-\alpha/2\), and \(F(U)\geq1-\alpha/2\). Discreteness can make its actual coverage exceed the nominal level.

For a Poisson variable with mean \(\lambda\), its probabilities satisfy \(p_0=\exp(-\lambda)\) and \(p_k=p_{k-1}\lambda/k\) for \(k\geq1\). Indeed, from \(p_k=\exp(-\lambda)\lambda^k/k!\), division by \(p_{k-1}\) gives \(p_k/p_{k-1}=\lambda/k\). Summing these recursively computed probabilities gives the cumulative probabilities used to locate the integer endpoints.

Worked Example: A 95% Plug-In Prediction Interval

Take \(\lambda=6.64\). Recursively computing the Poisson probabilities and summing them gives, to four decimal places,

F(1)\approx0.0100,\quad F(2)\approx0.0388,\quad F(11)\approx0.9613,\quad F(12)\approx0.9813.

For the lower endpoint, \(F(1)\approx0.0100<0.025\), whereas \(F(2)\approx0.0388\geq0.025\), so the lower equal-tail endpoint is \(L=2\). For the upper endpoint, \(F(11)\approx0.9613<0.975\), whereas \(F(12)\approx0.9813\geq0.975\), so \(U=12\). Thus the equal-tail 95% plug-in prediction interval is \([2,12]\).

Its actual coverage under this discrete plug-in distribution is

P(2\leq Y_{\mathrm{new}}\leq12) =F(12)-F(1) \approx0.9813-0.0100 =0.9713.

Using more precise cumulative probabilities gives approximately \(0.97134\). The coverage is above 0.95 because integer endpoints cannot generally allocate exactly 0.025 probability to each tail. The interval describes future-count variation conditional on the fitted mean being treated as fixed; it does not incorporate uncertainty in estimating the regression coefficients.

Putting the Quantities Together

The calculations answer distinct questions using the same fitted model. The interaction contrast is \(\exp(\hat\beta_{AB})\), with uncertainty determined by the variance of that coefficient contrast. The profile-specific mean comes from substituting covariates into the linear predictor and exponentiating. A plug-in prediction interval instead uses the entire Poisson distribution for a future observation at that profile.

Define the target contrast or profile.
State which rate ratio is being compared, or give the covariate and exposure values for the mean or future count.

Calculate on the log-mean scale.
Form the required coefficient contrast for an effect, or substitute the profile into the fitted linear predictor for a mean.

Transform and represent uncertainty appropriately.
Exponentiate a log-scale contrast or predictor; use a coefficient-based interval for a mean and a count distribution for a future observation.

A full predictive analysis can also account for parameter uncertainty, rather than fixing the mean at its fitted value. The interval calculated above is explicitly the plug-in version. Mixing over coefficient uncertainty generally changes predictive probabilities, so it should not be silently claimed that the plug-in interval incorporates that additional source of uncertainty.

Common Mistake

A common error is to report the confidence interval for the fitted mean, such as \((3.23,13.64)\), as a 95% prediction interval for a future count. The reasoning might be: “The true mean is likely to lie in this interval, so the next count should lie there too.” The conclusion does not follow. Even if the mean were known exactly to be 6.64, a Poisson count would still vary around that mean; moreover, a count is integer-valued. A prediction interval must use the predictive distribution of the count, as in the calculation yielding \([2,12]\). Conversely, that plug-in count interval does not quantify uncertainty in the fitted mean. Always name the target: a rate-ratio contrast, conditional mean, or future observation.

Key takeaway: A log-link model makes multiplicative rate comparisons and profile-specific means straightforward, but their uncertainty calculations answer different questions. Use the coefficient covariance for a contrast or fitted mean, and use the future-count distribution for a prediction interval.

Check Your Understanding

Try these questions using the model and fitted values above.

At \(B=0\), what is the fitted rate ratio for \(A=1\) versus \(A=0\)? What is the corresponding rate ratio at \(B=1\)?
Why does the interaction-specific ratio of rate ratios equal \(\exp(\beta_{AB})\), rather than \(\exp(\beta_A+\beta_{AB})\)?
What quantity is covered by the interval \((3.23,13.64)\), and why is it not a future-count prediction interval?
For a discrete distribution, why can an equal-tail interval with nominal 95% confidence have actual coverage greater than 95%?
What additional uncertainty is omitted when using \(\operatorname{Poisson}(6.64)\) as a plug-in predictive distribution?