Why Nonlinear Forecasts Need More Than a Best-Fitting Parameter
An individual’s observed concentrations inform their pharmacokinetic parameters, but usually do not determine them exactly. A Bayesian forecast carries that remaining uncertainty into predictions at future times. This matters particularly when the concentration model is nonlinear: averaging predictions across plausible parameter values need not give the same answer as predicting at an averaged parameter value.
The introductory material on posterior predictive distributions distinguishes uncertainty about an individual’s parameters from residual variability in a future measurement. Here the focus is on constructing the forecast when the concentration model is nonlinear and the individual parameter distribution, after conditioning on their data, may be skewed or have separated regions of support.
The first expression combines the individual’s population-based prior distribution with the likelihood of their observed concentrations. In a typical NLME model, individual parameters arise from random effects, for example through a transformation of \(\eta_i\sim N(0,\Omega)\); the induced distribution of \(\theta_i\) is the prior in this expression. The second expression averages the future-concentration distribution across the individual’s posterior parameter uncertainty. In this tutorial, population quantities \(\phi\) are treated as fixed. If their uncertainty is also to be propagated, it must be included in the conditioning and integration.
State the Concentration and Error Models First
Write the structural concentration as \(f(t;D_i,\theta_i)\), where the dosing history \(D_i\), time \(t\), and individual parameters determine the model-predicted concentration. An observation model specifies how measured concentrations vary around that structural prediction. For example, with lognormal residual error, \(\log Y_{ij}=\log f(t_{ij};D_i,\theta_i)+\epsilon_{ij}\), with \(\epsilon_{ij}\sim N(0,\sigma^2)\). This model requires positive structural predictions and positive observed concentrations. Under it, the future measurement distribution is lognormal conditional on \(\theta_i\); it is not a point mass at the structural concentration.
The forecast can target either the structural concentration \(f(t_*;D_i,\theta_i)\), which excludes future residual error, or a future observed concentration \(Y_{i*}\), which includes it. The target must be stated because their distributions and summaries differ. In either case, nonlinear forecasting averages the appropriate quantity over the posterior distribution of \(\theta_i\).
A Nonlinear Forecast Is an Average of Predictions
Suppose first that the forecast target is the structural concentration. Its posterior mean is \(E[f(t_*;D_i,\theta_i)\mid y_i,D_i,\phi]\), not generally \(f(t_*;D_i,E[\theta_i\mid y_i,D_i,\phi])\). The reason follows directly from conditioning and expectation:
This integral averages predictions, weighting each parameter value by its posterior probability. A local Taylor expansion gives a useful justification for why the plug-in alternative can differ. In a scalar case, let \(m=E[\theta_i\mid y_i]\). A second-order approximation requires suitable moment and remainder conditions and is most informative when posterior mass is concentrated near \(m\); it can be poor for broad or separated posterior distributions. For a sufficiently smooth \(f\),
Taking posterior expectations makes the linear term vanish because \(E[\theta_i-m\mid y_i]=0\). To second order, the difference between the posterior mean prediction and the plug-in prediction is therefore approximately \(\frac{1}{2}f''(m)\operatorname{Var}(\theta_i\mid y_i)\). Unless the model is linear in the uncertain parameter, or posterior uncertainty is negligible, there is no general reason for that difference to be zero.
Worked Example: Forecasting With Two Possible Clearances
Consider a hypothetical one-compartment intravenous model with volume \(V=10\) L and initial concentration \(10\) mg/L. Thus \(C(t;CL)=10\exp(-CL\,t/V)\) mg/L, with clearance \(CL\) in L/h and time in hours. Suppose an individual’s prior assigns probability 0.5 to each of \(CL=10\) L/h and \(CL=25\) L/h. At \(t=1\) h the observed concentration is \(1.779\) mg/L, and the observation model is lognormal with log-error standard deviation \(0.20\). These assumptions define a valid likelihood for the positive observation and a two-point prior for clearance.
For a candidate clearance, the log likelihood, apart from a common constant, is \(-[\log(y)-\log C(1;CL)]^2/(2(0.20)^2)\). The structural predictions at \(t=1\) are \(10e^{-1}=3.6788\) mg/L for \(CL=10\), and \(10e^{-2.5}=0.82085\) mg/L for \(CL=25\). The observed log concentration is \(\log(1.779)=0.576051\). Therefore the log likelihood ratio for \(CL=25\) versus \(CL=10\) is
The prior odds are 1, so the posterior odds are \(e^{-0.880}\). Consequently, \(P(CL=25\mid y)\approx 0.2932\) and \(P(CL=10\mid y)\approx 0.7068\). The posterior mean clearance is \(0.2932(25)+0.7068(10)=14.398\) L/h.
At \(t=2\) h, the structural concentrations for the two clearances are \(10e^{-4}=0.1832\) mg/L? That substitution would use \(CL=20\) L/h. Using the stated model and parameters instead gives \(10e^{-2(10)/10}=10e^{-2}=1.3534\) mg/L for \(CL=10\), and \(10e^{-2(25)/10}=10e^{-5}=0.06738\) mg/L for \(CL=25\). The posterior mean structural forecast is thus \(0.7068(1.35335)+0.2932(0.067379)=0.9763\) mg/L. In contrast, predicting at the posterior mean clearance gives \(10e^{-2(14.398)/10}=10e^{-2.8796}\approx0.5616\) mg/L. Both calculations use the same posterior, but answer different questions; \(0.9763\) mg/L is the posterior mean of the structural predictions, while \(0.5616\) mg/L is one prediction at the mean clearance.
Nonlinear Likelihoods Can Leave Separate Plausible Regions
A nonlinear likelihood need not yield a symmetric, single-peaked posterior. Different combinations of pharmacokinetic parameters can produce similar concentrations at the observed sampling times but make different predictions at other times. If the prior gives support to those distinct regions, the posterior may retain multiple modes. A mean parameter vector can then lie between plausible regions, and its corresponding prediction can represent neither region well.
Worked Example: Two Absorption-Rate Regimes Fit One Concentration
Consider a hypothetical one-compartment oral model with fixed elimination rate \(k=0.2\) h\(^{-1}\), an absorption rate \(k_a\), and a positive scale \(S\) in mg/L. Its structural concentration is \(C(t)=S\frac{k_a}{k_a-k}\{\exp(-kt)-\exp(-k_a t)\}\), for \(k_a\ne k\). Suppose the prior has equal probability on two parameter regimes: A, \(k_a=0.5\) h\(^{-1}\), \(S=14.1376\) mg/L; and B, \(k_a=1.0\) h\(^{-1}\), \(S=8.8721\) mg/L. Both rates exceed \(k\), so the displayed formula is defined and gives positive concentrations.
At \(t=1\) h, the rate-dependent factors are \(\frac{0.5}{0.3}(e^{-0.2}-e^{-0.5})=0.353667\) for A and \(\frac{1.0}{0.8}(e^{-0.2}-e^{-1})=0.563564\) for B. Multiplying by the respective scales gives approximately \(5.000\) mg/L in each regime. With an observed concentration of \(5.1\) mg/L and lognormal error with standard deviation \(0.10\), the predicted concentrations are both about \(5.000\) mg/L at the displayed precision, so the likelihoods and posterior probabilities are approximately equal (about \(0.5\) each).
At \(t=4\) h, regime A predicts \(14.1376\{\frac{0.5}{0.3}(e^{-0.8}-e^{-2})\}\approx7.399\) mg/L. Regime B predicts \(8.8721\{\frac{1.0}{0.8}(e^{-0.8}-e^{-4})\}\approx4.780\) mg/L. Thus the observation is compatible with two distinct parameter regimes whose later forecasts differ substantially. The two-state prior makes the posterior modes explicit; with continuous prior distributions concentrated around these regimes, the corresponding posterior could have two separated high-density regions. It is not valid to assume from the nonlinear model alone that the posterior must be multimodal: the actual posterior shape depends on the data, the sampling times, the prior, and the parameter constraints.
Posterior asymmetry can also arise without distinct modes. For instance, a positive parameter represented on its natural scale may have a long right tail, or a concentration likelihood may constrain one side of a parameter range more strongly than the other. A symmetric approximation on an unsuitable scale can obscure this shape. The relevant question is not whether the posterior looks Gaussian by default, but whether its shape is represented adequately for the intended forecast.
From Parameter Uncertainty to a Future Measurement
For a lognormal residual model, conditional on parameters the future measurement has median \(C_*=f(t_*;D_i,\theta_i)\) and mean \(C_*\exp(\sigma^2/2)\). After conditioning on past observations, the forecast is a mixture over the posterior parameter distribution. Even if each conditional distribution is simple, the mixture need not be symmetric or unimodal.
Worked Example: A Two-Component Predictive Mean
Use the two posterior regimes in the absorption example, each with posterior probability 0.5, and suppose the future measurement at \(t=4\) h has lognormal residual error with standard deviation \(\sigma=0.10\). The component structural concentrations are \(7.399\) and \(4.780\) mg/L. The conditional mean measurement in each regime is its structural concentration multiplied by \(\exp(0.10^2/2)=\exp(0.005)\approx1.00501\). These conditional means are approximately \(7.436\) and \(4.804\) mg/L.
Averaging over the two posterior probabilities gives the future measurement mean \(0.5(7.436)+0.5(4.804)=6.120\) mg/L. This is not the mean of the two structural predictions, \(0.5(7.399+4.780)=6.090\) mg/L, because the former also includes future residual variability. The forecast distribution is a mixture of two lognormal components, not one lognormal distribution centered at \(6.120\) mg/L. With sufficiently separated components and modest residual error, the mixture can retain two distinct peaks.
Checks, Alternatives, and Limitations
A defensible forecast begins by checking that the model and conditioning information match the question. Confirm the dose history, sampling times, units, residual-error scale, and parameter transformations. Posterior draws or other representations of the individual parameter distribution should reproduce the observed concentrations reasonably under the assumed likelihood. Strongly skewed weights, separated parameter regions, sensitivity to prior choices, or predictions concentrated near parameter bounds are reasons to inspect the model and data more carefully—not to replace the posterior with its mean automatically.
| Approach | What it computes | Main caution |
|---|---|---|
| Posterior-averaged nonlinear forecast | Averages future predictions over the conditioned individual parameter distribution | Requires representing the posterior adequately |
| Prediction at posterior mean parameters | Evaluates the model once at the posterior mean parameter vector | Generally differs from the posterior mean prediction in a nonlinear model |
| Gaussian random-effects benchmark | Uses a Gaussian updating structure under its conjugate assumptions | May not represent nonlinear likelihood shapes or separated posterior regions |
The Gaussian random-effects forecast is useful as a benchmark when its assumptions are appropriate, but it should not be treated as a general consequence of Bayesian updating. Conversely, nonlinear posterior averaging does not repair a misspecified structural model, an incorrect residual-error model, or poor dosing records. Forecasts can also be sensitive to fixed population parameters when those are treated as known. These issues limit interpretation: a posterior forecast is conditional on the model and information supplied, not a guarantee that future concentrations will follow the model.
Check Your Understanding
Try each question, then compare your reasoning with the solutions below.
- What two distributions are combined to obtain the posterior distribution of an individual’s parameters?
- Why can a forecast at the posterior mean parameter differ from the posterior mean forecast?
- In the absorption example, why did the two parameter regimes retain equal posterior probabilities after the \(5.1\) mg/L observation?
- Does a posterior predictive mean under lognormal residual error generally equal the mean structural concentration? Explain.
Solutions. (1) Bayes’ rule combines the likelihood of the individual’s observed concentrations with the population-based prior for their parameters. (2) The structural model is nonlinear, so in general \(E[f(\theta_i)\mid y_i]\ne f(E[\theta_i\mid y_i])\). (3) Their predicted concentrations were equal to the displayed precision, so their likelihoods and posterior probabilities were approximately equal (about 0.5 each). (4) No. Under the stated lognormal model, each conditional measurement mean is its structural concentration times \(\exp(\sigma^2/2)\), so future residual variability changes the predictive mean as well as its spread.