Posterior Predictive Checks Of Biomarker Trajectories

Why Replicate Biomarker Trajectories?

A fitted joint model may describe an average biomarker trend and account for subject-specific deviations, yet still fail to reproduce features that matter in the observed data. For example, it may smooth away unusually steep changes, generate too narrow a range of measurements, or underestimate short-term variation around a subject’s trajectory. Posterior predictive checks examine these features by asking: if the fitted model could generate another data set like this one, would that data set look like the observed biomarker data in the ways we care about?

The check is a model diagnostic, not a test that proves the model is correct. It depends on the chosen summaries, the fitted model, and what is held fixed when data are replicated. Here the focus is on longitudinal biomarker trajectories: their slopes, ranges, extrema, and within-subject variability. Replicated event histories are a separate diagnostic target.

Definition and Replication Setup

Let \(Y_{ij}\) be the biomarker measurement for subject \(i\) at time \(t_{ij}\), where \(i=1,\ldots,n\) and \(j=1,\ldots,n_i\). Write \(D\) for the observed joint data, including the biomarker measurements and event information, and let \(\theta\) denote the joint-model parameters. A typical longitudinal submodel has conditional distribution

where \(b_i\) is the subject-specific random effect and \(m_i(t;b_i,\theta)\) is the modelled biomarker trajectory. For example, a Gaussian marker model might set \(Y_{ij}=m_i(t_{ij};b_i,\theta)+\epsilon_{ij}\), with \(\epsilon_{ij}\) having the error distribution specified by the model. A joint model also has an event-time component; that component contributes information to the fitted posterior but does not mean that this particular check must simulate new events.

For a subject-conditional check, draw \((\theta^{(s)},b_1^{(s)},\ldots,b_n^{(s)})\) from the fitted posterior \(p(\theta,b_1,\ldots,b_n\mid D)\), then draw replicated values \(Y_{ij}^{\mathrm{rep},(s)}\) at the observed subjects’ measurement times and covariates. This asks whether the model reproduces the data for these subjects, using their posterior-informed latent effects. For a new-subject check, draw fresh effects \(b_i^{\mathrm{rep},(s)}\) from the population random-effects distribution implied by \(\theta^{(s)}\), then generate marker measurements. That asks whether the model can reproduce trajectories for subjects drawn from its population. The two checks target different aspects of fit and should not be treated as interchangeable.

Unless the visit or missingness process is itself part of the model and is being checked, the usual trajectory comparison holds the observed measurement schedule fixed. In particular, a replicate should not create a biomarker observation after a subject’s observed event and then compare it as if it had been measured. The intended prediction target must be stated when event-related observation or dropout affects which measurements are available.

Choosing Summaries and Comparing Replicates

A check is only informative about features captured by its summaries. For each subject with at least two distinct measurement times, one possible slope summary is the least-squares slope of the observed measurements against time. A range summary is the maximum minus the minimum measurement. The minimum and maximum can also be examined separately, since a similar range can arise from different extremes. Slopes and ranges have the biomarker’s units per unit time and the biomarker’s units, respectively.

Within-subject variability needs a precise definition. One useful summary is residual variation around a subject-specific fitted trajectory, rather than variation around that subject’s raw mean. For a linear trajectory fitted to one subject with \(n_i>2\) measurements, define the residual standard deviation as

where \(\widehat{a}_i\) and \(\widehat{c}_i\) are the fitted intercept and slope for that subject. This measures scatter around a linear trend; it is not a universal measure of within-subject variability. If the fitted trajectory is nonlinear, or if residuals are autocorrelated or heteroscedastic, the summary and its interpretation need to reflect that structure. In posterior predictive replication, calculate the same summary by the same rule for observed and replicated data.

A statistic \(T\) can summarize one trajectory or a feature across subjects—for example, the median subject-specific slope, the largest subject range, or the median residual standard deviation. For posterior draw \(s\), calculate \(T(Y^{\mathrm{rep},(s)})\) and compare it with \(T(Y)\). An upper-tail posterior predictive probability is estimated by

where \(B\) is the number of replicated data sets and \(I\{\cdot\}\) equals 1 when the condition is true and 0 otherwise. A value near 0 means the observed summary is larger than nearly all replicated summaries; a value near 1 means it is smaller than nearly all of them. A value near the middle indicates that the observed summary is not unusual relative to these replicates, but it does not establish adequate fit in other respects.

Worked Trajectory Calculations

Practical Workflow, Diagnostics, and Alternatives

A defensible check starts by stating the prediction target and the summary before looking at its result. Then use posterior draws that represent the fitted joint model, generate marker values under the longitudinal submodel, and calculate identical summaries for the observed and replicated data. The following sequence makes the comparison explicit.

Plots of observed summaries against replicated distributions can reveal features obscured by a single tail area. Check slopes at more than one time scale when scientifically relevant, and examine minima and maxima separately as well as ranges. A model could match a median slope while missing the spread of slopes, or match a range while placing its minimum and maximum incorrectly. For variability, compare the distribution across subjects and check whether discrepancies are concentrated at particular times or measurement levels.

The earlier material on longitudinal residual diagnostics provides complementary ways to look for systematic patterns in the observed data. Residual diagnostics focus on deviations from fitted values; posterior predictive checks ask whether whole replicated data sets reproduce selected observed features under parameter uncertainty. Neither replaces the other. A parametric bootstrap can also compare simulated data with observed summaries, but it generally uses a different treatment of parameter uncertainty than posterior predictive simulation. The appropriate comparison depends on the inferential framework and diagnostic question.

Several failures can make a check hard to interpret. If the model’s posterior simulation has not adequately explored its posterior distribution, apparent discrepancies may reflect computation rather than model fit. If observed and replicated summaries use different time grids, different rules for missing measurements, or different slope calculations, the comparison is not like-for-like. Very small numbers of measurements make individual slopes and extrema unstable. A maximum range is especially sensitive to sample size: larger groups have more opportunities to contain an extreme. Interpret it relative to replicates with the same design and group size.

Subject-conditional replications use posterior-informed random effects, so they can reproduce subject-specific patterns more closely than new-subject replications. Consequently, a good conditional check does not establish that the population random-effects distribution is adequate. Conversely, a new-subject check can expose poor population-level predictions but may differ because it draws fresh effects, not because the measurement-error model is wrong. These checks also cannot diagnose features they do not summarize. A model may pass a range check and still miss curvature, multimodality, or a time-dependent variance pattern.

Posterior predictive checks are most useful as a structured way to find consequential lack of fit, not as a score that certifies a model. Report what was held fixed, whether effects were conditioned on or redrawn, the summary and its units, the number of simulations, and the observed position in the replicated distribution. A discrepancy should prompt model investigation; it does not by itself identify which model component needs changing.