Worked Bayesian Forecasting Case with Concentration Data

From Observed Concentrations to a Future Prediction

An individual forecast is a distribution for a future concentration conditional on the information available now. Producing it requires more than substituting a single parameter estimate into a pharmacokinetic equation: uncertainty about the individual’s parameters and the future observation’s residual error both contribute. This case follows one hypothetical individual from observed concentration data through a Bayesian update to a future prediction interval.

The calculation uses a deliberately small, two-particle representation of the individual parameter distribution. Such a representation makes every step visible; in practical analyses, the posterior may instead be represented by many simulation draws or numerical integration. The introductory material on posterior predictive distributions establishes the general principle: average future-data distributions over the posterior distribution of the unknown parameters.

Definition: Let \(D_i\) denote the conditioning data for individual \(i\), including the observed concentrations and the doses and sampling times needed to interpret them. A posterior predictive distribution for a future concentration \(Y_{i,*}\) at a specified future time is the distribution of \(Y_{i,*}\) conditional on \(D_i\), after accounting for uncertainty in the individual parameters and the future residual error.

Case Model and Conditioning Data

Suppose an individual receives a hypothetical 200 mg intravenous bolus. Assume a one-compartment model with instantaneous distribution and first-order elimination:

C_i(t\mid CL_i,V_i)=\frac{\mathrm{Dose}}{V_i}\exp\left(-\frac{CL_i}{V_i}t\right).

Here \(C_i(t)\) is the model-predicted concentration in mg/L, \(V_i\) is apparent volume in L, and \(CL_i\) is clearance in L/h. For this worked case, \(V_i=10\) L is treated as known, while the individual clearance is uncertain. Thus the initial model concentration is \(200/10=20\) mg/L. The observed concentration at 1 h is 13.1411 mg/L.

Assume the observation model is additive normal error with known standard deviation \(\sigma=1.00\) mg/L:

Y_i(t_j)=C_i(t_j\mid CL_i,V_i)+\epsilon_{ij}, \qquad \epsilon_{ij}\sim N(0,\sigma^2).

For this calculation, errors at distinct observation times and the future time are conditionally independent given the individual parameters. The two possible clearances in the prior approximation are 6.00 L/h and 3.00 L/h, each with prior probability 0.5. The model, prior, error distribution, dose, volume, and observed concentration define the assumptions and information used in the update. The future prediction is for 2 h after the dose, using the same model and residual-error assumptions.

Worked Example: Updating the Individual Parameter Distribution

Let the two prior particles be \(CL_i=6.00\) L/h and \(CL_i=3.00\) L/h, with \(V_i=10\) L for both. Calculate their predictions at the observed time, then update their weights using the measured concentration \(y_i(1)=13.1411\) mg/L.

The elimination rates are \(6/10=0.600\) h\(^{-1}\) and \(3/10=0.300\) h\(^{-1}\). The model predictions at 1 h are \(20\exp(-0.600)=10.9762\) mg/L and \(20\exp(-0.300)=14.8164\) mg/L, respectively. The observed-minus-predicted residuals are therefore \(2.1649\) mg/L and \(-1.6753\) mg/L.

Under the stated normal error model, the likelihood for a particle is proportional to \(\exp[-(y_i(1)-C_i(1))^2/(2\sigma^2)]\). Since the prior weights are equal and the common normalizing factor cancels, the unnormalized posterior weights are approximately \(\exp[-(2.1649)^2/2]=0.0960\) and \(\exp[-(-1.6753)^2/2]=0.2458\).

Normalizing gives the posterior weight for the 6.00 L/h particle as \(0.0960/(0.0960+0.2458)\approx0.2809\); the weight for the 3.00 L/h particle is \(1-0.2809=0.7191\). The observed concentration is closer to the prediction from the 3.00 L/h particle, so that particle receives more posterior weight. The update changes the relative support; it does not make either particle certain.

Constructing the Future Predictive Distribution

At 2 h, the two particles predict concentrations of \(20\exp(-1.200)=6.0239\) mg/L and \(20\exp(-0.600)=10.9762\) mg/L. Rounded to three decimals, these are 6.024 and 10.976 mg/L. Conditional on either particle, the future measured concentration has a normal distribution centred on its model prediction, with standard deviation 1.00 mg/L. The posterior predictive distribution is therefore a weighted mixture of those two normal distributions.

