Quasi-Poisson Estimation and Inference

Why Quasi-Poisson Inference Is Different

Count data are often modeled with a mean that depends on predictors and a variance proportional to that mean. The Poisson model fixes the variance equal to the mean and specifies a full probability distribution. Quasi-Poisson regression instead requires the mean and variance relationship but does not assert a Poisson probability law. Its central practical consequence is precise: with a log link, estimating a common dispersion factor changes the estimated uncertainty, but not the fitted mean coefficients.

This separation is useful when counts vary more, or less, than the Poisson variance allows. The mean model can remain the same while the scale of uncertainty is adjusted. The calculations below derive both parts: the estimating equation for the coefficients and the covariance estimate used for inference.

Model and Estimating Equation

Definition: In a quasi-Poisson regression, the responses \(Y_1,\ldots,Y_n\) are independent, with \(E[Y_i]=\mu_i>0\) and \(\operatorname{Var}(Y_i)=\phi\mu_i\), where \(\phi>0\) is a common dispersion parameter. For a log-link regression with fixed covariate vectors \(x_i\in\mathbb{R}^p\), \(\log(\mu_i)=x_i^\mathsf{T}\beta\). With a known exposure \(t_i>0\), the model may instead be written \(\log(\mu_i)=\log(t_i)+x_i^\mathsf{T}\beta\).

The variance assumption specifies a first two moments model; it does not specify the distribution of \(Y_i\). In particular, the quasi-Poisson model is not a Poisson likelihood with an extra parameter. For inference based on the quasi-score, assume the mean model and variance relationship are appropriate, the observations are independent, and the design contains enough information to identify \(\beta\).

For a general differentiable mean \(\mu_i(\beta)\), the quasi-score contribution is the derivative of the mean multiplied by the residual and divided by its variance. Thus the estimating equation is

U(\beta) = \sum_{i=1}^n \frac{Y_i-\mu_i}{\phi\mu_i} \frac{\partial\mu_i}{\partial\beta} =0.

For the log link, \(\partial\mu_i/\partial\beta=\mu_i x_i\). The factor \(\mu_i\) cancels, leaving

U(\beta)=\frac{1}{\phi}\sum_{i=1}^n x_i(Y_i-\mu_i) =\frac{1}{\phi}X^\mathsf{T}(y-\mu)=0,

where \(X\) has rows \(x_i^\mathsf{T}\), and \(y\) and \(\mu\) are the vectors of observed and fitted counts. Since \(\phi\) is positive and common to every term, this equation is equivalent to \(X^\mathsf{T}(y-\mu)=0\). The same coefficient estimates therefore solve it for every value of \(\phi\). An exposure offset changes the \(\mu_i\) in the equation, but not this cancellation.

Worked Example: Estimating a Rate With Exposure

Suppose two independent count observations have exposures \(t_1=1\) and \(t_2=3\), and observed counts \(y_1=2\) and \(y_2=6\). Use an intercept-only log-rate model \(\mu_i=t_i\exp(\beta)\). Find the coefficient estimate and its model-based standard error if the dispersion is known to be \(\phi=2\).

The score equation is \(\sum_i(y_i-\mu_i)=0\). Here it becomes \(2+6-(1+3)\exp(\beta)=0\), so \(\exp(\hat\beta)=8/4=2\) and \(\hat\beta=\log(2)\). The fitted means are \(\hat\mu_1=2\) and \(\hat\mu_2=6\).

For an intercept-only log-link model, the information matrix for the coefficient is \(\sum_i\mu_i/\phi\). At the fitted values it is \(8/2=4\), so the estimated variance is \(1/4\) and the standard error is \(1/2\). The coefficient estimate comes from the mean equation; the specified dispersion determines its uncertainty.

Estimating Dispersion and the Covariance Matrix

A common estimate of \(\phi\) is the Pearson estimate. Given fitted means \(\hat\mu_i\), the Pearson residual contribution is \((y_i-\hat\mu_i)^2/\hat\mu_i\). If the fitted mean model has \(p\) coefficients, the residual-degrees-of-freedom estimate is

