Coverage Calibration Across Simulated Trials

Coverage Across Repeated Trials

A 90% posterior credible interval assigns 90% posterior probability to a parameter, conditional on the observed data and the fitted model. That statement alone does not guarantee that intervals from repeated studies will contain the true parameter 90% of the time. To assess that repeated-sampling property, simulate many trials under a specified pharmacokinetic model, refit each trial, and record whether each interval contains the generating value.

For a pharmacokinetic model, the parameters of interest might include clearance \(CL\), volume \(V\), and a diffusion scale \(\sigma\). The simulation must represent the full trial design being evaluated: for example, dose and sampling times, between-subject variation, stochastic paths, and observation error. Each replicate is then analyzed using the same fitting procedure and interval construction intended for the application.

Definition: Fix a true parameter value \(\theta_0\), a trial-generating process, a fitting procedure, and an interval rule. For a nominal \(100(1-\alpha)\%\) interval \(C(Y)\) computed from trial data \(Y\), its repeated-sampling coverage at \(\theta_0\) is \(P_{\theta_0}\{\theta_0\in C(Y)\}\). With \(M\) independently simulated trials, the estimated coverage is \(\hat p=M^{-1}\sum_{m=1}^{M} I_m\), where \(I_m=1\) if \(\theta_0\in C(Y_m)\), and \(I_m=0\) otherwise.

For marginal coverage, \(C(Y)\) is an interval for one parameter at a time. A joint region instead requires all specified parameter values to be covered simultaneously; its indicator is one only when the entire true parameter vector lies in the region. These are different events and need not have the same coverage. The earlier material on Marginal and Joint Calibration treats that distinction; here the emphasis is repeated-trial estimation and diagnosing departures from the nominal level.

Estimating Coverage and Its Simulation Error

If the simulated trials are independent and generated under the same fixed truth, the indicators \(I_1,\ldots,I_M\) are independent Bernoulli random variables with success probability \(p=P_{\theta_0}\{\theta_0\in C(Y)\}\). Therefore the number of covering intervals, \(K=\sum_{m=1}^{M}I_m\), has a binomial distribution with parameters \(M\) and \(p\).

Result: Under independent, identically generated trials, \(\hat p=K/M\) is unbiased for \(p\), and \(\operatorname{Var}(\hat p)=p(1-p)/M\). Since \(p\) is unknown, its estimated Monte Carlo standard error is \(\sqrt{\hat p(1-\hat p)/M}\).

To verify the result, \(E[I_m]=p\) and \(\operatorname{Var}(I_m)=p(1-p)\) for every replicate. Linearity of expectation gives \(E[\hat p]=M^{-1}\sum_{m=1}^{M}E[I_m]=p\). Independence makes the covariances between distinct indicators zero, so

\operatorname{Var}(\hat p) =\frac{1}{M^2}\sum_{m=1}^{M}\operatorname{Var}(I_m) =\frac{Mp(1-p)}{M^2} =\frac{p(1-p)}{M}.

Replacing the unknown \(p\) in this variance by \(\hat p\) gives the stated estimated standard error. This measures simulation noise in the estimated coverage; it does not measure uncertainty about whether the fitted pharmacokinetic model is scientifically appropriate. If trials share random numbers, are generated from different truths, or are otherwise dependent, the binomial variance calculation need not apply as stated. \(\square\)

Worked Example: Clearance Coverage in Ten Simulated Trials

Suppose the generating clearance is \(CL_0=5.0\) L/h, and each of ten trials produces a nominal 90% credible interval. The intervals are:

Trial	Interval for \(CL\), L/h	Width, L/h	Covers 5.0?
1	[4.5, 5.4]	0.9	Yes
2	[4.7, 5.3]	0.6	Yes
3	[5.1, 5.8]	0.7	No
4	[4.6, 5.2]	0.6	Yes
5	[4.9, 5.5]	0.6	Yes
6	[4.2, 4.9]	0.7	No
7	[4.8, 5.1]	0.3	Yes
8	[4.4, 5.6]	1.2	Yes
9	[5.0, 5.7]	0.7	Yes
10	[4.1, 4.8]	0.7	No

