Tutorials › Pharmacokinetic-Pharmacodynamic (PK/PD) Modeling › Simulation Checks for Baseline, Placebo, and Progression Models

Baseline, Placebo and Disease-Progression Components · Tutorial 67 of 598

Simulation Checks for Baseline, Placebo, and Progression Models

Use repeated simulated data sets to test parameter recovery, quantify uncertainty, and identify weaknesses in baseline, placebo, and progression models before relying on their estimates.

Intermediate 12 min read

What You'll Learn

  • Specify a data-generating model and simulation design for baseline, placebo, and untreated progression.
  • Calculate parameter bias, relative bias, RMSE, and confidence-interval coverage from simulation replicates.
  • Assess how visit timing affects recovery of placebo time-course and progression parameters.
  • Distinguish poor parameter recovery from structural model misspecification.
  • Compare simulation checks with standard errors, bootstrap resampling, and sensitivity analyses.
  • Recognize convergence, residual, and identifiability warnings in simulation results.

Why Check Recovery by Simulation?

A model can describe the response pattern in a data set and still estimate its components unreliably. Baseline, placebo response, and untreated disease progression may overlap in time: a changing response after treatment initiation can be attributed to more than one component, especially when measurements are sparse. Simulation checks ask a practical question before interpreting fitted parameters: if the model and sampling plan were used repeatedly under stated conditions, would the fitting procedure recover the parameters that generated the data?

Parameter recovery is not proof that a model is biologically correct. It evaluates the consequences of a specified model, error structure, participant design, and sampling schedule. The result is conditional on those choices. A useful check therefore makes the assumptions explicit, varies plausible conditions, and examines both numerical performance and the patterns of failure.

Definition: A simulation-based parameter-recovery check specifies parameter values and a data-generating model, repeatedly generates data under a defined design, fits the proposed model to each simulated data set, and summarizes how closely the fitted estimates and uncertainty statements recover the specified values.

Specify the Model and the Simulation Conditions

Consider an illustrative untreated response measured relative to a common time origin. Let \(Y_{ij}\) denote the response for participant \(i\) at time \(t_{ij}\), in months. Use the mean model

$$ E(Y_{ij}) = B + A\left(1-e^{-t_{ij}/\tau}\right) + \beta t_{ij}. $$

Here \(B\) is the mean response at baseline, \(A\) is the eventual magnitude of a time-dependent placebo response, \(\tau\) is its time constant in months, and \(\beta\) is the linear untreated progression rate in response units per month. In this illustration, positive \(A\) and positive \(\beta\) both increase the response. The placebo term is zero at \(t=0\), and approaches \(A\) as time increases. This is one particular model, not a universal decomposition of clinical response.

For a simple simulation, assume observations equal this mean plus independent normal error:

$$ Y_{ij} = B + A\left(1-e^{-t_{ij}/\tau}\right) + \beta t_{ij} + \epsilon_{ij}, \qquad \epsilon_{ij}\sim N(0,\sigma^2). $$

The error assumption means that, conditional on the specified mean, each observation has mean zero error and variance \(\sigma^2\); observations are independent in this basic version. A real design may require participant-level random effects, correlated residuals, dropout, or other features. If those features are relevant, they should be included in the simulation rather than ignored. The fitting procedure should match the proposed analysis, including any constraints and estimation method.

A recovery check needs more than parameter values. Specify the number and timing of visits, number of participants, error variability, participant variability if included, missingness assumptions, and the number of simulation replicates. Repeating the same design under different plausible parameter values or error conditions helps reveal whether performance depends on a particularly favorable scenario.

Summarize Recovery and Uncertainty

Suppose \(R\) simulated data sets are fitted successfully, giving estimates \(\hat{\theta}_1,\ldots,\hat{\theta}_R\) for a parameter whose generating value is \(\theta\). The average estimate is \(\bar{\hat{\theta}}=R^{-1}\sum_{r=1}^{R}\hat{\theta}_r\). Its difference from the generating value is the simulation bias. For a nonzero parameter, relative bias expresses that difference as a percentage of the generating value.

$$ \text{Bias} = \bar{\hat{\theta}}-\theta, \qquad \text{Relative bias (\%)}=100\frac{\bar{\hat{\theta}}-\theta}{\theta}. $$

The empirical standard deviation of the estimates describes their spread across simulated data sets. Root mean squared error (RMSE) combines spread and bias, so a method can have a small spread but still perform poorly if it is systematically biased. For each replicate, a nominal 95% confidence interval can also be computed using the fitting procedure. Empirical coverage is the fraction of those intervals that contain the generating value.

