Therapeutic Drug Monitoring With Bayesian Forecasts

Why Combine Population Information With TDM?

A measured drug concentration contains information about an individual’s pharmacokinetics, but a small number of measurements rarely determine an individual parameter precisely by themselves. Therapeutic drug monitoring (TDM) addresses this by combining the concentration data with prior information from a population pharmacokinetic (PK) model. The result is an individual-level posterior distribution: it represents which parameter values remain plausible after considering both sources of information. Forecasts under candidate regimens can then account for uncertainty in that distribution rather than treating one estimated parameter as known.

The calculation is only as credible as its PK and observation models. The example below uses a deliberately simple one-compartment model with known volume, linear elimination, intravenous bolus dosing, and lognormally distributed concentration error. It illustrates Bayesian forecasting; it is not a clinical dosing recommendation.

Model and Bayesian Updating

Definition: In Bayesian TDM, a prior distribution describes uncertainty about an individual’s PK parameter before using that individual’s concentration data. A likelihood describes how probable the observed data are under each parameter value. Bayes’ rule combines them to give the posterior distribution, which is the parameter distribution conditional on the observed data.

Let \(k\) be the individual elimination rate constant, in h\(^{-1}\), and let \(V\) be the apparent volume of distribution in L. After a single intravenous bolus dose \(D\), the one-compartment concentration at time \(t\) is

C(t\mid k)=\frac{D}{V}\exp(-kt).

Assume the observed concentration \(C_{\mathrm{obs},j}\) at time \(t_j\) has lognormal error with log-scale standard deviation \(\sigma\). Equivalently, conditional on \(k\), its logarithm is normally distributed around the logarithm of the model prediction:

\log C_{\mathrm{obs},j} \mid k \sim N\!\left( \log\!\left[\frac{D}{V}\exp(-kt_j)\right], \sigma^2 \right).

This model assumes the dose and sample time are correctly recorded, the dose was administered and absorbed as represented by the model, the assay is appropriate, and the residual errors are independent conditional on \(k\). It also assumes \(V\) is known for this calculation. If there are several observations, their conditional likelihoods multiply under the independence assumption.

For a finite set of plausible parameter values \(k_1,\ldots,k_m\), with prior probabilities \(\pi_i\), the posterior probability of \(k_i\) is

P(k_i\mid \mathbf{C}_{\mathrm{obs}}) = \frac{\pi_i L_i} {\sum_{r=1}^{m}\pi_r L_r}, \qquad L_i=P(\mathbf{C}_{\mathrm{obs}}\mid k_i).

For one lognormal observation, a factor common to all parameter values can be omitted when calculating the normalized probabilities, leaving the likelihood kernel

L_i \propto \exp\!\left[ -\frac{ \left\{ \log\!\left(C_{\mathrm{obs}}/C(t\mid k_i)\right) \right\}^{2}} {2\sigma^{2}} \right].

Worked Example: Updating Elimination-Rate Probabilities

Suppose a hypothetical patient receives a 100 mg intravenous bolus. Use \(V=20\) L, so the initial model concentration is \(100/20=5\) mg/L. The population prior assigns probabilities 0.25, 0.50, and 0.25 to \(k=0.08,\ 0.10,\) and \(0.12\) h\(^{-1}\), respectively. A concentration of 1.45 mg/L is measured 12 h after the dose. Assume log-scale residual standard deviation \(\sigma=0.15\). Calculate the posterior probabilities.

First calculate the model concentration at 12 h for each value:

C(12\mid k)=5\exp(-12k): \quad 1.9145,\quad 1.5060,\quad 1.1846\ \text{mg/L}.

The likelihood kernels are approximately 0.1798, 0.9686, and 0.4034. For example, for \(k=0.08\) h\(^{-1}\),

\exp\!\left[ -\frac{\{\log(1.45/1.9145)\}^{2}}{2(0.15)^2} \right] \approx \exp(-1.716) \approx 0.1798.

Here \(\log(1.45/1.9145)=-0.2779\) and \(2(0.15)^2=0.045\), so the exponent is \(-0.2779^2/0.045\approx-1.716\) and the kernel is approximately \(0.1798\). The prior-weighted kernels are therefore \(0.25(0.1798)=0.04495\), \(0.50(0.9686)=0.4843\), and \(0.25(0.4034)=0.1009\). Their sum is approximately 0.6301, giving

P(k_i\mid C_{\mathrm{obs}}) \approx (0.0713,\ 0.7686,\ 0.1600).

