Bootstrap Distributions and Consistency
A confidence interval requires information about the sampling error of an estimator. Sometimes an exact pivot or a randomization argument supplies that information. More often the error distribution depends on the unknown population. The bootstrap estimates this distribution by applying the sampling procedure to an estimated population. The issue is then to justify replacing the population while the sample size is also changing.
First fix the sampling model and the parameter. Let \(X^n = (X_1, \dots, X_n)\) denote a sample of \(n\) i.i.d. random variables taking values in a sample space \(\mathcal S\), and they have an unknown probability distribution \(P\). The distribution \(P\) belongs to a class of distributions, denoted \(\mathcal P\). We are interested in conducting inference on some parameter \(\theta = \theta(P)\). Because of the dual interpretation of a confidence interval (region), let us focus on discussing confidence region alone, and note that this discussion applies to hypothesis tests. Let
denote the collection of all parameters \(\theta(P)\).
To construct a confidence region around \(\theta(P)\), we use a statistic that measures discrepancy between the data and a proposed value of \(\theta(P)\). Its sampling distribution calibrates which discrepancies are compatible with \(\theta(P)\). A real-valued functional that depends on both the data and \(\theta(P)\) is called a root, denoted \(R_n (X^n, \theta(P))\). And, under \(P\), we may always want to use the plug-in estimator for estimating \(\theta(P)\), denoted \(\hat \theta_n = \theta(\hat P_n)\), where
is the empirical probability measure, assuming \(\theta(\hat P_n)\) is well-defined. We will return to this momentarily.
To study the distribution of the root, let \(J_n (P)\) denote the distribution of \(R_n (X^n, \theta(P))\) under the distribution \(P\). Let
denote the cumulative distribution function associated with \(J_n(P)\). Note that we need to know the appropriate quantiles or the sampling distribution \(J_n(P)\), to construct a confidence region for \(\theta(P)\) based on the root \(R_n(X^n, \theta(P))\). We begin with the case in which these quantiles do not depend on the unknown population.
Suppose the distribution \(J_n (P)\) of \(R_n (X^n, \theta(P))\), under \(P\), does not depend on \(P\). Then, we refer to the root \(R_n (X^n, \theta(P))\) as a pivot. Since the distribution \(J_n (P)\) is independent of \(P\), we can construct a confidence region at level \(1-\alpha\) for \(\theta(P)\) by choosing constants \(a,b\) that satisfy
Therefore, the confidence region
contains \(\theta(P)\) with probability at least \(1-\alpha\) under \(P\).
Remark 1. In general, if \(R_n (X^n, \theta(P))\) is a pivot, its distribution may not be explicitly computable or have a known tractable form. However, since there is only one distribution that needs to be known, the problem is much simpler than if the distribution depends on \(P\).
Now, suppose the distribution \(J_n (P)\) of \(R_n (X^n, \theta(P))\), under \(P\), does depend on \(P\). This is usually the case, because we do not know the underlying distribution \(P\) that generates the data. In a more tractable case, when we know \(J_n (P)\) converges weakly to a limiting distribution that is independent of \(P\). Then, we call \(R_n (X^n, \theta(P))\) an asymptotic pivot.
In the less tractable but even more common case, \(J_n (P)\) converges weakly to a limiting distribution that does depend on \(P\). Suppose that this limiting distribution has a known form, which depends on \(P\), but only through some unknown parameter. An example is for estimating the nonparametric mean, where the root \(\sqrt{n} (\overline X_n - \theta(P))\) converges to \(\mathcal{N}(0, \sigma^2(P))\). Therefore, the limiting distribution of \(J_n (P)\) depends on \(P\) through its dependence on \(\sigma^2(P)\). In this common situation, we need to then do a normal approximation based on estimated variance \(\sigma^2 (\hat P_n)\)m which converges in probability to \(\sigma^2(P)\), and the quantiles of \(J_n (P)\) may need to be approximated, too.
1. Bootstrap
The bootstrap takes the same plug-in idea one step further: it estimates the entire sampling distribution of the root. The use of the bootstrap to approximate a null distribution in the construction of hypothesis tests is discussed later. This is useful when the limiting distribution is difficult to calculate. It can also provide a different approximation when the limiting law is known but does not describe the finite-sample distribution accurately. To tackle the problem of estimating \(J_n (P)\), we note that an important aspect of this problem is that \(J_n(P)\) depends on \(n\). So, we can proceed (heuristically) as follows: for each \(n\), approximate \(J_n (P)\) by \(J_n (\hat P_n)\), i.e. obtain a plug-in estimate of \(J_n (P)\). Consistency requires this plug-in distribution and the true sampling distribution to approach the same limit \(J(P)\) under the stated assumptions.
The bootstrap law is random because it depends on the observed sample. Conditional on that sample it is an ordinary probability distribution. For any distribution function \(H\) and \(0<u<1\), the generalized quantile
is well defined and finite. Neither continuity nor strict increase is required. Right continuity gives \(H(H^{-1}(u))\ge u\), while \(H(H^{-1}(u)-)\le u\). In particular, with \(q=J_n^{-1}(1-\alpha,P)\),
Equality holds if the cdf is continuous at \(q\); an atom there may make the inequality strict. Thus the known-law confidence region is \(\{\theta\in\Theta:R_n(X^n,\theta)\le q\}\). The bootstrap replaces the unknown population quantile by the same generalized inverse of the conditional bootstrap cdf:
This defines the bootstrap confidence region
Remark 2. (1) is usually difficult to compute exactly. The definition here leads to involving \(n\)-fold convolutions in computing this exactly. The usual way we go about constructing bootstrap confidence region is to rely on Monte Carlo approximations of \(J_n (P)\). When \(P\) is the empirical distribution, this is equivalent to resampling the original observations \(X^n\) with replacement. The bootstrap estimator \(J_n (\hat P_n)\) is approximated by the empirical distribution of the values \(R_n ({X_j^n}^*, \hat \theta_n)\), where \({X_j^n}^*\) is the sample of \(n\) i.i.d. observations drawn with replacement from \(\hat P_n\), and \(\hat \theta_n\) is the estimator for \(\theta(P)\). We sample with replacement for \(B\) bootstraps, i.e. \(j\) ranges from \(1\) to \(B\). This procedure can be seen as a plug-in estimator of a distribution function.
1.1. Asymptotic Behaviour
We hope the bootstrap estimator \(J_n (\hat P_n)\) is consistent for the true sampling distribution \(J_n (P)\) for \(R_n (X^n, \theta(P))\), but this is not clear from the discussion above. We need to show \(J_n (P)\) converges weakly to a limit law \(J(P)\) (that is nondegenerate). The bootstrap replaces \(P\) by \(\hat P_n\) in \(J_n (\cdot)\), so we need to study \(J_n (P_n)\) under more general sequences \(\{P_n\}\). In particular, because \(\hat P_n\) is random, we should aim to first understand the limiting behaviour of \(J_n (P_n)\) for fixed sequences \(\{P_n\}\). From there, a preliminary observation we can make is, because we need to study the limiting behaviour of \(J_n (P)\), \(J_n (P)\) should be required to be smooth in \(P\) (in particular, we will replace \(P\) by \(\hat P_n\), so the smoothness requirement is natural). Under this requirement, if we have sequences of probability distributions \(\{P_n\}\) converging to \(P\) in some suitable sense, then we have some hope of obtaining a consistency result for fixed sequences \(\{P_n\}\). Then, it becomes straightforward to study the limiting behaviour for random sequences.
Proposition 3. Let \(\mathbf S_P\) denote the set of sequences \(\{P_n \in \mathcal P\}\) such that it contains the sequence \(\{P, P, \dots\}\). Assume for every sequence \(\{P_n\} \in \mathbf S_P\), \(J_n (P_n)\) converges weakly to a continuous limit law \(J(P)\), with distribution function \(J(x,P)\). Let \(X^n\) be a sample of size \(n\) drawn from \(P\). Assume \(\hat P_n\) is an estimator of \(P\) based on \(X^n\), such that \(\{\hat P_1, \hat P_2, \dots \} \in \mathbf S_P\) with probability one. Then, \(J_n (x, \hat P_n)\) converges uniformly almost surely to \(J_n (x, P)\) over \(x\), i.e.
If \(J(\cdot, P)\) is continuous and strictly increasing at \(J^{-1} (1-\alpha, P)\), then
Finally, (1) is pointwise consistent (in level \(\alpha\)), i.e. under \(P\),
Proof. On every sample path whose estimated-distribution sequence belongs to \(\mathbf S_P\), the assumed sequential convergence applies to \(J_n(\hat P_n)\). Weak convergence to a continuous cdf is uniform: choose a finite grid with small limiting cdf increments and small tail probabilities, then use monotonicity between grid points. Apply this both to the estimated sequence and to the constant sequence \(P_n=P\), and use the triangle inequality to obtain (2).
For completeness, if \(H_n\Rightarrow H\) and \(H\) is continuous, fix \(\eta>0\). Choose \(a<b\) with \(H(a)<\eta\) and \(1-H(b)<\eta\), and a finite grid \(a=x_0<\cdots<x_m=b\) on which every increment of \(H\) is at most \(\eta\). Once \(\max_j|H_n(x_j)-H(x_j)|<\eta\), monotonicity gives \(|H_n(x)-H(x)|<2\eta\) between grid points and in both tails. This proves uniform convergence without assuming a density. Apply it separately on each of the probability-one sample paths in the hypothesis.
Write \(q=J^{-1}(1-\alpha,P)\). The quantile condition means, precisely, that for every \(\varepsilon>0\),
Uniform cdf convergence traps the bootstrap quantile between these two endpoints, proving (3). Finally, the root converges in law to \(J(P)\) and its random critical value converges in probability to \(q\). Sandwiching the coverage event between the events with fixed thresholds \(q-\varepsilon\) and \(q+\varepsilon\), then letting \(\varepsilon\downarrow0\), proves (4). No independence between the root and its estimated critical value is needed. ◻
The last sandwich has the explicit form, for \(Q_n=J_n^{-1}(1-\alpha,\hat P_n)\) and \(R_n=R_n(X^n,\theta(P))\),
Its lower and upper limits are \(J(q-\varepsilon,P)\) and \(J(q+\varepsilon,P)\). Continuity and the quantile inequalities imply \(J(q,P)=1-\alpha\), so both bounds converge to the same desired coverage. The strict quantile condition is needed for this part of Proposition 3, in addition to continuity of the limiting cdf.
In a parametric model, the estimated population remains in the model. This gives the parametric bootstrap, whose consistency can be related directly to local asymptotic normality. Suppose \(X^n\) is a sample from a quadratic mean differentiable model \(\{P_\theta : \theta \in \mathbf T\}\), where \(\mathbf T \subset \mathbb R^d\). Assume the information matrix is positive definite, the parameter is an interior point, and let \(\hat \theta_n\) be an efficient likelihood estimator. Let \(g : \mathbf T \to \mathbb R\) be a differentiable map with nonzero gradient vector \(\dot g(\theta)\). Consider the root
Let \(J_n (x, \theta)\) denote the distribution function of \(R_n (X^n, \theta)\). By standard results on limiting distribution of likelihood estimators, we have
To estimate the distribution of \(\sqrt{n} (g(\hat \theta_n) - g(\theta))\), we can use the normal approximation \(\mathcal{N}(0, \sigma_\theta^2)\), provided we have access to a consistent estimator of \(\sigma_\theta^2\). For now, suppose we know the forms of \(\dot g(\cdot)\) and \(I(\cdot)\); and assume these maps are continuous at the true parameter, then a consistent estimator of \(\sigma_\theta^2\) is
If we are able to consistently estimate \(\sigma_\theta^2\), then we have the following result.
The local condition used below is stronger than a pointwise normal approximation: for every deterministic sequence \(h_n\to h\in\mathbb R^d\), with \(\theta_n=\theta+h_n/\sqrt n\) in the model, assume
Also assume \(\sqrt n(\hat\theta_n-\theta)=O_{P_\theta}(1)\). These properties hold for regular efficient estimators in the LAN setting. They state the precise local sampling-law input needed for the parametric bootstrap; knowledge of the variance formula alone is insufficient.
Proposition 4. In the parametric bootstrap setup above, under \(\theta\), \(J_n (x, \hat \theta_n)\) converges uniformly in probability to \(J(x,\theta)\), i.e.
and
Proof. Let \(d_n(h)=\sup_x|J_n(x,\theta+h/\sqrt n)-J(x,\theta)|\). The sequential local assumption and continuity of the normal cdf imply \(d_n(h_n)\to0\) whenever \(h_n\to h\). Consequently \(\sup_{\|h\|\le M}d_n(h)\to0\) for every fixed \(M\), restricting the supremum to parameters in the model: otherwise a violating sequence in this compact ball would have a convergent subsequence. Tightness of \(\sqrt n(\hat\theta_n-\theta)\) now gives \(d_n(\sqrt n(\hat\theta_n-\theta))\to0\) in probability. Taking \(h=0\) and applying the triangle inequality proves the second assertion. This argument works in every finite dimension. ◻
To spell out the compactness step in Proposition 4, a failure of uniform convergence would give \(\varepsilon>0\), indices \(n_j\), and \(\|h_j\|\le M\) with \(d_{n_j}(h_j)>\varepsilon\). Extract \(h_j\to h\) and extend these values to a full sequence by setting the remaining \(h_n=h\). The sequential assumption contradicts the displayed lower bound. With \(H_n=\sqrt n(\hat\theta_n-\theta)\),
First take \(n\to\infty\), then \(M\to\infty\). This is why tightness of the local estimator suffices and almost-sure convergence of \(H_n\) is unnecessary.
Gathering the discussion above, we have that a one-sided bootstrap lower confidence bound for \(g(\theta)\) is of the form \(g(\hat \theta_n) - \frac{1}{\sqrt{n}} J_n^{-1} (1-\alpha, \hat \theta_n)\). The quantile argument in Proposition 3, applied in probability to Proposition 4, implies that, under \(\theta\), we have \(J_n^{-1} (1-\alpha, \hat \theta_n) \xrightarrow{P} \sigma_\theta z_{1-\alpha}\). We are now ready to consider a testing problem of the following form: suppose we want to test \(H_0 : g(\theta)=0\) versus \(H_1 : g(\theta)>0\). Then, the ongoing discussion suggests (after appealing to the fact that confidence regions can be obtained from test inversions) that we can reject the null if \(\sqrt{n} g(\hat \theta_n) > J_n^{-1} (1-\alpha, \hat \theta_n)\). Note that under \(g(\theta)=0\), \(\sqrt{n} (g(\hat \theta_n) - g(\theta)) \xrightarrow{d} \mathcal{N}(0, \sigma_\theta^2)\) by Slutsky’s theorem. Therefore, this test is asymptotically level \(\alpha\). To characterize the limiting power of this test against a contiguous sequence of alternatives, we have the following result.
Proposition 5. In the parametric bootstrap setup above, with \(\theta\) satisfying \(g(\theta)=0\), the limiting power of the test that rejects when \(\sqrt{n} g(\hat \theta_n) > J_n^{-1} (1-\alpha, \hat \theta_n)\) against the sequence \(\theta_n = \theta + \frac{h}{\sqrt{n}}\) is characterized by
Proof. This follows from a result from the theory of efficient likelihood estimators. First, note that
Then, note that \(J_n^{-1} (1-\alpha, \hat \theta_n) \xrightarrow{P_\theta^n} \sigma_\theta z_{1-\alpha}\). Contiguity implies that under \(P_{\theta_n}^n\), \(J_n^{-1} (1-\alpha, \hat \theta_n) \xrightarrow{P_{\theta_n}^n} \sigma_\theta z_{1-\alpha}\). Therefore, because \(g(\theta)=0\) and \(g\) is differentiable, we have \(\sqrt{n} g(\theta_n) \to \langle \dot g(\theta)^{\top}, h \rangle\). Appealing to the limiting distribution results of efficient likelihood estimators, we have that \(\sqrt{n} (g(\hat \theta_n) - g(\theta_n)) \xrightarrow{d} \mathcal{N}(0,\sigma_\theta^2)\), hence, by Slutsky’s theorem,
which is the desired result. ◻
Remark 6. The bootstrap can produce an asymptotically optimal test, in the sense that the limiting power is optimal, as above. However, this can only happen if the initial estimator or test statistic is optimally chosen. If the root is based on a suboptimal estimator, then accurately approximating its sampling distribution does not remove that loss of efficiency.