Residual Plots for Binary GLMMs

Why Plot Event Rates for a Binary GLMM?

A binary GLMM can produce very small fitted event probabilities, especially when the outcome is rare. In that setting, individual residuals are mostly close to zero or one of a few possible values; a plot of each residual alone can be difficult to read. Grouping observations with similar fitted probabilities gives a more direct question: among observations assigned similar event probabilities, does the observed event rate agree with the model's average prediction?

This is a diagnostic for systematic calibration error, not a test that proves a model is correct. A persistent gap between observed and predicted rates may indicate that predictions are too high or too low in some range. But random variation, especially with few events, can also create large-looking gaps. Because responses in the same cluster can be dependent, the plot must be interpreted using the GLMM's marginal predictions and an understanding of cluster-level uncertainty.

Definition: For binary response $Y_{ij}$, let $p_{ij}^{M}=P(Y_{ij}=1\mid \mathbf{x}_{ij})$ denote the model-implied marginal fitted event probability, obtained by averaging the conditional event probability over the random-effects distribution. For a bin $B_g$ containing $n_g$ observations, define its observed event rate and mean fitted event probability as $$ O_g=\frac{1}{n_g}\sum_{(i,j)\in B_g}Y_{ij}, \qquad E_g=\frac{1}{n_g}\sum_{(i,j)\in B_g}p_{ij}^{M}. $$ The binned residual is $D_g=O_g-E_g$. A plot of $O_g$ against $E_g$ compares observed and predicted event rates; perfect agreement lies on the diagonal $O_g=E_g$.

The superscript $M$ emphasizes that these predictions average over random effects. A conditional fitted probability instead plugs in an estimated or specified cluster effect. Those are different prediction targets. For a diagnostic asking whether the model predicts population-average event rates across clusters, marginal fitted probabilities are the coherent comparison. The earlier material on Residuals Marginalized Over Random Effects distinguishes these targets.

Constructing and Reading the Plot

Calculate $p_{ij}^{M}$ for each observation, sort observations by that value, and form bins with similar fitted probabilities. Bins might have equal numbers of observations, or be defined by scientifically meaningful probability ranges. For each bin, plot $E_g$ on the horizontal axis and $O_g$ on the vertical axis. A companion plot of $D_g$ against $E_g$ centers the reference pattern at zero. Always report or display $n_g$ and, where useful, the number of observed events in each bin: the same rate difference has quite different evidential weight when it comes from two events than when it comes from twenty.

A point above the diagonal has $O_g>E_g$: the events occurred more often than the model predicted in that bin. A point below the diagonal has $O_g<E_g$: the model predicted more events than were observed. A sequence of points mostly above the line over a probability range suggests underprediction there; a sequence mostly below suggests overprediction. One isolated point far from the line is less persuasive, particularly if its bin contains few events.

Binning is a visualization choice, not a way to create additional information. Too few bins can hide a local pattern; too many create noisy rates, often zero in rare-event data. Bin boundaries based on fitted values can also put observations from the same cluster into several bins. Consequently, the points are not generally independent, and ordinary binomial standard errors that treat every observation as independent can be too small.

Worked Example: Comparing Three Bins of Rare-Event Predictions

Suppose a hypothetical binary GLMM is fitted to 60 observations from clustered units. Marginal fitted probabilities are grouped into three bins of 20 observations each. The observed event counts are 1, 5, and 7; the corresponding sums of fitted probabilities are 0.80, 1.20, and 6.00.

For the first bin, the observed rate is $O_1=1/20=0.050$, and the mean fitted probability is $E_1=0.80/20=0.040$. Thus $D_1=0.050-0.040=0.010$.

For the second bin, $O_2=5/20=0.250$, while $E_2=1.20/20=0.060$. Thus $D_2=0.250-0.060=0.190$.

For the third bin, $O_3=7/20=0.350$, and $E_3=6.00/20=0.300$. Thus $D_3=0.350-0.300=0.050$.