$$ \text{RMSE} = \sqrt{\frac{1}{R}\sum_{r=1}^{R}(\hat{\theta}_r-\theta)^2}, \qquad \text{Coverage} = \frac{1}{R}\sum_{r=1}^{R}I(\theta\text{ is in interval }r). $$

Here \(I(\cdot)\) equals 1 when the statement is true and 0 otherwise. Coverage near 95% is a long-run property under the specified simulation conditions, not a claim that a particular interval has a 95% probability of containing a fixed parameter. Report convergence failures as well as summaries among successful fits: excluding failures without disclosure can make a weak method look reassuring.

Worked Example: Calculating Bias and RMSE

Suppose 500 hypothetical data sets are generated with progression rate \(\beta=0.400\) response units per month and fitted using the generating model. Across the fits, the mean estimate is \(0.392\), and the empirical standard deviation is \(0.090\) response units per month. The mean squared deviation from the generating value, calculated across all 500 estimates, is \(0.008164\) squared response units per month squared.

Bias. The average estimate is below the generating value by

$$ 0.392-0.400=-0.008 \text{ response units per month}. $$

Relative bias. The percentage difference is

$$ 100\frac{-0.008}{0.400}=-2.00\%. $$

RMSE. Taking the square root of the stated mean squared deviation gives

$$ \sqrt{0.008164}=0.0904 \text{ response units per month}. $$

The average estimate has little bias in this hypothetical run, while its across-replicate spread is appreciable relative to the rate. The RMSE is close to the empirical standard deviation because the bias is small. These values alone do not establish acceptable performance; the acceptable uncertainty depends on the scientific use of the estimate.

Visit Timing Can Determine What Is Recoverable

The design must contain information about each component. A baseline measurement directly informs \(B\), while early post-baseline observations help distinguish the onset and curvature of the placebo term. Later observations may be needed to characterize progression over time and to separate a continuing trend from a placebo response approaching its limit. If all follow-up visits occur in a narrow time window, different combinations of \(A\), \(\tau\), and \(\beta\) can produce similar predictions there. That is a practical recovery problem even if an optimizer returns a numerical estimate.

For a nonlinear model, it is especially important to inspect joint behavior, not only each parameter’s average. Estimates of \(A\) and \(\tau\), for example, may trade off: a larger amplitude with a slower time course can resemble a smaller amplitude with a faster time course over a limited observation window. Examine scatterplots or correlations of estimates, interval widths, boundary estimates, and convergence rates alongside bias and coverage.

Worked Example: Comparing Two Sampling Schedules

Consider two hypothetical schedules for the same 80 participants, with independent observation error having standard deviation 3 response units. Schedule S measures at months 0, 1, 2, and 3; Schedule L measures at months 0, 1, 3, 6, and 12. In each of 500 simulated data sets, the generating values are \(B=100\) response units, \(A=6\) response units, \(\tau=2.5\) months, and \(\beta=0.400\) response units per month. Assume the same model is used for fitting. The following hypothetical simulation summaries illustrate how to compare recovery; they are not results from a clinical study.

ParameterGenerating valueSchedule S mean estimateSchedule S 95% coverageSchedule L mean estimateSchedule L 95% coverage
\(B\) (response units)100100.194%100.095%
\(A\) (response units)65.791%5.994%
\(\tau\) (months)2.53.482%2.694%
\(\beta\) (response units/month)0.4000.43088%0.40594%

For Schedule S, the mean estimate of \(\tau\) has bias \(3.4-2.5=0.9\) months, or \(100(0.9/2.5)=36\%\) relative bias. For Schedule L, its bias is \(2.6-2.5=0.1\) months, or \(4.00\%\). The late observations in Schedule L provide more opportunity to distinguish an approaching placebo plateau from ongoing progression in this example. The low coverage for \(\tau\) under Schedule S is also a warning: the stated interval procedure is not achieving its nominal rate under that design and scenario.

This comparison does not show that a later visit is always preferable. It demonstrates how simulation can evaluate a proposed schedule against a specific model and set of conditions. A design decision should also consider feasibility, missing visits, and whether the parameters are sufficiently precise for their intended interpretation.

Check Model Form as Well as Estimation

Recovery under a model that generates the data is an important baseline test, but it cannot establish that the model form is adequate for actual data. If the simulation generates data from the same equation used for fitting, it mainly evaluates estimation and design under that equation. To probe model-form risk, add scenarios in which data are generated from plausible alternatives, then fit the proposed model and examine bias, residual patterns, and coverage. Earlier material on disease-progression functional-form errors is relevant to choosing such alternatives; the simulation check quantifies their consequences under a design.

