Why Stratify a Predictive Check?
A model can reproduce the overall response distribution while missing systematic differences across a candidate covariate. For example, predictions may look adequate after pooling all participants even though responses are consistently too high in one age group and too low in another. A covariate-stratified visual predictive check (VPC) makes this kind of mismatch easier to see by comparing observed and simulated response distributions separately within covariate-defined groups.
A stratified VPC is a model diagnostic, not a stand-alone test that proves a covariate belongs in the model. Its message depends on the chosen strata, time bins, summaries, and simulation design. A visible mismatch can motivate further model investigation; agreement means only that the displayed summaries are compatible with the model in the data region examined.
Construction and Interpretation
Suppose participant \(i\) has covariate value \(C_i\), and response \(Y_{ij}\) is recorded at occasion \(j\). First, define strata using a categorical covariate's levels or, for a continuous covariate, clinically and statistically defensible ranges. If the response changes over time, also define time bins or other meaningful design-based groups. Within each covariate stratum and time bin, calculate observed summaries such as the median and selected quantiles.
Next, simulate \(B\) replicate data sets from the fitted model. Simulations should retain the observed dosing or treatment design, sampling times, covariate values, and relevant observation structure. In each replicate, calculate the same summaries within the same strata and bins. For a chosen summary, the central 90% simulation envelope is formed by its 5th and 95th percentiles across replicates. Plot the observed summary against these envelopes, often with time on the horizontal axis.
Here \(s\) indexes covariate strata, \(t\) indexes time bins, \(p\) is a response quantile such as 0.10, 0.50, or 0.90, and \(b=1,\ldots,B\). The displayed envelope for that summary is the interval from the 5th to the 95th percentile of \(\widehat{q}_{s,t,p}^{(b)}\) over replicates. It describes the range in which that statistic would fall in 90% of repeated data sets simulated under the fitted model and design, subject to the simulation procedure and parameter treatment. It is not a 90% confidence interval for a model parameter.
The model assumptions relevant to interpretation are those built into its simulation: the structural response model, covariate effects, random-effects distributions, residual error model, and any dependence across repeated observations. The design must also be represented faithfully. If those ingredients are wrong, the simulated envelope may be misleading even when the code runs successfully. A plug-in VPC typically treats estimated parameters as fixed; it therefore does not automatically include all uncertainty from estimating them.
Why the Simulation Envelope Is Informative
The justification is a comparison of like with like. Under a correctly specified model and the same fixed design, the observed data and a newly simulated data set are generated from the same distribution, if the model parameters are treated as known. Applying the same statistic to either data set therefore gives statistics with the same sampling distribution. The observed statistic should consequently fall within a central 90% simulation interval about 90% of the time, apart from finite simulation error and ties. This is a calibration argument, not a guarantee that any one observed statistic will lie inside its envelope.
For a given stratum, time bin, and response quantile, let \(T_{\mathrm{obs}}\) be the observed summary and \(T^{(1)},\ldots,T^{(B)}\) the replicate summaries. The simulation interval is obtained by ordering the replicate summaries and selecting the empirical 5th and 95th percentiles. If \(T_{\mathrm{obs}}\) lies outside that range, the observed summary is relatively unusual under the fitted model for that stratum and bin. The check does not identify which model component caused the discrepancy.
Worked Example: A Categorical Covariate
Consider hypothetical continuous PD responses measured at a specified post-dose time. A candidate categorical covariate defines two groups, A and B, with 40 participants in each. In group A, the 20th and 21st ordered observed responses are 8.0 and 8.4 response units, so the sample median is \((8.0+8.4)/2=8.2\) units. In group B, the corresponding values are 10.0 and 10.2 units, giving a median of \((10.0+10.2)/2=10.1\) units.
Suppose 500 data sets are simulated under the fitted model, preserving the observed group membership and sampling design. The 5th-to-95th percentile envelope of replicate medians is 7.7 to 8.7 units for group A and 8.5 to 9.5 units for group B.
The group A observed median, 8.2 units, is inside its simulated envelope. The group B median, 10.1 units, is above its upper limit of 9.5 units by \(10.1-9.5=0.6\) units. This indicates that the fitted model tends to predict group B medians lower than the observed median under this design. It is a diagnostic signal to investigate, not proof by itself that the candidate covariate has a particular causal effect or that a specific covariate model is correct.
Worked Example: Stratifying a Continuous Covariate
In a hypothetical study, body weight is a continuous candidate covariate. Before examining this VPC, the analyst defines three ranges: below 60 kg, 60 to below 80 kg, and at least 80 kg. At one post-dose time, the strata contain 25, 30, and 22 participants. Their observed response medians are 5.4, 6.1, and 7.0 mg/L, respectively. For the middle stratum, the two central ordered responses are 5.9 and 6.3 mg/L, giving a median of \((5.9+6.3)/2=6.1\) mg/L.
Across 500 simulated data sets, suppose the respective 90% envelopes for the median are 4.9–6.0, 5.6–6.6, and 6.3–7.5 mg/L. The first stratum's observed median exceeds its envelope by \(5.4-6.0=-0.6\) mg/L? It does not: \(5.4\) is within 4.9–6.0 mg/L. The second median, 6.1 mg/L, is also within its envelope, as is the third, 7.0 mg/L.
Thus these three medians do not show a conspicuous lack of fit against their respective envelopes. The rising observed medians across weight ranges are not, on their own, evidence of model failure: the model's simulated medians may also rise across the same strata. The comparison is between observed and simulated summaries within each stratum, not between raw medians across strata.
Choosing Strata and Reading the Display
Categorical covariates can generally be stratified by their modeled levels, provided each level has enough data to estimate the summaries meaningfully. For a continuous covariate, cut points should be chosen with the scientific question, data density, and expected covariate relationship in mind. Quantile-based ranges can balance sample sizes; scientifically motivated ranges can make interpretation easier. Either choice has a cost: broad ranges may hide a local trend, while many narrow ranges can produce unstable summaries.
When the response evolves over time, inspect whether bins preserve the important shape of the time course. Wide bins can average over a peak or transition; narrow bins can contain too few observations. Keep the same bin definitions for observed and simulated data. Report the number of participants and observations in each stratum and bin, since nominally large observation counts can still represent few independent participants.
A useful display often includes the observed median and lower and upper quantiles, along with their separate simulation envelopes. The median checks central tendency; the outer quantiles give information about spread and tails. A model may reproduce the median while missing dispersion, or match the overall range while misrepresenting the center. Compare all displayed summaries rather than treating one line as a complete assessment.
Worked Example: A Binary Response Within a Covariate Stratum
Suppose a hypothetical binary PD outcome is recorded once for each of 12 participants in a high-covariate stratum. Nine participants have outcome 1, so the observed proportion is \(\hat{p}=9/12=0.75\). In 500 simulations with the same 12 participants and stratum membership, the 5th and 95th percentiles of the simulated proportions are 0.42 and 0.83.
The observed proportion, 0.75, lies inside the 0.42–0.83 envelope. Under this fitted model, such a proportion is not unusual relative to the displayed central 90% of simulated proportions. The interval is wide because each simulated proportion changes in increments of \(1/12\), and the stratum is small. This example illustrates why sample size and outcome type matter: a proportion VPC is interpreted on its discrete scale and may have coarse envelopes.
Diagnostics, Alternatives, and Limitations
Investigate a VPC that looks unexpected by checking the data processing and simulation before changing the model. Confirm that covariate units and group labels are correct, observed and simulated records use identical strata and time bins, and simulations preserve each participant's design. Check that the simulation includes the intended random effects and residual variability. A mismatch confined to one time region, one tail, or one stratum is more informative than a pooled impression, but it still does not uniquely identify the source.
- Too many strata or bins: summaries become noisy and envelopes broad or jagged. Reduce the number of bins or report counts and uncertainty clearly.
- Unequal or sparse strata: apparent differences can reflect low information, especially for tail quantiles. Avoid strong conclusions from a small cell.
- Pooling away a pattern: an overall VPC may look adequate while stratum-specific summaries diverge. Inspect the covariate-specific displays rather than relying on pooled results.
- Unrepresentative simulation: omitting design features or dependence can create an envelope that is not comparable to the observed statistic. Simulate the observed design, not an idealized one.
- Monte Carlo noise: a small number of replicates makes envelope percentiles unstable. Increase \(B\) when the plotted boundary changes noticeably across repeated runs.
Compared with an unstratified VPC, a stratified VPC can reveal conditional misfit but has fewer observations per display and more opportunities for unstable summaries. A model comparison such as a likelihood-ratio test addresses a different question: whether a specified nested model improves fit under its assumptions. Covariate-effect plots or empirical Bayes plots can also help explore patterns, but they do not replace checking the model's replicated response distributions. These approaches are complementary, not interchangeable.
A VPC is limited by the summaries chosen, the number of observations, model-estimation uncertainty, and the covariate values represented in the data. It cannot establish reliable behavior outside that observed range. It is also not a formal multiplicity-adjusted test: inspecting many strata, bins, and quantiles increases the chance of seeing an unusual point by chance. Treat the visual pattern as evidence to assess alongside biological rationale and other diagnostics, not as a binary pass/fail rule.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- What information does a 90% simulation envelope for a stratum-specific median represent?
- A continuous covariate is divided into many narrow ranges, leaving only a few participants in each. What is the main diagnostic risk?
- An observed median is outside its envelope, while the 10th and 90th observed percentiles are inside theirs. What can you conclude?
- Why should simulations retain the observed participants' covariates and sampling design?
Solutions. (1) It is the central 90% of the replicate-specific medians generated under the fitted model and design; it is not a confidence interval for a parameter. (2) Stratum summaries become unstable and their envelopes may be broad or erratic, making apparent patterns hard to distinguish from sampling variation. (3) The central response summary is relatively unusual under the fitted model in that stratum and bin, while the displayed tails are not; this is a diagnostic signal, not proof of a particular cause or formal evidence by itself for adding a covariate. (4) Matching covariates and design makes observed and simulated summaries comparable under the same conditions; otherwise differences can arise from design mismatch rather than model fit.