\hat\phi = \frac{1}{n-p} \sum_{i=1}^n\frac{(y_i-\hat\mu_i)^2}{\hat\mu_i}, \qquad n>p.

This estimates the scale of variation relative to the Poisson variance \(\mu_i\). The divisor \(n-p\) is a degrees-of-freedom adjustment; it does not make the estimate exactly unbiased in every finite-sample model. When \(\hat\phi\) is used for inference, the log-link quasi-Poisson covariance estimate is

\widehat{\operatorname{Cov}}(\hat\beta) = \hat\phi\left(X^\mathsf{T}\widehat W X\right)^{-1}, \qquad \widehat W=\operatorname{diag}(\hat\mu_1,\ldots,\hat\mu_n).

This formula requires \(X^\mathsf{T}\widehat W X\) to be invertible. For large-sample inference, it also relies on regularity conditions such as a correctly specified mean and variance relationship, a fixed number \(p\) of coefficients, a design with information increasing in all coefficient directions, and conditions ensuring a central limit theorem for the estimating equation.

Result: Under the stated conditions, the log-link quasi-Poisson coefficient estimator has model-based covariance approximately \(\phi(X^\mathsf{T}WX)^{-1}\), evaluated at the fitted means in practice. If \(\phi\) is estimated, replacing it by \(\hat\phi\) multiplies every coefficient variance by \(\hat\phi\), and every coefficient standard error by \(\sqrt{\hat\phi}\).

To derive the covariance, let \(W=\operatorname{diag}(\mu_1,\ldots,\mu_n)\) and consider the estimating equation \(U(\beta)=X^\mathsf{T}(Y-\mu)/\phi\). Its negative expected derivative with respect to \(\beta\) is

A = -\operatorname{E}\left[\frac{\partial U}{\partial\beta^\mathsf{T}}\right] = \frac{1}{\phi}X^\mathsf{T}WX.

Indeed, the derivative of \(\mu\) with respect to \(\beta^\mathsf{T}\) is \(WX\), and the residual has expectation zero. Independence and \(\operatorname{Var}(Y_i)=\phi\mu_i\) give the variability of the score:

B=\operatorname{Var}(U) = \frac{1}{\phi^2}X^\mathsf{T}\operatorname{Var}(Y)X = \frac{1}{\phi}X^\mathsf{T}WX.

The first-order expansion of the estimating equation around the true coefficient gives the usual estimating-equation covariance \(A^{-1}BA^{-1}\). Substituting the two matrices yields

A^{-1}BA^{-1} = \left\{\phi(X^\mathsf{T}WX)^{-1}\right\} \left\{\frac{1}{\phi}X^\mathsf{T}WX\right\} \left\{\phi(X^\mathsf{T}WX)^{-1}\right\} = \phi(X^\mathsf{T}WX)^{-1}.

The plug-in estimate replaces \(\phi\) and the unknown means by their estimates. This derivation uses the assumed variance relationship; it is not a claim that quasi-Poisson responses have a Poisson distribution. \(\square\)

Worked Example: A Two-Group Log-Link Regression

There are two observations in each of two groups. The counts are \(1,3\) in group A and \(2,6\) in group B. Fit the log-link model \(\log(\mu_i)=\beta_0+\beta_1 x_i\), where \(x_i=0\) for A and \(x_i=1\) for B. Calculate the coefficient estimates, Pearson dispersion estimate, and covariance estimate.

The estimating equations make the fitted total in each group equal its observed total. The group means are therefore \(\hat\mu_A=(1+3)/2=2\) and \(\hat\mu_B=(2+6)/2=4\). It follows that \(\hat\beta_0=\log(2)\), and \(\hat\beta_0+\hat\beta_1=\log(4)\), giving \(\hat\beta_1=\log(2)\).

The Pearson sum is \( (1-2)^2/2+(3-2)^2/2+(2-4)^2/4+(6-4)^2/4=0.5+0.5+1+1=3\). There are \(n=4\) observations and \(p=2\) coefficients, so \(\hat\phi=3/(4-2)=1.5\).

