Tutorials › Population Pharmacokinetics › Parametric Bootstrap From a Fitted PopPK Model

Bootstrap and External Validation · Tutorial 406 of 599

Parametric Bootstrap From a Fitted PopPK Model

Build a parametric bootstrap that reproduces the study design and fitted PopPK model, and interpret what its refitted estimates can—and cannot—show.

Advanced 12 min read

What You'll Learn

  • Distinguish model-based parametric bootstrap samples from subject-level resampling.
  • Specify which covariates, random effects, residual errors, and sampling times are simulated.
  • Work through concentration generation under a hypothetical one-compartment model.
  • Explain why refitting each simulated dataset approximates the estimator’s sampling behavior.
  • Identify assumptions, diagnostics, and limitations of the resulting bootstrap distribution.

What a Parametric Bootstrap Repeats

A population PK analysis has two sources of structure: a model for how concentrations arise and a study design that determines who is observed, when, and under what conditions. A parametric bootstrap uses the fitted model to generate new datasets under that structure, then applies the estimator again to each dataset. It asks: if the fitted model and design were repeated, how would the estimation procedure behave?

This differs from resampling observed subjects. A subject-level bootstrap draws subjects from the empirical dataset, whereas a parametric bootstrap generates new observations from specified probability models. It can represent variation in random effects, residual errors, covariates, and sampling times, but its conclusions depend on the adequacy of those models.

Definition: Let \(D\) denote the observed dataset, \(\hat{\theta}\) the fitted population-model parameters, and \(T\) the estimation procedure. A parametric bootstrap replicate is generated as \(D^* \sim p_{\hat{\theta}}(\,\cdot\mid\text{design})\), and the model is refitted to obtain \(\hat{\theta}^*=T(D^*)\). Repeating this generates a conditional, model-based approximation to the sampling behavior of the estimator.

The parameter vector \(\theta\) commonly includes population parameters such as typical clearance and volume, covariate effects, random-effects covariance parameters, and residual-error parameters. In this tutorial, \(\Omega\) denotes the covariance matrix of the individual random effects. A fitted value such as \(\hat{\theta}\) is treated as the generator of bootstrap data; it is not assumed to equal the unknown true parameter.

Specify the Data-Generating Model

Consider a hypothetical one-compartment intravenous-bolus model with lognormal individual clearance and volume. For subject \(i\), let \(WT_i\) be weight in kilograms and let \(\eta_{CL,i}\) and \(\eta_{V,i}\) be individual random effects. One possible model is

$$ CL_i=\theta_{CL}\left(\frac{WT_i}{70}\right)^{\theta_{WT}}\exp(\eta_{CL,i}),\qquad V_i=\theta_V\left(\frac{WT_i}{70}\right)\exp(\eta_{V,i}), $$

with conditional concentration

$$ C_{ij}^{\mathrm{pred}}=\frac{D_i}{V_i}\exp\left(-\frac{CL_i}{V_i}t_{ij}\right),\qquad C_{ij}=C_{ij}^{\mathrm{pred}}\exp(\epsilon_{ij}). $$

Here \(D_i\) is dose in milligrams, \(t_{ij}\) is time in hours, \(CL_i\) is clearance in liters per hour, and \(V_i\) is volume in liters. The residual error \(\epsilon_{ij}\) is on the log scale. For example, \(\epsilon_{ij}\sim N(0,\sigma^2)\) gives a multiplicative error and positive simulated concentrations. This model treats \(C_{ij}^{\mathrm{pred}}\) as the conditional median; the arithmetic mean conditional on the random effects is \(C_{ij}^{\mathrm{pred}}\exp(\sigma^2/2)\).

A bootstrap generator must specify more than the concentration equation. It also needs a rule for covariates and the observation design. If the target is repeated sampling under a population covariate distribution, simulate covariates from an appropriate model or sampling scheme. If the target is uncertainty conditional on the enrolled subjects, retain their covariates instead. Likewise, either retain the observed sampling schedule or generate schedules from a stated design model. The choice changes the question being answered.

Typical assumptions include independent subjects conditional on the model, the stated distribution of random effects, the stated residual-error distribution, and a correctly represented covariate and sampling-time mechanism. Random effects are shared across observations from the same individual; residual errors are generated at the observation level. If random effects are correlated, draw their vector jointly from the fitted multivariate distribution rather than simulating each component independently. Informative dropout or sampling times related to a subject’s evolving condition require an appropriate joint model; treating them as independent when they are not can distort the bootstrap.

Worked Example: Generating One Subject’s Concentration

