Residual Reference Behavior Under a Correct GLM

What Does a Residual’s Reference Distribution Mean?

A residual is often compared with zero, with a unit scale, or with a standard normal distribution. These comparisons do not all make the same claim. Correctly specifying a GLM’s conditional mean and variance gives exact moment statements for response and Pearson residuals. It does not, by itself, make those residuals normal. Normal approximations need additional conditions, while randomized quantile residuals have an exact standard normal distribution when the full conditional response distribution is known.

To separate these claims, first imagine that the covariates and the true model parameters are fixed and known. For observation \(i\), let \(Y_i\) have conditional mean \(\mu_i\) and conditional variance \(\phi V(\mu_i)\), where \(\phi>0\) is the dispersion parameter and \(V\) is the variance function. In applications, residuals usually use fitted rather than true parameters; that distinction will matter at the end.

Exact Moment Statements

The first result needs only the correct conditional mean and variance. It does not require a particular exponential-family distribution, nor does it imply a normal distribution.

Proof. Subtracting the fixed number \(\mu_i\) from \(Y_i\) gives \(E[e_i\mid x_i]=E[Y_i\mid x_i]-\mu_i=0\), and \(\operatorname{Var}(e_i\mid x_i)=\operatorname{Var}(Y_i\mid x_i)=\phi V(\mu_i)\). The denominator of \(r_{P,i}\) is a fixed positive number conditional on \(x_i\). Therefore \(E[r_{P,i}\mid x_i]=0/\sqrt{\phi V(\mu_i)}=0\), and \(\operatorname{Var}(r_{P,i}\mid x_i)=\phi V(\mu_i)/\{\phi V(\mu_i)\}=1\). \(\square\)

The theorem gives an exact center and scale, but it does not determine the shape of the residual distribution. For example, a standardized Bernoulli variable has only two possible values, so it cannot have a normal distribution. Nor does the result assert independence across residuals; that requires a condition on the joint response distribution.

Why Deviance Residuals Often Resemble Pearson Residuals

The connection between Pearson and deviance residuals is local: when an observation is close to its fitted mean, the likelihood-ratio deviance has a quadratic approximation. The approximation is most useful when the data provide a regime in which the standardized response becomes close to normal and the deviance expansion’s higher-order terms become negligible.

For clarity, consider an exponential dispersion model with log-likelihood \(\ell(y;\theta,\phi)=\{y\theta-b(\theta)\}/\phi+c(y,\phi)\), where \(b\) is three times continuously differentiable near the relevant parameter values, \(b''(\theta)>0\), and \(\mu=b'(\theta)\). Then \(V(\mu)=b''(\theta)\). For a possible observed value \(y\) in the model’s mean domain, let \(\theta_y\) solve \(b'(\theta_y)=y\). The likelihood-ratio deviance contribution is \(D(y,\mu)=2\{\ell(y;\theta_y,\phi)-\ell(y;\theta,\phi)\}\); the common term \(c(y,\phi)\) cancels.

