Worked Residual Diagnosis for Poisson Counts

Why compare several residuals?

A residual is a way to express the gap between an observed count and the count predicted by a model. In a Poisson GLM, that gap can be measured on different scales: in counts, relative to the model’s standard deviation, through the likelihood, or on the linear-predictor scale. These scales answer related but not identical questions. Comparing them helps reveal whether a large-looking discrepancy is substantial for its fitted mean and whether a pattern remains across observations.

We will use the Poisson distribution with a log link. For observation \(i\), let \(Y_i\) be the count, \(\hat{\mu}_i\) its fitted mean, and \(\hat{\eta}_i=\log(\hat{\mu}_i)\). The Poisson variance equals its mean, so the dispersion parameter is taken to be 1. The five residuals considered here are response, Pearson, deviance, working, and randomized quantile residuals.

Definition: For a Poisson GLM with fitted mean \(\hat{\mu}_i>0\) and observed count \(y_i\), the response residual is \(r_i=y_i-\hat{\mu}_i\); the Pearson residual is \(r_{P,i}=(y_i-\hat{\mu}_i)/\sqrt{\hat{\mu}_i}\); the signed deviance residual is \(r_{D,i}=\operatorname{sign}(y_i-\hat{\mu}_i)\sqrt{d_i}\), where \(d_i=2\{y_i\log(y_i/\hat{\mu}_i)-(y_i-\hat{\mu}_i)\}\), with \(0\log(0/\hat{\mu}_i)=0\); and the working residual for the log link is \(r_{W,i}=(y_i-\hat{\mu}_i)/\hat{\mu}_i\). A randomized quantile residual is \(r_{Q,i}=\Phi^{-1}(U_i)\), where \(U_i\) is drawn uniformly on \([F_i(y_i-1),F_i(y_i)]\), \(F_i\) is the fitted Poisson CDF, and \(\Phi\) is the standard normal CDF.

Three residual scales and their derivations

The response residual is the direct difference in counts. It is easy to interpret, but its variability depends on the mean: a Poisson count with mean 20 is naturally more variable than one with mean 2. The Pearson residual accounts for this by dividing by the fitted standard deviation \(\sqrt{\hat{\mu}_i}\). The working residual expresses the discrepancy on the linear-predictor scale used by the log-link model.

To see the working-residual formula, the Poisson GLM’s local linearization uses the derivative of the mean with respect to the linear predictor. For the log link, \(\mu_i=\exp(\eta_i)\), so \(d\mu_i/d\eta_i=\mu_i\). The first-order change in the mean for a change in predictor is therefore \(d\mu_i\approx\mu_i\,d\eta_i\). A response discrepancy \(y_i-\hat{\mu}_i\), divided by this derivative, gives \((y_i-\hat{\mu}_i)/\hat{\mu}_i\). This is a local scale conversion, not a claim that the working residual is an exact observed difference in \(\eta_i\).

Derivation: The Poisson deviance contribution compares the fitted model’s log likelihood for observation \(i\) with the log likelihood at the saturated mean \(\mu_i=y_i\). For \(y_i>0\), the log likelihood, apart from terms not depending on \(\mu_i\), is \(y_i\log(\mu_i)-\mu_i\). Subtracting the fitted log likelihood from the saturated one gives \(y_i\log(y_i/\hat{\mu}_i)-(y_i-\hat{\mu}_i)\); twice this difference is \(d_i\). At \(y_i=0\), the saturated mean is zero and the limiting contribution is \(d_i=2\hat{\mu}_i\). The signed square root preserves whether the observed count is above or below its fitted mean.

The deviance residual is not simply a standardized count difference. It reflects the change in likelihood relative to the saturated model, and its sign records the direction of that change. The randomized quantile residual takes a different approach: because a count distribution is discrete, it randomizes within the probability interval assigned to the observed count, then maps that value to a normal scale.

Worked calculations for one fitted mean

Worked Example: Comparing Five Residuals at a Fitted Mean of 3.5

Suppose a Poisson model gives \(\hat{\mu}_i=3.5\) for an observation with \(y_i=1\). The response residual is \(1-3.5=-2.5\). Since the Poisson standard deviation is \(\sqrt{3.5}\), the Pearson residual is \(-2.5/\sqrt{3.5}\approx-1.336\). The working residual is \(-2.5/3.5\approx-0.714\).

For the deviance residual, \(d_i=2\{1\log(1/3.5)-(1-3.5)\}\). Since \(\log(1/3.5)\approx-1.2528\), this is \(2\{-1.2528+2.5\}\approx2.4945\). The count is below its fitted mean, so \(r_{D,i}=-\sqrt{2.4945}\approx-1.579\).

For a randomized quantile residual, the fitted Poisson probabilities give \(F_i(0)\approx0.0302\) and \(F_i(1)\approx0.1359\). One valid random draw is \(U_i=0.10\), which lies in this interval. Therefore \(r_{Q,i}=\Phi^{-1}(0.10)\approx-1.282\). A different valid draw produces a different randomized residual, so the draw should be recorded if exact reproducibility is needed.

Worked Example: Positive Discrepancies at the Same Mean

Keep \(\hat{\mu}_i=3.5\), but take \(y_i=7\). The response residual is \(7-3.5=3.5\), the Pearson residual is \(3.5/\sqrt{3.5}\approx1.871\), and the working residual is \(3.5/3.5=1\).