Y_{i,*}\mid D_i \sim 0.2809\,N(6.024,1.00^2) + 0.7191\,N(10.976,1.00^2).

This is a mixture distribution, not a single normal distribution. Its component weights express the posterior probabilities of the parameter particles. Its component variances represent future residual error. The separation between component means reflects uncertainty about the individual’s clearance, which persists even after observing the concentration.

Worked Example: Calculating the Predictive Mean and Variance

Use the displayed mixture to obtain a point forecast and describe the overall predictive spread. For squared deviations, use the displayed component means, whose difference is \(10.976-6.024=4.952\) mg/L.

The predictive mean is the posterior-weighted average of the component means: \(0.2809(6.024)+0.7191(10.976)=9.5849832\) mg/L, or approximately 9.585 mg/L. This is the mean of the future-concentration distribution conditional on the observed data and the stated model.

The variance of a mixture is the weighted average of the within-component variances plus the variance of the component means around the mixture mean. With a common within-component variance of \(1.00^2\), the between-component term for two components is \(w_1w_2(\mu_2-\mu_1)^2\). Thus \(0.2809(0.7191)(4.952)^2=4.953387\) (mg/L)\(^2\), and the total variance is \(1.000+4.953387=5.953387\) (mg/L)\(^2\). The predictive standard deviation is \(\sqrt{5.953387}\approx2.440\) mg/L.

The overall spread is substantially larger than the 1.00 mg/L residual standard deviation for either particle alone. That difference arises because the forecast also averages over uncertainty about which clearance particle is appropriate.

A Central Individual Prediction Interval

A central 95% posterior predictive interval has lower and upper endpoints at the 2.5th and 97.5th percentiles of the mixture. These are quantiles of future concentrations conditional on the observed data and the model assumptions. They are not obtained by taking 1.96 times the mixture standard deviation around the predictive mean: that shortcut would impose a single normal shape, whereas the posterior predictive distribution here is a mixture.

For a mixture with component means \(\mu_1,\mu_2\), common standard deviation \(\sigma\), and first-component weight \(w_1\), its cumulative distribution function is \(F(y)=w_1\Phi((y-\mu_1)/\sigma)+(1-w_1)\Phi((y-\mu_2)/\sigma)\), where \(\Phi\) is the standard normal cumulative distribution function. The interval endpoints solve \(F(y)=0.025\) and \(F(y)=0.975\). In this case, each endpoint lies near one component, allowing a transparent approximation that can be checked against the full mixture equation.

Worked Example: Finding the 95% Predictive Interval

For the lower endpoint, the second component’s contribution is negligible near the 2.5th percentile because its mean is nearly 5 mg/L above the first component. Solving \(0.2809\Phi((q_{0.025}-6.024)/1.00)\approx0.025\) gives a standard-normal probability \(0.025/0.2809\approx0.0890\). Its quantile is approximately \(-1.347\), so \(q_{0.025}\approx6.024-1.347(1.00)=4.677\) mg/L.

For the upper endpoint, the first component is essentially entirely below the endpoint. The remaining cumulative probability required from the second component is \((0.975-0.2809)/0.7191\approx0.965234\). The corresponding standard-normal quantile is approximately 1.815, giving \(q_{0.975}\approx10.976+1.815(1.00)=12.791\) mg/L.

The resulting central 95% posterior predictive interval is approximately 4.677 to 12.791 mg/L. The endpoint approximations neglect the opposite component’s small tail contribution; numerical inversion of the full mixture CDF is the appropriate way to obtain more precise endpoints. The interval describes a range containing 95% of the predictive probability under this posterior and model, conditional on the stated data and assumptions. It does not mean there is a 95% probability that a fixed unknown future value is inside an already calculated interval.

Why the Full Forecast Matters

The sequence of calculations distinguishes three quantities that are easy to conflate: the observed data update the individual parameter distribution; the model maps each parameter value to a future concentration; and residual error adds variation between a model-predicted concentration and a future measured concentration. The predictive distribution combines all three. A forecast of the model-predicted concentration alone would omit the final observation-error contribution and would not answer the same question as a prediction for a future measured concentration.

