Summarizing Simulation Uncertainty for Design Selection

Why a Point Estimate Is Not Enough

A simulation can estimate how a MAMS design behaves under a specified scenario, but the estimate itself varies across simulation runs. For example, if a design rejects a null hypothesis in 760 of 1,000 simulated trials, the estimated rejection probability is 0.760. A different set of 1,000 simulated trials would usually produce a somewhat different estimate. Reporting only 0.760 hides that Monte Carlo variability.

This tutorial concerns uncertainty caused by using a finite number of simulation replicates. It is distinct from uncertainty about the true treatment effect, nuisance parameters, or which scenario is most relevant. Those are represented by the scenario choices and assumptions supplied to the simulation. A precise estimate within one scenario does not establish that the scenario is realistic, and averaging over scenarios requires a separately justified weighting scheme.

Definition: For a fixed design, simulation scenario, and analysis rule, a Monte Carlo estimate is a summary calculated from a finite number \(R\) of simulated trials. Monte Carlo uncertainty is the variability in that estimate across repeated sets of \(R\) simulation replicates, conditional on the specified design, scenario, and correctly implemented simulation procedure.

Uncertainty for a Simulated Operating Characteristic

Many operating characteristics are probabilities: familywise type I error, power for a specified treatment effect, or the probability of stopping at an interim analysis. For one such outcome, define \(X_r=1\) if replicate \(r\) has the event of interest and \(X_r=0\) otherwise. For example, \(X_r=1\) could mean that at least one experimental arm is declared effective in replicate \(r\). If replicates are independently generated under the same scenario and the event has probability \(p\), then \(X_r\) are independent Bernoulli random variables with probability \(p\).

The estimate is the observed proportion of replicates with the event. Under these assumptions, its Monte Carlo standard error estimates the standard deviation of the estimate over repeated simulation runs.

\hat{p}=\frac{1}{R}\sum_{r=1}^{R}X_r, \qquad \widehat{SE}_{MC}(\hat{p}) =\sqrt{\frac{\hat{p}(1-\hat{p})}{R}}.

For sufficiently large \(R\), a commonly used approximate 95% Monte Carlo interval is \(\hat{p}\pm1.96\widehat{SE}_{MC}(\hat{p})\). It describes the precision of the simulated estimate for this fixed scenario; it is not a confidence interval for an unknown treatment effect and does not include uncertainty about scenario specification. Because a probability lies between 0 and 1, the approximate interval can be unsatisfactory near either boundary. The Wilson interval or an exact binomial interval is generally preferable when the estimated probability is near 0 or 1, or when the event count is small. In particular, a Wald interval based on the displayed standard error collapses to a point when there are zero simulated events, despite substantial uncertainty about how small the probability may be.

Interpretation: A Monte Carlo interval describes the precision expected from repeating the simulation with the same number of independent replicates under the same conditions. It does not assign a probability to the fixed value \(p\), and it does not account for model misspecification, coding errors, or uncertainty about which scenario applies.

Worked Example: Monte Carlo Precision for Estimated Power

Suppose a hypothetical MAMS design is simulated under a specified alternative scenario with a treatment effect of 0.30 log-odds units. In \(R=2{,}500\) independently generated trials, the prespecified power event occurs 1,900 times. Assume the simulation correctly implements the design and analysis, and that the replicates are independent under this fixed scenario.

The estimated power is \(1{,}900/2{,}500=0.760\). Its estimated Monte Carlo standard error is

\widehat{SE}_{MC} =\sqrt{\frac{0.760(1-0.760)}{2{,}500}} =\sqrt{0.00007296} =0.008542.

The approximate 95% Monte Carlo interval is

0.760\pm1.96(0.008542) =0.760\pm0.01674 =(0.7433,\;0.7767).

Thus the estimated power is 76.0%, with a Monte Carlo margin of error of about 1.67 percentage points by this approximation. The interval concerns repeated simulation estimates under the stated effect scenario; it does not imply that the design’s power under every plausible effect is between 74.3% and 77.7%.

Rare Events and Estimates at the Boundary

When an operating characteristic is a rare error probability, the number of simulated events matters as much as the total number of replicates. A zero count is not evidence that the event is impossible. The ordinary Wald calculation is especially unsuitable here: if \(\hat{p}=0\), it reports a standard error of zero.

