Multiplier Bootstrap for High-Dimensional Maxima

Why Resample a Maximum?

In high-dimensional inference, a statistic may report the largest standardized deviation across many coordinates. For example, a researcher might ask whether any of $p$ treatment effects differs substantially from zero. The uncertainty of the largest deviation depends not only on the uncertainty of each coordinate but also on their dependence. A multiplier bootstrap approximates that joint behavior by attaching random weights to the observed, centered data vectors.

The method is particularly useful when the statistic is nonsmooth: the maximum can switch from one coordinate to another after a small change in the data. Rather than resampling individual coordinates independently, the multiplier method applies one random weight to each whole observation vector. That shared weight preserves the observed cross-coordinate structure.

Target and Bootstrap Statistics

Definition: Let $X_1,\ldots,X_n$ be independent random vectors in $\mathbb{R}^p$, with $E[X_i]=0$ and common covariance matrix $\Sigma$. Write $\sigma_j^2=\operatorname{Var}(X_{ij})>0$. The maximum of the standardized coordinates of the normalized sum is $$ T_n=\max_{1\le j\le p}\left|\frac{1}{\sqrt{n}}\sum_{i=1}^n\frac{X_{ij}}{\sigma_j}\right|. $$ Given observations, let $\bar X_j=n^{-1}\sum_{i=1}^n X_{ij}$ and let $\hat\sigma_j>0$ be an estimate of $\sigma_j$. For independent multipliers $\xi_1,\ldots,\xi_n$, independent of the data, with mean zero and variance one, define $$ T_n^*=\max_{1\le j\le p}\left|\frac{1}{\sqrt{n}}\sum_{i=1}^n \xi_i\frac{X_{ij}-\bar X_j}{\hat\sigma_j}\right|. $$ The multiplier-bootstrap distribution is the conditional distribution of $T_n^*$ given the observations.

Two common choices are Gaussian multipliers, $\xi_i\sim N(0,1)$, and Rademacher multipliers, for which $P(\xi_i=1)=P(\xi_i=-1)=1/2$. Conditional on the data, Gaussian multipliers make the vector inside the maximum exactly Gaussian. Rademacher multipliers instead give a sum of randomly signed observed vectors; that conditional distribution is not generally Gaussian in a finite sample.

Let $U_{ij}=X_{ij}/\sigma_j$, and let $R$ be the correlation matrix with entries $R_{jk}=E[U_{ij}U_{ik}]$. The Gaussian reference statistic is $T_G=\max_j|G_j|$, where $G\sim N(0,R)$. Its covariance structure records the dependence relevant to the maximum.

Worked Example: Two Independent Coordinates

Suppose $n=4$, $p=2$, and the observed vectors are $(1,0),(0,1),(-1,0),(0,-1)$. Their sample mean is $(0,0)$. Using divisor $n$, the sample variances are $$ \hat\sigma_1^2=\frac{1^2+0^2+(-1)^2+0^2}{4}=\frac12,\qquad \hat\sigma_2^2=\frac{0^2+1^2+0^2+(-1)^2}{4}=\frac12. $$ The sample covariance is zero, so the sample correlation matrix is the $2$-by-$2$ identity matrix.

For Gaussian multipliers, the two bootstrap coordinates are conditionally independent $N(0,1)$ variables after standardization. Thus, for $t\ge0$, $$ P(T_n^*\le t\mid X_1,\ldots,X_4)=\{2\Phi(t)-1\}^2, $$ where $\Phi$ is the standard normal distribution function. The conditional $0.95$ quantile solves $\{2\Phi(t)-1\}^2=0.95$, so $\Phi(t)=(1+\sqrt{0.95})/2$, giving $t$ approximately $2.24$. This calculation illustrates how a joint maximum distribution, rather than two separate coordinate intervals, determines the cutoff.

Worked Example: A Correlation the Bootstrap Must Preserve

Consider four centered observations with first coordinate $a=(1,-1,1,-1)$. Define $b=(1,1,-1,-1)$ and set the second coordinate to $y=0.5a+\sqrt{0.75}\,b$. Both $a$ and $b$ have sample variance $1$ using divisor $4$, and their sample covariance is zero: $$ \frac{1}{4}\sum_{i=1}^4 a_i b_i=0. $$ Consequently, the second coordinate has sample variance $0.5^2+0.75=1$, and its covariance with the first coordinate is $0.5$. The sample correlation matrix is therefore $$ \hat R=\begin{pmatrix}1&0.5\\0.5&1\end{pmatrix}. $$