Worked Example: Detecting a Progression-Form Mismatch

Suppose a hypothetical true progression component is saturating, \(D(t)=12t/(6+t)\), in response units, while the fitted model assumes linear progression \(D(t)=\beta t\). To isolate the functional-form issue, suppose the baseline and placebo parameters are known and correctly specified. At months 1, 3, 6, and 12, the true progression values are respectively \(12/7=1.714\), \(36/9=4.000\), \(72/12=6.000\), and \(144/18=8.000\) response units.

For these follow-up times, a zero-intercept least-squares linear fit to the progression values has slope

$$ \hat{\beta} = \frac{\sum tD(t)}{\sum t^2} = \frac{1(1.714)+3(4)+6(6)+12(8)} {1^2+3^2+6^2+12^2} = \frac{145.714}{190} = 0.767 \text{ response units per month}. $$

The fitted linear values are 0.767, 2.301, 4.602, and 9.204 response units. Thus, observed-minus-fitted progression residuals are approximately 0.947, 1.699, 1.398, and \(-1.204\) response units. The residuals first rise and then turn negative rather than fluctuating without a time pattern. This systematic shape warns that a linear progression term does not represent the saturating generating process well, even though a best-fitting slope exists.

In a repeated simulation with measurement error, that pattern may be less obvious in any one data set. Across replicates, the fitted slope may also be precise but systematically different from a scientifically meaningful rate over the full time range. The remedy is not to declare the parameter recovered because the optimizer converged; evaluate the residual structure and test plausible alternative generating scenarios.

Practical Diagnostics and Alternatives

A useful simulation report should include the number of attempted and successful fits, convergence warnings, estimates at parameter boundaries, bias, empirical spread, RMSE, interval coverage, and interval width. Inspect representative fitted curves and residuals against time. Unexpected parameter correlations, broad or asymmetric estimate distributions, and a few extreme fits can be hidden by reporting only the mean and standard deviation. Increase the number of replicates when simulation summaries themselves are noisy; for example, a coverage estimate based on 500 replicates has finite Monte Carlo uncertainty and should not be treated as exact.

Simulation checks complement other approaches rather than replacing them. Model-based standard errors use the fitted model’s local uncertainty approximation and can be quick to compute, but may not reveal finite-sample bias or convergence problems. Bootstrap resampling can assess variability using observed data, but it generally cannot test a structural generating scenario absent from those data. Sensitivity analysis compares conclusions under alternative assumptions; simulation makes the repeated-sampling consequences of selected assumptions measurable. None removes the need for scientific judgment about plausible models and designs.

Common pitfall: A researcher simulates data from a proposed model, fits that same model, sees small average bias, and concludes that the model is validated. This only shows that the fitting procedure can recover parameters under the chosen model and simulation conditions. It does not show that the model describes the real response process. Include plausible misspecification scenarios and evaluate residual patterns, uncertainty, and failures as well as mean bias.

Check Your Understanding

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

  1. A parameter has generating value 2.0 months and mean estimate 2.3 months across simulations. Calculate its bias and relative bias.
  2. Why can a schedule with only early follow-up visits make placebo time course and progression difficult to distinguish?
  3. In 400 simulations, 360 nominal 95% confidence intervals contain the generating value. What is the empirical coverage, and what does it describe?
  4. A model converges in every replicate and has low average bias, but residuals repeatedly curve with time under an alternative generating scenario. What conclusion is justified?
  5. Name two useful checks beyond average parameter bias and explain what each can reveal.

Solutions. (1) Bias is \(2.3-2.0=0.300\) months; relative bias is \(100(0.300/2.0)=15.0\%\). (2) Over a short time range, the placebo component may not have approached its plateau, so changes in its amplitude or time constant can resemble a progression trend. This can produce correlated estimates and weak recovery. (3) Coverage is \(360/400=0.900\), or 90.0%. It is the fraction of intervals containing the fixed generating value across these simulated data sets; it is below the nominal 95% rate under the stated conditions. (4) Convergence and low average bias in the matched scenario do not validate the model. Curved residuals under a plausible alternative indicate potential functional-form failure and motivate evaluating that scenario and alternative models. (5) Convergence and boundary rates reveal numerical or estimation failures; empirical coverage checks uncertainty calibration. Residuals over time can reveal systematic model mismatch, while estimate correlations can reveal parameter trade-offs.

Clinical Biostats