Tutorials › Bayesian Clinical Trials › Diagnosing Monte Carlo Error in Predictive Success Estimates

Posterior and Predictive Probability Rules · Tutorial 167 of 696

Diagnosing Monte Carlo Error in Predictive Success Estimates

Use the number of simulated successes to quantify simulation error, judge whether an estimate is precise enough for a decision, and plan additional simulations when needed.

Intermediate 10 min read

What You'll Learn

  • Define the predictive success probability and its replicate-based Monte Carlo estimate.
  • Derive the binomial Monte Carlo standard error from independent simulated outcomes.
  • Calculate and interpret an approximate interval for Monte Carlo error.
  • Plan a simulation size to meet a stated precision target.
  • Recognize why zero simulated successes do not imply zero predictive success.
  • Diagnose dependence, coding errors, and other failures of the binomial error calculation.

Why a Predictive Success Estimate Needs an Error Check

A predictive success probability summarizes the chance that a future trial will meet its prespecified success criterion, given the current data and the specified model. It is often estimated by simulating many possible future trials and counting how many meet that criterion. The resulting estimate is itself subject to simulation error: a different set of random draws could produce a somewhat different count, even when the model and current data are unchanged.

This Monte Carlo error is distinct from the uncertainty about future trial outcomes represented by the predictive probability. The predictive probability describes uncertainty about what may happen in the future under the model. Monte Carlo error describes how accurately a finite simulation estimates that probability. A small Monte Carlo error does not establish that the model or success criterion is appropriate; it only indicates that the simulation has estimated the probability precisely under the stated setup.

Definition: Let \(\theta\) be the predictive success probability under a specified model, current data, and success rule. For each of \(M\) independent simulated future trials, define \(I_m=1\) if replicate \(m\) meets the success rule and \(I_m=0\) otherwise. The number of simulated successes is \(K=\sum_{m=1}^{M} I_m\), and the Monte Carlo estimate is \(\hat{\theta}_{MC}=K/M\). If the replicates are independent and each has success probability \(\theta\), then \(K\sim\operatorname{Binomial}(M,\theta)\).

The binomial description assumes that each replicate is generated from the intended predictive distribution, uses the same success rule, and is independent of the other replicates. In practice, a replicate may involve drawing unknown parameters from their posterior distribution, generating future trial data, and applying the final decision rule. The indicator \(I_m\) records only whether that complete replicate met the rule. A single simulated trial can have complicated internal structure; the binomial calculation concerns the resulting success or failure indicator across replicates.

Deriving the Monte Carlo Standard Error

Because an indicator with success probability \(\theta\) has variance \(\theta(1-\theta)\), independence gives the variance of the replicate average directly:

$$ \operatorname{Var}(\hat{\theta}_{MC}) =\operatorname{Var}\left(\frac{1}{M}\sum_{m=1}^{M} I_m\right) =\frac{1}{M^2}\sum_{m=1}^{M}\operatorname{Var}(I_m) =\frac{\theta(1-\theta)}{M}. $$

The standard deviation of this simulation estimate is the Monte Carlo standard error. Since \(\theta\) is unknown, it is estimated by substituting \(\hat{\theta}_{MC}\):

$$ \widehat{SE}_{MC} =\sqrt{\frac{\hat{\theta}_{MC}(1-\hat{\theta}_{MC})}{M}}. $$

This estimated standard error measures the typical fluctuation in \(\hat{\theta}_{MC}\) across repeated sets of \(M\) independent replicates under the same predictive distribution. It decreases approximately in proportion to \(1/\sqrt{M}\), not \(1/M\). Consequently, reducing the standard error by half requires about four times as many replicates.

Worked Example: Checking Precision for a High Predictive Success Estimate

Suppose \(M=2{,}000\) independent future-trial simulations are run under a specified predictive model. A replicate is counted as successful if it meets the prespecified final posterior go criterion. Of the simulations, \(K=1{,}604\) meet that criterion. Estimate the predictive success probability and its Monte Carlo error.

The estimated predictive success probability is \(\hat{\theta}_{MC}=1{,}604/2{,}000=0.8020\), or 80.20%. The estimated Monte Carlo standard error is

$$ \widehat{SE}_{MC} =\sqrt{\frac{0.8020(1-0.8020)}{2{,}000}} =\sqrt{\frac{0.158796}{2{,}000}} =\sqrt{0.000079398} \approx 0.0089106. $$

Thus, the estimated standard error is about 0.00891, or 0.891 percentage points. Under the normal approximation, an approximate 95% Monte Carlo interval is the estimate plus or minus 1.96 estimated standard errors:

$$ 0.8020 \pm 1.96(0.0089106) =0.8020\pm0.01746 =(0.7845,\;0.8195). $$

This interval describes approximate simulation precision for \(\theta\), conditional on the specified predictive distribution and success rule. It is not a 95% probability statement about the unknown model parameters, nor does it quantify uncertainty about whether the predictive model is scientifically correct.

If the goal is an approximate 95% half-width of at most 0.01 near \(\theta=0.80\), the normal approximation suggests planning for \(M\geq 1.96^2(0.80)(0.20)/(0.01)^2=6{,}146.56\). Rounding up gives 6,147 replicates. This is a planning approximation based on the assumed success probability; if that probability is uncertain, a conservative calculation can use \(\theta(1-\theta)\leq0.25\), which is largest at \(\theta=0.5\).

Using an Approximate Interval and a Precision Target

When \(M\) is large enough for a normal approximation to be reasonable and the estimated success and failure counts are not too small, a common approximate interval is \(\hat{\theta}_{MC}\pm z_{1-\alpha/2}\widehat{SE}_{MC}\). For a 95% interval, \(z_{1-\alpha/2}\) is approximately 1.96. This is an interval for the success probability governing the independent simulation indicators, not for a treatment effect or a parameter of the clinical model.

The normal approximation can behave poorly when the simulated success count or failure count is small, particularly when \(\hat{\theta}_{MC}\) is close to 0 or 1. A Wald interval can then be too narrow or extend outside the possible range from 0 to 1. Wilson or exact binomial intervals are alternatives in such cases. The underlying diagnostic question remains the same: how much could the estimated success probability vary simply because only a finite number of replicates were generated?

A useful way to decide whether more simulation is needed is to define an acceptable Monte Carlo precision before interpreting a borderline estimate. If a two-sided interval half-width no greater than \(\delta\) is desired, an approximate planning formula is \(M\geq z_{1-\alpha/2}^{\,2}\theta(1-\theta)/\delta^2\). Because \(\theta\) is unknown before simulation, an initial estimate can be used for a first calculation, or the conservative bound \(\theta(1-\theta)\leq0.25\) can be used to plan for any possible success probability. The required precision should be chosen in relation to the decision: a tiny numerical error may be unnecessary when the estimate is far from a decision threshold, while a more precise estimate can matter when simulation error could change the interpretation.

Worked Example: A Low Estimate Near a Decision Boundary

In a hypothetical analysis, \(500\) independent predictive replicates produce \(38\) successes. Suppose the trial team uses \(0.10\) as a reference value when interpreting predictive success. Assess the Monte Carlo precision of the estimate and whether the simulation clearly places it on one side of that reference.

The estimate is \(\hat{\theta}_{MC}=38/500=0.0760\). Its estimated standard error is

$$ \widehat{SE}_{MC} =\sqrt{\frac{0.0760(1-0.0760)}{500}} =\sqrt{\frac{0.070224}{500}} =\sqrt{0.000140448} \approx 0.01185. $$

The approximate 95% Monte Carlo interval is \(0.0760\pm1.96(0.01185)=0.0760\pm0.02323\), or approximately \((0.0528,0.0992)\). This interval is just below 0.10 at its upper endpoint, but the estimate is close to the reference and the normal approximation should be treated cautiously because there are only 38 simulated successes. A Wilson or exact binomial interval is preferable for a more reliable boundary assessment. If the decision depends on distinguishing values near 0.10, the simulation should be planned or extended to achieve adequate precision, and the chosen interval method should be reported.

For illustration, using the current estimate to plan for a normal-approximation half-width of 0.01 gives \(M\geq1.96^2(0.0760)(0.9240)/(0.01)^2\approx2{,}698\) total independent replicates. This is a planning estimate, not a guarantee: the true \(\theta\) is unknown, and the final precision should be checked using the actual replicate count and an appropriate interval.

Zero or Very Few Simulated Successes

A plug-in standard error of zero occurs if all simulations have the same outcome. For example, if \(K=0\), then \(\hat{\theta}_{MC}=0\) and the plug-in formula returns zero. That does not demonstrate that \(\theta=0\); it demonstrates only that no success was observed among the simulated replicates. A simulation can miss a rare event.

Worked Example: Interpreting Zero Successes

Suppose a hypothetical predictive simulation generates \(M=1{,}000\) independent future trials and observes no successes. The estimate is \(0/1{,}000=0\), but the plug-in standard error is not a useful indication of uncertainty. Find an exact one-sided 95% binomial upper confidence bound for the success probability.