The shared multiplier $\xi_i$ is attached to the entire $i$th pair, so the conditional covariance of the standardized multiplier sum is $\hat R$. Replacing the two coordinates by independently resampled versions would instead impose zero covariance. That would change the distribution of the maximum and generally its quantiles. The observed correlation is part of the uncertainty calculation, not a nuisance that can be discarded.

A Validity Argument

The central reasoning is a comparison through Gaussian distributions. It separates two tasks: approximating the sampling distribution by a Gaussian maximum, and estimating the covariance that determines that Gaussian maximum. The following result makes the required approximations explicit.

Theorem: Consider a sequence of problems in which $n$ and $p=p_n$ may grow. Let $G\sim N(0,R)$, and let $\hat R$ be the sample correlation matrix formed from the centered observations and positive sample standard deviations. Let $\hat G$, conditionally on the observations, have distribution $N(0,\hat R)$. Suppose $$ \delta_n=\sup_t|P(T_n\le t)-P(\max_j|G_j|\le t)|\longrightarrow0, $$ and, for the chosen multipliers, $$ \delta_n^*=\sup_t|P(T_n^*\le t\mid X_1,\ldots,X_n)-P(\max_j|\hat G_j|\le t)|\xrightarrow{p}0. $$ If $$ \eta_n=\sup_t|P(\max_j|G_j|\le t)-P(\max_j|\hat G_j|\le t)|\xrightarrow{p}0, $$ then $$ \sup_t|P(T_n^*\le t\mid X_1,\ldots,X_n)-P(T_n\le t)|\xrightarrow{p}0. $$

Proof. Fix any $t$. Add and subtract the two Gaussian distribution functions. The triangle inequality gives $$ \begin{aligned} |P(T_n^*\le t\mid X_1,\ldots,X_n)-P(T_n\le t)| &\le |P(T_n^*\le t\mid X_1,\ldots,X_n)-P(\max_j|\hat G_j|\le t)|\\ &\quad+|P(\max_j|\hat G_j|\le t)-P(\max_j|G_j|\le t)|\\ &\quad+|P(\max_j|G_j|\le t)-P(T_n\le t)|. \end{aligned} $$ Taking the supremum over $t$ bounds the left-hand side by $\delta_n^*+\eta_n+\delta_n$. The first two terms tend to zero in probability, and the deterministic term $\delta_n$ tends to zero. Their sum therefore tends to zero in probability, proving the claim. $\square$

For Gaussian multipliers, the conditional multiplier sum is Gaussian exactly, with covariance $\hat R$; hence the theorem's $\delta_n^*$ is zero when the bootstrap uses the same coordinate standardization as $\hat R$. For Rademacher multipliers, $\delta_n^*$ requires a conditional high-dimensional normal approximation. This is an additional condition, not a finite-sample identity.

Covariance Accuracy and Dimension Growth

The covariance comparison condition $\eta_n\to0$ says precisely what covariance estimation must accomplish for the maximum. A useful sufficient condition is $$ \Delta_n\{\log(2p)\}^2\longrightarrow0\quad\text{in probability},\qquad \Delta_n=\max_{j,k}|\hat R_{jk}-R_{jk}|. $$ For correlation matrices with unit diagonal, a Gaussian comparison bound gives, up to a universal constant, a difference between the two maximum distribution functions bounded by $$ \{\Delta_n\log^2(2p/\Delta_n)\}^{1/3} $$ when $\Delta_n$ is small. Thus the displayed condition makes the Gaussian maximum distributions close. The logarithmic factor matters: an entrywise covariance error that is small at fixed $p$ may be insufficient when $p$ increases.

One familiar sufficient regime makes the connection to sample size more concrete. Suppose the standardized coordinates are uniformly bounded in absolute value by $B_n$, their variances are bounded away from zero, and the observations are independent. Standard high-dimensional Gaussian approximation and multiplier central limit results for bounded vectors give sufficient growth conditions of the form $$ \frac{B_n^2\{\log(pn)\}^7}{n}\longrightarrow0. $$ In the same setting, concentration of sample covariances typically gives $$ \Delta_n=O_p\left(B_n^2\sqrt{\frac{\log p}{n}}\right) $$ (up to lower-order centering and variance-estimation terms). For bounded $B_n$, this rate implies the covariance comparison condition if $(\log p)^5/n\to0$. These are conservative sufficient conditions, not necessary thresholds. With heavier tails, dependence across observations, or inaccurate variance estimates, different assumptions and methods are needed.