The middle value, \(0.10\) h\(^{-1}\), receives the greatest posterior probability because its predicted concentration is closest to the observation. The posterior does not prove that this is the patient’s exact elimination rate; it quantifies uncertainty under the stated model and prior.

Forecasting Concentrations Under Candidate Regimens

For repeated intravenous boluses of dose \(D\) every \(\tau\) hours, with linear one-compartment elimination, the steady-state trough immediately before the next dose is the sum of the residual concentrations from all earlier doses:

C_{\min,\mathrm{ss}}(k;D,\tau) = \frac{D}{V} \left[ \exp(-k\tau)+\exp(-2k\tau)+\cdots \right] = \frac{D/V}{\exp(k\tau)-1}.

The equality follows from the geometric series: for \(r=\exp(-k\tau)\), \(0<r<1\), the sum \(r+r^2+\cdots\) is \(r/(1-r)=1/(\exp(k\tau)-1)\). The formula presumes the repeated-dose regimen has reached steady state. For a finite set of posterior parameter values, the forecast mean of the latent trough concentration is

E(C_{\min,\mathrm{ss}}\mid\mathbf{C}_{\mathrm{obs}},D,\tau) = \sum_{i=1}^{m} P(k_i\mid\mathbf{C}_{\mathrm{obs}}) C_{\min,\mathrm{ss}}(k_i;D,\tau).

This is an average of the regimen-specific predictions over the posterior. It is generally not equal to the prediction at the posterior mean of \(k\), because the concentration function is nonlinear. The posterior distribution of latent troughs can also be reported directly; adding future assay error would produce a predictive distribution for a future measured concentration.

Worked Example: Forecasting Two Schedules With the Posterior

Use the posterior probabilities from the preceding example to forecast steady-state troughs for two hypothetical regimens: 100 mg every 24 h and 50 mg every 12 h. Keep \(V=20\) L and the three posterior values of \(k\). These calculations compare model-based exposure forecasts; they do not establish which regimen is appropriate for a patient.

For 100 mg every 24 h, the trough predictions for \(k=0.08,\ 0.10,\) and \(0.12\) h\(^{-1}\) are approximately 0.8590, 0.4988, and 0.2974 mg/L. For example, at \(k=0.10\) h\(^{-1}\),

C_{\min,\mathrm{ss}} = \frac{100/20}{\exp(0.10\times24)-1} = \frac{5}{\exp(2.4)-1} \approx 0.4988\ \text{mg/L}.

Weighting all three troughs by their posterior probabilities gives

(0.0713)(0.8590)+(0.7686)(0.4988)+(0.1600)(0.2974) \approx 0.4922\ \text{mg/L}.

For 50 mg every 12 h, the corresponding trough predictions are approximately 1.551, 1.078, and 0.7762 mg/L. Their posterior-weighted mean is

(0.0713)(1.551)+(0.7686)(1.078)+(0.1600)(0.7762) \approx 1.063\ \text{mg/L}.

Thus, under this bolus model, the more frequent schedule has a higher predicted trough despite the same total daily dose. The discrete posterior assigns probabilities 0.0713, 0.7686, and 0.1600 to the three trough values for either schedule. For the 100 mg every 24 h schedule, the posterior median trough is 0.4988 mg/L, while the central 95% interval spans approximately 0.2974 to 0.8590 mg/L. This wide interval reflects the coarse three-value prior and substantial remaining parameter uncertainty; it is not a 95% probability statement about a fixed parameter.

Additional Measurements and Model Checking

A further concentration can update the posterior if its sampling time is informative and the observation model is credible. Under conditional independence given \(k\), the new posterior is proportional to the current posterior probability multiplied by the new observation’s likelihood. The next example shows the arithmetic and also illustrates how a measurement can concentrate the posterior.

Worked Example: Updating With a Second Sample

In the same hypothetical single-dose case, suppose a second concentration is 0.48 mg/L at 24 h after the dose. Assume the same \(\sigma=0.15\) and independent log-scale residual error conditional on \(k\). The model predictions at 24 h are approximately 0.7330, 0.4536, and 0.2807 mg/L. The corresponding likelihood kernels are approximately 0.0186, 0.9313, and 0.00166.

Multiply these by the previous posterior probabilities, not by the original prior probabilities: the unnormalized weights are approximately \(0.0713(0.0186)=0.00133\), \(0.7686(0.9313)=0.7158\), and \(0.1600(0.00166)=0.000266\). The sum is approximately 0.7174, so the updated probabilities are

