Monte Carlo Error in Bootstrap Confidence Limits

Why the Number of Bootstrap Replicates Matters

A bootstrap confidence limit is usually calculated from a percentile of bootstrap parameter estimates. A bootstrap standard error is calculated from their sample standard deviation. Both are computed from a finite number of refits, so they vary slightly if the bootstrap is repeated with a different set of replicates. That variation is Monte Carlo error: simulation noise in the numerical approximation, not uncertainty in the original parameter estimate itself.

For example, increasing the number of bootstrap replicates can make an estimated 2.5th percentile more stable. It does not make the original study larger, remove uncertainty in its parameter estimate, or correct a misspecified model. It improves the precision with which the chosen bootstrap procedure is calculated.

Let $\hat{\theta}^{*(b)}$ be the estimate from bootstrap replicate $b$, for $b=1,\ldots,B$, where the replicates are independently generated and fitted using the same bootstrap procedure. The bootstrap sample standard deviation is

s_{\mathrm{boot}} = \sqrt{\frac{1}{B-1}\sum_{b=1}^{B} \left(\hat{\theta}^{*(b)}-\bar{\theta}^{*}\right)^2}, \qquad \bar{\theta}^{*}=\frac{1}{B}\sum_{b=1}^{B}\hat{\theta}^{*(b)}.

Here $B$ is the number of valid, independent bootstrap estimates used in the calculation. The value $s_{\mathrm{boot}}$ estimates the bootstrap standard error; it is itself subject to Monte Carlo error. Likewise, an interval endpoint calculated as a bootstrap quantile is an estimate of the corresponding quantile of the bootstrap distribution.

Definition: Monte Carlo error of a bootstrap summary is the variability in that summary across repeated finite sets of bootstrap replicates, when the underlying data and bootstrap procedure are held fixed. It is distinct from the sampling uncertainty that the bootstrap is intended to describe.

Monte Carlo Error in a Bootstrap Standard Error

If the bootstrap estimates are independent draws from a distribution with finite variance and a reasonably well-behaved fourth moment, the Monte Carlo standard error of their sample standard deviation can be approximated using the distribution's kurtosis. Under an approximately normal bootstrap distribution, a commonly used approximation is

\operatorname{MCSE}(s_{\mathrm{boot}}) \approx \frac{s_{\mathrm{boot}}}{\sqrt{2(B-1)}}.

The approximation is useful as a quick precision check. It is not a statement that the bootstrap estimates themselves must be normally distributed for the bootstrap procedure to be valid. Rather, it uses approximate normality to describe the sampling variation of the calculated standard deviation across repeated simulation runs. Strong skewness or heavy tails can make this approximation inaccurate.

Justification: If the bootstrap estimates are normally distributed with standard deviation $\sigma$, their sample variance $s_{\mathrm{boot}}^2$ has variance $2\sigma^4/(B-1)$. Applying the delta method to the square-root function, whose derivative at $\sigma^2$ is $1/(2\sigma)$, gives $$ \operatorname{Var}(s_{\mathrm{boot}}) \approx \left(\frac{1}{2\sigma}\right)^2\frac{2\sigma^4}{B-1} =\frac{\sigma^2}{2(B-1)}. $$ Replacing $\sigma$ with $s_{\mathrm{boot}}$ and taking the square root gives the stated approximation.

Worked Example: Monte Carlo Error in a Bootstrap Standard Error

Suppose 500 hypothetical subject-level bootstrap refits produce a bootstrap standard error for clearance of $12.0\ \mathrm{L/h}$. Assume the replicate estimates are independent and their distribution is not markedly nonnormal. Estimate the Monte Carlo standard error of this bootstrap standard error.

Substitute $s_{\mathrm{boot}}=12.0\ \mathrm{L/h}$ and $B=500$:

\operatorname{MCSE}(s_{\mathrm{boot}}) \approx \frac{12.0}{\sqrt{2(500-1)}} =\frac{12.0}{\sqrt{998}} =0.380\ \mathrm{L/h}.

The estimated bootstrap standard error is $12.0\ \mathrm{L/h}$, with an estimated Monte Carlo standard error of about $0.380\ \mathrm{L/h}$, or $0.380/12.0=3.17\%$ of the estimate. This describes the approximate run-to-run simulation variation in the reported standard error; it is not an additional $0.380\ \mathrm{L/h}$ of clearance uncertainty.

If the same calculation used $B=2{,}000$ replicates and gave a similar $s_{\mathrm{boot}}$, the Monte Carlo standard error would be approximately $12.0/\sqrt{2(1{,}999)}=0.190\ \mathrm{L/h}$. Quadrupling $B$ roughly halves this error because the error decreases at a rate proportional to $1/\sqrt{B}$.