The distinction also clarifies the main alternatives. A plug-in forecast uses one estimated clearance and returns one model prediction; it is simple but suppresses parameter uncertainty unless an uncertainty calculation is added. A posterior mean parameter inserted into a nonlinear concentration equation generally need not equal the posterior mean concentration, because averaging and nonlinear transformation do not in general commute. A posterior predictive calculation instead transforms each posterior parameter value and averages the resulting future-data distributions. A normal approximation based on the predictive mean and variance is another possible summary, but can obscure mixture shape or skewness; here the mixture CDF directly supplies the interval quantiles.

A two-particle posterior is useful pedagogically but coarse. It can only represent two clearance values and may miss substantial posterior mass between or beyond them. In an applied analysis, the forecast should use a sufficiently representative posterior sample or a suitable integration method, and numerical accuracy should be checked by increasing the number of draws or refining the integration. The present case also treats volume, residual error, and model form as known. If those quantities are uncertain, or if the one-compartment model and additive error assumption are inadequate, the displayed interval will not account for that uncertainty or misspecification.

Diagnostics and Common Pitfalls

Check the conditioning set. Record which concentrations, doses, and sampling times informed the posterior. Accidentally including a future concentration in \(D_i\) would make the forecast no longer a prediction made before that observation.
Check the likelihood and units. The residuals and \(\sigma\) must use the same concentration units. Here both are in mg/L, so the standardized residuals are dimensionless.
Check posterior weights. Each weight must be nonnegative and the weights must sum to one. A particle with a poorer data fit can still retain positive posterior weight when its prior weight is positive and its likelihood is nonzero.
Check predictive quantiles using the full CDF. Component-wise normal intervals are not, by themselves, mixture quantiles. The approximations above are justified by negligible opposite-component tails; for overlapping components, solve the full mixture equation numerically.
Do not omit future residual error when predicting a measured concentration. Using only the parameter-based model predictions produces a distribution for the latent model concentration, not for the future measurement specified here.
Do not infer calibration from one case. A model-based 95% predictive interval has nominal conditional probability under its assumptions. Whether intervals achieve appropriate coverage across individuals requires evaluation across suitable cases, not inspection of this one forecast alone.

Common mistake: “The posterior mean clearance is enough, so I will calculate one concentration and attach the residual standard deviation as the prediction uncertainty.” This neglects posterior uncertainty in clearance and can also misrepresent a nonlinear transformation. In this case, the mixture standard deviation is about 2.440 mg/L, not 1.00 mg/L. The correction is to propagate the parameter posterior through the concentration model and then include future residual error before calculating predictive summaries.

Check Your Understanding

Use the case assumptions and results above. Each solution is included after the questions.

Why does the 3.00 L/h clearance particle receive greater posterior weight than the 6.00 L/h particle?
What would the forecast represent if the future residual error were omitted?
Why is the 95% interval not calculated as the predictive mean plus or minus 1.96 predictive standard deviations?
Suppose the observed concentration had been measured in micrograms per litre but was entered as though it were in mg/L. What part of the analysis would be affected?
State a correct contextual interpretation of the 4.677 to 12.791 mg/L interval.

Solutions. (1) At 1 h, the prediction from 3.00 L/h is 14.816 mg/L, closer to the observed 13.1411 mg/L than the 10.976 mg/L prediction from 6.00 L/h; with equal prior weights and the stated normal likelihood, it therefore receives the larger posterior weight. (2) The result would describe uncertainty in the model-predicted concentration due to the parameter posterior, not the distribution of a future measured concentration under the stated observation model. (3) The predictive distribution is a mixture of normals and need not itself be normal; its quantiles should be found from its CDF. (4) The likelihood residuals would be calculated on incompatible scales, giving an invalid update and forecast; convert the measurement to mg/L before analysis. (5) Conditional on the observed concentration and the stated prior, one-compartment model, and error assumptions, 95% of the posterior predictive probability for the future measured concentration at 2 h lies approximately between 4.677 and 12.791 mg/L.