From a Resampling Result to a Reproducible Analysis
A bootstrap interval or permutation \(p\)-value is not determined by a data set and a method name alone. The result also depends on how observations are sampled or reassigned, how many replicates are generated, how randomness is managed, and how numerical failures are handled. If those choices are not recorded, a later run may produce a different result with no clear way to determine why.
A reproducible resampling pipeline makes these choices explicit. It identifies the analysis inputs, assigns each replicate a stable identity, records the random-number configuration and software environment, and retains enough diagnostics to distinguish Monte Carlo variation from numerical error. Reproducibility does not mean that every computer or software version must produce identical floating-point values. It means that the computation is specified and auditable, and that any differences can be investigated.
Specify the Computation Before Running It
For a bootstrap, state what constitutes one sampling unit and how each replicate is drawn. If observations are clustered, paired, stratified, or otherwise dependent, the resampling mechanism must preserve the relevant design rather than silently sampling individual rows. For a permutation procedure, state which labels or assignments are exchangeable under the null hypothesis and which features of the design must remain fixed. A reproducible seed cannot repair a resampling scheme that does not match the study design.
Also define the statistic precisely, including any preprocessing, missing-data rule, optimization criterion, and treatment of undefined results. Identify the data by an immutable file or version and record the order of observations when that order can affect computation. A useful run record includes the code revision, software and random-number generator versions, master seed or stream key, replicate range, worker configuration, and numerical settings. These details support exact repetition within a fixed environment; bit-for-bit agreement across different libraries or hardware is not guaranteed.
Worked Example: A Named Bootstrap Replicate
Suppose the observed values are \(2,4,7\), and the statistic is the sample mean. A nonparametric bootstrap replicate draws three observations with replacement from these values. Imagine replicate 17 draws \(2,2,7\). Its recorded result is
A useful record for this replicate identifies the input data version, bootstrap scheme, statistic, replicate ID 17, and the sampled indices or equivalent random draws. Storing the sampled values or indices makes this particular result directly auditable; storing only the final mean does not reveal whether the replicate was generated as intended. In a large analysis, the full draw record may be costly, so the pipeline should at least retain enough metadata and deterministic stream information to regenerate it.
Random Streams That Survive Parallel Execution
A single seed can reproduce a serial run if the generator, software, draw order, and environment are unchanged. Parallel execution complicates that guarantee: task scheduling can change which worker consumes which random values. Giving every worker the same seed is worse, because workers may generate identical or strongly overlapping sequences. Giving workers different integer seeds is not, by itself, proof that their streams are independent.
A robust design assigns randomness to replicate identity rather than to worker identity. For example, a documented counter-based generator can generate the random values for replicate \(b\) from a master key and the index \(b\). Another valid design uses a generator that provides documented independent substreams, assigning one substream to each replicate or fixed task. The required property is that a replicate's random draws do not depend on which worker happens to execute it. Store the generator family and version, master key or seed, stream-allocation rule, and replicate IDs. If the method consumes a variable number of random draws, use replicate-local streams or another documented scheme that prevents one replicate's draw count from shifting later replicates.
Parallel-safe assignment also makes output assembly predictable: save each result with its replicate ID, then sort or join by that ID rather than relying on completion order. With replicate-indexed streams, fixed inputs, and deterministic replicate code, changing the worker count changes scheduling but not the random input assigned to any replicate.
Worked Example: Parallelizing 1,000 Replicates
Suppose a bootstrap requires \(B=1000\) replicates and four workers are available. Assign IDs 1 through 1000 before launching any work. One possible schedule gives worker 1 IDs \(1,5,9,\ldots\), worker 2 IDs \(2,6,10,\ldots\), and so on. A different schedule may give each worker a contiguous block. Either schedule is reproducible if each ID uses its own documented random input, such as \(R_b=g(K,b)\), and results are stored with \(b\).
For instance, if replicate 17 finishes on worker 3 in the first schedule and worker 1 in the second, both runs use \(g(K,17)\) for that replicate. Thus its resample and statistic agree, assuming the replicate calculation itself is deterministic in the fixed environment. By contrast, if each worker simply starts from the same seed, the schedules can repeat random sequences across workers and do not implement 1,000 safely separated replicates.
Convergence Logs and Monte Carlo Precision
A convergence log records what the computation has produced at prespecified checkpoints, not just its final answer. Depending on the target, record the running estimate, the number of completed and failed replicates, and a Monte Carlo standard error (MCSE). For an interval endpoint, a quantile, or a multivariate summary, state the estimator and any method used to assess its simulation error. A running estimate changing little is not, by itself, evidence that it is accurate.
For example, if \(B\) independent resampling replicates estimate a tail probability \(p\), let \(I_b\) equal 1 when replicate \(b\) meets the tail event and 0 otherwise. The Monte Carlo estimate is \(\hat p=B^{-1}\sum_{b=1}^{B}I_b\). If the replicate indicators are independent Bernoulli variables with common probability \(p\), the MCSE is estimated by \(\sqrt{\hat p(1-\hat p)/B}\).
Worked Example: Reporting Tail-Probability Precision
In \(B=2000\) independently generated resampling replicates, suppose \(K=74\) statistics are at least as extreme as the observed statistic. The estimated tail probability is
Its estimated Monte Carlo standard error is
A normal-approximation summary of simulation precision is \(0.037\pm1.96(0.00422)\), or approximately \((0.0287,0.0453)\). This interval describes Monte Carlo error under the independent Bernoulli approximation; it does not replace a confidence interval for a population parameter. The report should state the replicate count and tail-count rule, since both are needed to interpret the estimate.
Choose the replicate budget and any stopping rule before examining the final result. Fixed budgets are straightforward to reproduce. If computation stops adaptively when an estimate appears stable, the stopping criterion and its statistical justification must also be specified; repeatedly checking an ordinary fixed-sample error formula does not automatically make an adaptive procedure valid. Keep checkpoint logs so that a later reviewer can see how the estimate and its uncertainty changed as replicates accumulated.
Numerical Diagnostics Are Not Monte Carlo Diagnostics
A resampling statistic may require an eigendecomposition, a matrix solve, or an optimization. Its numerical error concerns the accuracy of the computation for a given replicate; Monte Carlo error concerns variation across the random replicates. More resampling cannot fix an unstable solve, and a small numerical residual does not imply that the Monte Carlo estimate is precise.
For each replicate, define explicit rules for convergence and failure. Record optimizer status, iteration limits, tolerances, matrix dimensions, rank decisions where relevant, and warnings or exceptions. For matrix computations, useful diagnostics can include a residual relative to the scale of the input and a condition estimate. State thresholds in advance where possible, and report how many replicates were omitted, retried, or treated as failures. Silently dropping difficult replicates can change the distribution being estimated.
Worked Example: A Rank Threshold in a Bootstrap Replicate
Suppose a replicate produces the diagonal matrix \(A=\operatorname{diag}(1,10^{-10})\), and the recorded numerical-rank rule counts singular values greater than \(10^{-8}\) times the largest singular value. The singular values are \(1\) and \(10^{-10}\); the threshold is
Because \(1>10^{-8}\) but \(10^{-10}<10^{-8}\), this rule reports numerical rank 1. The diagnostic record should include the threshold, singular values, and rank outcome, not just the resulting statistic. If another replicate has a small singular value just above the threshold, tiny floating-point changes could alter the rank and downstream result. That is a numerical sensitivity to flag and investigate, not evidence that the bootstrap itself has converged.
Common Pitfall: Treating a Seed as a Full Specification
The reasoning “I recorded the seed, so the analysis is reproducible” is incomplete. A seed does not identify the data version, resampling unit, generator implementation, stream assignment, replicate count, software, numerical tolerances, or rules for failed fits. Even with all of those fixed, changing software or hardware can change floating-point results. The correction is to record the complete pipeline specification and distinguish exact reruns in a fixed environment from reproducibility of the scientific conclusion across environments.
A concise report should state the resampling design and statistic, number of replicates, random-number generator and stream strategy, software environment, key numerical settings, failures and their treatment, and Monte Carlo precision for the reported quantities. Preserve the code and run log alongside the result. These practices make it possible to reproduce the computation, diagnose discrepancies, and assess whether simulation or numerical limitations matter for the conclusion.
Check Your Understanding
Use these questions to check whether a resampling result is both reproducible and interpretable.
- Why can assigning the same seed to every worker fail to provide safe parallel random streams?
- What information should be recorded so that replicate 17 can be regenerated independently of which worker ran it?
- For \(B\) independent Bernoulli tail indicators, derive the variance of their sample mean and state what uncertainty that variance measures.
- Why is a small Monte Carlo standard error not evidence that an ill-conditioned matrix calculation is numerically stable?
- A pipeline drops every replicate whose optimizer does not converge. What should be documented, and why can silent dropping change the result?