Derivation. Define \(g(u)=u\theta_u-b(\theta_u)\), where \(b'(\theta_u)=u\). Differentiation gives \(g'(u)=\theta_u\), because the terms involving \(d\theta_u/du\) cancel. Differentiating \(b'(\theta_u)=u\) gives \(d\theta_u/du=1/b''(\theta_u)=1/V(u)\). Hence \(g''(u)=1/V(u)\). The deviance is \(D(y,\mu)=2\{g(y)-g(\mu)-g'(\mu)(y-\mu)\}/\phi\). Taylor’s theorem at \(\mu\) therefore gives \(D(y,\mu)=(y-\mu)^2/\{\phi V(\mu)\}+O(|y-\mu|^3)\), since \(g'''\) is bounded on a sufficiently small closed neighborhood of \(\mu\).

To obtain the signed-square-root approximation, write \(h=y-\mu\) and \(a=\phi V(\mu)>0\). The expansion has the form \(D=h^2\{1/a+O(|h|)\}\). By differentiability of the square-root function at \(1/a>0\), \(\sqrt{D}=|h|\{1/\sqrt a+O(|h|)\}\). Multiplying by \(\operatorname{sign}(h)\) gives \(h/\sqrt a+O(h|h|)\), the Pearson residual plus a smaller local term. At \(y=\mu\), both residuals are zero. This proves a local approximation, not an exact identity for arbitrary observations. \(\square\)

When Is a Normal Reference Reasonable?

A normal reference for Pearson residuals requires more than correct first two moments. It is reasonable when the conditional response itself has a suitable normal approximation. One precise way to state the needed condition for a sequence of models is that, conditional on the covariates, the standardized variables \((Y_i-\mu_i)/\sqrt{\phi V(\mu_i)}\) converge in distribution to \(N(0,1)\) in the regime being considered. For independent responses, familiar examples include Poisson means tending to infinity and binomial counts with both expected successes and expected failures growing. The details of a normal approximation depend on the response family and the sequence of parameters.

In such a regime, a Pearson residual using the true mean and variance inherits that approximate normal reference. The local deviance expansion can then make the signed deviance residual approximately equivalent as well, provided its remainder is negligible on the scale of the fluctuations. For instance, a Poisson observation with large mean \(\mu\) fluctuates on the order of \(\sqrt{\mu}\), while the cubic remainder in its deviance is of smaller order after standardization. At a fixed small mean, that argument need not apply.

These statements concern residuals formed with true parameters. Replacing them with estimates changes the distribution: fitted residuals share estimated quantities, and fitting can impose constraints on their aggregate values. Consequently, fitted residuals generally are not exactly independent standard normal variables, even when the model is correct. Their approximation may also depend on sample size, the number of estimated parameters, and the design. A unit-variance calculation for a true-parameter Pearson residual cannot simply be transferred to every fitted residual.

Randomized Quantile Residuals: An Exact Distributional Reference

Unlike Pearson and deviance residuals, randomized quantile residuals use the entire conditional response distribution. For a discrete response, define \(F_i(y)=P(Y_i\le y\mid x_i)\), \(F_i(y-)=P(Y_i<y\mid x_i)\), and let \(V_i\sim\operatorname{Uniform}(0,1)\) be independent of the response. Set \(U_i=F_i(Y_i-)+V_iP(Y_i\mid x_i)\), and \(R_{Q,i}=\Phi^{-1}(U_i)\), where \(\Phi\) is the standard normal cumulative distribution function. The randomization fills the probability jump at the observed response value.

Proof for a discrete response. Fix \(x_i\), and write \(p_k=P(Y_i=k\mid x_i)\), \(F_k=P(Y_i\le k\mid x_i)\), and \(F_{k-}=P(Y_i<k\mid x_i)\). Conditional on \(Y_i=k\), \(U_i\) is uniform on the interval from \(F_{k-}\) to \(F_k\). For any \(t\in[0,1]\), choose a support value \(k\) such that \(F_{k-}\le t\le F_k\). The probability contributed by outcomes below \(k\) is \(F_{k-}\). The probability contributed by \(Y_i=k\) and \(U_i\le t\) is \(p_k(t-F_{k-})/p_k=t-F_{k-}\) when \(p_k>0\). Outcomes above \(k\) cannot contribute. Thus \(P(U_i\le t\mid x_i)=F_{k-}+t-F_{k-}=t\). Values at interval endpoints give the same result by continuity of the uniform distribution. Therefore \(U_i\) is uniform. Since \(\Phi\) is continuous and strictly increasing, \(P(\Phi^{-1}(U_i)\le z\mid x_i)=P(U_i\le\Phi(z)\mid x_i)=\Phi(z)\), proving the normal claim. \(\square\)

Common Mistakes and Practical Meaning

• Confusing unit variance with normality. The Bernoulli example has Pearson mean zero and variance one exactly, but its distribution has only two points. Moments alone do not specify a distribution.
• Treating a deviance residual as exactly Pearson. The quadratic expansion is local. In general, the likelihood-ratio deviance and Pearson statistic can differ, as the Poisson calculation shows.
• Assuming every correct GLM gives normal residuals. A correct conditional mean and variance establish exact Pearson moments, not a normal approximation. The response distribution and its asymptotic regime matter.
• Calling fitted residuals exact standard normals. Exact reference statements above use true conditional distributions and parameters. Estimation changes residual dependence and distribution.
• Ignoring what randomized quantile residuals require. Their exact normal reference uses the full correct conditional distribution and fresh independent randomization, not merely a correct mean and variance.

These distinctions matter because diagnostic comparisons are meaningful only when the reference is appropriate. A standard normal reference is particularly well justified for randomized quantile residuals under a known correct distribution, and can be an approximation for Pearson or deviance residuals in suitable large-sample regimes. For sparse counts, small binomial denominators, or fitted models with substantial estimation effects, treating normality as exact can make ordinary discreteness or fitting constraints look like model failure.