Worked Example: Checking the Covariance-Error Requirement

Suppose $p=100$, the true variances are one, and the largest absolute entrywise error in the estimated correlation matrix is $\Delta_n=0.01$. The covariance comparison scale is $$ \Delta_n\{\log(2p)\}^2 =0.01\{\log(200)\}^2 \approx0.01(5.30)^2 \approx0.281. $$ This number is not close to zero, so this calculation alone does not justify treating the estimated-covariance Gaussian maximum as indistinguishable from the target Gaussian maximum. If, at a larger sample size, the error were $0.001$ with the same $p$, the scale would be approximately $0.0281$, which is smaller. Asymptotic validity requires the scale to tend to zero along the sequence of problems; one small numerical value is not a proof of that limit.

Rademacher Multipliers and Finite-Sample Discreteness

The two multiplier choices agree in mean and variance, but their conditional distributions differ. Gaussian multipliers make the conditional law a Gaussian maximum once the covariance estimate is fixed. Rademacher multipliers produce a discrete distribution for every finite data set. Under suitable moment, growth, and covariance conditions, that discreteness becomes negligible for the distributional approximation of interest.

Worked Example: A Discrete Rademacher Bootstrap

Take one coordinate and two centered observations $X_1=1$, $X_2=-1$, with $\hat\sigma=1$. The Rademacher bootstrap statistic is $$ \left|\frac{\xi_1-\xi_2}{\sqrt{2}}\right|. $$ The four equally likely multiplier pairs are $(1,1),(1,-1),(-1,1),(-1,-1)$. The statistic is zero for two pairs and $\sqrt{2}$ for two pairs, so its conditional distribution assigns probability $1/2$ to each value. With Gaussian multipliers, $(\xi_1-\xi_2)/\sqrt{2}\sim N(0,1)$, and the absolute value has a continuous distribution. The difference is real at this sample size; asymptotic approximation requires conditions that control it as $n$ and possibly $p$ grow.

Common Pitfalls and Practical Use

A multiplier maximum can provide simultaneous critical values or approximate a familywise tail probability, but its validity depends on matching the target's dependence and scale. In applications, the data vectors should correspond to the independent sampling units. If observations are clustered or dependent, the multiplier construction must operate at an appropriate unit or use a dependence-aware method; treating dependent rows as independent does not meet the independent-vector conditions above.

Using uncentered observations. Multipliers have mean zero, but that does not make it safe to substitute $X_{ij}$ for $X_{ij}-\bar X_j$. The conditional covariance of the uncentered weighted sum estimates a second moment about zero, not the centered covariance required for the sampling fluctuation.
Resampling each coordinate separately. Independent coordinate-wise weights erase the empirical cross-coordinate dependence. Use the same $\xi_i$ across all coordinates of observation $i$.
Assuming Gaussian multipliers fix every problem. They remove the conditional central-limit approximation step, but do not fix a poor covariance estimate, incorrect standardization, or invalid independence assumptions.
Confusing marginal and simultaneous cutoffs. A $95\%$ cutoff for each coordinate separately is not automatically a $95\%$ cutoff for the maximum. The maximum's joint distribution is the relevant one.

The practical purpose of the method is to estimate the distribution of the largest standardized fluctuation while retaining estimated dependence. Its high-dimensional validity is a joint statement: the normalized sum must be approximable by a Gaussian vector, the multiplier law must approximate its Gaussian counterpart, and the estimated covariance must be accurate enough at the scale induced by $p$.

Check Your Understanding

Use the definitions and validity argument to answer these questions.

Why does a single multiplier $\xi_i$ need to be shared across all coordinates of observation $i$?
What is the conditional distribution of the multiplier vector when the multipliers are Gaussian?
State the three approximation errors controlled in the validity theorem.
Why can a small entrywise covariance error still matter when the number of coordinates grows?
What finite-sample feature distinguishes the Rademacher multiplier distribution from the Gaussian multiplier distribution?