set.seed(2450)
n1 <- 12
x1 <- round(rlnorm(n1, log(20), 0.5), 1)
B <- 4000
boot_mean <- replicate(B, mean(sample(x1, n1, replace = TRUE)))
boot_med <- replicate(B, median(sample(x1, n1, replace = TRUE)))
d_mean <- length(unique(round(boot_mean, 6))) # distinct resampled means
d_med <- length(unique(boot_med)) # distinct resampled medians
atom_med <- max(table(boot_med)) / B # largest single spikeChecking a bootstrap: when resampling fails
Corrected 14 October 2026: the coverage section said the bias correction helped for the skewed ratio in the first post, as it does here; in that post’s simulation the bias-corrected interval covered less often than the percentile interval, and the text now says so.
Corrected 7 October 2026: the check in the maximum section asked whether a small change in one observation moves the statistic smoothly, which the maximum passes, and the introduction, a heading, the take-away and the description still blamed statistics that are not smooth; they now point to statistics set by an edge of the data or of the parameter space. The median section also said the resampled medians sit on order statistics; most sit on midpoints between pairs of them, and the text and caption now say so.
Corrected 30 September 2026: the maximum section said the failure extends to any statistic that is not smooth in the data; the median is not smooth and its bootstrap still works, so the text now names statistics set by an edge of the data or parameter space.
The bootstrap is easy to run and easy to over-trust. Across this series it kept undercovering a little, even when it worked; and there are cases where it does not work at all, quietly returning a confident answer that is wrong. This closing post collects three ways resampling fails: when the sample is too small for the statistic to resolve, when the statistic is set by an edge of the data, and when the interval simply does not cover at the rate it claims. Each comes with the check that would have caught it.
Small samples cannot resolve a rough statistic
With few observations the bootstrap distribution is built from a limited set of resamples, and how well it fills in depends on how smooth the statistic is. Take twelve values and bootstrap two summaries: the mean and the median.
From the same twelve numbers the bootstrap mean takes 1016 distinct values, forming a smooth distribution, while the bootstrap median takes only 41, landing on a coarse lattice of order statistics and midpoints between pairs of them, with a single value carrying 9% of the resamples. The median is not smooth in the data, so at this sample size its bootstrap distribution is a handful of spikes, and any quantile read from it, including the ends of a percentile interval, can only sit on that lattice. The check is direct: look at the bootstrap distribution before trusting its tails. If it is a few tall bars rather than a filled-in shape, the interval is coarser than its decimals suggest, and you need more data or a smoother statistic.
A statistic set by an edge breaks it outright
Small samples only sharpen a deeper problem: for some statistics the bootstrap is inconsistent, and no amount of data repairs it. The textbook case is the maximum of a sample from a uniform distribution, the natural estimator when the parameter is an upper bound. We draw forty values from a uniform on the interval up to one, and estimate that upper bound by the sample maximum.
set.seed(2451)
n2 <- 40; theta <- 1
x2 <- runif(n2, 0, theta)
mx <- max(x2)
boot_max <- replicate(B, max(sample(x2, n2, replace = TRUE)))
atom <- mean(boot_max == mx) # share of resamples pinned at the max
atom_theory <- 1 - (1 - 1/n2)^n2 # exact for n = 40; the limit is 1 - 1/e = 0.632
ci_boot <- quantile(boot_max, c(.025, .975)) # bootstrap percentile interval
ci_exact <- c(mx, mx / (0.05)^(1/n2)) # exact 95% interval from max/theta ~ Beta(n,1)The observed maximum is 0.973. A resample includes it unless every draw misses it, so a fraction 0.65 of bootstrap maxima are exactly the observed value, matching the 0.64 that this sample size implies, a shade above the large-sample limit of 1 - 1/e. The bootstrap distribution is a huge spike at the maximum with a short tail below, and it never exceeds the observed maximum, because resampling cannot invent a larger value. So the bootstrap interval [0.888, 0.973] sits entirely at or below the sample maximum, while the truth is above it. The exact interval, from the fact that the maximum divided by the parameter follows a known distribution, runs [0.973, 1.049] and reaches past the data as it must. Its whole 5% sits in the upper tail: the sample maximum is below the true bound with probability one, so a lower limit below it would buy nothing.
covmax <- function(R = 2000, Bb = 500) {
hb <- 0; he <- 0
for (r in seq_len(R)) {
xx <- runif(n2, 0, 1); m <- max(xx)
bb <- replicate(Bb, max(sample(xx, n2, TRUE)))
hb <- hb + (1 >= quantile(bb, .025) & 1 <= quantile(bb, .975))
ce <- c(m, m / (0.05)^(1/n2)); he <- he + (1 >= ce[1] & 1 <= ce[2])
}
c(bootstrap = hb / R, exact = he / R)
}
set.seed(22)
cmax <- covmax()Repeated over many samples, the bootstrap interval covers the true bound 0% of the time, against 95% for the exact one. Not merely poor: a flat failure. The lesson extends past the maximum to other statistics set by an edge of the data or of the parameter space, including minima, ranges, and parameters tested at the edge of their space, the boundary case behind the parametric bootstrap. Before bootstrapping, ask whether the statistic is set by such an edge, as the maximum is. If it is, the bootstrap may be inconsistent, and simulating from a model, or an exact result, is the safer route.
Coverage is a promise, so check it
Even for a smooth statistic at a reasonable sample size, a bootstrap interval labelled 95% rarely covers exactly 95%. We saw it throughout this series. Here it is head on: the variance of twenty normal observations, whose sampling distribution is right-skewed. We compare the plain percentile interval with the bias-corrected and accelerated interval from the first post.
nn <- 20; sig2 <- 1; zc <- qnorm(c(.025, .975))
one_var <- function(v) mean((v - mean(v))^2)
cov_pb <- function(R = 1500, Bb = 600) {
hp <- 0; hb <- 0
for (r in seq_len(R)) {
xx <- rnorm(nn, 0, 1); th <- one_var(xx)
bv <- replicate(Bb, one_var(sample(xx, nn, TRUE)))
hp <- hp + (sig2 >= quantile(bv, .025) & sig2 <= quantile(bv, .975))
z0 <- qnorm(mean(bv < th))
jk <- sapply(seq_len(nn), function(i) one_var(xx[-i])); jm <- mean(jk)
a <- sum((jm - jk)^3) / (6 * (sum((jm - jk)^2))^1.5)
ad <- pnorm(z0 + (z0 + zc) / (1 - a * (z0 + zc)))
cb <- quantile(bv, ad); hb <- hb + (sig2 >= cb[1] & sig2 <= cb[2])
}
c(percentile = hp / R, bca = hb / R)
}
set.seed(50)
cvar <- cov_pb()The percentile interval covers 82% of the time and the bias-corrected interval 88%, both short of 95%. The correction helps here, but it is not a guarantee: for the skewed ratio in the first post, the bias-corrected interval covered less often than the percentile interval in that post’s simulation. At this sample size neither reaches the promised level. The point is not which interval wins; it is that the nominal level is optimistic and you cannot know by how much without checking. When you can simulate, as we have throughout, do it and read the actual coverage. When you cannot, prefer the bias-corrected interval, report the resampling scheme and the number of resamples, and treat 95% as a hopeful label. For an automatic correction, the double bootstrap re-calibrates the nominal level using an inner resample, at a squared cost in computation.
What to take away
Resampling is a tool, not a guarantee, and this series has tried to show both its reach and its edges. It rescues honest uncertainty from awkward statistics like a ratio, it builds null distributions that asymptotics get wrong, and it adapts to dependence through blocks and clusters. But it needs a statistic that is not set by an edge of the data or of the parameter space, enough data to resolve it, and the right resampling unit; and even then its intervals lean a little narrow. The habit that ties the series together is to check rather than assume: inspect the bootstrap distribution, ask whether the statistic is set by such an edge and whether the observations are independent, and verify coverage by simulation whenever you can. A bootstrap you have checked is a strong instrument. A bootstrap you have merely run is a guess with decimals.
References
Bickel and Freedman 1981 Annals of Statistics 9(6):1196-1217 (10.1214/aos/1176345637)
Andrews 2000 Econometrica 68(2):399-405 (10.1111/1468-0262.00114)
Canty, Davison, Hinkley and Ventura 2006 Canadian Journal of Statistics 34(1):5-27 (10.1002/cjs.5550340103)
Puth, Neuhauser and Ruxton 2015 Journal of Animal Ecology 84(4):892-897 (10.1111/1365-2656.12382)
Davison and Hinkley 1997 Bootstrap Methods and their Application. ISBN 978-0-521-57391-7