The second bin has the clearest apparent discrepancy: its observed event rate is 19 percentage points above its mean prediction. The first and third bins are closer to their predictions. This pattern warrants investigation, but it is not by itself a significance test. In particular, five observed events in the second bin and dependence among clustered observations affect how surprising the difference is.

Why Dependence Changes the Uncertainty

The point estimate $O_g-E_g$ is easy to calculate, but its uncertainty cannot generally be calculated as if all binary responses were independent Bernoulli trials. Let $S_g=\sum_{(i,j)\in B_g}Y_{ij}$. Under a correctly specified marginal model, $E(S_g)=\sum_{B_g}p_{ij}^{M}$, so $E(D_g)=0$ when the fitted probabilities are treated as known. The variance, however, includes covariance between observations in the same cluster.

The variance identity follows by expanding the variance of a sum. Dividing by $n_g^2$ gives the variance of the observed rate, since the fitted probabilities are fixed in this calculation:

\operatorname{Var}(D_g) = \operatorname{Var}(O_g) = \frac{1}{n_g^2} \left[ \sum_{(i,j)\in B_g}p_{ij}^{M}(1-p_{ij}^{M}) + 2\sum_{\substack{a<b\\a,b\in B_g}} \operatorname{Cov}(Y_a,Y_b) \right].

Here $a$ and $b$ index observations in the bin. The first sum is the Bernoulli variance contribution; the second adds the covariance contributions for pairs. For observations in different independent clusters, those covariances may be zero under the model. Within a cluster, shared random effects commonly induce positive covariance. Parameter estimation adds further uncertainty, which this formula does not include. Thus the formula explains why a simple independent-binomial error bar is not generally adequate; it does not, by itself, provide a complete uncertainty interval for a fitted-model diagnostic.

Worked Example: Independent and Cluster-Aware Rate Uncertainty

Consider a hypothetical bin of 20 observations, each with marginal fitted probability $0.05$. Suppose they form 10 clusters of two observations, and the two responses in each cluster have marginal correlation $0.20$. The expected event rate is $E_g=0.05$. For one response, the variance is $0.05(0.95)=0.0475$. Within a cluster the covariance is $0.20\sqrt{0.0475(0.0475)}=0.0095$.

There are 10 within-cluster pairs, so the variance of the total event count is $20(0.0475)+2(10)(0.0095)=1.14$. The approximate standard error of the observed rate is therefore $\sqrt{1.14}/20=0.0534$.

If the 20 observations were incorrectly treated as independent, the standard error would be $\sqrt{20(0.0475)}/20=\sqrt{0.95}/20=0.0487$, which is smaller. If this bin happened to contain 3 events, then $O_g=3/20=0.15$, and $D_g=0.15-0.05=0.10$. The difference is about $0.10/0.0534=1.87$ cluster-aware standard errors using the simplified calculation. This is a scale comparison, not a formal $p$-value or a complete uncertainty assessment; it illustrates why treating correlated responses as independent can overstate the apparent precision.

Model Assumptions and Prediction Target

A meaningful plot assumes that the binary outcome, covariates, link, random-effects structure, and random-effects distribution used to obtain the marginal predictions are appropriate for the intended prediction target. The observations must be assigned to the correct clusters, and the response must be coded consistently so that $Y_{ij}=1$ denotes the event of interest. If the model omits an important predictor or uses an unsuitable random-effects distribution, the resulting predictions may be systematically wrong even if the plotting arithmetic is correct.

The plotting procedure also assumes that each marginal prediction is computed for the covariates of its observation and averaged over the random-effects distribution, rather than evaluated only at a typical random effect. In a nonlinear model, averaging conditional probabilities is not generally the same as applying the inverse link to an average linear predictor. The observed rate is a sample quantity; it should be compared with the mean of the individual marginal probabilities, not with an average linear predictor on the link scale.

Worked Example: Conditional and Marginal Predictions Can Differ

Suppose a hypothetical random-intercept logistic GLMM gives an average conditional fitted probability of $0.020$ for ten observations when the estimated cluster effects are plugged in. Integrating over the random-intercept distribution instead gives a mean marginal fitted probability of $0.050$. One event is observed among the ten responses, so the observed rate is $1/10=0.10$.