Suppose a fitted hypothetical model has \(\hat{\theta}_{CL}=5.0\) L/h, \(\hat{\theta}_{WT}=0.75\), and \(\hat{\theta}_{V}=50\) L. Generate a subject weighing 56 kg with simulated random effects \(\eta_{CL}=-0.20\) and \(\eta_V=0.10\). The subject receives 100 mg. For illustration, the design assigns each sample time independently as 0.5, 2, or 6 hours with probabilities 0.20, 0.50, and 0.30. A uniform draw of 0.63 selects 2 hours because it falls in the interval from 0.20 to 0.70.

The covariate ratio is \(56/70=0.80\). Thus,

$$ CL_i=5.0(0.80)^{0.75}\exp(-0.20) =5.0(0.8459)(0.8187) \approx 3.463\text{ L/h}. $$

Similarly,

$$ V_i=50(0.80)\exp(0.10) =40(1.1052) \approx 44.21\text{ L}. $$

At 2 hours, the conditional prediction is

$$ C_{ij}^{\mathrm{pred}} =\frac{100}{44.21}\exp\left(-\frac{3.463}{44.21}\times 2\right) \approx 2.262\exp(-0.1567) \approx 1.934\text{ mg/L}. $$

If the simulated residual draw is \(\epsilon_{ij}=0.10\), then \(C_{ij}=1.934\exp(0.10)\approx 2.137\) mg/L. The random effect for clearance is used throughout this subject’s profile, while a new residual draw is used for each observation. The numerical values are hypothetical and illustrate one data-generation step, not a study result.

From Simulated Dataset to Bootstrap Distribution

A complete replicate applies the same estimation procedure to the generated dataset as to the original analysis. The fit should use the same structural model, parameterization, likelihood or estimation criterion, constraints, and relevant data-handling rules. Otherwise, variation in the bootstrap estimates may reflect a changed analysis procedure rather than repeated use of the estimator being evaluated.

1
Choose the target design.
Decide whether covariates and sampling times are held at their observed values or generated. State the target population and the design assumptions.
2
Generate subjects and observations.
Draw covariates when required, draw one random-effects vector per subject, generate the observation schedule, and then draw residual errors at the observations.
3
Refit the model.
Apply the original estimation procedure to the replicate and save the estimates and relevant fit information.
4
Repeat and summarize.
Across replicates, examine the distribution of each estimate and record failures or boundary behavior as well as successful fits.

The justification is a plug-in approximation. If the true data-generating parameter is \(\theta_0\), repeated datasets would produce estimates \(T(D)\) under \(D\sim p_{\theta_0}\). Since \(\theta_0\) is unknown, the bootstrap substitutes the fitted value \(\hat{\theta}\), generates \(D^*\sim p_{\hat{\theta}}\), and calculates \(T(D^*)\). When the model is adequate and the estimator behaves regularly, the distribution of \(\hat{\theta}^*-\hat{\theta}\), conditional on the observed data, approximates the sampling variation of \(\hat{\theta}-\theta_0\). This is an approximation, not a guarantee: model misspecification and weak identification are not repaired by generating more replicates.

A simple special case shows why refitting matters. Suppose a scalar estimate is based on \(n\) independent subject-level values \(Y_i=\theta+\eta_i\), where \(\eta_i\sim N(0,\omega^2)\), and the estimator is \(\bar{Y}\). Under the generating model,

$$ E(\bar{Y})=\theta,\qquad \operatorname{Var}(\bar{Y})=\operatorname{Var}\left(\frac{1}{n}\sum_{i=1}^n\eta_i\right) =\frac{1}{n^2}\sum_{i=1}^n\omega^2 =\frac{\omega^2}{n}. $$

A parametric bootstrap that draws new subject effects and recalculates \(\bar{Y}\) reproduces this estimator’s variation under the assumed model. In a nonlinear PK model, the refitting step is more complex, but the principle is the same: generate new data, then let the estimator respond to them.

Worked Example: Summarizing Hypothetical Refit Results

Suppose five hypothetical bootstrap datasets were generated from a fitted model with \(\hat{\theta}_{CL}=5.00\) L/h, and the same model was refitted successfully to each. The clearance estimates were 4.80, 5.30, 4.90, 5.10, and 4.70 L/h. Their mean is

$$ \bar{\theta}_{CL}^* =\frac{4.80+5.30+4.90+5.10+4.70}{5} =\frac{24.80}{5} =4.960\text{ L/h}. $$

The average difference from the generating value is \(4.960-5.00=-0.0400\) L/h for these five replicates. The sample standard deviation of the five estimates is calculated from deviations \(-0.16,\ 0.34,\ -0.06,\ 0.14,\ -0.26\) L/h:

$$ s^*=\sqrt{\frac{(-0.16)^2+(0.34)^2+(-0.06)^2+(0.14)^2+(-0.26)^2}{5-1}} =\sqrt{\frac{0.232}{4}} \approx 0.2408\text{ L/h}. $$

These results illustrate how replicate estimates can be summarized to assess behavior under the fitted generator. Five replicates are far too few for a stable characterization of estimator behavior; the numbers are for arithmetic illustration only. The mean difference here is a Monte Carlo summary relative to the generating value, not evidence that the fitted model is unbiased in reality.

Worked Example: Simulating a Sampling-Time Design

Suppose a hypothetical design has 10 subjects with two samples per subject, and each sampling time is independently selected from 0.5, 2, and 6 hours with probabilities 0.20, 0.50, and 0.30. There are \(10\times2=20\) observations. The expected numbers at those times are \(20(0.20)=4\), \(20(0.50)=10\), and \(20(0.30)=6\), respectively. These are expected counts, not fixed counts: one generated dataset may have fewer or more than each expectation.

For example, a uniform draw of 0.15 selects 0.5 hours, while a draw of 0.82 selects 6 hours. In this example, independent draws for every observation are part of the stated design assumption. If the real design instead assigns a fixed pair of times to each subject, or uses correlated visit schedules, the generator should preserve that structure rather than use independent time draws.

Diagnostics, Alternatives, and Limits

Check that simulated datasets follow the intended design before interpreting refits. Useful checks include covariate distributions and stratum counts, the distribution of random effects, numbers and timing of observations per subject, concentration ranges, and the frequency of values near assay limits. Compare these features with the assumptions used to build the generator and with the intended study design. A mismatch may reveal a coding error or an inappropriate simulation rule.

Refitting also tests how the estimation procedure responds to data generated by its own fitted model. Track convergence, estimates at parameter bounds, implausible estimates, and warnings; do not silently omit problematic replicates. A high failure rate is itself important evidence about estimator behavior in the simulated setting. It may also indicate a numerical or implementation problem, so inspect individual failures rather than assuming a single cause.

ApproachWhat is generatedMain strengthMain limitation
Parametric bootstrapNew data from specified fitted distributions and design rulesCan represent model structure and generate covariate or sampling patterns beyond those observedDepends on the assumed data-generating model
Subject-level nonparametric bootstrapResampled observed subjects, often with their records intactUses empirical subject-level variation with fewer distributional assumptionsCannot readily create covariate patterns absent from the observed sample
Asymptotic covariance approximationNo new datasets; uses local uncertainty from the fitted estimation procedureComputationally direct when its approximation is reliableMay be inaccurate with weak information, strong nonlinearity, or boundary estimates

The parametric approach is useful when the model and design are explicit and the question concerns the estimator under that model. It cannot establish that the fitted model describes the real population, and it does not automatically represent uncertainty about the model form or the generating assumptions. If the covariate distribution, residual structure, or sampling process is misspecified, the bootstrap can give a precise account of the wrong repeated-sampling scenario.

A common mistake is to treat the observed individual random effects as fixed subject characteristics and reuse them in every replicate. That fails to generate between-subject variation: the intended bootstrap draws new random effects from the fitted population distribution, unless the question explicitly conditions on particular individuals. A related mistake is to resample observed concentrations after also generating new residual errors; that mixes resampling and parametric generation and generally no longer corresponds to the stated model.

Key takeaway: A parametric bootstrap is defined jointly by its fitted model, covariate and sampling-time rules, random-effects and residual distributions, and refitting procedure. State each element, check that generated data obey it, and interpret the result as model-conditional estimator behavior.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

  1. In a bootstrap replicate, should one random effect be drawn per observation or per subject? Why?
  2. A fitted model uses a lognormal residual error. What is the consequence of drawing a new residual error for each subject rather than for each observation?
  3. In the sampling-time example, how many observations are expected at 2 hours, and does every replicate have to contain exactly that number?
  4. What does a parametric bootstrap estimate if covariates are simulated from a population distribution rather than held at their observed values?

Solutions. (1) Draw the random-effects vector once per subject and use it for that subject’s repeated observations, because it represents persistent between-subject differences. (2) It imposes a shared residual perturbation across a subject’s observations instead of the specified observation-level residual variation, changing the within-subject error structure. (3) The expected number is \(20(0.50)=10\); a replicate need not contain exactly 10 because the times are randomly drawn. (4) The bootstrap represents repeated sampling that includes variation in covariate composition under the chosen covariate-generation model, rather than uncertainty conditional on the enrolled subjects’ covariates.

Clinical Biostats