Why Use Posterior Predictive Checks?
Sparse outcomes can make familiar diagnostic summaries unstable: a small number of successes, events, or nonempty cells may dominate the evidence about fit. A posterior predictive check asks a focused question: if the fitted model could generate new outcomes under the uncertainty supported by the data, how often would a chosen summary look as unusual as the one observed? The answer depends both on which summary is targeted and on how parameter and random-effects uncertainty enter the replication.
For a GLMM, outcomes within a cluster are linked through shared random effects. A replication that simply inserts estimated coefficients and estimated random effects as if they were known can understate predictive variation. A posterior predictive check instead averages over a posterior distribution, while keeping the intended prediction target clear: the observed clusters or newly generated clusters.
Posterior Predictive Distribution
Let \(Y_{ij}\) be the outcome for observation \(i\) in cluster \(j\), let \(x_{ij}\) denote its covariates, and let \(b_j\) denote the cluster random effect. Write \(\theta\) for the model parameters, including fixed-effect coefficients and random-effects distribution parameters. A fitted Bayesian GLMM supplies a posterior distribution for the unknown quantities given the observed data \(y\).
The expression makes two sources of variation explicit. Conditional on a particular draw of \(\theta\) and \(b\), outcomes vary according to the response distribution in the GLMM. Across posterior draws, the fitted parameters and latent cluster effects also vary. For any discrepancy statistic \(T\), the law of total variance gives
The first term is outcome variation conditional on a draw of the unknown quantities; the second is variation in the statistic’s conditional expectation across posterior draws. Holding estimates fixed omits the second contribution and can make a predictive distribution too narrow.
A discrepancy statistic is a function chosen to emphasize a feature relevant to the model and data. Examples for sparse outcomes include the number of all-zero clusters, the number of clusters with exactly one event, the total number of events, or the count above a specified threshold. The statistic should be chosen before inspecting which choice gives the most reassuring result; otherwise, a reported check can conceal the feature that motivated it.
Choosing What Is Replicated
The random-effects part of the replication must match the question. For a check of outcomes in the same observed clusters, one can draw from the posterior distribution of the existing \(b_j\) values and then simulate new outcomes conditional on those draws. This assesses the response model given the cluster effects inferred from the data. For a check of the model’s behavior across new clusters, draw new random effects from their distribution under each posterior draw of \(\theta\), then generate outcomes. This assesses population-level cluster variation as well as outcome variation.
These targets are not interchangeable. Reusing posterior random effects is not the same as drawing fresh effects from the population distribution. Likewise, a marginal check that averages over random effects answers a different question from a check conditional on cluster effects. State the target alongside the statistic so readers know what sources of variation the replication includes.
Worked Example: All-Zero Probability for a Rare Binary Outcome
Suppose a hypothetical logistic random-intercept GLMM describes five binary outcomes in one cluster. At one posterior draw, the intercept is \(\beta_0=-2.5\) and the cluster effect is \(b_j=0\). The conditional event probability for each observation is therefore
Under conditional independence given the draw, the probability that all five outcomes are zero is
Thus, at this particular draw, a zero-event cluster is not surprising: its conditional probability is about 0.674. That is not yet a posterior predictive probability, because it uses just one parameter-and-random-effect draw.
For a simplified illustration of posterior averaging, suppose three equally weighted posterior draws give cluster event probabilities \(p=0.05,\ 0.10,\ 0.20\). The posterior predictive probability of five zeros is
The averaged probability, about 0.564, incorporates variation across the illustrative draws. In an actual check, the draws would come from the joint posterior and would be numerous enough to characterize the predictive distribution. If the target were new clusters, the random effect for each replicated cluster would also be generated according to the model’s random-effects distribution.
Estimating a Predictive Tail Area
A common summary compares the observed statistic \(T(y)\) with replicated values \(T(y^{\mathrm{rep}})\). For an upper-tail check, the posterior predictive tail area is
With \(M\) posterior predictive replicates, estimate this probability by the fraction at least as large as the observed statistic. The direction of comparison should reflect the scientific concern: too many empty clusters calls for an upper-tail check of the empty-cluster count, whereas too few calls for a lower-tail check. A two-sided check may be useful for a clearly specified discrepancy, but there is no universal rule that converts a discrete statistic into a unique two-sided tail area.
Worked Example: Checking the Number of Empty Clusters
In a hypothetical study, 20 clusters each contain a fixed set of observations, and the observed number of all-zero clusters is \(T(y)=4\). The target is replication for these same clusters, including posterior uncertainty in the relevant parameters and cluster effects. Across \(M=1000\) posterior predictive replications, the simulated numbers of empty clusters are summarized below.
| Replicated empty clusters | Number of replications |
|---|---|
| 0 | 20 |
| 1 | 180 |
| 2 | 280 |
| 3 | 270 |
| 4 | 150 |
| 5 or more | 100 |
For the upper-tail check, replications with four or more empty clusters number \(150+100=250\). The estimated posterior predictive tail area is
The simulation-based estimate has approximate Monte Carlo standard error
So the observed four empty clusters are not in an extreme upper tail under this check. This is not proof that the model is correct: it says that this chosen statistic, under this replication target, is not unusually large relative to the posterior predictive distribution. A lower-tail check would address the distinct concern that the model produces too many empty clusters compared with the observed data.
Accounting for Posterior Uncertainty
The following calculation shows how posterior mixing changes a predictive tail probability. Consider a hypothetical statistic \(T\) that counts, among ten independent clusters, how many meet a specified sparse-outcome criterion. Suppose that, conditional on each of three equally weighted posterior draws, each cluster meets the criterion with probability \(q\), where the three values of \(q\) are \(0.05,\ 0.15,\ 0.30\). This simplified setup treats the cluster criterion indicators as independent conditional on the draw. The posterior draws stand in for uncertainty in the model quantities relevant to the check.
Worked Example: Mixing Tail Probabilities Across Posterior Draws
Suppose the observed statistic is \(T(y)=2\). Conditional on a draw with probability \(q\), \(T(Y^{\mathrm{rep}})\) has a binomial distribution with 10 trials. For each draw,
For \(q=0.05\), the upper-tail probability is \(1-0.95^{10}-10(0.05)0.95^9\approx 0.08614\). For \(q=0.15\), it is \(1-0.85^{10}-10(0.15)0.85^9\approx 0.4557\). For \(q=0.30\), it is \(1-0.70^{10}-10(0.30)0.70^9\approx 0.8507\). Averaging the three equally weighted draws gives
The posterior predictive upper-tail probability is therefore about 0.464 in this illustration. It is the average of the conditional probabilities over the posterior draws, not the probability evaluated at an average value of \(q\). In general, nonlinear tail probabilities cannot be obtained by simply substituting a posterior mean parameter.
Interpreting the Check and Recognizing Failure Modes
A posterior predictive tail area is a diagnostic summary, not a classical \(p\)-value. The observed data contribute to the posterior used to generate the replicates, so the observed statistic and its predictive reference distribution are not independent in the way required for a standard null-distribution test. Consequently, posterior predictive tail areas are not generally uniform under a correctly specified model, and values near one half do not certify model adequacy. Very small or very large values can flag a mismatch for the chosen statistic, but the check should be read alongside other summaries and substantive knowledge.
Sparse data make the choice of statistic especially consequential. A total event count can look adequate even if the model allocates events to the wrong clusters; an empty-cluster count can reveal that feature but miss tail behavior among nonempty counts. Use statistics that target plausible failure patterns, and inspect their replicated distributions rather than reporting only a single number.
Common warning signs include a replication that treats fitted random effects as fixed despite a population-level target, simulated data that do not preserve the observed covariates or cluster sizes, or a statistic selected only after trying many alternatives. For each replication, keep the design features fixed when they are conditioned on in the analysis, and reproduce the model’s sampling structure. A check also cannot diagnose features it does not summarize: a satisfactory empty-cell count, for example, does not guarantee satisfactory calibration or dependence structure.
Posterior predictive checks complement rather than replace other approaches. A parametric bootstrap predictive envelope can assess replicated behavior under fitted frequentist estimates; a posterior predictive distribution additionally averages over a posterior distribution for unknown quantities. Conditional checks based on fitted cluster effects focus on within-cluster response behavior, whereas checks that generate fresh cluster effects assess population-level cluster variation. Neither target is automatically preferable; the question determines the relevant one. In all cases, results depend on the assumed GLMM and the chosen discrepancy statistic.
Check Your Understanding
Answer each question, then compare with the solutions below.
- Two equally weighted posterior draws give event probabilities \(p=0.10\) and \(p=0.20\). For six conditionally independent outcomes, calculate the posterior predictive probability that all six are zero.
- In 500 posterior predictive replications, 80 statistics are at least as large as the observed statistic. Calculate the estimated upper-tail area and its approximate Monte Carlo standard error.
- For a check about the behavior of new clusters, should the replication reuse each observed cluster’s posterior random effect or generate a new effect from the fitted population distribution? Explain.
- Why is a posterior predictive tail area of 0.50 not proof that the model fits well?
Solutions. (1) Average the two conditional probabilities: \(\{0.90^6+0.80^6\}/2=(0.531441+0.262144)/2=0.3968\). (2) The estimate is \(80/500=0.160\); its approximate Monte Carlo standard error is \(\sqrt{0.16(0.84)/500}\approx 0.01640\). (3) Generate a new random effect from the fitted population distribution for each replicated new cluster; reusing observed effects targets the existing clusters instead. (4) The data are used to form the posterior that generates the replicates, and the check examines only its chosen statistic. A value near 0.50 is not a classical test result and cannot rule out mismatches in other features.