For zero events in \(R\) independent replicates, an exact one-sided 95% binomial upper confidence bound \(u\) solves \((1-u)^R=0.05\). This follows because, if the true event probability were \(u\), the probability of seeing zero events in \(R\) trials would be \(0.05\). Solving gives \(u=1-0.05^{1/R}\). This is a bound on the event probability based on the binomial model, not proof that it is below the bound with certainty.

Worked Example: Zero Simulated Type I Errors

Under a hypothetical global-null scenario, suppose a design produces zero familywise type I error events in \(R=1{,}000\) independent simulation replicates. The estimated error rate is 0/1,000, or 0%, but its Monte Carlo uncertainty is not zero. The exact one-sided 95% upper bound is

u=1-0.05^{1/1{,}000} =1-\exp\left(\frac{\ln(0.05)}{1{,}000}\right) \approx 1-\exp(-0.002996) \approx0.002991.

The simulation is therefore consistent, at this one-sided confidence level, with a familywise error probability as high as about 0.299% under this scenario. The result does not establish that the error rate is zero. If a design requirement is, for example, a probability below 0.5%, these simulations provide stronger evidence of adequate precision than if the requirement is below 0.1%; the latter threshold is not supported by this bound.

Comparing Designs Without Hiding Simulation Error

A design choice should be informed by the scenario and the relevant trade-offs, not just by a ranking of point estimates. For each design, report the operating characteristic, its Monte Carlo uncertainty, the number of replicates, and the scenario that generated it. Across a scenario grid, present results scenario by scenario or state a defensible way of combining them. Do not treat one scenario’s precision as evidence that the design performs similarly in other scenarios.

When candidate designs are simulated using common random numbers, outcomes from the same replicate index are paired. For a probability difference, let \(X_{Ar}\) and \(X_{Br}\) be the event indicators for designs A and B in paired replicate \(r\), and let \(D_r=X_{Ar}-X_{Br}\). The estimated difference is the mean of the paired differences. Its Monte Carlo standard error is the sample standard deviation of the \(D_r\) values divided by \(\sqrt{R}\). This accounts for the pairing; calculating two separate standard errors and treating the estimates as independent discards the covariance induced by common random numbers.

Worked Example: Uncertainty in a Paired Difference in Power

In a hypothetical treatment-effect scenario, two candidate designs are simulated using the same \(R=2{,}000\) replicate-level random inputs. Design A has a power event in 1,500 replicates, and design B has one in 1,460. Both designs have the event in 1,400 paired replicates. Thus A has the event but B does not in 100 replicates; B has it but A does not in 60. The remaining 1,840 pairs have the same event status.

The estimated difference in power, A minus B, is \((1{,}500-1{,}460)/2{,}000=0.0200\), or 2.00 percentage points. Each \(D_r\) is 1 for the 100 A-only events, \(-1\) for the 60 B-only events, and 0 otherwise. Therefore the sample mean of \(D_r^2\) is \(160/2{,}000=0.0800\), and the sample variance is

s_D^2 =\frac{2{,}000}{1{,}999}\left(0.0800-0.0200^2\right) =0.07964.

The paired Monte Carlo standard error and approximate 95% interval for the difference are

\widehat{SE}_{MC}(\bar{D}) =\sqrt{\frac{0.07964}{2{,}000}} =0.006310, \qquad 0.0200\pm1.96(0.006310) =(0.00763,\;0.03237).

Under this scenario, A’s estimated power is 2.00 percentage points higher than B’s, and the approximate paired interval for the difference ranges from 0.763 to 3.24 percentage points. This quantifies Monte Carlo uncertainty in the comparison, not uncertainty about the scenario’s clinical relevance. The paired approach is appropriate here because the stated simulation uses aligned replicate pairs.

Choosing the Number of Replicates

A simulation size should be chosen in relation to the precision needed for the decision. For a probability estimate, the approximate standard error is largest at \(p=0.5\). If a planning value for \(p\) is available and the desired approximate two-sided margin of error is \(h\), the normal approximation gives

R\approx\frac{z_{0.975}^2p(1-p)}{h^2},

where \(z_{0.975}\approx1.96\). This is a planning calculation, not a guarantee: the achieved standard error depends on the event rate, and rare-event requirements may call for methods and replicate counts that do not rely on a normal approximation.

Worked Example: Planning Replicates for a Target Margin