For a nonnormal bootstrap distribution, the approximation depends on more than $s_{\mathrm{boot}}$ and $B$. A highly skewed distribution or a few very large estimates can make the standard deviation especially variable. In such cases, compare estimates from independent batches of replicates or use a suitable simulation-based assessment as a sensitivity check. Batch comparisons are diagnostic, not a substitute for checking whether all replicates represent the intended bootstrap procedure.

Monte Carlo Error in a Quantile Endpoint

For a percentile endpoint, the main issue is how precisely a finite set of bootstrap estimates locates the desired quantile. Let $q_p$ denote the $p$th quantile of the bootstrap distribution, and let $f(q_p)$ be its density at that point. For independent, identically distributed bootstrap estimates and a smooth distribution with positive density near $q_p$, the large-sample approximation is

\operatorname{MCSE}(\hat q_p) \approx \frac{\sqrt{p(1-p)/B}}{f(q_p)}.

This formula shows why extreme endpoints can be hard to stabilize. The numerator reflects uncertainty in the empirical fraction of estimates below the quantile. The density in the denominator converts uncertainty in percentile rank into uncertainty in the endpoint's units. If the distribution is sparse near the target quantile, a small change in rank can correspond to a large change in the parameter value.

The density is not usually known. It can be estimated from nearby empirical quantiles or with a local density estimate, but the resulting endpoint MCSE is itself approximate. For example, if $q_{0.01}$ and $q_{0.04}$ are available, a local density estimate near $q_{0.025}$ is

\hat f(q_{0.025}) \approx \frac{0.04-0.01}{\hat q_{0.04}-\hat q_{0.01}}.

This estimate is only useful if those nearby quantiles are themselves reasonably well resolved. Quantile endpoints based on adjusted tail probabilities, as in BCa intervals, are subject to the same finite-replicate quantile error at the adjusted probability. Uncertainty in estimated adjustment quantities can add further error.

Worked Example: Assessing a Lower Percentile Endpoint

In a hypothetical analysis, $B=1{,}600$ bootstrap refits give a 95% percentile interval for clearance. The estimated lower endpoint is $\hat q_{0.025}=4.20\ \mathrm{L/h}$. Nearby quantiles are $\hat q_{0.01}=3.80\ \mathrm{L/h}$ and $\hat q_{0.04}=5.80\ \mathrm{L/h}$. Estimate the endpoint's Monte Carlo standard error.

First estimate the local density:

\hat f(q_{0.025}) \approx \frac{0.04-0.01}{5.80-3.80} =\frac{0.03}{2.00} =0.0150\ (\mathrm{L/h})^{-1}.

Then substitute $p=0.025$, $B=1{,}600$, and this density estimate:

\operatorname{MCSE}(\hat q_{0.025}) \approx \frac{\sqrt{0.025(0.975)/1{,}600}}{0.0150} =\frac{0.00390}{0.0150} =0.260\ \mathrm{L/h}.

The estimated lower endpoint is $4.20\ \mathrm{L/h}$, but its approximate Monte Carlo standard error is $0.260\ \mathrm{L/h}$. Under the approximation, a roughly 95% range for the endpoint's Monte Carlo variation is $4.20\pm1.96(0.260)$, or $3.69$ to $4.71\ \mathrm{L/h}$. This is a range for simulation variation in the endpoint estimate, not a confidence interval for the true clearance parameter.

The calculation depends on a local density estimate based on only the 1st-to-4th percentile span. It should therefore be treated as a precision diagnostic, not as an exact uncertainty bound. A sparse or irregular lower tail can make the estimate unstable.

Planning More Refits and Checking Stability

Once a tolerable Monte Carlo error is chosen, the approximate $1/\sqrt{B}$ rate can guide the number of additional refits. If the current endpoint MCSE is $m_{\mathrm{current}}$ at $B_{\mathrm{current}}$, and the target MCSE is $m_{\mathrm{target}}$, use

B_{\mathrm{target}} \approx B_{\mathrm{current}} \left(\frac{m_{\mathrm{current}}}{m_{\mathrm{target}}}\right)^2.

This calculation assumes the same bootstrap distribution and a similarly reliable density estimate as $B$ increases. It provides a planning approximation; it does not guarantee that the endpoint will move smoothly or that the chosen precision target is scientifically adequate.

Worked Example: Deciding Whether to Add Bootstrap Refits

Suppose a hypothetical analysis has $B_{\mathrm{current}}=800$ independent bootstrap refits. A local density calculation gives an estimated MCSE of $0.20\ \mathrm{L/h}$ for the 2.5th-percentile endpoint. The analysis team wants to reduce this endpoint MCSE to $0.10\ \mathrm{L/h}$. Estimate the required total number of refits.

Using the scaling formula:

B_{\mathrm{target}} \approx 800\left(\frac{0.20}{0.10}\right)^2 =800(2)^2 =3{,}200.

