Why the Selection Rule Must Be Simulated
In a multi-arm multi-stage trial, the final estimate is affected not only by sampling variation but also by the rule that determines which arms continue and which arm is reported. An arm can be selected partly because its interim estimate was unusually high. The final estimate may therefore be biased upward relative to that arm’s treatment effect, even when the stagewise estimates are unbiased before selection.
A parametric bootstrap makes this mechanism explicit. It simulates stagewise outcomes from a fitted probability model, applies the trial’s actual interim rules to each simulated data set, and records the final estimator for the arm selected in that replicate. The key is to reproduce the whole path, not merely to generate a new estimate for the arm selected in the observed trial.
Model and Bootstrap Procedure
For clarity, consider \(J\) treatment arms compared with a control, with \(K\) possible stages. Let \(Y_k=(Y_{1k},\ldots,Y_{Jk})^\mathsf{T}\) be the vector of stage-\(k\) treatment-effect estimates, and let \(\theta=(\theta_1,\ldots,\theta_J)^\mathsf{T}\) be the vector of true effects. A standard working model is
where \(V_k\) is the covariance matrix for the stage-\(k\) estimates. The model assumes that stage increments are independent across stages, conditional on \(\theta\), and that their covariance matrices are correctly specified. Estimates for different treatment arms at the same stage need not be independent: a shared control group, for example, can induce covariance. If stages have different information sizes, the stage estimates and the final estimator must use the trial’s prespecified information weights rather than an unweighted average.
Fit the model to the observed data to obtain \(\hat\theta\) and any required covariance or nuisance-parameter estimates. In each replicate, generate data under that fitted model, run the exact interim decision rules—including futility, efficacy, arm ranking, and any sample-size changes—and calculate the final estimate for the selected arm. The same rules and analysis estimator used by the trial must be used in the bootstrap.
Here \(B\) is the number of bootstrap replicates, \(S_b^*\) is the arm selected in replicate \(b\), and \(\hat\theta^{*(b)}_{S_b^*}\) is its simulated final estimate. The reference value \(\hat\theta_{S_b^*}\) is the fitted effect for the arm selected in that replicate. Thus, the errors are measured against the corresponding arm’s fitted effect, not against the observed winner’s effect in every replicate.
Why the Replicate Average Estimates Plug-In Bias
Conditional on the fitted model, define the replicate error as \(D^*=\hat\theta^*_{S^*}-\hat\theta_{S^*}\). If the bootstrap draws are independent and use the same simulation and selection procedure, each replicate error has the same plug-in distribution. Therefore,
The sample average is consequently a Monte Carlo estimate of the expected selection error under the fitted model. Its Monte Carlo standard error is estimated by \(s_D/\sqrt{B}\), where \(s_D\) is the sample standard deviation of the replicate errors. This quantifies simulation noise only; it does not account for uncertainty in the fitted model or guarantee that the model represents the trial correctly.
Subtracting the estimated bias from the observed estimate is a first-order bias correction. It does not automatically produce an unbiased estimate, a confidence interval, or a correction conditional on the observed selection path. The ordinary procedure above targets an average over simulated selection paths under the plug-in model. A conditional target requires a procedure that represents the relevant conditioning event.
Worked Example: A Two-Stage Selection Rule
Suppose two arms, A and B, have independent stage estimates with unit variance, and the fitted effects are \(\hat\theta_A=0.60\) and \(\hat\theta_B=0.40\), in outcome units. At stage 1, discard any arm with an estimate at or below 0. The arm with the larger stage-1 estimate among those remaining continues; if neither is above 0, the trial stops. The continuing arm’s final estimate is the equal-weight average of its two stage estimates. Assume the illustrative bootstrap replicates below all continue to a final estimate.
| Replicate | Stage 1 A | Stage 1 B | Stage 2 for continuing arm |
|---|---|---|---|
| 1 | 0.90 | 0.70 | A: 0.50 |
| 2 | 0.20 | 0.80 | B: 0.60 |
| 3 | -0.10 | 0.30 | B: 0.10 |
| 4 | 0.60 | 0.50 | A: 0.40 |
In replicate 1, both arms pass the interim threshold and A has the larger interim estimate. Its final estimate is \((0.90+0.50)/2=0.70\), so the error is \(0.70-0.60=0.10\). In replicate 2, B is selected; its estimate is \((0.80+0.60)/2=0.70\), and its error is \(0.70-0.40=0.30\). In replicate 3, A is discarded and B continues; the estimate is \((0.30+0.10)/2=0.20\), giving error \(0.20-0.40=-0.20\). In replicate 4, A is selected, giving estimate \((0.60+0.40)/2=0.50\) and error \(0.50-0.60=-0.10\).
The four-replicate estimate is \( \hat b_{\mathrm{PB}}=(0.10+0.30-0.20-0.10)/4=0.025 \) outcome units. If the observed selected-arm estimate were 0.80 outcome units, its illustrative corrected estimate would be \(0.80-0.025=0.775\) outcome units. Four replicates are far too few for a stable practical correction; the calculation demonstrates the procedure, not a recommended simulation size.
Representing Dependence and Trial Operations
The covariance model matters because selection compares arms. For two arms whose estimates share a control, simulating independent arm estimates when their true covariance is positive changes the chance that one appears much larger than the other. It can therefore alter both the selected-arm frequencies and the estimated bias.
Worked Example: Simulating Correlated Stage Estimates
Suppose the fitted effects are \(\hat\theta_A=0.50\) and \(\hat\theta_B=0.30\), in outcome units, and the stage covariance matrix is
A valid simulation uses independent standard normal draws \(U_1,U_2\) and sets \(Y_A^*=0.50+U_1\) and \(Y_B^*=0.30+0.50U_1+\sqrt{0.75}U_2\). The variances are both 1, while \( \operatorname{Cov}(Y_A^*,Y_B^*)=\operatorname{Cov}(U_1,0.50U_1+\sqrt{0.75}U_2)=0.50 \). Thus this construction reproduces the stated covariance.
For one illustrative draw, let \(U_1=0.40\) and \(U_2=-0.20\). Then \(Y_A^*=0.90\), and \(Y_B^*=0.30+0.50(0.40)+\sqrt{0.75}(-0.20)\approx0.327\). If the interim rule retains the larger estimate when both exceed 0, this replicate selects A. The remaining stages must then be simulated for the arms allowed by the actual design, and the full rule must be applied again. Generating correlated estimates without replaying those later decisions is not a simulation of the MAMS selection process.
In an actual trial, the final estimator may be an information-weighted combination, and the decision rule may include multiple interim looks, stopping for efficacy, or varying follow-up after an arm is dropped. These details belong in the bootstrap code. If an arm stops early, its simulated final estimator must use the same data and stopping convention as the trial’s reported estimator. An arm that is not observed after dropping cannot be treated as though later outcomes had been collected unless the fitted model and target explicitly justify that construction.
Selection Frequencies and a Frozen-Arm Comparison
A useful check is to inspect how often each arm is selected and how often each stopping path occurs. The correction is an average across those paths, so a coding error that freezes the observed winner changes the quantity being estimated. The following calculation uses the same illustrative draws as the first example to show the distinction.
Worked Example: Replaying Selection Versus Freezing the Winner
In the four replicates above, the selection rule chooses A in replicates 1 and 4 and B in replicates 2 and 3. Replaying the rule gives errors \(0.10,0.30,-0.20,-0.10\), with average \(0.025\) outcome units.
Suppose instead that A is the observed winner and the bootstrap code always analyzes A, ignoring the replicate-specific selection rule. Using the illustrative A stage estimates from the table and hypothetical stage-2 A estimates \(0.50,0.40,0.30,0.40\), the fixed-A final estimates are \(0.70,0.30,0.10,0.50\). Relative to \(\hat\theta_A=0.60\), their errors are \(0.10,-0.30,-0.50,-0.10\), averaging \(-0.20\) outcome units.
The two averages differ because they answer different questions. The first estimates error after the selection rule has operated; the second estimates error for A under a simulation that does not select among arms. This small constructed comparison is not evidence that either correction has a particular direction in general. It demonstrates why the arm-selection step cannot be omitted.
Diagnostics, Alternatives, and Limitations
Before interpreting a correction, verify the model and the simulated paths. At a minimum, check that stage estimates have the intended means and covariance, the code reproduces the prespecified decision rules, and simulated selection and stopping frequencies are plausible. Increase \(B\) until Monte Carlo error is small relative to the correction being reported. Also examine sensitivity to reasonable changes in fitted effects and covariance assumptions: the plug-in procedure treats those estimates as fixed, so model uncertainty can be consequential.
A parametric bootstrap is attractive when the stagewise distribution and covariance structure are credibly specified and the complete design can be simulated. Compared with a nonparametric bootstrap, it can generate outcomes for a structured multi-stage design without requiring resampling units in a way that may fail to represent stage increments or shared-control dependence. But it relies more directly on its parametric assumptions. If the distribution, covariance, or selection rule is misspecified, a large number of simulations will estimate the wrong plug-in quantity precisely.
The estimate is also sensitive to how nuisance parameters are fitted. Fitting effects only from data available after selection can make the plug-in model reflect the selection event itself. State the fitting procedure and its data inputs, and use a likelihood or other fitting method appropriate to the observed stages and stopping process. This bootstrap is a point-estimation tool; it does not by itself quantify the sampling uncertainty of the corrected estimate. Nor should an unconditional correction be described as removing bias conditional on a particular observed path.
State the stagewise distribution, covariance structure, parameter-fitting method, and any information weights used by the final estimator.
Generate outcomes under fitted parameters and apply the actual interim rules, including stopping, dropping, and selection.
Average replicate-specific selected-arm errors, subtract that average from the observed selected estimate, and report selection-path diagnostics and Monte Carlo error.
Check Your Understanding
Answer each question, then compare with the solutions.
- Why must the bootstrap select an arm anew in each replicate rather than always analyze the observed winner?
- For replicate errors \(0.12,-0.04,0.08\) outcome units, what is the bootstrap bias estimate?
- If two stage estimates have covariance 0.50, what feature of the simulation must reflect that covariance?
- Does the ordinary plug-in bias estimate target bias conditional on the observed selection path?
- Name one diagnostic that addresses simulation noise and one that addresses model or implementation problems.
Solutions. (1) The selected arm is random; freezing the observed winner omits the selection process and estimates a different quantity. (2) The average is \((0.12-0.04+0.08)/3=0.0533\) outcome units. (3) The joint simulation must generate correlated stage estimates, using the specified covariance rather than independent draws. (4) No. It estimates an average over paths under the fitted model; a path-conditional target requires appropriate conditioning. (5) Increasing \(B\) and checking the Monte Carlo standard error addresses simulation noise. Checking simulated covariance, selection frequencies, and rule implementation addresses model and coding problems.