Suppose a hypothetical design is expected to have power near 0.80, and the goal is an approximate 95% Monte Carlo margin of error no greater than 0.010, or one percentage point, for that scenario. Substituting \(p=0.80\) and \(h=0.010\) gives

R\approx \frac{1.96^2(0.80)(0.20)}{0.010^2} =\frac{0.614656}{0.0001} =6{,}146.56.

Rounding up gives a planning total of at least 6,147 replicates. If the event probability is unknown and a conservative approximation is wanted, use \(p(1-p)\leq0.25\); the corresponding calculation is \(R\approx1.96^2(0.25)/0.010^2=9{,}604\). These totals address precision for one probability estimate under the specified simulation conditions; several scenarios or a rare error-rate target may need separate planning.

Diagnostics, Comparisons, and Limitations

The formulas rely on a well-defined event or outcome, valid replicates, and independent replicates for a single-design binomial estimate. For paired comparisons, the pairs must be correctly aligned and independent across replicate indices. Check that the number of valid replicates is the denominator used in each estimate, and that exclusions or invalid replicates are handled transparently. A Monte Carlo interval cannot repair an incorrect simulation or a scenario that does not answer the design question.

Useful diagnostics include checking event counts, comparing estimates across independent batches, and tracking whether Monte Carlo standard errors are small relative to decision thresholds. Large changes across batches can indicate too few replicates, unusual simulation behavior, or a problem requiring investigation. If the output is an expected sample size rather than a probability, estimate it by the replicate mean and use the sample standard deviation divided by \(\sqrt{R}\) as its Monte Carlo standard error, when the replicates are independent and the variance is finite. The precision should be reported in the same units as the quantity, such as participants per trial.

Approach	Useful for	Important limitation
Point estimate only	A quick descriptive summary	Does not show whether apparent differences may be simulation noise
Separate intervals for each design	Summarizing each design’s operating characteristic	Does not directly quantify a paired difference; overlap is not a test of equality
Replicate-level paired difference	Comparing designs simulated with common random numbers	Requires correctly matched pairs and a standard error based on their differences

A narrow interval is not, by itself, a reason to select a design: the estimated operating characteristic must also be acceptable under the relevant scenarios, and practical trade-offs matter. Conversely, overlapping separate intervals do not establish that two designs perform equally. If a choice depends on a small difference relative to its Monte Carlo error, increase the simulation effort or describe the ranking as uncertain rather than treating the larger point estimate as decisive. The general principle from the earlier material on comparing candidate designs under common random numbers is to preserve the paired structure when estimating a difference.

Common pitfall: “Design A has a higher estimated power than design B, so A is better.” This ignores Monte Carlo uncertainty, the scenario in which power was estimated, and other design criteria. Report scenario-specific estimates and uncertainty, quantify the difference using the correct paired or unpaired calculation, and make the selection in light of the prespecified performance requirements and trade-offs.

Check Your Understanding

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

A hypothetical design has 420 power events in 600 independent replicates. What is \(\hat{p}\), and what is the approximate Monte Carlo standard error?
Why does zero simulated familywise error events not establish that the true simulated-scenario error probability is zero?
Two designs are simulated with common random numbers. Why is the standard error of their estimated difference based on replicate-level differences rather than just the separate estimates?
In a single scenario, what does a 95% Monte Carlo interval describe, and what sources of uncertainty does it not include?
Using \(p=0.5\), estimate the number of replicates needed for an approximate 95% margin of error of 0.02 for one probability estimate.

Solutions. (1) \(\hat{p}=420/600=0.700\). The standard error is \(\sqrt{0.700(0.300)/600}=\sqrt{0.000350}=0.01871\), or about 1.87 percentage points. (2) With finite \(R\), zero events can occur even when the event probability is positive; for example, the exact one-sided 95% upper bound is \(1-0.05^{1/R}\), which is greater than zero. (3) Pairing induces dependence within each replicate pair. The variance of the difference depends on that pairing, so the standard error must use the observed differences \(D_r\), not an independence calculation that discards their covariance. (4) It describes the repeated-simulation precision of the estimate under the fixed scenario and simulation procedure. It does not include uncertainty about the clinical effect, the scenario assumptions, or the correctness of the code. (5) \(R\approx1.96^2(0.5)(0.5)/0.02^2=0.9604/0.0004=2{,}401\) replicates.