The approximate plan is therefore 3,200 total refits, or $3{,}200-800=2{,}400$ additional refits. This targets an estimated MCSE of $0.10\ \mathrm{L/h}$; it does not promise that a particular endpoint will change by less than $0.10\ \mathrm{L/h}$. After adding replicates, the endpoint and its MCSE should be recalculated, since the estimated density and endpoint may have changed.

A practical stability check is to calculate the endpoint at prespecified cumulative replicate counts, such as $B=1{,}000$, $2{,}000$, and $4{,}000$, and record its movement alongside an MCSE estimate. If independent batches are used, compare their endpoint estimates as well. Large batch-to-batch differences or a density estimate that changes substantially with the neighborhood used are warnings that the tail is not yet well resolved. A stable-looking endpoint alone is not proof of adequate precision: endpoints can move irregularly, and two estimates can happen to be close.

Assumptions, Diagnostics, and Common Pitfalls

The formulas above require the bootstrap estimates being summarized to be independent draws from the same intended bootstrap procedure. They also rely on enough replicates for the relevant approximation to be useful. For quantile MCSE, the distribution should have a positive, reasonably smooth density near the target quantile. For the standard-error approximation, a normal-like distribution is a useful condition; strong skewness or heavy tails can require more care.

Check the effective replicate count. Use the number of independent valid estimates actually included, not the number of attempted refits. Excluding failed fits must follow a defensible procedure; silently discarding failures can alter the distribution being summarized.
Check the tail resolution. A 2.5th percentile lies near the lower tail. With $B=400$, only about $0.025(400)=10$ estimates are expected below that percentile. This makes the endpoint sensitive to individual order statistics.
Check local density sensitivity. If modest changes to the quantile neighborhood produce very different density estimates, the endpoint MCSE is not reliable enough to justify a precise claim.
Check for dependent replicates. Reusing or correlating simulated replicate data can reduce the effective information relative to $B$, making formulas based on independent replicates too optimistic.
Check the reported units and scale. An endpoint MCSE for clearance is in clearance units; the MCSE for its logarithm is on the log scale. Transformations change the scale and interpretation of endpoint variation.

Common mistake: “The bootstrap standard error is $12\ \mathrm{L/h}$, so the uncertainty in the standard error is also $12\ \mathrm{L/h}$.” This confuses the estimated bootstrap standard error with the Monte Carlo error of that estimate. In the 500-refit example, the corresponding approximate Monte Carlo standard error was only $0.380\ \mathrm{L/h}$. Another mistake is to treat a more precise Monte Carlo calculation as evidence that the model or bootstrap design is correct. Increasing $B$ reduces simulation noise; it does not repair model misspecification, an inappropriate resampling unit, or biased handling of failed fits.

The principal alternatives are complementary rather than interchangeable. The analytic standard-error approximation is quick but relies on distributional assumptions. Independent-batch comparisons reveal instability without requiring a density estimate, but can be noisy and depend on how batches are formed. Local-density quantile MCSEs provide an interpretable endpoint-scale measure but may be unreliable in sparse tails. Increasing $B$ is the direct way to reduce Monte Carlo error, though NLME refits can be computationally expensive. There is no universal replicate count that guarantees stable limits for every parameter and tail shape.

Key takeaway: Evaluate bootstrap summaries at the precision that matters for the analysis. Estimate standard-error MC error with a replicate-count approximation, assess endpoint error using tail rank and local density, and add refits when the resulting Monte Carlo uncertainty is too large for the intended interpretation.

Check Your Understanding

Answer each question, then compare with the solutions below.

A bootstrap standard error is $8.0\ \mathrm{mg/L}$ from $B=799$ independent, approximately normal bootstrap estimates. Estimate its Monte Carlo standard error.
Why can a 2.5th-percentile endpoint have greater Monte Carlo error than a median calculated from the same number of replicates?
An endpoint MCSE is estimated as $0.30\ \mathrm{L/h}$ at $B=500$. Approximately how many total refits are needed to target an MCSE of $0.15\ \mathrm{L/h}$?
Does doubling $B$ halve the Monte Carlo error? State the approximate relationship and one condition needed for it to apply.

Solutions. (1) $8.0/\sqrt{2(799-1)}=8.0/\sqrt{1{,}596}=0.200\ \mathrm{mg/L}$, approximately. (2) The tail may have lower density than the center, so a similar uncertainty in percentile rank translates to a larger change in endpoint value; there are also fewer observations in the tail. (3) $B_{\mathrm{target}}\approx500(0.30/0.15)^2=2{,}000$ total refits. (4) No. Under the usual independent-replicate approximation, error scales as $1/\sqrt{B}$; doubling $B$ multiplies error by $1/\sqrt{2}$, about $0.707$. The replicates must be independent draws from the same procedure, and the relevant approximation must be reasonable.