P(k_i\mid C_{12}=1.45,\ C_{24}=0.48) \approx (0.0018,\ 0.9978,\ 0.0004).

The second observation is close to the prediction for \(k=0.10\) h\(^{-1}\), so the posterior becomes concentrated there. This apparent precision is conditional on the assumed model, known \(V\), and independent residuals. If the samples share a systematic assay bias, or if the actual dose time is wrong, treating their errors as independent can overstate the information.

Assumptions, Diagnostics, and Alternatives

A useful Bayesian forecast requires more than applying Bayes’ rule. Check the dose history, actual sample times, assay units, and model structure before interpreting posterior summaries. Plot observed concentrations against time alongside predictions from posterior parameter values. Systematic residual trends, such as later concentrations consistently above predictions, may indicate that the one-compartment elimination model, dose history, or assumed error model is inadequate. A single concentration can also be weakly informative when multiple parameters, such as \(V\) and \(k\), are unknown and can compensate for one another.

The calculation here treats the population information as a discrete prior over \(k\). In applied analyses, population PK models often provide continuous distributions and covariate-dependent priors. Numerical integration, simulation, or Bayesian sampling can then approximate posterior forecasts. Regardless of computation method, the forecast must propagate the posterior parameter uncertainty. If population parameters or covariates are uncertain, that uncertainty may also need to be propagated.

Approach	What it does	Main limitation
Population prediction	Uses the population model without conditioning on this patient’s TDM data	Does not incorporate individual concentration evidence
Plug-in individual estimate	Predicts using one fitted parameter value, such as a posterior mean or maximum a posteriori estimate	Can conceal parameter uncertainty and nonlinear forecast variation
Bayesian posterior forecast	Averages or otherwise summarizes predictions over parameter values weighted by their posterior probabilities	Depends on the prior, PK model, and observation assumptions being appropriate

A plug-in forecast is convenient, but for a nonlinear model \(E\{C(k)\mid\text{data}\}\) need not equal \(C\{E(k\mid\text{data})\}\). Bayesian forecasts address this by carrying the parameter distribution through the prediction. They do not automatically account for model misspecification, nor do they establish that a forecasted concentration will produce a particular clinical outcome.

Common Pitfalls and Practical Use

Using the prior as if it were the posterior. Population predictions do not incorporate the individual’s measured concentrations; update the prior with the likelihood.
Confusing a likelihood with a probability of a parameter. A likelihood compares how well parameter values explain the data. Posterior probabilities require prior weights and normalization.
Reporting only a point forecast. A posterior mean alone can hide a broad or multimodal forecast distribution; report uncertainty and, where useful, posterior quantiles.
Ignoring the observation process. Sampling-time error, assay error, incorrect dose records, or correlated residuals can distort the update.
Extrapolating beyond the model. A one-compartment bolus model does not represent every route, absorption process, time-varying clearance, or nonlinear disposition.

In practice, posterior forecasts help organize the evidence available from a population model and a patient’s concentrations. Their value is greatest when the sampling schedule is informative, the prior applies to the patient, and predictive checks support the PK model. A forecast should be interpreted as conditional on those choices, not as a guarantee of future concentration or response.

Check Your Understanding

Answer each question, then compare with the solutions below.

In the one-compartment model, what does \(k\) represent, and what units must it have if time is measured in hours?
Why are prior-weighted likelihoods divided by their sum when calculating discrete posterior probabilities?
For a fixed \(k\), how does the steady-state trough formula change if the dose is halved while \(V\) and \(\tau\) remain fixed?
Why can a posterior-weighted mean forecast differ from the forecast calculated at the posterior mean of \(k\)?
Name two diagnostics or data checks that could reveal a poor basis for a TDM forecast.

Solutions. (1) \(k\) is the first-order elimination rate constant, with units h\(^{-1}\), so that \(kt\) is dimensionless. (2) Bayes’ rule requires posterior probabilities to sum to one; dividing each prior-times-likelihood weight by their total performs this normalization. (3) The trough is proportional to \(D\) in the stated linear model, so halving \(D\) halves the predicted trough. (4) The concentration prediction is nonlinear in \(k\), so averaging predictions over parameter values generally differs from evaluating the function at the average parameter. (5) Examples include checking dose and sample-time records, plotting observed concentrations against posterior predictions, and examining residuals for systematic trends or unexpectedly large discrepancies.