The fitted weighted design matrix is

X^\mathsf{T}\widehat W X = \begin{pmatrix} 12&8\\ 8&8 \end{pmatrix}, \qquad (X^\mathsf{T}\widehat W X)^{-1} = \begin{pmatrix} 0.25&-0.25\\ -0.25&0.375 \end{pmatrix}.

Multiplying by \(\hat\phi=1.5\) gives \(\widehat{\operatorname{Cov}}(\hat\beta)= \begin{pmatrix}0.375&-0.375\\-0.375&0.5625\end{pmatrix}\). Thus the standard error of the group effect \(\hat\beta_1\) is \(\sqrt{0.5625}=0.75\). The fitted group means remain 2 and 4: dispersion scales the covariance, not these fitted means.

What Changes, and What Does Not

The coefficient estimates depend on the mean model and the observed counts, not on a common positive \(\phi\). The fitted means, and consequently fitted rates or relative rates derived from those means, do not change when the common dispersion estimate changes. The covariance matrix does change: if \(\hat\phi=2\), standard errors are multiplied by \(\sqrt{2}\) relative to the \(\phi=1\) calculation; if \(\hat\phi=1/4\), they are multiplied by \(1/2\). Coefficient contrasts inherit the same scale change in their estimated variance.

This result is specific to a common dispersion factor in the log-link quasi-score. It does not mean that every modification to a count model leaves its fitted means unchanged. Changing the link, the mean structure, the offset, or the variance weights can change the estimating equation itself.

Worked Example: Overdispersion Changes a Standard Error

Consider an intercept-only log-link model for four independent counts \(0,0,4,4\). Calculate the fitted mean, Pearson dispersion estimate, and standard error for the log mean, then compare with the \(\phi=1\) standard error.

The estimating equation sets the fitted total equal to the observed total, so \(\hat\mu=(0+0+4+4)/4=2\). The Pearson sum is \(4(2-0)^2/2=8\), because each of the four squared residual contributions is \(4/2=2\). With \(n-p=4-1=3\), \(\hat\phi=8/3\).

For the intercept, \(X^\mathsf{T}\widehat W X=4(2)=8\), so the estimated variance of the log mean is \((8/3)/8=1/3\), and its standard error is \(1/\sqrt{3}\), approximately \(0.577\). With \(\phi=1\), the variance would be \(1/8\), and the standard error approximately \(0.354\). Both calculations use the same fitted mean, 2; the larger Pearson dispersion increases the standard error by \(\sqrt{8/3}\).

Common Pitfall

Wrong reasoning: “The dispersion estimate is larger than one, so the quasi-Poisson fit must have a larger fitted mean.” This confuses the mean equation with its covariance. In a log-link quasi-Poisson model, \(\phi\) is a common multiplier of the score, so it cancels from the equation for \(\hat\beta\). It enters the covariance after the coefficients and fitted means have been obtained.

Correction: First solve \(X^\mathsf{T}(y-\mu)=0\) for the mean coefficients. Then estimate or specify \(\phi\), and use it to scale the coefficient covariance. Also, do not interpret quasi-Poisson fitting as proof that the counts follow a Poisson distribution: only the specified mean and variance relationship is being used.

Check Your Understanding

Use the definitions and calculations above to check your understanding.

For a log-link quasi-Poisson regression, show why the common factor \(1/\phi\) does not change the coefficient estimating equation.
What assumptions about \(Y_i\) are made by the quasi-Poisson model, and what distributional claim is not made?
If \(\hat\phi=4\), by what factor are coefficient standard errors multiplied relative to the calculation with \(\phi=1\), assuming the fitted means are unchanged?
In a model with \(n=20\) observations and \(p=3\) coefficients, what divisor is used in the residual-degrees-of-freedom Pearson estimate?
Why does changing the common dispersion estimate not change fitted means, while changing the link or mean model may change them?