If the success probability is \(\theta\), the probability of observing zero successes is \((1-\theta)^{1{,}000}\). Set this probability equal to 0.05 and solve for the upper bound:

$$ (1-\theta)^{1{,}000}=0.05 \quad\Longrightarrow\quad \theta=1-0.05^{1/1{,}000} \approx0.002991. $$

Thus, observing zero successes in 1,000 replicates is compatible with success probabilities as large as about 0.003 at this one-sided confidence level. The familiar approximation \(3/M=0.003\) is close here, but the bound does not prove that the predictive success probability is zero. More simulations may be needed if probabilities of this magnitude would affect the decision.

Diagnostics, Alternatives, and Limitations

Before using the binomial standard error, check how the simulation was produced. Each replicate should represent one independently generated future trial from the intended predictive distribution. The success rule must be applied consistently, including any posterior threshold, endpoint definition, and direction of effect. A coding error that changes the rule or fails to propagate uncertainty in a required model component is not repaired by increasing \(M\).

The binomial formula is not appropriate without adjustment when replicate indicators are dependent, as can occur when draws are produced by a correlated Markov chain or when simulation outputs share random components. Correlation means that \(K\) does not have the independent-binomial variance used above, and the naive standard error can be too small. Such simulations require an error assessment that accounts for dependence. Likewise, if replicates have different success probabilities because they are generated under different scenarios, pooling them does not estimate one common \(\theta\) without a clearly defined weighting target.

When an analytic predictive calculation is available and can be evaluated reliably, it can avoid Monte Carlo error for that calculation. Simulation is valuable when the future-trial and decision rules are too complex for convenient analytic evaluation. The practical alternatives are therefore not interchangeable: analytic evaluation addresses the same target without finite-replicate sampling error when feasible, whereas simulation approximates the target and must report its Monte Carlo precision.

A simulation estimate may appear stable across repeated runs, but visual stability alone is not a quantitative precision guarantee. Report at least \(K\), \(M\), \(\hat{\theta}_{MC}\), and an estimated Monte Carlo standard error or interval. When a decision boundary is relevant, compare the interval with that boundary and state whether the estimated simulation error is small enough for the intended interpretation.

Common Mistakes

  • Calling Monte Carlo error predictive uncertainty. The standard error quantifies finite-simulation variability in the estimated probability; it is not the variation in future clinical outcomes themselves.
  • Interpreting a Monte Carlo interval as a probability that the estimate is correct. A frequentist binomial interval has its stated coverage over repeated sets of simulations under the assumptions. It is not a posterior probability for \(\theta\).
  • Assuming zero simulated successes means zero predictive success. A finite simulation may fail to generate a rare success. Use a suitable binomial bound or run more replicates.
  • Assuming more simulations fix a wrong model or rule. They reduce numerical Monte Carlo error only when the simulated distribution and event definition are correct.
  • Using the independent-binomial formula for correlated replicates. Dependence changes the variance; the ordinary formula may understate error.
Key takeaway: For \(M\) independent predictive replicates with \(K\) successes, report \(K/M\) together with its Monte Carlo error. Increase \(M\) when that error is material to the intended decision, while separately checking that the predictive model, success rule, and independence assumptions are appropriate.

Check Your Understanding

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

  1. A simulation has \(M=1{,}000\) independent replicates and \(K=250\) successes. Calculate the estimated predictive success probability and its estimated Monte Carlo standard error.
  2. Why does quadrupling the number of independent replicates approximately halve the Monte Carlo standard error?
  3. If no successes occur in a finite simulation, what can and cannot be concluded about the predictive success probability?
  4. What assumption fails if the replicate indicators are correlated, and what is the consequence of using the ordinary binomial standard error anyway?
  5. Give one reason to run more simulations and one problem that more simulations cannot solve.

Solutions. (1) \(\hat{\theta}_{MC}=250/1{,}000=0.25\). The standard error is \(\sqrt{0.25(0.75)/1{,}000}=\sqrt{0.0001875}\approx0.01369\), or about 1.369 percentage points. (2) The standard error is proportional to \(1/\sqrt{M}\); replacing \(M\) by \(4M\) multiplies it by \(1/\sqrt{4}=1/2\). (3) One can conclude only that no success was observed in the simulated replicates. A finite simulation does not establish \(\theta=0\); a suitable upper bound can quantify probabilities still compatible with the result. (4) Independence fails. The binomial variance formula no longer applies, and the resulting standard error may be too small. (5) More simulations can reduce Monte Carlo error when an estimate is too imprecise for a decision. They cannot correct a misspecified predictive model, an incorrectly coded success rule, or a simulation that targets the wrong quantity.

Clinical Biostats