For the deviance, \(d_i=2\{7\log(7/3.5)-(7-3.5)\}=2\{7\log(2)-3.5\}\). Using \(\log(2)\approx0.6931\), \(d_i\approx2(4.8520-3.5)=2.7041\), so \(r_{D,i}\approx\sqrt{2.7041}=1.644\).

Here \(F_i(6)\approx0.9347\) and \(F_i(7)\approx0.9733\). A valid randomized value is \(U_i=0.96\), yielding \(r_{Q,i}=\Phi^{-1}(0.96)\approx1.751\). The five values all indicate a count above expectation, but their magnitudes differ because the scales and constructions differ.

Worked diagnosis of a count dataset

Worked Example: Residuals Versus Predictor and Fitted Mean

Consider eight observations with predictor \(x=1,\ldots,8\) and counts \(y=(1,1,2,2,4,4,7,7)\). Fit an intercept-only Poisson model. Its log likelihood, ignoring terms not depending on the mean, is \(\ell(\mu)=\sum_i\{y_i\log(\mu)-\mu\}\). Differentiating gives \(\ell'(\mu)=\sum_i y_i/\mu-8\), so the maximum occurs at \(\hat{\mu}=\sum_i y_i/8=28/8=3.5\). Thus every observation has fitted mean 3.5.

The following table gives four residuals, rounded to three decimals. Rows with the same count have the same values because every fitted mean is identical.

x	y	Response	Pearson	Deviance	Working
1	1	-2.500	-1.336	-1.579	-0.714
2	1	-2.500	-1.336	-1.579	-0.714
3	2	-1.500	-0.802	-0.873	-0.429
4	2	-1.500	-0.802	-0.873	-0.429
5	4	0.500	0.267	0.261	0.143
6	4	0.500	0.267	0.261	0.143
7	7	3.500	1.871	1.644	1.000
8	7	3.500	1.871	1.644	1.000

For example, at \(y=2\), the deviance contribution is \(2\{2\log(2/3.5)-(2-3.5)\}\approx0.762\), giving the signed residual \(-\sqrt{0.762}\approx-0.873\). The negative sign follows because 2 is below 3.5.

Against fitted means, the residual plot is a vertical stack at 3.5. This happens because the model fits the same mean to every observation; that plot cannot show how residuals vary with \(x\). Against the predictor, however, residuals move from negative at small \(x\) to positive at large \(x\). The model leaves a clear increasing pattern unexplained. A model allowing the mean to change with \(x\) is worth investigating, but this plot alone does not establish which alternative model is appropriate.

The randomized quantile residuals also vary with \(x\), but each is random within its count-specific probability interval. That randomization can make the plotted points less tied and more comparable to a continuous normal scale; it does not remove a systematic relationship between the distribution of counts and the predictor.

What the residual plots can and cannot show

A residual-versus-fitted plot checks whether residual behavior changes across the model’s predicted means. Look for a systematic curve, changing spread, or clusters that suggest the mean structure or distributional assumptions may be inadequate. A residual-versus-predictor plot checks whether a covariate retains a relationship with the response after fitting the model. It can reveal structure that a fitted-value plot hides, especially when many observations share similar fitted means.

For Poisson data, the response residuals are discrete and their spread naturally varies with the mean. Pearson residuals account for the fitted Poisson standard deviation, while deviance residuals express discrepancy through likelihood. Neither transformation guarantees independent, exactly normal residuals in a fitted model. Randomized quantile residuals have a standard normal reference distribution when the observations are independent and the fitted distributions are the true distributions; with parameters estimated from the data, the residuals are not generally independent standard normal variables.

Result: If \(Y\) has a discrete CDF \(F\), define \(U=F(Y-1)+V\{F(Y)-F(Y-1)\}\), where \(V\) is independent Uniform\((0,1)\). Then \(U\) is Uniform\((0,1)\), and \(\Phi^{-1}(U)\) is standard normal. To verify the first claim, condition on \(Y=y\): \(U\) is uniform on the interval \([F(y-1),F(y)]\), whose length is \(P(Y=y)\). These intervals partition \([0,1]\), and each interval receives probability equal to its length. Thus for any subinterval \(A\subseteq[0,1]\), summing the conditional probabilities over the intervals gives \(P(U\in A)=\operatorname{length}(A)\). Applying the standard normal quantile function to a uniform variable gives a standard normal variable. \(\square\)

Common mistake

Wrong reasoning: “A Pearson residual of 1.5 and a deviance residual of 1.5 are the same amount of evidence, because the numbers match.” Correction: The residual types use different transformations. Pearson residuals divide the count difference by the fitted standard deviation; deviance residuals take a signed square root of a likelihood contribution. Compare patterns and investigate observations in context, rather than treating equal numerical values as interchangeable measures of evidence.

Another practical caution is that randomized quantile residuals depend on random draws. A changed random seed can change individual points, even when the data and fitted model are unchanged. Record the seed or repeat the diagnostic to check that a conclusion does not depend on one realization.

Check Your Understanding

Use the formulas and plotting logic above to answer these questions.

For \(y=0\) and fitted mean \(\hat{\mu}=2\), calculate the response, Pearson, working, and deviance residuals.
Why can a residual-versus-fitted plot fail to reveal a pattern that appears against a predictor?
For a Poisson observation with \(y=3\) and fitted mean 3, what are its response, Pearson, and working residuals? What is the deviance residual?
Why does a randomized quantile residual for a count need a random draw within a CDF interval?
In the eight-observation example, what feature of the residual-versus-predictor plot suggests that the intercept-only model leaves structure unexplained?