Seven intervals contain 5.0, so \(\hat p=7/10=0.70\). The estimated Monte Carlo standard error is \[ \sqrt{\frac{0.70(1-0.70)}{10}}=\sqrt{0.021}\approx 0.145. \] This estimate is too imprecise to conclude much about whether the procedure truly misses 90% coverage. The widths sum to \(7.0\) L/h, giving mean width \(7.0/10=0.70\) L/h. Both quantities should be reported: coverage describes containment, while width describes interval precision.

Comparing Parameters and Reading Simulation Uncertainty

Coverage should be calculated separately for clearance, volume, and diffusion parameters when those are the scientific targets. A parameter can be poorly calibrated even when the other parameters perform adequately. For each parameter, report the number of trials, the number of covering intervals, the estimated coverage, and its Monte Carlo standard error. Also summarize interval widths, since a procedure can attain coverage by producing intervals too broad to be useful.

Worked Example: Comparing Clearance, Volume, and Diffusion

In 100 independently simulated trials, suppose nominal 90% intervals cover the generating values in 91 trials for clearance, 87 for volume, and 72 for the diffusion scale \(\sigma\). The estimated coverages are \(91/100=0.91\), \(87/100=0.87\), and \(72/100=0.72\).

The plug-in Monte Carlo standard errors are, respectively, \[ \sqrt{\frac{0.91(0.09)}{100}}\approx 0.0286,\qquad \sqrt{\frac{0.87(0.13)}{100}}\approx 0.0336,\qquad \sqrt{\frac{0.72(0.28)}{100}}\approx 0.0449. \] For the diffusion parameter, the difference from the nominal target is \(0.72-0.90=-0.18\). Under the nominal value, the binomial standard error would be \(\sqrt{0.90(0.10)/100}=0.03\), so this discrepancy is six nominal standard errors. It is a strong signal to investigate the diffusion interval procedure, although the result alone does not identify the cause. The clearance and volume estimates are much closer to 0.90 relative to their simulation uncertainty.

These calculations assume the 100 trials are independent and use the same true parameter values and trial design. If the goal is to assess performance over a range of truths or designs, organize and report coverage by those conditions rather than pooling unlike trials into one percentage.

How Bias and Width Affect Coverage

Coverage depends on where intervals lie relative to the truth and how wide they are. A useful derivation isolates those effects under a simple approximation. Suppose an interval is centered on an estimator \(\hat\theta\), with fixed half-width \(h>0\), and across repeated trials \(\hat\theta-\theta_0\sim N(b,s^2)\). Here \(b=E[\hat\theta]-\theta_0\) is the estimator's bias and \(s\) is its repeated-sampling standard deviation. The interval is \([\hat\theta-h,\hat\theta+h]\).

Result: Under the normal approximation just stated, repeated-sampling coverage is \[ P\{\theta_0\in[\hat\theta-h,\hat\theta+h]\} =\Phi\left(\frac{h-b}{s}\right)-\Phi\left(\frac{-h-b}{s}\right), \] where \(\Phi\) is the standard normal cumulative distribution function.

To derive the expression, containment is equivalent to \(\hat\theta-h\leq\theta_0\leq\hat\theta+h\), or \(-h\leq\hat\theta-\theta_0\leq h\). Standardizing the normally distributed error \(Z=\hat\theta-\theta_0\) gives

P(-h\leq Z\leq h) =P\left(\frac{-h-b}{s}\leq\frac{Z-b}{s}\leq\frac{h-b}{s}\right) =\Phi\left(\frac{h-b}{s}\right)-\Phi\left(\frac{-h-b}{s}\right).

Thus, increasing \(h\) increases coverage in this fixed-center, fixed-bias setting, while nonzero bias shifts the error distribution away from zero and can reduce coverage. The result is an explanatory approximation, not a general formula for data-dependent posterior intervals: in an actual simulation, both interval location and width may vary with the trial.

Worked Example: A Biased Diffusion Estimate With Narrow Intervals

Consider a positive diffusion parameter \(\sigma\), and suppose the interval is constructed on the log scale. Let \(\theta_0=\log(\sigma_0)\). Assume the repeated-sampling log-scale estimate satisfies \(\hat\theta-\theta_0\sim N(0.10,0.20^2)\), and each nominal interval has log-scale half-width \(h=1.96(0.12)=0.2352\). Exponentiating its endpoints preserves containment because the exponential function is strictly increasing.