Against the conditional prediction, the rate difference is $0.10-0.02=0.08$. Against the marginal prediction, it is $0.10-0.05=0.05$. These are different comparisons because the predictions answer different questions. For assessing population-average event rates across clusters, the marginal comparison is the relevant one. The numerical predictions here are hypothetical outputs; they illustrate the distinction rather than asserting that marginal predictions must always exceed conditional predictions.

Ten observations and one event provide little information about a stable pattern. The calculation identifies the direction and size of the discrepancy, but it cannot establish that the model is miscalibrated without considering uncertainty and the broader pattern across bins.

Diagnostics, Alternatives, and Limitations

Use the plot as one part of a diagnostic workflow, alongside checks of model specification and predictive behavior. A few practical warning signs are:

Many empty bins or bins with one event. Rates then change in large steps, so apparent swings may reflect discreteness rather than a systematic error.
Unequal or very small bin sizes. A rate from a small bin is much more variable than the same rate from a large bin; show bin sizes and event counts.
Dependence ignored in uncertainty displays. Independent-binomial error bars can be too narrow when observations share clusters.
Misleading prediction scale. Comparing observed rates to conditional plug-in predictions can diagnose a different target from comparison with marginal predictions.
Overinterpreting a single gap. A visual difference is a prompt to investigate, not proof of a specific omitted effect or a formal rejection of the model.

A binned plot compresses individual information and depends on bin boundaries. Changing the number or placement of bins can change the appearance. With rare events, bins with no events are common even under a well-calibrated model; zero observed events do not imply zero underlying event probability. When uncertainty matters, simulation-based residual checks or predictive comparisons can assess whether patterns like those observed are plausible under the fitted GLMM, with the simulation preserving the relevant clustering and refitting choices where appropriate.

The main alternative displays answer related but distinct questions. Individual response-residual plots retain observation-level detail but are difficult to interpret for binary outcomes. Randomized quantile residuals provide a more continuous display for discrete responses, while simulation-based residuals compare observed patterns with replicated data under the fitted model. Binned observed-versus-fitted rates are especially transparent about event-rate calibration, but they can conceal within-bin structure and do not diagnose every form of model misspecification.

The short derivation above also clarifies a limitation: even if the model is correct, the binned differences vary from sample to sample, and their variation depends on the cluster structure. Conversely, a plot close to the diagonal does not establish that the random-effects distribution, covariate effects, or individual predictions are correct. Different errors can cancel when averaged within a bin.

Key takeaway: Plot marginal observed event rates against the corresponding mean marginal fitted probabilities, display bin sizes and event counts, and interpret coherent patterns rather than isolated points. In clustered binary data, dependence and rare-event discreteness make visual gaps uncertain; the plot is a diagnostic, not a stand-alone test.

Check Your Understanding

Answer each question, then compare your reasoning with the brief solutions.

A bin has 25 observations, 2 events, and a sum of marginal fitted probabilities of 1.25. Calculate $O_g$, $E_g$, and $D_g$, and interpret the sign.
Why can an independent-binomial standard error be too small for a bin containing observations from the same clusters?
A plot has one point well above the diagonal, but that bin contains only one event and few observations. What conclusion is justified?
For population-average calibration across clusters, should the plot use conditional plug-in predictions or marginal fitted probabilities? Explain briefly.

Solutions. (1) $O_g=2/25=0.080$, $E_g=1.25/25=0.050$, and $D_g=0.080-0.050=0.030$. The observed rate is 3 percentage points above the mean prediction. (2) Positive within-cluster covariance contributes to the variance of the summed responses; an independence calculation omits those covariance terms. (3) The bin's observed rate exceeds its prediction, but the point is highly uncertain and does not establish systematic underprediction. Examine the broader pattern and uncertainty. (4) Use marginal fitted probabilities, because they average over the random-effects distribution and target average event rates across clusters; conditional plug-in predictions target particular estimated cluster effects.