The coverage formula gives \[ \Phi\left(\frac{0.2352-0.10}{0.20}\right) -\Phi\left(\frac{-0.2352-0.10}{0.20}\right) =\Phi(0.676)-\Phi(-1.676) \approx 0.7505-0.0469=0.7036. \] So coverage is about 70%, despite the intervals being labelled nominal 95% under the posterior construction. The positive bias moves the interval centers above the truth, and the relatively small half-width does not compensate.

For comparison, if the same half-width were unbiased (\(b=0\)), coverage would be \(2\Phi(0.2352/0.20)-1=2\Phi(1.176)-1\approx 0.760\). If the biased procedure instead had \(h=0.392\), then coverage would be \(\Phi(1.46)-\Phi(-2.46)\approx0.9279-0.0069=0.921\). Increasing width can restore coverage in this example, but it does not remove bias and may make inference less informative.

Simulation Design, Diagnostics, and Limitations

Specify the target.
Choose fixed generating values for \(CL\), \(V\), and the diffusion parameter, along with the dose, sampling schedule, population structure, and observation model to be evaluated.

Generate independent trials.
For each replicate, simulate the intended trial under the specified model, including relevant random effects, stochastic paths, and measurement error.

Refit and form intervals.
Apply the same estimation, convergence checks, and credible-interval rule to every replicate. Record endpoints, widths, and coverage indicators.

Summarize and investigate.
Report coverage and Monte Carlo error by parameter, then inspect width and estimation error across replicates to identify plausible failure mechanisms.

Before interpreting a coverage deficit as a property of the inferential method, check whether every replicate completed successfully and whether posterior computation was reliable. Failed fits, divergent sampling, or inaccurate approximations can distort both endpoints and coverage. Record failures rather than silently discarding them: a coverage estimate based only on successful trials may describe a selected subset of the simulated studies.

For each parameter, useful diagnostics include the distribution of interval widths, the distribution of point-estimate errors \(\hat\theta-\theta_0\), and plots of interval endpoints relative to the truth. Many intervals entirely on one side of the truth suggest a location or bias problem; intervals frequently close to the truth but just too narrow suggest a width or uncertainty problem. These patterns are clues, not proofs of a cause. Dependence between parameters, skewed posterior distributions, and boundary constraints can make a simple bias-and-width explanation incomplete.

Fixed-truth repeated-sampling coverage answers a conditional question: how often would this procedure cover this specified truth under this specified trial design? It is distinct from averaging performance over a prior distribution of parameter values. It is also distinct from posterior probability conditional on one observed dataset. Results may change with dose, sampling schedule, sample size, parameter value, or model misspecification; a calibration result at one setting does not establish universal performance. For positive parameters, intervals on a log scale and intervals on the original scale have the same coverage when one is exactly the monotone transformation of the other, but their reported widths are on different scales.

Common Pitfall

Common mistake: “The interval is a 90% credible interval, so 90% of repeated intervals must contain the true value.” The first statement is a posterior probability conditional on the data and model; the second is a repeated-sampling property of the complete procedure. It must be assessed by simulation or justified by a separate result. Another mistake is to interpret a low estimated percentage without its Monte Carlo error: a small number of replicates can produce a noisy estimate even when the underlying coverage is nominal.

Coverage is only one performance measure. A method can attain nominal coverage with very wide intervals, while a narrow method can miss the truth too often. Report coverage together with interval width, and investigate systematic displacement as well as computational and design-related causes.

Check Your Understanding

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

In 200 independent trials, 174 intervals contain the fixed true clearance. What is estimated coverage and its plug-in Monte Carlo standard error?
For a nominal 90% procedure, what is the expected number of covering intervals in 500 independent trials if its true coverage is 0.90?
Under the normal approximation, let \(b=0\), \(s=0.25\), and \(h=0.49\). What is coverage, approximately?
Why is a 90% posterior probability statement not by itself a guarantee of 90% repeated-trial coverage?

Solutions. (1) The estimate is \(174/200=0.87\). Its plug-in standard error is \(\sqrt{0.87(0.13)/200}=\sqrt{0.0005655}\approx0.0238\). (2) The expected count is \(500(0.90)=450\), since each indicator has expectation 0.90 and expectations add. (3) With zero bias, coverage is \(2\Phi(h/s)-1=2\Phi(0.49/0.25)-1=2\Phi(1.96)-1\approx0.95\). (4) Posterior probability conditions on the observed data and fitted model, whereas repeated-trial coverage averages the interval's containment indicator over new datasets generated at a fixed truth. Those are different probabilities.