We propose a method for constructing a confidence region for the solution to a
conditional moment equation. The method is built around a class of algorithms
for nonparametric regression based on subsampled kernels. This class includes
random forest regression. We bound the error in the confidence region’s
nominal coverage probability, under the restriction that the conditional
moment equation of interest satisfies a local orthogonality condition. The
method is applicable to the construction of confidence regions for conditional
average treatment effects in randomized experiments, among many other similar
problems encountered in applied economics and causal inference. As a by-
product, we obtain several new order-explicit results on the concentration and
normal approximation of high-dimensional $U$-statistics.

Keywords: Random Forest Regression, Half-Sample Bootstrap, $U$-statistics
JEL: C01, C14, C12

##  1\. Introduction

Consider an independent and identically distributed sample
$\mathbf{D}_{n}=(D_{i})_{i=1}^{n}$, where each observation $D_{i}$ can be
partitioned $D_{i}=(A_{i},X_{i})$. We study a method for constructing a
uniform confidence region for the parameter vector
$\theta_{0}(\bm{x}^{(d)})=(\theta_{0}(x^{(j)}))_{j=1}^{d}$, where each
$\theta_{0}(x)$ is the unique scalar solution to the conditional moment
equation

$M(x;\theta,g_{0})=\mathbb{E}\left[m(D_{i};\theta,g_{0})\mid X_{i}=x\right]=0$
(1.1)

in $\theta$ and $\bm{x}^{(d)}=(x^{(j)})_{j=1}^{d}$ is a specified $d$-vector
in the domain of $X_{i}$. Here, $g_{0}$ is an unknown nuisance parameter,
identified via an auxiliary statistical problem, and $m(\cdot;\theta,g)$ is a
known moment function. Many problems in applied economics and causal inference
can be formulated as instances of (1.1), including nonparametric regression,
quantile regression, and estimation of conditional average treatment effects.

To fix ideas, consider Banerjee et al., (2015), who study the effects of a
poverty alleviation program implemented in Ghana.111Banerjee et al., (2015)
study data collected from several similar graduation programs. We focus on the
data from their evaluation of the program implemented in Ghana. Appendix E
gives further details on our treatment of the Banerjee et al., (2015) data.
For each individual in their sample, they observe the data
$D_{i}=(Y_{i},W_{i},Z_{i})$, where $Y_{i}$ is a measurement of total assets
taken two years after the implementation of the program, $W_{i}$ is an
indicator denoting assignment to the program, and $Z_{i}$ is a vector of
covariates. A broad aim of the study is to determine the conditions under
which recipients of aid experience lasting improvements in welfare. One
quantity that can inform this determination is the conditional average
treatment effect (CATE)

$\theta_{0}(x)=\mathbb{E}_{P}\left[Y_{i}(1)-Y_{i}(0)\mid
X_{i}=x\right]\leavevmode\nobreak\ ,$ (1.2)

where $Y_{i}(1)$ and $Y_{i}(0)$ are the potential outcomes generated by the
intervention $W_{i}$ and $X_{i}$ is some chosen subvector of $Z_{i}$. A
canonical approach to estimating (1.2) is premised on the observation that
$\theta_{0}(x)$ is the solution to the moment equation

$M(x;\theta,g_{0})=\mathbb{E}\left[(\mu_{1}(Z_{i})-\mu_{0}(Z_{i}))+\beta(W_{i},Z_{i})(Y_{i}-\mu_{W_{i}}(Z_{i}))-\theta\mid
X_{i}=x\right]=0\leavevmode\nobreak\ ,$ (1.3)

of Robins et al., (1994) and Hahn, (1998), where the nuisance parameter
$g_{0}$ collects the conditional outcome regression and Horvitz-Thompson
weight

$\mu_{w}(z)=\mathbb{E}_{P}\left[Y_{i}\mid
W_{i}=w,Z_{i}=z\right]\quad\text{and}\quad\beta(w,z)=\frac{w}{\pi(z)}-\frac{1-w}{1-\pi(z)}\leavevmode\nobreak\
,$ (1.4)

for the propensity score $\pi(z)=P\\{W_{i}=1\mid Z_{i}=z\\}$. See e.g., Nie
and Wager, (2021), Foster and Syrgkanis, (2023), and references therein for
further discussion.

Often, estimates of solutions to conditional moment equations of the form
(1.1) or (1.3) are obtained by solving the empirical conditional moment
equation

$M_{n}(x;\theta,\hat{g}_{n},\mathbf{D}_{n})=\sum_{i=1}^{n}K(x,X_{i})m(D_{i};\theta,\hat{g}_{n})=0$
(1.5)

in $\theta$, where $\hat{g}_{n}$ is some first-stage estimator of the nuisance
parameter $g_{0}$ and $K(x,x^{\prime})$ is some, potentially random and data-
dependent, kernel function measuring the distance between $x$ and
$x^{\prime}$. Subsampled or bagged kernels, introduced by Breiman, (1996), and
popularized by Wager and Athey, (2018) and Athey et al., (2019), are
particularly convenient, due in part to their computational efficiency and
robustness to tuning parameter choice. Popular examples of subsampled kernel
estimators include $k$-NN regression (Fix and Hodges,, 1989) and random forest
regression (Breiman,, 2001). Solutions to conditional moment equations (1.5)
constructed with random forest regression are referred to as Generalized
Random Forests (GRF) (Athey et al.,, 2019).

Figure 1 displays GRF estimates of the CATE (1.2) for the experiment studied
in Banerjee et al., (2015), where the chosen conditioning covariates $X_{i}$
are pretreatment measurements of monthly consumption and total assets.222These
estimates are constructed by approximating the solution to the conditional
moment (1.3) with a subsampled kernel estimator of the form (1.5) at each
point $x$ on a grid. Both the nuisance parameter estimator $\hat{g}_{n}$ and
the kernel $K(x,x^{\prime})$ are constructed with the implementation of random
forest regression made available through the “GRF” R package (Athey et al.,,
2019). The graduation program appears to be most effective for individuals
with high level of baseline consumption and a low level of baseline
assets.333The quartiles of baseline log consumption are 3.33, 3.76, and 4.20.
The quartiles of baseline assets are -0.45, -0.71, and 0.03. Panel A of Figure
E.4 displays of scatter plot of the joint distribution of baseline log
consumption and assets. That is, individuals with an opportunity to increase
their assets are able to do so only if they have a high level of baseline
consumption. The effect of the program for individuals with low baseline
consumption or high baseline assets appears more muted. These results are
suggestive of a poverty trap: individuals without a stable source of
consumption may be incentivized to sell productive assets (Kraay and
McKenzie,, 2014; Balboni et al.,, 2022).

Figure 1. CATE Estimates

---

Notes: Figure 1 displays a heat map giving CATE estimates for the intervention
studied in Banerjee et al., (2015) on post-treatment assets. The color of each
rectangle indicates the estimate of the CATE queried at the rectangle’s
central point. The horizontal and vertical axes display the baseline monthly
consumption, normalized to dollars and measured in logs base 10, and the
baseline value of an index for total assets. CATE estimates are obtained by
solving the empirical moment equation (1.5) for each value $x$ on an evenly
spaced grid on both axes. The nuisance parameter estimate $\hat{g}_{n}$ and
the kernel $K(x,x^{\prime})$ are constructed with the implementation of random
forest regression made available through the “GRF” R package (Athey et al.,,
2019). See Appendix E for further details.

The information communicated by Figure 1 is rich and granular. This contrasts
with the more widely encountered practice of reporting regression coefficients
on linear interactions of pretreatment covariates with treatment indicators.
In fact, as we will see, the latter approach gives a substantively different
picture of the heterogeneity in the effect of the Banerjee et al., (2015)
graduation program.

We contribute a method for assessing the statistical significance of estimates
typified by Figure 1. In particular, we propose a computationally simple
procedure for constructing uniform upper and lower confidence bounds for
solutions to conditional moment equations (1.1) centered around subsampled
kernel estimators of the form (1.5). Formally, we construct a family of random
intervals

$\hat{\mathcal{C}}(\bm{x}^{(d)})=\left\\{\hat{\mathcal{C}}(x^{(j)})=[c_{L}(x^{(j)}),c_{U}(x^{(j)})]:j\in[d]\right\\}\leavevmode\nobreak\
,$ (1.6)

on the basis of the observed data, such that

$\sup_{P\in\mathbf{P}}\big{|}P\left\\{\theta_{0}(\bm{x}^{(d)})\in\hat{\mathcal{C}}(\bm{x}^{(d)})\right\\}-(1-\alpha)\big{|}\leq
r_{n,d}$ (1.7)

for some sequence $r_{n,d}$, where $\mathbf{P}$ is some statistical family
that contains the distribution $P$ of the data $D_{i}$. We say that a region
(1.6) satisfying (1.7) is uniformly asymptotically valid at the rate
$r_{n,d}$. Here, uniformity operates over both the $d$-dimensional query-
vector $\bm{x}^{(d)}$ and the statistical family $\mathbf{P}$.444Following Li,
(1989), uniform validity of a confidence region for nonparametric regression
is often referred to as “Honesty” (see., e.g., Chernozhukov et al., (2014);
Armstrong and Kolesár, (2020); Kuchibhotla et al., (2023)). We use the label
“uniform validity,” following e.g., Romano and Shaikh, (2012), to distinguish
this property from “Honest” construction of subsampled kernels (Athey and
Imbens,, 2016; Wager and Athey,, 2018), which will be studied in detail in
Section 3. The main theoretical contribution of this paper is the construction
of uniform asymptotically valid confidence regions whose error rate $r_{n,d}$
converges to zero in asymptotic regimes where the number of points $d$ in the
query-vector $\bm{x}^{(d)}$ may increase much more quickly than the sample
size $n$.555The leading examples for choices of the query-vector
$\bm{x}^{(d)}$ are cases where $\bm{x}^{(d)}$ is taken to be the observed
values of the covariates $X_{1},\ldots,X_{n}$ or where $\bm{x}^{(d)}$ gives a
fine grid over the domain of $X_{i}$. In both cases the dimension $d$ of the
query-vector $\bm{x}^{(d)}$ is either equal to or potentially large relative
to the sample size $n$.

We begin, in Section 2, by defining the proposed confidence region and
illustrating its application to the Banerjee et al., (2015) experiment. Our
construction can be seen as an instance of subsampling (Politis et al.,, 1999;
Politis and Romano,, 1994), although our formal analysis is more directly
connected to the exchangeably weighted bootstrap (Præstgaard and Wellner,,
1993; Chernozhuokov et al.,, 2022).

In Section 3, we give a bound on the accuracy of the proposed confidence
region. We require that the confidence region be built around the solution to
a Neyman orthogonal moment. This restriction mitigates the error induced by
estimation of nuisance parameters. Our result is comparable to the generic
bounds on the accuracy of Gaussian multiplier bootstrap confidence regions for
nonparametric regression and $Z$-estimation given in Chernozhukov et al.,
(2014) and Belloni et al., (2018), respectively. We generalize these results,
in the sense that we treat inference for conditional $Z$-estimators whose
score function is potentially unknown to the researcher. We document that the
proposed confidence region is accurate and informative at empirically relevant
sample sizes with a simulation calibrated to the Banerjee et al., (2015) data.

As a by-product of this analysis, we give several new results on the
concentration and normal approximation of large order, high-dimensional,
$U$-statistics. In particular, we give a concentration inequality and central
limit theorem for high-dimensional $U$-statistics with explicit order-
dependence. These bounds are applicable to non-degenerate $U$-statistics whose
order $b$ satisfies $b=o(n)$, up to a dimension dependent logarithmic factor.
This generality represents a substantial improvement over existing results
(Song et al.,, 2019; Minsker,, 2023), that apply to the regime $b=o(n^{1/3})$,
and is essential for our application. Our results hinge on a new concentration
inequality for the difference between a $U$-statistic and its Hájek projection
(Hájek,, 1968), obtained through a Hoffman-Jørgensen type argument enabled by
a symmetrization inequality due to Sherman, (1994). We collect these results
in Section 4. Section 5 concludes.

###  1.1. Related Literature

There is an extensive literature on estimation of solutions to conditional
moment equations. See, for example, Newey, (1993), Ai and Chen, (2003), Chen
and Pouzo, (2012), and Chernozhukov et al., (2023). Chen and Christensen,
(2018) and Chen et al., (2024) propose related methods for constructing
uniform confidence bands for parameters identified by conditional moments,
emphasizing achieving minimax rates in Hölder classes by building confidence
regions around carefully constructed sieve estimators with Lepski’s method
(Chernozhukov et al.,, 2014). See also Singh and Vijaykumar, (2023) for an
analysis of inference for kernel ridge regression. By contrast, our aim is to
provide a simple procedure for uniform inference based on estimators whose
precise structure may be unknown to the user.

We contribute to a large literature on the role of Neyman orthogonality in
estimation of solutions to moment equations with nuisance parameters.
Chernozhukov et al., (2018), Chernozhukov et al., (2022), and Ichimura and
Newey, (2022) provide extensive discussion and guidance on the derivation of
orthogonal moments. Nie and Wager, (2021), Foster and Syrgkanis, (2023), and
Kennedy, (2023) apply aspects of this analysis to conditional moment
estimation. The closest paper, in this literature, is Belloni et al., (2018),
who give a related, general, treatment of $Z$-estimation. We build on this
analysis by studying conditional $Z$-estimators with unknown score functions.

Our paper is motivated by Athey and Imbens, (2016), Wager and Athey, (2018),
and Athey et al., (2019), who popularized estimation of solutions to
conditional moment equations with subsampled kernel regression. The estimators
considered in Section 3 are closely related to the “Orthogonal Random Forests”
estimator proposed in Oprescu et al., (2019). Methods for constructing
confidence regions for random forests regression are studied in Sexton and
Laake, (2009), Wager et al., (2014), Mentch and Hooker, (2016) and Athey et
al., (2019). We contribute to this literature by providing methods for
constructing uniformly asymptotically valid confidence regions.

Our formal analysis builds on a groundbreaking sequence of papers on central
limit theorems for maxima of sums initiated by Chernozhukov et al., (2013).
Extensions and refinements of these results are given in Chernozhukov et al.,
2017a  and Chernozhuokov et al., (2022). Similar approaches for applying these
results to the construction of uniform confidence regions are given in
Chernozhukov et al., (2014) and Belloni et al., (2018). Other aspects of our
analysis draw on the consideration of exchangeably weighted bootstrap
approximations to asymptotically linear statistics given in Præstgaard and
Wellner, (1993), Chung and Romano, (2013), and Yadlowsky et al., (2023).

The asymptotic analysis of $U$-statistics has a long and involved history (see
e.g., Lee,, 1990 for a textbook introduction). We provide a more detailed
literature review in Section 4. Recently, several papers have demonstrated
that central limit theorems for non-degenerate, real-valued, $U$-statistics of
order $b$ can hold, even if $b$ is increasing at some rate that satisfies
$b=o(n)$ as $n$ increases to infinity (Wager and Athey,, 2018; DiCiccio and
Romano,, 2022; Peng et al.,, 2022; Minsker,, 2023). However, to the best of
our knowledge, there are no general, order-explicit, exponential moment
inequalities or high-dimensional central limit theorems that obtain in this
regime. In particular, Song et al., (2019) and Minsker, (2023) establish high-
dimensional moment inequalities and central limit theorems that apply to the
regime $b=o(n^{1/3})$.666Minsker, (2023) additionally gives an exponential
moment inequality that holds in the regime $b=o(n)$ by placing stringent
restrictions on the smoothness of the kernel of the $U$-statistic of interest.
These smoothness restrictions will not hold in our application to subsampled
kernel regression. As we will see, this is insufficient for our main
application to subsampled kernel regression. We generalize these results,
obtaining exponential moment inequalities and non-asymptotic, high-
dimensional, central limit theorems, with explicit order dependence, that hold
in the regime $b=o(n)$.

Finally, we contribute to a large literature on the statistical analysis of
subsampled kernel regression and random forest regression. A wide variety of
consistency results are given in, e.g., Bühlmann and Yu, (2002), Lin and Jeon,
(2006), Biau et al., (2008), Mentch and Hooker, (2014), Scornet et al.,
(2015), and Cattaneo et al., (2024). High dimensional consistency results are
given in Syrgkanis and Zampetakis, (2020), Chi et al., (2022), and Huo et al.,
(2023).

###  1.2. Notation

The data $Z_{i}$ and $X_{i}$ take values in the spaces $\mathcal{Z}$ and
$\mathcal{X}$, respectively. Define the generic norm $\|\cdot\|_{\infty}$ on
$\mathcal{X}$. The nuisance parameter $g_{0}$ is a finite collection of $h$
real-valued functions $g_{0}=(g_{0}^{(k)})_{k=1}^{h}$, each having domain
$\mathcal{Z}$. Define the norm

$\|g-g_{0}\|_{2,\infty}=\sup_{k\in[h]}\sup_{j\in[d]}\left(\mathbb{E}\left[(g^{(k)}(Z_{i})-g_{0}^{(k)}(Z_{i}))^{2}\mid
X_{i}=x^{(j)}\right]\right)^{1/2}$ (1.8)

for any $g=(g^{(k)})_{k=1}^{h}$. The nuisance parameter $g_{0}$ takes values
in the space $\mathcal{G}$.

The quantities $c$ and $C$ denote universal positive constants, whose values
are allowed to depend only on the family of distributions $\mathbf{P}$. For
two real-valued functions $f$ and $g$ on a domain $\mathcal{X}$, we say
$g(x)\lesssim f(x)$ if $g(x)\leq Cf(x)$ for each $x$ in $\mathcal{X}$. The set
$\mathcal{S}_{n,b}$ collects all of the subsets of $[n]=\\{1,\ldots,n\\}$ of
size $b$ and $D_{\mathsf{s}}$ denotes the subset of the observed data
$\mathbf{D}_{n}$ with indices in the set $\mathsf{s}$. For a functional $F$ on
$\mathcal{F}$, we use the notation

$\partial_{f}F(f)[h]=\frac{\text{d}}{\text{d}t}F(f+th)\big{|}_{t=0}\quad\text{and}\quad\partial_{f,f}F(f)[h]=\frac{\text{d}^{2}}{\text{d}t^{2}}F(f+th)\big{|}_{t=0}$

to denote first and second order directional derivatives, respectively.
Throughout, for any function $f(x)$ and vector
$\bm{x}^{(d)}=(x^{(j)})_{j=1}^{d}$, we let $f(\bm{x}^{(d)})$ denote the vector
$(f(x^{(1)}),\cdots,f(x^{(d)}))$.

##  2\. Implementation

We build confidence regions around solutions to the empirical conditional
moment equation

$M_{n}(x;\theta,\hat{g}_{n},\mathbf{D}_{n})=0$ (2.1)

in the variable $\theta$, evaluated at each $x^{(j)}$ in $\bm{x}^{(d)}$. Let
$\hat{\theta}_{n}(x)$ denote the solution to (2.1) evaluated at $x$. Our
construction is agnostic to the structure of the empirical moment (2.1). This
property is essential for the application to subsampled kernel regression
considered in Section 3.

###  2.1. Construction

For any $\mathsf{s}$ in $[n]$, let $\hat{\theta}_{\mathsf{s}}(\bm{x}^{(d)})$
denote the vector of solutions to (2.1) evaluated at each $x^{(j)}$ in
$\bm{x}^{(d)}$ with the data $\mathbf{D}_{n}$ replaced by the subsample
$D_{\mathsf{s}}$. Our proposal is premised on approximating the sampling
distribution of the root

$R_{n}(\bm{x}^{(d)})=\hat{\theta}_{n}(\bm{x}^{(d)})-\theta_{0}(\bm{x}^{(d)})$
(2.2)

with the conditional distribution of the half-sample bootstrap root

$R^{*}_{n}(\bm{x}^{(d)})=\hat{\theta}_{\mathsf{h}}(\bm{x}^{(d)})-\hat{\theta}_{n}(\bm{x}^{(d)})\leavevmode\nobreak\
,$ (2.3)

where $\mathsf{h}$ denotes a random element of $\mathcal{S}_{n,n/2}$, i.e., a
random half-sample of $[n]$. The nuisance parameter estimator $\hat{g}_{n}$
does not need to be re-estimated when computing (2.3). A version of the half-
sample bootstrap is implemented by default in the GRF R package (Athey et
al.,, 2019).777The version of the half-sample bootstrap considered in Athey et
al., (2019) is based on combining estimates of the variances of components of
a linearization of the moment $M(x;\theta,g)$ with a Delta method type
argument. By contrast, the bootstrap root Equation 2.3 is agnostic to the
structure of the conditional moment under consideration. The half-sample
bootstrap is an instance of subsampling (Politis and Romano,, 1994; Politis et
al.,, 1999).

Let $\hat{\lambda}^{2}_{n,j}$ denote the variance of
$\sqrt{n}R^{*}_{n}(x^{(j)})$, conditioned on the data $\mathbf{D}_{n}$, and
let $\mathsf{cv}_{n}(\alpha)$ denote the $1-\alpha$ quantile of the
distribution of the studentized process

$\hat{S}^{*}_{n}(\bm{x}^{(d)})=\sqrt{n}\|\hat{\Lambda}_{n}^{-1/2}R^{*}_{n}(\bm{x}^{(d)})\|_{\infty}\leavevmode\nobreak\
,$ (2.4)

again conditioned on the data $\mathbf{D}_{n}$. Here, $\hat{\Lambda}_{n}$
denotes the diagonal matrix with elements $\hat{\lambda}^{2}_{n,j}$.888The
quantities $\hat{\lambda}^{2}_{n,j}$ and $\mathsf{cv}_{n}(\alpha)$ are easily
approximated by resampling the bootstrap root (2.3). To simplify exposition,
we omit explicit consideration of residual randomness induced by this
approximation. The confidence region considered in this paper has the
following structure.

######  Definition 2.1 (Uniform Confidence Region).

Define the intervals

$\hat{\mathcal{C}}(x^{(j)})=\hat{\theta}_{n}(x^{(j)})\pm
n^{-1/2}\hat{\lambda}_{n,j}\mathsf{cv}_{n}(\alpha)\quad\text{for each}\quad
j\text{ in }[d]\leavevmode\nobreak\ .$ (2.5)

The level-$\alpha$ uniform confidence region for $\theta_{0}(\bm{x}^{(d)})$ is
given by $\hat{\mathcal{C}}(\bm{x}^{(d)})$.

Confidence regions with the same structure, based on different choices of
bootstrap root, are studied in, e.g., Chernozhukov et al., (2014) and Belloni
et al., (2018). The essential feature of the bootstrap root (2.3) is that it
can be computed without knowing anything about the structure of the estimator
$\hat{\theta}_{n}(\bm{x}^{(d)})$. In particular, approaches based on the
Rademacher or Gaussian multiplier bootstrap rely on knowledge of a linear
approximation to $\hat{\theta}_{n}(\bm{x}^{(d)})$.

To gain intuition, suppose that the estimator $\hat{\theta}_{n}(\bm{x}^{(d)})$
satisfies a linear representation

$\hat{\theta}_{n}(\bm{x}^{(d)})=\frac{1}{n}\sum_{i=1}^{n}\overline{u}(\bm{x}^{(d)},D_{i})$
(2.6)

for some function $\overline{u}(\cdot,\cdot)$. Let $V_{1},\ldots,V_{n}$ denote
a collection of random variables, where $V_{i}$ takes the value $1$ if the
index $i$ is an element of the random set $\mathsf{h}$ used to define the
half-sample bootstrap, and takes the value $-1$ otherwise. Observe that

$\displaystyle
R^{*}_{n}(\bm{x}^{(d)})=\hat{\theta}_{\mathsf{h}}(\bm{x}^{(d)})-\hat{\theta}_{n}(\bm{x}^{(d)})$
$\displaystyle=\frac{2}{n}\sum_{i=1}^{n}\mathbb{I}\\{i\in\mathsf{h}\\}\overline{u}(\bm{x}^{(d)},D_{i})-\frac{1}{n}\sum_{i=1}^{n}\overline{u}(\bm{x}^{(d)},D_{i})$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}V_{i}\left(\overline{u}(\bm{x}^{(d)},D_{i})-\theta_{0}(\bm{x}^{(d)})\right)\leavevmode\nobreak\
.$ (2.7)

The representation (2.7) is due to Yadlowsky et al., (2023), who draw on a
similar observation made in the context of two-sample testing in Chung and
Romano, (2013). The weights $V_{i}$ are exchangeable Rademacher random
variables, i.e., they are uniformly distributed on $\\{1,-1\\}$. If the
weights were fully independent, the representation (2.7) reduces to the
Rademacher bootstrap.999In Section D.2, we consider a variant of the bootstrap
root (2.3) based on re-estimating and re-scaling
$\hat{\theta}_{n}(\bm{x}^{(d)})$ on a subsample of a random size
$\mathsf{Bin}(n,1/2)$. For this construction, the equivalent objects to the
weights $V_{i}$ are fully independent. That is, this bootstrap root is
equivalent to the Rademacher bootstrap root for linear statistics. We show
that the resulting confidence regions obtain the same error rates on coverage
accuracy. Chernozhuokov et al., (2022) give a central limit theorem for
statistics of the form (2.6). We obtain a central limit theorem for the
statistic (2.7) by adapting a coupling argument due to Yadlowsky et al.,
(2023). In particular, we show (in a suitable sense) that limiting conditional
distribution of the bootstrap root (2.3) is the same as the limiting
distribution of the root (2.2).

In practice, many widely applied estimators are not perfectly linearly
decomposable. Often, however, estimators do satisfy an approximate linear
decomposition, in the sense that the equality (2.6) holds with a remainder
term of order, say, $o_{p}(n^{-\gamma})$ for some positive constant $\gamma$.
As we will see, estimators constructed with subsampled kernels are
approximately linear around some unknown function $\overline{u}(\cdot,\cdot)$.
If an estimator $\hat{\theta}_{n}(\bm{x}^{(d)})$ is approximately linear, then
the subsampled estimate $\hat{\theta}_{\mathsf{h}}(\bm{x}^{(d)})$ is
_immediately_ also approximately linear.101010Approximate linearity does not
immediately imply a representation analogous to (2.7) for a root constructed
with an empirical bootstrap, as approximate linearity would not necessarily
hold under the empirical distribution. Thus, for approximately linear
statistics, the representation (2.7) continues to hold, now with a remainder
term of order $o_{p}(n^{-\gamma})$. As a consequence, the validity of the
confidence region formulated in 2.1, for approximately linear estimators,
follows from the central limit theorems discussed in the preceding paragraph
and an appropriate generalization of Slutsky’s theorem.

###  2.2. Application

We now return to the application to the data studied in Banerjee et al.,
(2015), considered in Section 1. Figure 2 displays upper and lower confidence
bounds for the CATE (1.2) on post-treatment assets. These bounds are built
with the confidence region formulated in 2.1.

Consistent with the qualitative features of the estimates displayed in Figure
1, the null hypothesis that the CATE is equal to zero is only rejected for
individuals with low baseline assets and high baseline consumption. That is,
the graduation program has a positive impact on individuals who do not have
many assets to begin with, but who do have access to a stable source of
consumption. On the other hand, the confidence regions contain zero for
individuals who have either low baseline consumption or high baseline assets.

Figure 2. Half-Sample Confidence Region

Panel A: Upper Bound
---
Panel B: Lower Bound

Notes: Figure 2 displays heat maps giving half-sample upper and lower
confidence bounds for the CATE of the intervention studied in Banerjee et al.,
(2015) on post-treatment total assets. The confidence bounds are constructed
at level $\alpha=0.1$. The upper and lower bounds are displayed with different
color palettes to emphasize the use of different scales. A contour line has
been superimposed over the lower bound to demarcate where the bound crosses
zero. The axes and estimator are the same as in Figure 1.

It is illustrative to contrast the estimates and confidence bounds displayed
in Figures 1 and 2 with a more frequently encountered method for assessing
treatment effect heterogeneity—interacted linear regression. Table 1 reports
estimates and standard errors associated with several linear regression
specifications, constructed with the same data. The first column reports the
coefficient from a regression of post-treatment assets on a treatment
indicator. Consistent with results reported in Banerjee et al., (2015), the
average effect of the program is positive and significant. The second and
third specifications interact the treatment indicator with pre-treatment
assets and consumption, respectively. In both cases, the estimate of the
coefficient on the interaction is statistically insignificant. The fourth
specification interacts both pre-treatment assets and consumption with a
treatment indicator. Here, all coefficients lose statistical significance.

The implicit view of much of applied economics appears to be that the
flexibility afforded by nonparametric methods is not worth sacrificing the
statistical precision of more parsimonious, linear, alternatives.111111We
conduct a small survey of papers published in the American Economic Review in
the first six months of 2023. Of 38 empirical papers, 30 assess treatment
effect heterogeneity in some way. As best as we can tell, only two papers
display nonparametric estimates of effect heterogeneity. By contrast, 10
display the results of an interacted linear regression typified by Table 1.
The rest are either structural papers, or are only interested in interactions
with binary covariates. See Section D.1 for more details. The exercise here
suggests otherwise. The linearity imposed by interacted regression masks the
structure in the effect heterogeneity recovered by the GRF estimator. The
resulting bias is so substantial that statistical significance is lost. The
half-sample confidence regions developed in this paper enable the recovery of
statistically significant measurements of effect heterogeneity. The remainder
of the paper is devoted to developing theoretical guarantees on the accuracy
of confidence bounds typified by Figure 2.

Table 1. Interacted Linear Regression

| Dependent Variable: Post-Treatment Assets
---|---
| (1) | (2) | (3) | (4)
Treatment | 0.218 (0.031) | 0.194 (0.032) | 0.102 (0.172) | 0.111 (0.201)
Assets |  | 0.770 (0.051) |  | 1.229 (0.288)
Consumption |  |  | 0.06 (0.022) | -0.026 (0.026)
Assets $\times$ Consumption |  |  |  | -0.119 (0.076)
Treatment $\times$ Assets |  | 0.047 (0.093) |  | -0.065 (0.552)
Treatment $\times$ Consumption |  |  | 0.031 (0.045) | 0.021 (0.050)
Treatment $\times$ Assets $\times$ Consumption |  |  |  | 0.027 (0.141)
Observations: | 2,438 | 2,438 | 2,438 | 2,438

Notes: Table 1 reports estimates of the coefficients of four linear regression
specifications, constructed with the Banerjee et al., (2015) data. Robust
standard errors are displayed in parentheses. Coefficients on intercepts are
not displayed.

##  3\. Subsampled Kernel Regression

We establish a bound on the accuracy of the nominal coverage probability for
the confidence region introduced in 2.1. We begin in Section 3.1 by discussing
subsampled kernel regression and introducing several quantities that take a
prominent role in our analysis. Our results apply to conditional moments that
satisfy a Neyman orthogonality condition, in addition to several simple
regularity conditions. We overview these restrictions in Section 3.2. Our main
result is stated in Section 3.3. The results of a simulation calibrated to the
Banerjee et al., (2015) data are reported in Section 3.4.

###  3.1. Subsampled Kernel Regression

Subsampled kernel regression is a broad class of algorithms for solving
regression problems of the form (1.5), based on constructing a data-driven
kernel function $K(x,x^{\prime})$ with subsampling. Formally, fix some
positive integer $r$ and let $(\mathsf{s}_{q})_{q=1}^{r}$ collect a sequence
of subsets of $[n]$ drawn independently and uniformly from
$\mathcal{S}_{n,b}$. Let $\xi$ denote some auxiliary source of randomness and
let $(\xi_{\mathsf{s}_{q}})_{q=1}^{r}$ collect a set of independent random
variables with the same distribution as $\xi$. We study conditional empirical
moment estimators of the form (1.5), where the kernel function
$K(x,x^{\prime})$ admits the decomposition

$K(x,X_{i})=\sum_{q=1}^{r}\mathbb{I}\\{i\in\mathsf{s}_{q}\\}\kappa(x,X_{i},D_{\mathsf{s}_{q}},\xi_{\mathsf{s}_{q}})$
(3.1)

for some known kernel $\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi_{\mathsf{s}})$.
Several widely applied instances of subsampled kernels are as follows.

######  Example 3.1 (Subsampled $k$-Nearest Neighbor Regression).

Nearest-neighbors regression is a simple example of a kernel with the
structure (3.1) (Fix and Hodges,, 1989). Here, the kernel
$\kappa(x,X_{i^{\prime}},D_{\mathsf{s}},\xi_{\mathsf{s}})$ is non-zero if and
only if $X_{i^{\prime}}$ is one of the $k$ closest points to $x$ among the
points in the subsample $D_{\mathsf{s}}$. $\blacksquare$

######  Example 3.2 (Random Forest Regression).

Random forest regression, introduced by Breiman, (2001), is another example of
a kernel with the structure (3.1). In this case, each pair
$(D_{\mathsf{s}},\xi_{\mathsf{s}})$ generates some partition of the domain of
$X_{i}$. The kernel $\kappa(x,x^{\prime},D_{\mathsf{s}},\xi_{\mathsf{s}})$ is
non-zero if and only if $x$ and $x^{\prime}$ are in the same element of the
partition generated by $(D_{\mathsf{s}},\xi_{\mathsf{s}})$. Often, such
partitions are constructed with recursive algorithms, e.g., the “CART”
algorithm of Breiman et al., (1984). $\blacksquare$

We impose the following restrictions on the kernels under consideration.

######  Assumption 3.1 (Honesty and Symmetry).

(i) The kernel $\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi)$ is Honest in the sense
that

$\kappa(x,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})\perp\\!\\!\\!\perp
m(D_{i};\theta,g)\mid X_{i},D_{\mathsf{s}_{-i}}\leavevmode\nobreak\ ,$ (3.2)

where $\perp\\!\\!\\!\perp$ denotes conditional independence and
$\mathsf{s}_{-i}$ denotes the set $\mathsf{s}\setminus\\{i\\}$.

(ii) The kernel $\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi)$ is positive and
satisfies the restriction
$\sum_{i\in{s}}\kappa(\cdot,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})=1$ almost
surely. Moreover, the conditional expectation
$\mathbb{E}\left[\kappa(\cdot,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})\mid
D_{\mathsf{s}}\right]$ is invariant to permutations of the data
$D_{\mathsf{s}}$.

The “Honesty” condition stipulated in Part (i) of Assumption 3.1 imposes the
restriction that any part of the data $D_{i}$ that can affect the value of the
moment $m(D_{i};\theta,g)$ cannot affect the value of the kernel
$\kappa(x,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})$. This condition was
introduced in Athey and Imbens, (2016). Honesty is often achieved through
kernel construction schemes based on sample-splitting; see Wager and Athey,
(2018) and Athey et al., (2019) for further discussion. Part (i) of Assumption
3.1 imposes several weak regularity conditions.

The following two quantities restrict the “size” and “variability” of the
chosen kernel.

######  Definition 3.1 (Shrinkage and Incrementality).

(i) We say that the kernel
$\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi_{\mathsf{s}})$ has a uniform shrinkage
rate $\varepsilon_{b}$ if

$\sup_{P\in\mathbf{P}}\sup_{j\in[d]}\mathbb{E}\left[\max\left\\{\|X_{i}-x^{(j)}\|_{\infty}:\kappa(x^{(j)},X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})>0\right\\}\right]\leq\varepsilon_{b}\leavevmode\nobreak\
.$ (3.3)

(ii) We say that a kernel
$\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi_{\mathsf{s}})$ is uniformly incremental
if

$\inf_{P\in\mathbf{P}}\sup_{j\in[d]}\operatorname{Var}\left(\mathbb{E}\left[\sum_{i\in\mathsf{s}}\kappa(x,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})m(D_{i};\theta,g)\mid
l\in\mathsf{s},D_{l}=D\right]\right)\gtrsim b^{-1}$ (3.4)

where $D$ is an independent random variable with distribution $P$.

The shrinkage rate of a kernel
$\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi_{\mathsf{s}})$ is analogous to the
bandwidth of a classical, deterministic, kernel. The incrementality
restriction ensures that the chosen kernel is not overly dependent on a single
data point. Both notions were introduced by Wager and Athey, (2018) and have
been characterized explicitly for various widely applied subsampled kernel
estimators.121212The terminology “shrinkage” was introduced in Oprescu et al.,
(2019), and is not intended to connote (explicit) regularization.

######  Example 3.3 (continues=ex: knn).

In the case of honest subsampled $k$-NN regression, Khosravi et al., (2019)
show that $\varepsilon_{b}\lesssim b^{-1/p}$, where $p$ is the “intrinsic
dimension” of the measure of the covariates $X_{i}$. Roughly speaking, a
distribution has an intrinsic dimension of $p$ if it is (locally) well
approximated by a measure supported on a subspace of $\mathcal{X}$ of
dimension $p$. See e.g., Kpotufe, (2011) for further discussion. In turn,
Khosravi et al., (2019) and Peng et al., (2022) show that the kernels
associated with both honest and non-honest variants of $k$-NN regression are
incremental, up to logarithmic factors that depend on the dimension of the
covariates. $\blacksquare$

######  Example 3.4 (continues=ex: random forest).

Analogously, for honest random forest regression, Wager and Athey, (2018)
establish that $\varepsilon_{b}\lesssim b^{-c/p}$, where $p$ is the dimension
of the domain of $X_{i}$. See e.g., Wager and Athey, (2018) and Oprescu et
al., (2019) for further discussion. Bounds adaptive to the intrinsic dimension
of the measure of $X_{i}$ are given in Huo et al., (2023) under further
restrictions. Wager and Athey, (2018) and Peng et al., (2022) give simple
conditions under which the kernel associated with subsampled, honest, random
forest regression is uniformly incremental, again up to dimension dependent
logarithmic factors. $\blacksquare$

###  3.2. Moment Restrictions

The uniform confidence region introduced in 2.1 is based on undersmoothing, in
the sense that the construction makes no explicit correction for bias. In
other words, the confidence regions that we consider are reliant on the use of
estimators $\hat{\theta}_{n}(\bm{x}^{(d)})$ whose bias is of a smaller
stochastic order than the sampling variance. To this end, we emphasize the use
of conditional moments $M(\cdot;\theta_{0},g_{0})$, that satisfy a local
Neyman Orthogonality condition.

######  Definition 3.2 (Local Neyman Orthogonality).

We say that a conditional moment $M(\cdot;\theta_{0},g_{0})$ is uniformly
locally Neyman orthogonal if

$\partial_{g}M(x^{(j)};\theta_{0},g_{0})[g-g_{0}]=0$ (3.5)

for all $P$ in $\mathbf{P}$ and $x^{(j)}$ in $\bm{x}^{(d)}$.

The use of Neyman orthogonal moments ensures that the bias induced by the
estimation of the nuisance parameter $g_{0}$ with $\hat{g}_{n}$ is small
(Newey,, 1994).

In addition to Neyman orthogonality, we require several smoothness
restrictions on the function $m(\cdot;\theta,g)$. In the main text, to ease
exposition, we impose the following linearity and boundedness restriction.

######  Assumption 3.2 (Moment Linearity and Boundedness).

The moment function $m(\cdot;\theta,g)$ satisfies the linear representation

$m(D_{i};\theta,g)=m^{(1)}(D_{i};\theta,g)\cdot\theta+m^{(2)}(D_{i};g)$ (3.6)

for some known functions $m^{(1)}(\cdot;\theta,g)$ and $m^{(2)}(\cdot;g)$.
Moreover, the absolute value of the function $m(\cdot;\theta,g)$ is bounded by
the constant $(\theta+1)\phi$ almost surely.

The linearity restriction entailed in Assumption 3.2 is inessential and is
imposed for the sake of simplicity.131313In Appendix A, we show that moment
linearity can be replaced by the high-level assumption that
$\hat{\theta}_{n}(\bm{x}^{(d)})$ is consistent for $\theta_{0}(\bm{x}^{(d)})$.
This state of affairs is standard in $M$-estimation problems (see e.g., Newey
and McFadden,, 1994). As our running examples use linear moments, and
sufficient conditions for the consistency of $\hat{\theta}_{n}(\bm{x}^{(d)})$
have been established (see e.g., Assumption 4.1 and Theorem 4.3 of Oprescu et
al.,, 2019), we omit a detailed consideration of nonlinear moments. The
boundedness restriction is easily weakened to a slightly more involved
assumption stated in terms of the sub-exponential norm. Again, we impose
boundedness to ease exposition.

We maintain the following mild smoothness restrictions on the moment function
$M(\cdot;\theta_{0},g_{0})$.

######  Assumption 3.3 (Moment Smoothness).

(i) The moment $M(\cdot;\theta,g_{0})$ is second order smooth, in the sense
that

$\sup_{P\in\mathbf{P}}\sup_{j\in[d]}\big{|}\partial_{g,g}M(x^{(j)};\theta_{0},g_{0})[g-g_{0}]\big{|}\lesssim\|g-g_{0}\|^{2}_{2,\infty}$
(3.7)

for each $g$ in $\mathcal{G}$.

(ii) The variogram

$V(x;g)=\mathbb{E}\left[(m(Z_{i};\theta_{0}(x),g)-m(Z_{i};\theta_{0}(x),g_{0}))^{2}\mid
X=x\right]$

is uniformly Lipschitz in both of its components, in the sense that

$\sup_{P\in\mathbf{P}}\sup_{g\in\mathcal{G}}\big{|}V(x;g)-V(x^{\prime};g)\big{|}\lesssim\|x-x^{\prime}\|_{\infty}$
(3.8)

holds for all $x$ and $x^{\prime}$ in $\mathcal{X}$ and

$\sup_{P\in\mathbf{P}}\sup_{j\in[d]}|V(x^{(j)},g)-V(x^{(j)},g^{\prime})|\lesssim\|g-g^{\prime}\|^{2}_{2,\infty}$
(3.9)

holds for each $g$ and $g^{\prime}$ in $\mathcal{G}$.

(iii) Define the moments

$M^{(1)}(x;\theta,g)=\mathbb{E}\left[m^{(1)}(D_{i};\theta,g)\mid
X_{i}=x\right]\quad\text{and}\quad
M^{(2)}(x;g)=\mathbb{E}\left[m^{(2)}(D_{i};g)\mid
X_{i}=x\right]\leavevmode\nobreak\ ,$

associated with the functions $m^{(1)}(\cdot;\theta,g)$ and $m^{(2)}(\cdot;g)$
introduced in Assumption 3.2. Both moments are uniformly Lipschitz in their
first component, in the sense that

$\sup_{P\in\mathbf{P}}\sup_{g\in\mathcal{G}}\big{|}M^{(j)}(x;g)-M^{(j)}(x^{\prime};g)\big{|}\lesssim\|x-x^{\prime}\|_{\infty}$
(3.10)

for each $j$ in $\\{1,2\\}$ and all $x$ and $x^{\prime}$ in $\mathcal{X}$.
Moreover, the first moment is uniformly Lipschitz in its third component and
bounded from below in the sense that

$\displaystyle\sup_{P\in\mathbf{P}}\sup_{j\in[d]}\big{|}M^{(1)}(x^{(j)};g)-M^{(1)}(x^{(j)};g_{0})\big{|}$
$\displaystyle\lesssim\|g-g_{0}\|_{2,\infty}\quad\text{and}$ (3.11)
$\displaystyle\inf_{P\in\mathbf{P}}\inf_{j\in[d]}\big{|}M^{(1)}(x^{(j)};g)\big{|}$
$\displaystyle\geq c$ (3.12)

for each $g$ in $\mathcal{G}$ and some positive constant $c$.

Neyman orthogonal moments satisfying the smoothness restrictions specified in
Assumptions 3.2 and 3.3 are available for many widely considered statistical
problems.

######  Example 3.5 (Nonparametric Regression).

Suppose that $D_{i}=(Y_{i},Z_{i})$ contains a measurement of an outcome
$Y_{i}$ and a vector of covariates $Z_{i}$ and that the data $X_{i}$ are some
function, e.g., a sub-vector, of $Z_{i}$. In this setting, we are often
interested in estimating the conditional expectation function

$\theta_{0}(x)=\mathbb{E}_{P}\left[Y_{i}\mid
X_{i}=x\right]\leavevmode\nobreak\ .$

Here, the parameter $\theta_{0}(x)$ is identified via the linear moment
function $m(D_{i},\theta)=Y_{i}-\theta$. As there are no nuisance parameters,
local Neyman orthogonality is immediate. $\blacksquare$

######  Example 3.6 (Conditional Average Treatment Effects).

Now, suppose that $D_{i}$ additionally contains a binary valued variable
$W_{i}$ that indicates whether unit $i$ has been randomly assigned to an
intervention. In this case, interest may be in estimating a CATE (1.2). Under
strong ignorability, the CATE is identified by several moment functions, e.g.,
those implied by inverse propensity weighting or outcome regression (Imbens
and Rubin,, 2015). The moment (1.3) is the unique Neyman orthogonal
identifying moment for the parameter (1.2) (see e.g., Hahn,, 1998;
Chernozhukov et al.,, 2018). Note that, in this case, Part (i) of Assumption
3.3 is implied by the more refined bound

$\sup_{P\in\mathbf{P}}\sup_{j\in[d]}\big{|}\partial_{g,g}M(x^{(j)};\theta_{0},g_{0})[g-g_{0}]\big{|}\lesssim\|\mu-\mu_{0}\|_{2,\infty}\|\beta-\beta_{0}\|_{2,\infty}\leavevmode\nobreak\
,$ (3.13)

where $\mu$ and $\beta$ are defined in (1.4). This structure yields the
celebrated “Double Robustness” result for estimation of average treatment
effects and conditional average treatment effects; see Chernozhukov et al.,
(2018) and Kennedy, (2023) for further discussion. We impose the more general
condition (3.7), as this bound exhibits many of the same features and will
hold for a wider variety of problems. $\blacksquare$

######  Example 3.7 (Conditional Local Average Treatment Effects).

In turn, suppose that $D_{i}$ additionally contains a measurement of a binary
instrument $B_{i}$ and that interest is now in estimating the conditional
local average treatment effect

$\theta_{0}(x)=\mathbb{E}_{P}\left[Y_{i}(1)-Y_{i}(0)\mid
W_{i}(1)>W_{i}(0),X_{i}=x\right]\leavevmode\nobreak\ ,$ (3.14)

where $W_{i}(1)$ and $W_{i}(0)$ are the potential outcomes for the treatment
$W_{i}$ generated by the instrument $B_{i}$. Under standard assumptions
(Angrist et al.,, 1996), Tan, (2006) and Frölich, (2007) show that the unique
Neyman orthogonal moment function for the parameter (3.14) is given by

$m(D_{i};\theta,g)={1\choose-\theta}^{\top}\left((\gamma_{1}(Z_{i})-\gamma_{0}(Z_{i}))+\beta^{\varrho}(B_{i},Z_{i})\left({Y_{i}\choose
W_{i}}-\gamma_{B_{i}}(Z_{i})\right)\right)\leavevmode\nobreak\ ,$ (3.15)

where the nuisance parameter $g_{0}$ now collects the nuisance functions

$\gamma_{b}(z)=\mathbb{E}_{P}\left[{Y_{i}\choose W_{i}}\mid
B_{i}=b,Z_{i}=z\right]\quad\text{and}\quad\beta^{\varrho}(b,z)=\frac{b}{\varrho(z)}-\frac{1-b}{\varrho(z)}$
(3.16)

and $\varrho(z)=P\\{B_{i}=1\mid Z_{i}=z\\}$ denotes the instrumental
propensity score. $\blacksquare$

Additional examples of problems where smooth Neyman orthogonal identifying
moments are available include estimation of partially linear regression and
partially linear instrumental variable regression (Chernozhukov et al.,,
2018), dynamic treatment effects (Lewis and Syrgkanis,, 2021), and long-term
treatment effects identified by surrogate outcomes (Athey et al.,, 2020; Chen
and Ritzwoller,, 2023). Ichimura and Newey, (2022) and Belloni et al., (2018)
provide extensive discussion on the derivation of Neyman orthogonal moments.
Further discussion of the role of Neyman orthogonality in semiparametric
estimation is given in Chernozhukov et al., (2018) and Foster and Syrgkanis,
(2023).

###  3.3. Coverage

The following theorem gives a bound on the error in the nominal coverage
probability of the confidence regions introduced in 2.1. We emphasize that the
result is applicable to asymptotic regimes where the dimension of the query-
vector $\bm{x}^{(d)}$ can be exponentially larger than the sample size $n$.

######  Theorem 3.1 (Coverage Error Bound).

Suppose that the kernel $\kappa(\cdot,\cdot,D_{\mathsf{s}},\xi_{\mathsf{s}})$
satisfies Assumption 3.1, has uniform shrinkage rate $\varepsilon_{b}$, and is
uniformly incremental and that the Neyman orthogonal moment function
$M(\cdot;\theta_{0},g_{0})$ satisfies Assumptions 3.2 and 3.3. Moreover,
suppose that quantity $\|\theta_{0}(\bm{x}^{(d)})\|_{\infty}$ is uniformly
bounded as $P$ varies over $\mathbf{P}$ and that $r$ has been chosen to
satisfy $n\leq b\sqrt{r}$. If the nuisance parameter estimator $\hat{g}_{n}$
satisfies the probability bound

$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\|\hat{g}_{n}-g_{0}\|^{2}_{2,\infty}\geq\sqrt{\frac{b}{n}}\delta^{2}_{n,g}\right\\}\lesssim\frac{1}{n}$
(3.17)

for some sequence $\delta_{n,g}$, then the confidence region formulated in 2.1
satisfies the bound

$\displaystyle\sup_{P\in{\mathbf{P}}}\big{|}P\left\\{\theta_{0}(\bm{x}^{(d)})\in\hat{\mathcal{C}}(\bm{x}^{(d)})\right\\}-(1-\alpha)\big{|}$
$\displaystyle\quad\quad\lesssim\left(\frac{b\log^{5}(dn)}{n}\right)^{1/4}+\left(\delta^{2}_{n,g}+\sqrt{\frac{n}{b}}\varepsilon_{b}\right)\sqrt{\log(d)}\leavevmode\nobreak\
,$ (3.18)

for all sufficiently large $n$ and $b$.141414The statistical family
$\mathbf{P}$ is defined implicitly by the omitted constants in the uniform
bounds stated in 3.1 and Assumption 3.3, in addition to the restriction that
$\|\theta_{0}(\bm{x}^{(d)})\|_{\infty}$ is bounded as $P$ varies over
$\mathbf{P}$.

######  Remark 3.1.

Theorem 3.1 follows from an application of a more general result stated in
Appendix A. This result applies to conditional moment estimators that are not
necessarily constructed with subsampled kernels. Rather, the result holds
under the high-level assumption that the estimator
$\hat{\theta}_{n}(\bm{x}^{(d)})$ is approximately linear, with a sufficiently
small remainder term, and has sufficiently small bias and stochastic
equicontinuity. These conditions are verified for subsampled kernel regression
in Appendix B.

Several aspects of this argument are new. In particular, through a standard
series of expansions (see e.g., Chernozhukov et al.,, 2018), we show that the
root $R_{n}(\bm{x}^{(d)})$ is approximated by

$\displaystyle\overline{U}_{n,b}(\bm{x}^{(d)})$
$\displaystyle=\frac{1}{N_{b}}\sum_{\mathsf{s}\in\mathcal{S}_{n,b}}\mathbb{E}\left[u(\bm{x}^{(d)};D_{\mathsf{s}},\xi_{\mathsf{s}},\theta_{0},g_{0})\mid
D_{\mathsf{s}}\right]\leavevmode\nobreak\ ,\quad\text{where}$ (3.19)
$\displaystyle u(x;D_{\mathsf{s}},\xi_{\mathsf{s}},\theta,g)$
$\displaystyle=M^{(1)}(x;g_{0})^{-1}\sum_{i\in\mathsf{s}}\big{(}\kappa(x,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})m(D_{i};\theta,g)-\mathbb{E}\big{[}\kappa(x,X_{i},D_{\mathsf{s}},\xi_{\mathsf{s}})m(D_{i};\theta,g)\big{]}\big{)}\leavevmode\nobreak\
,$

with a remainder term given by the second term in (3.18). The quantity (3.19)
can be recognized as a complete, deterministic, $U$-statistic of order
$b$.151515The incrementality condition specified in 3.1 ensures that this
$U$-statistic is non-degenerate. We then apply a new result that demonstrates
that $U$-statistics of order $b$ are approximately linear with a remainder
term of order $(b/n)^{b/2}$, up to logarithmic factors. This is a dramatic
improvement over analogous results given in Song et al., (2019) and Minsker,
(2023), whose remainder terms exhibit polynomial decay as $b$ and $n$ grow and
only apply to regime $b=o(n^{1/3})$. This regime is unsuitable for our
application. This result has other applications and is discussed in detail in
Section 4.

In other words, we establish that the root $R_{n}(\bm{x}^{(d)})$ satisfies a
linear representation of the form (2.6) up to a small remainder term. It
immediately follows that the bootstrap root $R^{*}_{n}(\bm{x}^{(d)})$
satisfies the linear representation (2.7), up to a small remainder term, as it
is constructed with subsampling. We conclude by applying suitable high-
dimensional central limit theorems (Chernozhuokov et al.,, 2022) to the linear
terms (2.6) and (2.7). To handle the half-sample bootstrap, we apply an
argument, based on a coupling proposed in Yadlowsky et al., (2023), similar to
the Poissonization trick used in Præstgaard and Wellner, (1993).
$\blacksquare$

######  Remark 3.2.

The bound (3.18) can be interpreted as a bias-variance decomposition. The
first term in (3.18) results from a bound on the accuracy of a normal
approximation to (3.19). The term involving the kernel shrinkage
$\varepsilon_{b}$ is a remnant of a bound on the supremum of the bias of the
estimator $\hat{\theta}_{n}(\bm{x}^{(d)})$. $\blacksquare$

######  Remark 3.3.

It is worth emphasizing that the coverage bound given in Theorem 3.1 is
achieved without assuming that the estimator $\hat{\theta}_{n}(\bm{x}^{(d)})$
is constructed with sample-splitting. That is, the nuisance parameter
estimator $\hat{g}_{n}$ can be computed using the same data used to evaluate
the conditional moment $M_{n}(\cdot;\theta,\hat{g}_{n},\bm{D}_{n})$. Often,
sample-splitting is necessary to ensure that stochastic equicontinuity is
sufficiently small (Chernozhukov et al.,, 2018). Avoiding sample-splitting can
be practically important; Ritzwoller and Romano, (2023) show that the
randomness induced by sample-splitting can be large. $\blacksquare$

Theorem 3.1 states a bound on coverage error in terms of two generic
sequences: $\delta_{n,g}$, expressing the rate of convergence of the nuisance
parameter estimator $\hat{g}_{n}$, and $\varepsilon_{b}$, measuring the
effective bandwidth of the kernel. The bound (3.18) exhibits an interesting
tradeoff between these objects and the choice of subsample size $b$. To see
this, suppose that

$b=n^{\gamma_{b}},\quad\varepsilon_{b}\lesssim
b^{-\gamma_{\varepsilon}},\quad\text{and}\quad\|\hat{g}_{n}-g_{0}\|_{2,\infty}\lesssim
n^{-\gamma_{g}}\leavevmode\nobreak\ ,$ (3.20)

with probability greater than $1-1/n$, for some constants $\gamma_{b}$,
$\gamma_{\varepsilon}$, and $\gamma_{g}$ between $0$ and $1$. In this case,
ignoring logarithmic factors and other constants, the bound (3.18) can be re-
expressed as

$n^{\frac{\gamma_{b}-1}{4}}+n^{\frac{1-\gamma_{b}-4\gamma_{g}}{2}}+n^{\frac{1-\gamma_{b}(1+2\gamma_{\varepsilon})}{2}}\leavevmode\nobreak\
.$ (3.21)

In other words, the confidence region defined in 2.1 is consistent if

$1\leq\gamma_{b}+4\gamma_{g}\quad\text{and}\quad
1\leq\gamma_{b}(1+2\gamma_{\varepsilon})\leavevmode\nobreak\ ,$ (3.22)

respectively. That is, we are able to accommodate larger values of the
shrinkage rate $\varepsilon_{b}$ and nuisance parameter estimation error
$\delta_{n,g}$ if the subsample size $b$ is larger, relative to the sample
size $n$. However, as the subsample size $b$ increases, the normal
approximation error (i.e., the first term in (3.18)) increases.

Recall from Section 3.1 that, for many popular, honest, subsampled kernel
estimators, the shrinkage rate $\varepsilon_{b}$ satisfies a bound
$\varepsilon_{b}\lesssim b^{-c/p}$ for some small constant $c$ and some
integer $p$ measuring the (potentially, intrinsic) dimension of the covariates
$X_{i}$. Thus, in order to ensure that the second inequality in consistency
condition (3.22) is satisfied, it is essential to accommodate subsample sizes
$\gamma_{b}$ close to one, i.e., the regime $b=o(n)$. This is enabled by the
general results on the asymptotic linearity of $U$-statistics given in Section
4.

On the other hand, when the subsample size scaling factor $\gamma_{b}$ is
close to one, the restriction imposed by the consistency condition (3.22) on
the rate of convergence of the nuisance parameter estimator $\hat{g}_{n}$ is
very weak. Observe that the case $\gamma_{b}=0$ implies the familiar condition
that nuisance parameters can be estimated at the rate $n^{-1/4}$ in root mean
squared error (Chernozhukov et al.,, 2018). If the nuisance parameter
estimator $\hat{g}_{n}$ is itself estimated with random forest regression,
Syrgkanis and Zampetakis, (2020), Chi et al., (2022), and Huo et al., (2023),
among others, give conditions under which sufficient rates of convergence can
be achieved.

###  3.4. Performance

We now measure the performance of the confidence region formulated in 2.1. We
apply a method for simulation design proposed by Athey et al., (2021). In
particular, we calibrate a simulation to the Banerjee et al., (2015) data
using a Generative Adversarial Network (GAN) (Goodfellow et al.,, 2014).
Further details on this calibration are given in Appendix E. In effect, we
construct a data generating process that approximates the Banerjee et al.,
(2015) data, where we know the true value of the CATE $\theta_{0}(x)$ queried
at each value $x$ used to construct the grids displayed in Figures 1 and 2.

Measurements are taken as two parameters vary. First, we consider several
values of the sample size $n$. In particular, we consider settings with
$n=h\cdot n_{0}$, for $h$ in $\\{1,2.5,5,7.5,10\\}$, where $n_{0}$ is the
sample size of the Banerjee et al., (2015) data. Second, we vary the
proportion $b/n$. We consider three regimes: $b/n=0.05$, $b/n=(2/(h+1))0.05$,
and $b/n=(1/h)0.05$. Observe that $b$ increases in proportion to $n$ in the
first regime and that $b$ is constant as $n$ varies in the third
regime.161616We choose the same value of $b$ for constructing the empirical
moment $M_{n}(\cdot;\theta,g,\mathbf{D}_{n})$ and the nuisance parameter
estimator $\hat{g}_{n}$. The second regime resides between these two extremes.

Figure 3 displays measurements of the coverage and width of the confidence
region formulated in 2.1, in addition to measurements of the bias of the
estimator $\hat{\theta}_{n}(\bm{x}^{(d)})$. The first row displays
measurements of the coverage of the confidence region, the coverage of the
lower bound (i.e., Panel B of Figure 2), and the coverage of the upper bound
(i.e., Panel A of Figure 2). The nominal level is $\alpha=0.1$. At the
observed sample size $n_{0}=2438$, i.e., $h=1$, the confidence region is
somewhat anti-conservative. Consistent with the structure of the bound (3.18),
coverage does not improve as the sample size increases unless the proportion
$b/n$ also decreases. In the regimes where the proportion $b/n$ decreases as
the sample multiplier $h$ increases, coverage becomes moderately conservative.

Figure 3. Performance

---

Notes: Figure 3 displays several measurements of the performance of the
confidence intervals formulated in 2.1 in a simulation calibrated to the
Banerjee et al., (2015) data. The nominal level is $\alpha=0.1$; a horizontal
dotted line is displayed at the nominal coverage $1-\alpha$ in the first
panel. The confidence bounds considered are constructed analogously to the
confidence bounds displayed in Figure 2. The $x$-axis of each panel is the
sample multiplier $h$. The color of each measurement varies with the choice of
$b/n$. Further details on this design and implementation of this simulation
are given in Appendix E.

The second row of Figure 3 illustrates a bias-variance trade-off with the
subsample size $b$. The first panel displays measurements of the average width
of the confidence region. Here, the average is taken over both simulation
draws and the query-vector $\mathbf{x}^{(d)}$. The width of the confidence
region is increasing in the proportion $b/n$ and is essentially constant if
$b/n$ is constant as $n$ increases. By contrast, the second two panels display
the maximum and average bias of the estimator
$\hat{\theta}_{n}(\bm{x}^{(d)})$, again taken over the query-vector
$\mathbf{x}^{(d)}$. The bias is decreasing in the proportion $b/n$ and is
essentially constant if $b$ is constant as $n$ varies.171717The default value
for the proportion $b/n$ in the GRF R package is $0.5$. The results of this
simulation suggest that this proportion should be reduced, at least in
settings where the dimension of the covariate vector $X_{i}$ is small.

##  4\. General Results for High-Dimensional $U$-Statistics

An essential step in the proof of Theorem 3.1 follows from a new order-
explicit bound on the remainder in a linear approximation to a high-
dimensional $U$-statistic. In this section, we present this result and state
several corollaries. In particular, we give new order-explicit results on the
concentration and normal approximation of high-dimensional $U$-statistics.
That is, we consider the asymptotic behavior of the $b$ order $U$-statistic

$\overline{U}_{n,b}(\bm{x}^{(d)})=\frac{1}{N_{b}}\sum_{i\in\mathcal{S}_{n,b}}u(\bm{x}^{(d)};D_{\mathsf{s}})\leavevmode\nobreak\
,$ (4.1)

where the vector $u(\bm{x}^{(d)};\cdot)$ collects the deterministic,
symmetric, real-valued kernel function $u(x;\cdot),$ evaluated at the
$d$-vector of points $\bm{x}^{(d)}=(x^{(j)})_{j=1}^{d}$ in the space
$\mathcal{X}$. We assume that each component of the kernel function
$u(\bm{x}^{(d)};D_{\mathsf{s}})$ has mean zero. Proofs for results stated in
this section are given in Appendix C.

###  4.1. Context

The asymptotic analysis of $U$-statistics was initiated by Hoeffding, (1948),
who established a central limit theorem in the regime where the order $b$ is
fixed and the sample size $n$ is increasing. The Hoeffding central limit
theorem has been extended only recently to the regime where the order $b$
increases with the sample size $n$. DiCiccio and Romano, (2022) give a result
with this flavor in the regime where $b=o(n^{1/2})$. Wager and Athey, (2018),
Peng et al., (2022), and Minsker, (2023) are able to strengthen this result to
the regime where $b=o(n)$. We state and prove this more general result for the
sake of completeness, and because its main ideas will serve as useful touch
points in the more involved analysis to follow. We use an argument similar to
the proof of Theorem 3.1 of Minsker, (2023). Define the kernel variance
$\nu^{2}_{j}=\operatorname{Var}(u(x^{(j)};D_{[b]}))$, the Hájek projection

$u^{(1)}(x^{(j)};D)=\mathbb{E}\left[u(x^{(j)};D_{[b]})\mid
D_{1}=D\right]\leavevmode\nobreak\ ,$ (4.2)

and the Hájek projection variance
$\sigma^{2}_{b,j}=\operatorname{Var}(\tilde{u}^{(1)}(x^{(j)};D_{i}))$.

######  Theorem 4.1 (Hoeffding Central Limit Theorem).

For any sequence of kernel orders $b=b_{n}$, where

$\frac{1}{n}\frac{\nu^{2}_{j}}{\sigma^{2}_{b,j}}\to 0$ (4.3)

as $n\to 0$, we have that

$\sqrt{\frac{n}{\sigma^{2}_{b,j}b^{2}}}\frac{1}{N_{b}}\sum_{\mathsf{s}\in\mathsf{S}_{n,b}}u(x^{(j)};D_{\mathsf{s}})\overset{d}{\to}\mathsf{N}(0,1)\leavevmode\nobreak\
,$ (4.4)

as $n\to\infty$, where $\overset{d}{\to}$ denotes convergence in distribution.

######  Remark 4.1.

Theorem 4.1 is established by considering the decomposition

$\displaystyle\sqrt{\frac{1}{\sigma^{2}_{b,j}n}}\sum_{i=1}^{n}u^{(1)}(x^{(j)};D_{\mathsf{s}})-\sqrt{\frac{n}{\sigma^{2}_{b,j}b^{2}}}\left(\frac{1}{N_{b}}\sum_{\mathsf{s}\in\mathsf{S}_{n,b}}u(x^{(j)};D_{\mathsf{s}})-\frac{b}{n}\sum_{i=1}^{n}u^{(1)}(x^{(j)};D_{\mathsf{s}})\right)\leavevmode\nobreak\
.$ (4.5)

In particular, we give a high-probability bound for the second term and apply
a standard central limit theorem to the first term. The condition (4.3) is
needed to bound the second term. As a by-product of the proof, we show that
$b\sigma^{2}_{b,j}\leq\nu^{2}_{j}$. Thus, the normalization (4.3) implies that
$b=o(n)$. $\blacksquare$

Large deviation bounds for high-dimensional $U$-statistics, i.e.,
$U$-statistics with vector-valued kernels, were not given until Hoeffding,
(1963). This result is now more standard; a modern version is stated as
follows. A proof is given in Song et al., (2019). The norm
$\|\cdot\|_{\psi_{1}}$ denotes the $\psi_{1}$-Orlicz norm.181818Random
variables are sub-Exponential if and only if they have a finite
$\psi_{1}$-Orlicz norm (Section 2.7, Vershynin,, 2018).

######  Lemma 4.1 (Lemma A.5, Song et al., (2019)).

If the $\psi_{1}$-Orlicz norm bound

$\|u(x^{(j)};D_{\mathsf{s}})\|_{\psi_{1}}\leq\phi$ (4.6)

is satisfied for each $j$ in $[d]$, then

$\|\overline{U}_{n,b}(\bm{x}^{(d)})\|_{\infty}\lesssim\sqrt{\frac{b\overline{\nu}^{2}\log(dn)}{n}}+\frac{b\phi\log^{2}(dn)}{n}$
(4.7)

with probability greater than $1-C/n$, where
$\overline{\nu}^{2}=\max_{j\in[d]}\nu^{2}_{j}$.

Again, Lemma 4.1 demonstrates that $\overline{U}_{n,b}(\bm{x}^{(d)})$
concentrates in the regime that $b=o(n)$, up to a logarithmic factor that
depends on the dimension $d$. Here, however, concentration is expressed in
terms of the quantity $n^{-1}b\nu_{j}^{2}$, rather than the more appropriate,
and potentially substantively smaller, normalizing quantity
$n^{-1}b^{2}\sigma^{2}_{b,j}$ used in Theorem 4.1. In part motivated by this
incongruity, Arcones and Giné, (1993), Arcones, (1995), Giné et al., (2000),
establish a series of refined large deviation bounds for high-dimensional
$U$-statistics that use the appropriate normalizing factor (among many other
related results). See De la Pena and Giné, (1999) for a textbook treatment.
However, the constants used to express these bounds depend implicitly on the
order $b$, and so are not applicable to asymptotic regimes where $b$ may be
growing with the sample size $n$.

More recently, Chen, (2018), Chen and Kato, (2019), and Song et al., (2019)
have studied central limit theorems for high-dimensional $U$-statistics. Of
these papers, only Song et al., (2019) gives results with explicit dependence
on the order $b$. Their results are only applicable to the regime where
$b=o(n^{1/3})$. Minsker, (2023) gives a large deviation bound with the correct
normalizing factor and explicit order dependence, but this result is again
only applicable to the regime $b=o(n^{1/3})$.

###  4.2. Concentration of the Hájek Residual

We obtain a large deviation bound on the difference

$\overline{U}_{n,b}(\bm{x}^{(d)})-\frac{b}{n}\sum_{i=1}^{n}u^{(1)}(\bm{x}^{(d)};D_{i})\leavevmode\nobreak\
.$ (4.8)

We refer to the quantity (4.8) as the Hájek residual (Hájek,, 1968). This
bound is used in the proof of Theorem 3.1 and implies a new large deviation
bound and central limit theorem for high-dimensional $U$-statistics, stated in
the following subsection.

######  Theorem 4.2.

Define the terms

$\displaystyle\overline{\psi}^{2}_{b}=\max_{j\in[d]}\left\\{\nu^{2}_{j}-b\sigma^{2}_{b,j}\right\\}\quad\text{and}\quad\underline{\sigma}^{2}_{b}=\min_{j\in[d]}\sigma^{2}_{b,j}\leavevmode\nobreak\
.$ (4.9)

If the kernel function $u(x^{(j)};D_{\mathsf{s}})$ satisfies the bound (4.6)
for each $j$ in $[d]$, then

$\displaystyle\sqrt{\frac{n}{b^{2}\underline{\sigma}^{2}_{b}}}\Big{\|}\overline{U}_{n,b}(\bm{x}^{(d)})-\frac{b}{n}\sum_{i=1}^{n}u^{(1)}(\bm{x}^{(d)};D_{i})\Big{\|}_{\infty}\lesssim\xi_{n,b}\leavevmode\nobreak\
,\quad\text{where}$ (4.10)
$\displaystyle\xi_{n,b}=\left(\frac{Cb\log(dn)}{n}\right)^{b/2}\left(\left(\frac{n\overline{\psi}^{2}_{b}}{b^{2}\underline{\sigma}^{2}_{b}}\right)^{1/2}+\left(\frac{\phi^{2}b\log^{4}(dn)}{\underline{\sigma}^{2}_{b}}\right)^{1/2}\right)\leavevmode\nobreak\
,$

with probability greater than $1-C/n$.

######  Remark 4.2.

Roughly speaking, the bound (4.10) follows by first demonstrating that the
Hájek residual can be expressed as a degenerate $U$-statistic of order $b$.
This allows us to derive Hoffman-Jørgensen type bounds on higher moments of
(4.8) with a symmetrization argument. Here, we make essential use of a
symmetrization inequality for completely degenerate kernels, with explicit
dependence on the order, due to Sherman, (1994). This symmetrization
inequality was also used in Song et al., (2019) and Minsker, (2023).
$\blacksquare$

######  Remark 4.3.

The bound (4.10) implies that if $b^{-C}\leq\underline{\sigma}^{2}_{b}$ for
some positive constant $C$, then the Hájek residual is (roughly) of stochastic
order $(b/n)^{b/2}$ for all sufficiently large $b$. This is a dramatic
improvement over existing bounds (Song et al.,, 2019; Minsker,, 2023), which
decay polynomially as $b$ and $n$ increase and only converge to zero in the
regime $b=o(n^{1/3})$. In our view, this result hints at an explanation for
the widespread success of subsampling in machine learning and statistical
inference. Subsampled statistics, of a large order, are _essentially_ linear.
$\blacksquare$

###  4.3. Concentration and Normal Approximation

We now state a large deviation bound and central limit theorem for the high-
dimensional $U$-statistic (4.1). Both results are corollaries of Theorem 4.2,
apply to the regime $b=o(n)$, and depend on the correct normalizing factor.

######  Corollary 4.1.

Let $\Sigma$ be the diagonal matrix with components $\sigma_{b,j}^{2}$.

(i) Under the same conditions as Theorem 4.2, we have that

$\sqrt{\frac{n}{b^{2}}}\Sigma^{-1/2}\overline{U}_{n,b}(\bm{x}^{(d)})\lesssim\log^{1/2}(dn)+\frac{\phi\log^{2}(dn)}{\underline{\sigma}_{b}n^{1/2}}+\xi_{n,b}$
(4.11)

with probability greater than $1-C/n$.

(ii) Let $Z$ denote a centered Gaussian random vector with covariance matrix
$\operatorname{Var}(\tilde{u}^{(1)}(\bm{x}^{(d)},D_{i}))$. Under the same
conditions as Theorem 4.2, we have that

$\displaystyle\sup_{\mathsf{R}\in\mathcal{R}}\bigg{|}P\left\\{\sqrt{\frac{n}{b^{2}}}\Sigma^{-1/2}\overline{U}_{n,b}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{\Sigma^{-1/2}Z\in\mathsf{R}\right\\}\bigg{|}$
$\displaystyle\quad\quad\lesssim\left(\frac{\phi^{2}\log^{5}(dn)}{\underline{\sigma}^{2}_{b}n}\right)^{1/4}+\xi_{n,b}\sqrt{\log(d)}\leavevmode\nobreak\
,$ (4.12)

where $\mathcal{R}$ denotes the set of hyper-rectangles in $\mathbb{R}^{d}$.

######  Remark 4.4.

De la Pena and Giné, (1999) state a result analogous to Part (i) of 4.1, in
the sense that the $U$-statistic $\overline{U}_{n,b}(\bm{x}^{(d)})$ is
normalized by the correct quantity $n^{1/2}b^{-2}\Sigma^{-1/2}$. The constants
used in their bound depend implicitly on $b$. Part (i) of 4.1 improves
substantially on Lemma 4.1 in contexts where the Hájek projection variances
$\sigma_{b,j}^{2}$ are smaller than $b^{-1}$.191919Section 4 of Song et al.,
(2019) gives several examples of statistics where the Hájek projection
variance $\sigma_{b,j}^{2}$ is smaller than $b^{-1}$. Under the restriction
$\sigma_{b,j}^{2}\gtrsim b^{-1}$, imposed in the application to subsampled
kernel regression considered in Section 3, Part (i) of 4.1 only improves on
Lemma 4.1 by a constant factor. Part (ii) of 4.1 gives a high-dimensional
equivalent to Theorem 4.1. An analogous half-sample bootstrap central limit
theorem follows from arguments very similar to parts of the proof of Theorem
3.1. $\blacksquare$

##  5\. Conclusion

We propose a confidence region for solutions to conditional moment equations.
The confidence region is built around an estimator based on subsampled kernel
regression. As a running example, we consider the construction of confidence
regions for conditional average treatment effects around a Generalized Random
Forest (Athey et al.,, 2019). Empirically, we document that the proposed
confidence region is able to recover treatment effect heterogeneity undetected
by interacted linear regression. Theoretically, we establish a bound on
coverage accuracy that illustrates a bias-variance tradeoff in the user-chosen
subsample size. In order to do this, we obtain several new results on the
asymptotics of high-dimensional $U$-statistics.

We give conditions sufficient for the asymptotic validity of the proposed
confidence regions. However, the confidence region is not necessarily optimal,
in any particular sense, under the maintained assumptions. It is likely to be
the case an optimal confidence region would need to incorporate a bias
estimate and procedure for choosing tuning parameters to balance bias and
variance (see e.g., Chernozhukov et al., (2014), Chen et al., (2024) for
instantiations of these ideas). Adapting this approach to subsampled moment
regression is an interesting direction for further research.

## References

  * Adler and Taylor, (2009) Adler, R. J. and Taylor, J. E. (2009).  Random fields and geometry.  Springer Science & Business Media.
  * Ai and Chen, (2003) Ai, C. and Chen, X. (2003).  Efficient estimation of models with conditional moment restrictions containing unknown functions.  Econometrica, 71(6):1795–1843.
  * Angrist et al., (1996) Angrist, J. D., Imbens, G. W., and Rubin, D. B. (1996).  Identification of causal effects using instrumental variables.  Journal of the American statistical Association, 91(434):444–455.
  * Arcones, (1995) Arcones, M. A. (1995).  A bernstein-type inequality for u-statistics and u-processes.  Statistics & probability letters, 22(3):239–247.
  * Arcones and Giné, (1993) Arcones, M. A. and Giné, E. (1993).  Limit theorems for u-processes.  The Annals of Probability, pages 1494–1542.
  * Armstrong and Kolesár, (2020) Armstrong, T. B. and Kolesár, M. (2020).  Simple and honest confidence intervals in nonparametric regression.  Quantitative Economics, 11(1):1–39.
  * Athey et al., (2020) Athey, S., Chetty, R., Imbens, G., and Kang, H. (2020).  Estimating treatment effects using multiple surrogates: The role of the surrogate score and the surrogate index.
  * Athey and Imbens, (2016) Athey, S. and Imbens, G. (2016).  Recursive partitioning for heterogeneous causal effects.  Proceedings of the National Academy of Sciences, 113(27):7353–7360.
  * Athey et al., (2021) Athey, S., Imbens, G. W., Metzger, J., and Munro, E. (2021).  Using wasserstein generative adversarial networks for the design of monte carlo simulations.  Journal of Econometrics.
  * Athey et al., (2019) Athey, S., Tibshirani, J., and Wager, S. (2019).  Generalized random forests.  The Annals of Statistics, 47(2):1148–1178.
  * Balboni et al., (2022) Balboni, C., Bandiera, O., Burgess, R., Ghatak, M., and Heil, A. (2022).  Why do people stay poor?  The Quarterly Journal of Economics, 137(2):785–844.
  * Banerjee et al., (2015) Banerjee, A., Duflo, E., Goldberg, N., Karlan, D., Osei, R., Parienté, W., Shapiro, J., Thuysbaert, B., and Udry, C. (2015).  A multifaceted program causes lasting progress for the very poor: Evidence from six countries.  Science, 348(6236):1260799.
  * Belloni et al., (2018) Belloni, A., Chernozhukov, V., Chetverikov, D., and Wei, Y. (2018).  Uniformly valid post-regularization confidence regions for many functional parameters in z-estimation framework.  Annals of statistics, 46(6B):3643.
  * Biau et al., (2008) Biau, G., Devroye, L., and Lugosi, G. (2008).  Consistency of random forests and other averaging classifiers.  Journal of Machine Learning Research, 9(9).
  * Boucheron et al., (2013) Boucheron, S., Lugosi, G., and Massart, P. (2013).  Concentration inequalities: A nonasymptotic theory of independence,(2013).
  * Breiman, (1996) Breiman, L. (1996).  Bagging predictors.  Machine learning, 24:123–140.
  * Breiman, (2001) Breiman, L. (2001).  Random forests.  Machine learning, 45:5–32.
  * Breiman et al., (1984) Breiman, L., Friedman, J., Olshen, R., and Stone, C. J. (1984).  Classification and regression trees.  Routledge.
  * Bühlmann and Yu, (2002) Bühlmann, P. and Yu, B. (2002).  Analyzing bagging.  The annals of Statistics, 30(4):927–961.
  * Cattaneo et al., (2024) Cattaneo, M. D., Klusowski, J. M., and Tian, P. M. (2024).  On the pointwise behavior of recursive partitioning and its implications for heterogeneous causal effect estimation.
  * Chen and Ritzwoller, (2023) Chen, J. and Ritzwoller, D. M. (2023).  Semiparametric estimation of long-term treatment effects.  Journal of Econometrics, 237(2):105545.
  * Chen, (2018) Chen, X. (2018).  Gaussian and bootstrap approximations for high-dimensional u-statistics and their applications.
  * Chen et al., (2024) Chen, X., Christensen, T., and Kankanala, S. (2024).  Adaptive estimation and uniform confidence bands for nonparametric structural functions and elasticities.  The Review of Economic Studies.
  * Chen and Christensen, (2018) Chen, X. and Christensen, T. M. (2018).  Optimal sup-norm rates and uniform inference on nonlinear functionals of nonparametric iv regression: Nonlinear functionals of nonparametric iv.  Quantitative Economics, 9(1):39–84.
  * Chen and Kato, (2019) Chen, X. and Kato, K. (2019).  Randomized incomplete $U$-statistics in high dimensions.  The Annals of Statistics, 47(6):3127 – 3156.
  * Chen and Pouzo, (2012) Chen, X. and Pouzo, D. (2012).  Estimation of nonparametric conditional moment models with possibly nonsmooth generalized residuals.  Econometrica, 80(1):277–321.
  * Chernozhukov et al., (2018) Chernozhukov, V., Chetverikov, D., Demirer, M., Duflo, E., Hansen, C., Newey, W., and Robins, J. (2018).  Double/debiased machine learning for treatment and structural parameters: Double/debiased machine learning.  The Econometrics Journal, 21(1).
  * Chernozhukov et al., (2013) Chernozhukov, V., Chetverikov, D., and Kato, K. (2013).  Gaussian approximations and multiplier bootstrap for maxima of sums of high-dimensional random vectors.
  * Chernozhukov et al., (2014) Chernozhukov, V., Chetverikov, D., and Kato, K. (2014).  Anti-concentration and honest, adaptive confidence bands.  The Annals of Statistics, 42(5):1787–1818.
  * (30) Chernozhukov, V., Chetverikov, D., and Kato, K. (2017a).  Central limit theorems and bootstrap in high dimensions.  The Annals of Probability, 45(4):2309 – 2352.
  * (31) Chernozhukov, V., Chetverikov, D., and Kato, K. (2017b).  Detailed proof of nazarov’s inequality.  arXiv preprint arXiv:1711.10696.
  * Chernozhukov et al., (2022) Chernozhukov, V., Escanciano, J. C., Ichimura, H., Newey, W. K., and Robins, J. M. (2022).  Locally robust semiparametric estimation.  Econometrica, 90(4):1501–1535.
  * Chernozhukov et al., (2023) Chernozhukov, V., Newey, W. K., and Santos, A. (2023).  Constrained conditional moment restriction models.  Econometrica, 91(2):709–736.
  * Chernozhuokov et al., (2022) Chernozhuokov, V., Chetverikov, D., Kato, K., and Koike, Y. (2022).  Improved central limit theorem and bootstrap approximations in high dimensions.  The Annals of Statistics, 50(5):2562–2586.
  * Chi et al., (2022) Chi, C.-M., Vossler, P., Fan, Y., and Lv, J. (2022).  Asymptotic properties of high-dimensional random forests.  The Annals of Statistics, 50(6):3415–3438.
  * Chung and Romano, (2013) Chung, E. and Romano, J. P. (2013).  Exact and asymptotically robust permutation tests.  The Annals of Statistics, 41(2):484–507.
  * De la Pena and Giné, (1999) De la Pena, V. and Giné, E. (1999).  Decoupling: from dependence to independence.  Springer.
  * DiCiccio and Romano, (2022) DiCiccio, C. and Romano, J. (2022).  Clt for u-statistics with growing dimension.  Statistica Sinica, 32(1).
  * Efron and Stein, (1981) Efron, B. and Stein, C. (1981).  The jackknife estimate of variance.  The Annals of Statistics, pages 586–596.
  * Fix and Hodges, (1989) Fix, E. and Hodges, J. L. (1989).  Discriminatory analysis. nonparametric discrimination: Consistency properties.  International Statistical Review/Revue Internationale de Statistique, 57(3):238–247.
  * Foster and Syrgkanis, (2023) Foster, D. J. and Syrgkanis, V. (2023).  Orthogonal statistical learning.  The Annals of Statistics, 51(3):879–908.
  * Frölich, (2007) Frölich, M. (2007).  Nonparametric iv estimation of local average treatment effects with covariates.  Journal of Econometrics, 139(1):35–75.
  * Giné et al., (2000) Giné, E., Latała, R., and Zinn, J. (2000).  Exponential and moment inequalities for u-statistics.  In High Dimensional Probability II, pages 13–38. Springer.
  * Goodfellow et al., (2014) Goodfellow, I., Pouget-Abadie, J., Mirza, M., Xu, B., Warde-Farley, D., Ozair, S., Courville, A., and Bengio, Y. (2014).  Generative adversarial nets.  Advances in neural information processing systems, 27.
  * Götze et al., (2021) Götze, F., Sambale, H., and Sinulis, A. (2021).  Concentration inequalities for polynomials in $\alpha$-sub-exponential random variables.  Electronic Journal of Probability, 26(none):1 – 22.
  * Hahn, (1998) Hahn, J. (1998).  On the role of the propensity score in efficient semiparametric estimation of average treatment effects.  Econometrica, pages 315–331.
  * Hájek, (1968) Hájek, J. (1968).  Asymptotic normality of simple linear rank statistics under alternatives.  The Annals of Mathematical Statistics, pages 325–346.
  * Hanson and Wright, (1971) Hanson, D. L. and Wright, F. T. (1971).  A bound on tail probabilities for quadratic forms in independent random variables.  The Annals of Mathematical Statistics, 42(3):1079–1083.
  * Hoeffding, (1948) Hoeffding, W. (1948).  A class of statistics with asymptotically normal distribution.  Ann. Math. Statist., 19(4):293–325.
  * Hoeffding, (1963) Hoeffding, W. (1963).  Probability inequalities for sums of bounded random variables.  Journal of the American Statistical Association, pages 13–30.
  * Huo et al., (2023) Huo, Y., Fan, Y., and Han, F. (2023).  On the adaptation of causal forests to manifold data.  arXiv preprint arXiv:2311.16486.
  * Ichimura and Newey, (2022) Ichimura, H. and Newey, W. K. (2022).  The influence function of semiparametric estimators.  Quantitative Economics, 13(1):29–61.
  * Imbens and Rubin, (2015) Imbens, G. W. and Rubin, D. B. (2015).  Causal inference in statistics, social, and biomedical sciences.  Cambridge University Press.
  * Kennedy, (2023) Kennedy, E. H. (2023).  Towards optimal doubly robust estimation of heterogeneous causal effects.  Electronic Journal of Statistics, 17(2):3008–3049.
  * Khosravi et al., (2019) Khosravi, K., Lewis, G., and Syrgkanis, V. (2019).  Non-parametric inference adaptive to intrinsic dimension.  arXiv preprint arXiv:1901.03719.
  * Kpotufe, (2011) Kpotufe, S. (2011).  k-nn regression adapts to local intrinsic dimension.  Advances in neural information processing systems, 24.
  * Kraay and McKenzie, (2014) Kraay, A. and McKenzie, D. (2014).  Do poverty traps exist? assessing the evidence.  Journal of Economic Perspectives, 28(3):127–148.
  * Kuchibhotla et al., (2023) Kuchibhotla, A. K., Balakrishnan, S., and Wasserman, L. (2023).  Median regularity and honest inference.  Biometrika, 110(3):831–838.
  * Lee, (1990) Lee, A. (1990).  U-statistics: theory and practice.  Taylor & Francis.
  * Lewis and Syrgkanis, (2021) Lewis, G. and Syrgkanis, V. (2021).  Double/debiased machine learning for dynamic treatment effects.  In NeurIPS, pages 22695–22707.
  * Li, (1989) Li, K.-C. (1989).  Honest confidence regions for nonparametric regression.  The Annals of Statistics, 17(3):1001–1008.
  * Lin and Jeon, (2006) Lin, Y. and Jeon, Y. (2006).  Random forests and adaptive nearest neighbors.  Journal of the American Statistical Association, 101(474):578–590.
  * Mentch and Hooker, (2014) Mentch, L. and Hooker, G. (2014).  Ensemble trees and clts: Statistical inference for supervised learning.  stat, 1050:25.
  * Mentch and Hooker, (2016) Mentch, L. and Hooker, G. (2016).  Quantifying uncertainty in random forests via confidence intervals and hypothesis tests.  The Journal of Machine Learning Research, 17(1):841–881.
  * Minsker, (2023) Minsker, S. (2023).  U-statistics of growing order and sub-gaussian mean estimators with sharp constants.  Mathematical Statistics and Learning.
  * Montgomery-Smith, (1993) Montgomery-Smith, S. J. (1993).  Comparison of sums of independent identically distributed random variables.  Probability and Mathematical Statistics, 14:281–285.
  * Newey, (1993) Newey, W. K. (1993).  Efficient estimation of models with conditional moment restrictions.  In Econometrics, volume 11 of Handbook of Statistics, pages 419–454. Elsevier.
  * Newey, (1994) Newey, W. K. (1994).  The asymptotic variance of semiparametric estimators.  Econometrica: Journal of the Econometric Society, pages 1349–1382.
  * Newey and McFadden, (1994) Newey, W. K. and McFadden, D. (1994).  Large sample estimation and hypothesis testing.  Handbook of econometrics, 4:2111–2245.
  * Nie and Wager, (2021) Nie, X. and Wager, S. (2021).  Quasi-oracle estimation of heterogeneous treatment effects.  Biometrika, 108(2):299–319.
  * Oprescu et al., (2019) Oprescu, M., Syrgkanis, V., and Wu, Z. S. (2019).  Orthogonal random forest for causal inference.  In International Conference on Machine Learning, pages 4932–4941. PMLR.
  * Peng et al., (2022) Peng, W., Coleman, T., and Mentch, L. (2022).  Rates of convergence for random forests via generalized u-statistics.  Electronic Journal of Statistics, 16(1):232–292.
  * Politis and Romano, (1994) Politis, D. N. and Romano, J. P. (1994).  Large sample confidence regions based on subsamples under minimal assumptions.  The Annals of Statistics, pages 2031–2050.
  * Politis et al., (1999) Politis, D. N., Romano, J. P., and Wolf, M. (1999).  Subsampling.  Springer Science & Business Media.
  * Præstgaard and Wellner, (1993) Præstgaard, J. and Wellner, J. A. (1993).  Exchangeably weighted bootstraps of the general empirical process.  The Annals of Probability, pages 2053–2086.
  * Ritzwoller and Romano, (2023) Ritzwoller, D. M. and Romano, J. P. (2023).  Reproducible aggregation of sample-split statistics.  arXiv preprint arXiv:2311.14204.
  * Robins et al., (1994) Robins, J. M., Rotnitzky, A., and Zhao, L. P. (1994).  Estimation of regression coefficients when some regressors are not always observed.  Journal of the American statistical Association, 89(427):846–866.
  * Romano and Shaikh, (2012) Romano, J. P. and Shaikh, A. M. (2012).  On the uniform asymptotic validity of subsampling and the bootstrap.  The Annals of Statistics, 40(6).
  * Scornet et al., (2015) Scornet, E., Biau, G., and Vert, J.-P. (2015).  Consistency of random forests.  Annals of Statistics, 43(4):1716–1741.
  * Sexton and Laake, (2009) Sexton, J. and Laake, P. (2009).  Standard errors for bagged and random forest estimators.  Computational Statistics & Data Analysis, 53(3):801–811.
  * Sherman, (1994) Sherman, R. P. (1994).  Maximal inequalities for degenerate u-processes with applications to optimization estimators.  The Annals of Statistics, 22(1):439–459.
  * Singh and Vijaykumar, (2023) Singh, R. and Vijaykumar, S. (2023).  Kernel ridge regression inference.  arXiv preprint arXiv:2302.06578.
  * Song et al., (2019) Song, Y., Chen, X., and Kato, K. (2019).  Approximating high-dimensional infinite-order $u$-statistics: Statistical and computational guarantees.  Electronic Journal of Statistics, 13(2).
  * Syrgkanis and Zampetakis, (2020) Syrgkanis, V. and Zampetakis, M. (2020).  Estimation and inference with trees and forests in high dimensions.  In Conference on learning theory, pages 3453–3454. PMLR.
  * Tan, (2006) Tan, Z. (2006).  Regression and weighting methods for causal inference using instrumental variables.  Journal of the American Statistical Association, 101(476):1607–1618.
  * van der Vaart and Wellner, (2013) van der Vaart, A. W. and Wellner, J. A. (2013).  Weak convergence and empirical processes: with applications to statistics.  Springer.
  * Vershynin, (2018) Vershynin, R. (2018).  High-dimensional probability: An introduction with applications in data science, volume 47.  Cambridge university press.
  * Wager and Athey, (2018) Wager, S. and Athey, S. (2018).  Estimation and inference of heterogeneous treatment effects using random forests.  Journal of the American Statistical Association, 113(523):1228–1242.
  * Wager et al., (2014) Wager, S., Hastie, T., and Efron, B. (2014).  Confidence intervals for random forests: The jackknife and the infinitesimal jackknife.  The Journal of Machine Learning Research, 15(1):1625–1651.
  * Yadlowsky et al., (2023) Yadlowsky, S., Fleming, S., Shah, N., Brunskill, E., and Wager, S. (2023).  Evaluating treatment prioritization rules via rank-weighted average treatment effects.

Supplemental Appendix to: Uniform Inference for Subsampled Moment
Regression111Date:
David M. Ritzwoller Vasilis Syrgkanis Stanford University Stanford University

Contents
appendix.Asubsection.A.1subsection.A.2subsection.A.3subsection.A.4subsubsection.A.4.1subsubsection.A.4.2subsubsection.A.4.3appendix.Bsubsection.B.1subsubsection.B.1.1subsubsection.B.1.2subsubsection.B.1.3subsubsection.B.1.4appendix.Csubsection.C.1subsection.C.2subsection.C.3subsection.C.4appendix.Dsubsection.D.1subsection.D.2subsubsection.D.2.1subsubsection.D.2.2appendix.Esubsection.E.1subsection.E.2subsection.E.3

##  Appendix A An Abstract Bound on Coverage Error

In this appendix, we give an abstract bound on the accuracy of the nominal
coverage probability for the confidence region introduced in 2.1. That is, we
make no use of the kernel structure expressed in (1.5). We specify a set of
high-level assumptions in Section A.1. The abstract bound is stated in Section
A.2 and proved in Section A.3. Proofs for supporting Lemmas are given in
Section A.4.

###  A.1. Assumptions

We require several mild smoothness restrictions on the moment function
$M(\cdot;\theta,g)$. In contrast to the set of assumptions specified in
Section 3.2, we do not require moment linearity. Instead, we impose the
following generalization of Part (iii) of Assumption 3.3.

######  Assumption A.1 (Moment Restrictions).

The moment function $M(x;\theta_{0},g_{0})$ is twice continuously
differentiable in its second argument. Let

$M^{(1)}(x;\theta,g)=\frac{\partial}{\partial\theta^{\prime}}M(x;\theta^{\prime},g_{0})|_{\theta^{\prime}=\theta}\quad\text{and}\quad
H(x;\theta,g)=\frac{\partial^{2}}{\partial^{2}\theta^{\prime}}M(x;\theta^{\prime},g_{0})|_{\theta^{\prime}=\theta}$
(A.1)

denote the Jacobian and Hessian of $M(\cdot;\theta_{0},g_{0})$ in $\theta$,
respectively. The Jacobian $M^{(1)}(x;\theta,g)$ is uniformly Lipschitz in its
second argument and bounded from below in the sense that

$\displaystyle\sup_{P\in\mathbf{P}}\sup_{j\in[d]}\big{|}M^{(1)}(x^{(j)};g)-M^{(1)}(x^{(j)};g_{0})\big{|}$
$\displaystyle\lesssim\|g-g_{0}\|_{2,\infty}\quad\text{and}$ (A.2)
$\displaystyle\inf_{P\in\mathbf{P}}\inf_{j\in[d]}\big{|}M^{(1)}(x^{(j)};g)\big{|}$
$\displaystyle\geq c$ (A.3)

for each $g$ in $\mathcal{G}$ and some positive constant $c$. The Hessian
$H(x;\theta,g)$ is uniformly bounded as $x$, $\theta$, and $g$ vary over their
respective domains.

We additionally require Neyman orthogonality and second order smoothness,
i.e., Part (i) of Assumption 3.3.

Moreover, we impose an analogous smoothness restriction on the centered
empirical moment

$\overline{M}_{n}(x;\theta,g)=M_{n}(x;\theta,g,\mathbf{D}_{n})-\mathbb{E}\left[M_{n}(x;\theta,g,\mathbf{D}_{n})\right]\leavevmode\nobreak\
,$ (A.4)

where we have made the dependence on $\mathbf{D}_{n}$ implicit to ease
notation.

######  Assumption A.2 (Empirical Smoothness).

The centered empirical moment $\overline{M}_{n}(x;\theta,g)$ is twice
continuously differentiable in its second argument. Let

$\displaystyle\frac{\partial}{\partial\theta^{\prime}}\overline{M}_{n}(x;\theta^{\prime},g_{0})|_{\theta^{\prime}=\theta}=\overline{M}^{(1)}_{n}(x;\theta,g)\quad\text{and}\quad\frac{\partial^{2}}{\partial^{2}\theta^{\prime}}\overline{M}_{n}(x;\theta^{\prime},g_{0})|_{\theta^{\prime}=\theta}=\overline{H}_{n}(x;\theta,g)$
(A.5)

denote the Jacobian and Hessian of $\overline{M}_{n}(\cdot;\theta,g)$ in
$\theta$, respectively. The Hessian $H_{n}(x;\theta,g)$ is uniformly bounded
almost surely as $x$, $\theta$, and $g$ vary over their respective domains.

Next, we impose a set of high-level restrictions on the structure of the
empirical conditional moment (2.1). At times we refer to the normalized
statistic

$U_{n}(x)=-(M^{(1)}(x;g_{0}))^{-1}\overline{M}_{n}(x;\theta_{0},g_{0})\leavevmode\nobreak\
.$ (A.6)

First, we impose a condition that ensures that (A.6) is approximately linear.
Recall that the norm $\|\cdot\|_{\psi_{1}}$ denotes the $\psi_{1}$-Orlicz
norm.

######  Assumption A.3 (Approximate Linearity).

There exists a function $\overline{u}(\cdot,\cdot)$, a constant $\varphi\geq
1$, and real-valued sequences $\delta_{n,u}$ and $\rho_{n,u}$ such that
$\mathbb{E}\left[\overline{u}(x^{(j)},D_{i})\right]=0$,

$\|\overline{u}(x^{(j)},D_{i})\|_{\psi_{1}}\leq\varphi\leavevmode\nobreak\
,\quad\text{and}\quad\mathbb{E}\left[\overline{u}^{4}(x^{(j)},D_{i})\right]\leq\operatorname{Var}(\overline{u}(x^{(j)},D_{i}))\varphi^{2}$
(A.7)

hold for all $j$ in $[d]$ and $P$ in $\mathbf{P}$. Moreover, if
$\lambda^{2}_{j}$ denotes $\operatorname{Var}(\overline{u}(x^{(j)},D_{i}))$
and $\underline{\lambda}^{2}=\min_{j\in[d]}\lambda^{2}_{j}$, then

$\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\big{\|}U_{n}(\bm{x}^{(d)})-\frac{1}{n}\sum_{i=1}^{n}\overline{u}(\bm{x}^{(d)},D_{i})\big{\|}_{\infty}\geq\delta_{n,u}\right\\}\leq\rho_{n,u}\leavevmode\nobreak\
.$ (A.8)

######  Remark A.1.

Condition (A.8) enables the application of the representation (2.7). The
bounds (A.7) hold, for example, if the quantities
$\overline{u}_{n}(x^{(j)},D_{i})$ are bounded by $\varphi$. $\blacksquare$

Second, we impose several restrictions relating to the estimators
$\hat{\theta}_{n}(\bm{x}^{(d)})$ and $\hat{g}_{n}$. Throughout, we measure the
error in the estimator of $\hat{\theta}_{n}(\bm{x}^{(d)})$ in terms of the
norm

$\|\hat{\theta}_{n}(\bm{x}^{(d)})-\theta_{0}(\bm{x}^{(d)})\|_{\infty}=\sup_{j\in[d]}|\hat{\theta}_{n}(x^{(j)})-\theta_{0}(x^{(j)})|\leavevmode\nobreak\
.$ (A.9)

A closely related collection of conditions, in the context of estimation of
the solution to unconditional orthogonal moments, is stated as Assumption 3.2
of Chernozhukov et al., (2018).

######  Assumption A.4 (Bias, Consistency, and Stochastic Equicontinuity).

Recall the definition of the object $\underline{\lambda}^{2}$ introduced in
Assumption A.3.
(i) Define the quantity

$\mathsf{Bias}_{n}(x;\theta,g)=M(x;\theta,g)-\mathbb{E}\left[M_{n}(x;\theta,g,\bm{D}_{n})\right]\leavevmode\nobreak\
.$ (A.10)

There exists a sequence $\delta_{n,B}$ such that

$\sup_{P\in\mathbf{P}}\sup_{g\in\mathcal{G}}\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\mathsf{Bias}_{n}(\bm{x}^{(d)};\theta(\bm{x}^{(d)}),g)\|_{\infty}\lesssim(1+\|\theta(\bm{x}^{(d)})\|_{\infty})\delta_{n,B}$
(A.11)

uniformly over any vector
$\theta(\bm{x}^{(d)})=\\{\theta(x^{(j)})\\}_{j=1}^{d}$.

(ii) There exist sequences $\delta_{n,m}$, $\delta_{n,g}$,
$\delta_{n,\theta}$, $\rho_{n,m}$, $\rho_{n,g}$, and $\rho_{n,\theta}$ such
that

$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\big{\|}\overline{M}^{(1)}_{n}(\bm{x}^{(d)};\theta_{0}(\bm{x}^{(d)}),g_{0})\big{\|}^{2}_{\infty}\geq\delta_{n,m}\right\\}$
$\displaystyle\leq\rho_{n,m}\leavevmode\nobreak\ ,$ (A.12)
$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\hat{g}_{n}-g_{0}\|_{2,\infty}^{2}\geq\delta^{2}_{n,g}\right\\}\leq\rho_{n,g}\leavevmode\nobreak\
,\quad\text{and}$ (A.13)
$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\hat{\theta}_{n}(\bm{x}^{(d)})-\theta_{0}(\bm{x}^{(d)})\|_{\infty}^{2}\geq\delta^{2}_{n,\theta}\right\\}\leq\rho_{n,\theta}\leavevmode\nobreak\
.$ (A.14)

(iii) There exist sequences $\delta_{n,S}$, $\delta_{n,J}$, $\rho_{n,S}$, and
$\rho_{n,J}$ such that

$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\big{\|}\overline{M}_{n}(\bm{x}^{(d)};\theta_{0}(\bm{x}^{(d)}),\hat{g}_{n})-\overline{M}_{n}(\bm{x}^{(d)};\theta_{0}(\bm{x}^{(d)}),g_{0})\big{\|}_{\infty}\geq\delta_{n,S}\right\\}\leq\rho_{n,S}\leavevmode\nobreak\
.$ (A.15)

and

$\displaystyle\sup_{P\in\mathbf{P}}P\left\\{\sqrt{\frac{n}{\underline{\lambda}^{2}}}\big{\|}\overline{M}^{(1)}_{n}(\bm{x}^{(d)};\theta_{0}(\bm{x}^{(d)}),\hat{g}_{n})-\overline{M}^{(1)}_{n}(\bm{x}^{(d)};\theta_{0}(\bm{x}^{(d)}),g_{0})\big{\|}_{\infty}\geq\delta_{n,J}\right\\}$
$\displaystyle\leq\rho_{n,J}\leavevmode\nobreak\ .$ (A.16)

respectively.

###  A.2. Coverage

The following theorem gives a non-asymptotic bound on the error in the nominal
coverage probability of the confidence regions introduced in 2.1.

######  Theorem A.1 (Generic Coverage Error Decomposition).

Collect the error sequences

$\displaystyle\delta_{n}$
$\displaystyle=\delta^{2}_{n,g}+\delta^{2}_{n,\theta}+\delta_{n,B}+\delta_{n,S}+\delta_{n,u}+\delta_{n,\theta}\left(\delta_{n,m}+\delta_{n,g}+\underline{\lambda}^{1/2}n^{-1/4}\left(\delta_{n,B}+\delta_{n,J}\right)\right)\quad\text{and}$
$\displaystyle\rho_{n}$
$\displaystyle=\rho_{n,m}+\rho_{n,g}+\rho_{n,\theta}+\rho_{n,S}+\rho_{n,J}+\rho_{n,u}$

and assume that $\delta_{\varepsilon n}\leq C_{\varepsilon}\delta_{n}$ and
$\rho_{\varepsilon n}\leq C_{\varepsilon}\rho_{n}$ for any $0<\varepsilon<1$.
Suppose that the Neyman orthogonal moment function $M(x;\theta_{0},g_{0})$
satisfies Assumption A.1 and Part (i) of Assumption 3.3 and that the centered
empirical moment function $\overline{M}_{n}(x;\theta_{0},g_{0})$ satisfies
Assumption A.2. If Assumptions A.3 and A.4 hold, then the confidence region
defined in 2.1 satisfies

$\sup_{P\in\mathbf{P}}\big{|}P\left\\{\theta_{0}(\bm{x}^{(d)})\in\hat{\mathcal{C}}(\bm{x}^{(d)})\right\\}-(1-\alpha)\big{|}\lesssim\left(\frac{\varphi^{2}\log^{5}\left(dn\right)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}+\rho_{n}\leavevmode\nobreak\
.$ (A.17)

######  Remark A.2.

Theorem A.1 is verified by considering the decomposition

$\sqrt{n}R_{n}(\bm{x}^{(d)})=\sqrt{n}\left(R_{n}(\bm{x}^{(d)})-\frac{1}{n}\sum_{i=1}^{n}\overline{u}(\bm{x}^{(d)},D_{i})\right)+\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\overline{u}(\bm{x}^{(d)},D_{i})\leavevmode\nobreak\
.$ (A.18)

First, through a standard series of expansions (see e.g., Chernozhukov et
al.,, 2018), we show that the first term in (A.18) is bounded above by
$\delta_{n}$ with probability greater than $1-\rho_{n}$. An analogous bound
holds for the bootstrap root $R^{*}_{n}(\bm{x}^{(d)})$, as the half-sample
bootstrap root is obtained with subsampling. These bounds produce the latter
terms on the right-hand side of (A.17).

The first term on the right-hand side of (A.17) is obtained through the
application of suitable of high-dimensional central limit theorems for the
second term in (A.18) (Chernozhuokov et al.,, 2022). In the case of the
bootstrap root, the central limit theorem that we apply leverages the
representations (2.7). In particular, we apply an argument based on a coupling
proposed in Yadlowsky et al., (2023). Basically, the weights in the sums (2.7)
are coupled with a sequence of independent and identically distributed weights
and the difference between the sums using the two types of weights are
bounded. This bound uses a Lévy-type generalization of an appropriate
Bernstein-type bound. $\blacksquare$

###  A.3. Proof of Theorem A.1

Throughout, we let $\mathcal{R}$ denote the set of hyper-rectangles in
$\mathbb{R}^{d}$. Fix a rectangle $\mathsf{R}=[a_{l},a_{u}]$ in $\mathcal{R}$,
where $a_{l}$ and $a_{u}$ are vectors in $\mathbb{R}^{d}$ with $a_{l}\leq
a_{u}$, interpreted componentwise. For each $t>0$, define the enlarged
rectangle $\mathsf{R}_{t}=[a_{l}-t\bm{1}_{d},a_{u}+t\bm{1}_{d}]$. Let $Z$
denote a centered Gaussian random vector with covariance matrix
$\operatorname{Var}(\overline{u}(\bm{x}^{(d)},D_{i}))$. Let $\Lambda$ be the
diagonal matrix with components
$\lambda_{j}^{2}=\operatorname{Var}(\overline{u}(x^{(j)},D_{i}))$. We have
that

$\displaystyle\big{|}P\left\\{\sqrt{n}\hat{\Lambda}_{n}^{-1/2}R_{n}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}\big{|}$
$\displaystyle\quad\quad\leq\big{|}P\left\\{\sqrt{n}\Lambda^{-1/2}R_{n}(\bm{x}^{(d)})\in\mathsf{R}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}\big{|}$
(A.19)
$\displaystyle\quad\quad\quad+\big{|}P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}\big{|}$
(A.20)
$\displaystyle\quad\quad\quad+P\left\\{\sqrt{n}\|(\Lambda^{-1/2}-\hat{\Lambda}_{n}^{-1/2})R_{n}(\bm{x}^{(d)})\|_{\infty}\geq
t\right\\},$ (A.21)

and that similarly

$\displaystyle\big{|}P\left\\{\sqrt{n}\hat{\Lambda}_{n}^{-1/2}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}\big{|}$
$\displaystyle\quad\quad\leq|P\left\\{\sqrt{n}\Lambda^{-1/2}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}_{t}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}|$
(A.22)
$\displaystyle\quad\quad\quad+|P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}|$
(A.23)
$\displaystyle\quad\quad\quad+P\left\\{\sqrt{n}\|(\Lambda^{-1/2}-\hat{\Lambda}_{n}^{-1/2})R_{n}^{*}(\bm{x}^{(d)})\|_{\infty}\geq
t\mid\mathbf{D}_{n}\right\\},$ (A.24)

for any $t>0$. We provide bounds for each term, (A.19) through (A.24).

To bound the Gaussian approximation errors (A.19) and (A.22), we apply the
following Theorem, which establishes a generic quantitative central limit for
the statistic $R_{n}$ in addition to generic quantitative conditional central
limit theorems for both bootstrap procedures.

######  Theorem A.2.

Suppose that the moment function $M(x;\theta_{0},g_{0})$ satisfies Assumption
A.1 and Part (i) of Assumption 3.3 and that Assumptions A.3 and A.4 hold.

(i) The inequality

$\displaystyle\sup_{R\in\mathcal{R}}\sup_{P\in\mathbf{P}}\Big{|}P\left\\{\sqrt{n}R_{n}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}\Big{|}$
$\displaystyle\lesssim\frac{\varphi^{1/2}}{\underline{\lambda}^{1/2}}\left(\frac{\log^{5}(dn)}{n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}+\rho_{n}$
(A.25)

holds.

(ii) If the bootstrap root is constructed with the Half-Sample bootstrap, then
the inequality

$\displaystyle\sup_{R\in\mathcal{R}}\sup_{P\in\mathbf{P}}\Big{|}P\left\\{\sqrt{n}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}\mid\mathbf{D}_{n}\right\\}-P\left\\{Z_{n}\in\mathsf{R}\right\\}\Big{|}$
$\displaystyle\lesssim\frac{\varphi^{1/2}}{\underline{\lambda}^{1/2}}\left(\frac{\log^{5}\left(dn\right)}{n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}$
(A.26)

holds with probability greater than
$1-Cn^{-1/2}\varphi\underline{\lambda}^{-1}\log^{3/2}(dn)-\rho_{n}$.

To bound the differences in the Gaussian probabilities (A.20) and (A.23), we
apply the following anti-concentration inequality, stated in Chernozhukov et
al., 2017b  and often referred to as Nazarov’s inequality.

######  Lemma A.1 (Chernozhukov et al., 2017b, , Theorem 1, ).

Let $Z=(Z_{j})_{j=1}^{d}$ be a centered Gaussian random vector in
$\mathbb{R}^{d}$ such that $\mathbb{E}[Z_{j}^{2}]\geq c$ for all $j$ in $[d]$
and some constant $c$. For every $z\in\mathbb{R}^{d}$ and $t>0$, the
inequality

$P\left\\{Z\leq z+t\right\\}-P\left\\{Z\leq
z\right\\}\lesssim\frac{t}{c}\sqrt{\log d}$

holds.

In particular, we have that

$|P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}|\leq
t\sqrt{\log d}$

for all $t>0$.

Finally, to bound the terms (A.21) and (A.24) resulting from variance
estimation, we apply the following bound on the accuracy of bootstrap variance
estimate $\hat{\lambda}^{2}_{n,j}$.

######  Lemma A.2.

Suppose that Assumption A.3 holds. If the bootstrap root is constructed with
the Half-Sample bootstrap, then

$P\left\\{\sup_{j\in[d]}\bigg{|}\frac{\hat{\lambda}^{2}_{n,j}}{\lambda^{2}_{j}}-1\bigg{|}\geq
C\frac{\varphi^{2}}{\underline{\lambda}^{2}n}\log^{2}(dn)+C\frac{1}{n}\delta_{n}^{2}\right\\}\lesssim
1-C(\rho_{n}+n^{-1})\leavevmode\nobreak\ .$ (A.27)

To apply Lemma A.2, observe that the Borell-TIS inequality (e.g., Theorem
2.1.1 of Adler and Taylor,, 2009) implies that

$P\left\\{\|\Lambda^{-1/2}Z\|_{\infty}\geq C\sqrt{\log{dn}}\right\\}\leq
n^{-1}\leavevmode\nobreak\ .$

Thus, Theorem A.2 implies that

$\displaystyle P\left\\{\|\Lambda^{-1/2}R_{n}(\bm{x}^{(dn)})\|_{\infty}\geq
C\sqrt{\log{dn}}\right\\}$
$\displaystyle\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}+\rho_{n}\leavevmode\nobreak\
.$

and that

$\displaystyle P\left\\{\|\Lambda^{-1/2}R^{*}_{n}(\bm{x}^{(d)})\|_{\infty}\geq
C\sqrt{\log{dn}}\mid\mathbf{D}_{n}\right\\}$
$\displaystyle\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}$

with probability greater than
$1-Cn^{-1/2}\varphi\underline{\lambda}^{-1}\log^{3/2}(dn)-\rho_{n}$. Hence,
Lemma A.2 implies that

$\sup_{j\in[d]}\big{|}{\hat{\lambda}_{n,j}}/{\lambda_{j}}-1\big{|}\leq\sup_{j\in[d]}\big{|}{\hat{\lambda}_{n,j}^{2}}/{\lambda_{j}^{2}}-1\big{|}\lesssim\frac{\varphi^{2}}{\underline{\lambda}^{2}n}\log^{2}(dn)+\frac{1}{n}\delta_{n}^{2}$

with probability at least $1-C(\rho_{n}+n^{-1})$. Thus, we have that

$\displaystyle
P\left\\{\|(\Lambda^{-1/2}-\hat{\Lambda}_{n}^{-1/2})R_{n}(\bm{x}^{(d)})\|_{\infty}\gtrsim\left(\frac{\varphi^{2}}{\underline{\lambda}^{2}n}\log^{2}(dn)+\frac{1}{n}\delta_{n}^{2}\right)\sqrt{\log{dn}}\right\\}$
$\displaystyle\quad\quad\quad\quad\quad\quad\quad\quad\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}+\rho_{n}$

and that

$\displaystyle
P\left\\{\|(\Lambda^{-1/2}-\hat{\Lambda}_{n}^{-1/2})R_{n}^{*}(\bm{x}^{(d)})\|_{\infty}\gtrsim\left(\frac{\varphi^{2}}{\underline{\lambda}^{2}n}\log^{2}(dn)+\frac{1}{n}\delta_{n}^{2}\right)\sqrt{\log{dn}}\mid\mathbf{D}_{n}\right\\}$
$\displaystyle\quad\quad\quad\quad\quad\quad\quad\quad\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{d}}$

with probability greater than
$1-Cn^{-1/2}\varphi\underline{\lambda}^{-1}\log^{3/2}(dn)-\rho_{n}$.

Putting the pieces together, by setting

$t=\left(\frac{\varphi^{2}}{\underline{\lambda}^{2}n}\log^{2}(dn)+\frac{1}{n}\delta_{n}^{2}\right)\sqrt{\log{dn}}\leavevmode\nobreak\
,$

we have that

$\displaystyle\sup_{R\in\mathcal{R}}\sup_{P\in\mathbf{P}}|P\left\\{\hat{\Lambda}_{n}^{-1/2}R_{n}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}|$
$\displaystyle\quad\quad\quad\quad\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log(d)}+\left(\frac{\varphi^{2}\log^{2}(dn}{\underline{\lambda}^{2}n})+\frac{1}{n}\delta_{n}^{2}\right)\log(dn)+\rho_{n}$
(A.28)

and that

$\displaystyle\sup_{R\in\mathcal{R}}\sup_{P\in\mathbf{P}}|P\left\\{\hat{\Lambda}_{n}^{-1/2}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}\mid\mathbf{D}_{n}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}|$
$\displaystyle\quad\quad\quad\quad\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\left(\frac{\varphi^{2}\log^{2}(dn)}{\underline{\lambda}^{2}n}+\frac{1}{n}\delta_{n}^{2}\right)\log(dn)$
(A.29)

with probability greater than
$1-Cn^{-1/2}\varphi\underline{\lambda}^{-1}\log^{3/2}(dn)-\rho_{n}$. To
complete the proof, observe that

$\displaystyle\big{|}P\left\\{\bm{\theta}_{0}(\bm{x}^{(d)})\in\hat{\mathcal{C}}(\bm{x}^{(d)},\mathbf{D}_{n})\right\\}-(1-\alpha)\big{|}$
$\displaystyle\quad=\big{|}P\left\\{\sqrt{n}\|\hat{\Lambda}_{n}^{-1/2}R_{n}(\bm{x}^{(d)})\|_{\infty}\leq\mathsf{cv}_{n}(\alpha)\right\\}-P\left\\{\sqrt{n}\|\hat{\Lambda}_{n}^{-1/2}R^{*}_{n}(\bm{x}^{(d)})\|_{\infty}\leq\mathsf{cv}_{n}(\alpha)\right\\}\big{|}$
$\displaystyle\quad\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log{dn}}+\rho_{n}$

by (A.28) and (A.29), as required.

###  A.4. Proofs for Supporting Lemmas

####  A.4.1. Proof of Theorem A.2, Part (i)

Throughout, we take $\theta_{0}(x)=0$ for all $x$ without loss of generality.
We begin by showing that the root $R_{n}(\bm{x}^{(d)})$ is well-approximated
by the statistic $U_{n}(\bm{x}^{(d)})$. Take $x$ to be any component of the
vector $\bm{x}^{(d)}$. By a Taylor expansion about $\theta_{0}(x)$, we have
that

$\displaystyle
M(x;\hat{\theta}_{n}(x),\hat{g}_{n})-M(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle=(\hat{\theta}_{n}(x)-\theta_{0}(x))M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle+(\hat{\theta}_{n}(x)-\theta_{0}(x))^{2}H(x;\tilde{\theta}_{0}(x),\hat{g}_{n})$
(A.30)

for some $\tilde{\theta}_{0}(x)$ between $\hat{\theta}_{n}(x)$ and
$\theta_{0}(x)$. Moreover, we can write

$\displaystyle(\hat{\theta}_{n}(x)-\theta_{0}(x))M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle=(\hat{\theta}_{n}(x)-\theta_{0}(x))M^{(1)}(x;\theta_{0}(x),g_{0})$
$\displaystyle+(\hat{\theta}_{n}(x)-\theta_{0}(x))\left(M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})-M^{(1)}(x;\theta_{0}(x),g_{0})\right)$
(A.31)

and

$\displaystyle
M(x;\hat{\theta}_{n}(x),\hat{g}_{n})-M(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle=(M(x;\hat{\theta}_{n}(x),\hat{g}_{n})-M_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n},\mathbf{D}_{n}))$
$\displaystyle\quad+(M(x;\theta_{0}(x),g_{0})-M(x;\theta_{0}(x),\hat{g}_{n}))$
$\displaystyle=\left(M(x;\hat{\theta}_{n}(x),\hat{g}_{n})-\mathbb{E}\left[M_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n},\mathbf{D}_{n})\right]\right)$
$\displaystyle\quad+\left(\mathbb{E}\left[M_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n},\mathbf{D}_{n})\right]-M_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n},\mathbf{D}_{n})\right)$
$\displaystyle\quad+(M(x;\theta_{0}(x),g_{0})-M(x;\theta_{0}(x),\hat{g}_{n}))\leavevmode\nobreak\
.$ (A.32)

Thus, by the identity

$\displaystyle\overline{M}_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n})$
$\displaystyle=\overline{M}_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n})-\overline{M}_{n}(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle+\overline{M}_{n}(x;\theta_{0}(x),\hat{g}_{n})-\overline{M}_{n}(x;\theta_{0}(x),g_{0})+\overline{M}_{n}(x;\theta_{0}(x),g_{0})\leavevmode\nobreak\
,$ (A.33)

the equalities (A.30), (A.31), and (A.32) imply that

$\displaystyle
M^{(1)}(x;\theta_{0}(x),g_{0})(\hat{\theta}_{n}(x)-\theta_{0}(x))$
$\displaystyle=M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})(\hat{\theta}_{n}(x)-\theta_{0}(x))$
$\displaystyle-\left(M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})-M^{(1)}(x;\theta_{0},g_{0})\right)(\hat{\theta}_{n}(x)-\theta_{0}(x))$
$\displaystyle=-\overline{M}_{n}(x;\theta_{0}(x),g_{0})$ (A.34)
$\displaystyle\quad+\mathsf{Bias}(x;\hat{\theta}_{n}(x),\hat{g}_{n})+\mathsf{Nuis}(x;\theta_{0}(x),\hat{g}_{n})$
(A.35)
$\displaystyle\quad+\mathsf{Stoch}^{(1)}(x;\hat{\theta}_{n}(x),\hat{g}_{n})+\mathsf{Stoch}^{(2)}(x;\hat{g}_{n})$
(A.36)
$\displaystyle\quad-(\hat{\theta}_{n}(x)-\theta_{0}(x))^{2}H(x;\tilde{\theta}_{0}(x),\hat{g}_{n})$
(A.37)
$\displaystyle\quad-(\hat{\theta}_{n}(x)-\theta_{0}(x))\left(M^{(1)}(x;\theta_{0}(x),\hat{g}_{n})-M^{(1)}(x;\theta_{0}(x),g_{0})\right)\leavevmode\nobreak\
$ (A.38)

where

$\displaystyle\mathsf{Bias}(x;\hat{\theta}_{n}(x),\hat{g}_{n})$
$\displaystyle=M(x;\hat{\theta}_{n}(x),\hat{g}_{n})-\mathbb{E}\left[M_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n})\right]\leavevmode\nobreak\
,$ (A.39) $\displaystyle\mathsf{Nuis}(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle=M(x;\theta_{0}(x),g_{0})-M(x;\theta_{0}(x),\hat{g}_{n})\leavevmode\nobreak\
,$ (A.40)
$\displaystyle\mathsf{Stoch}^{(1)}(x;\hat{\theta}_{n}(x),\hat{g}_{n})$
$\displaystyle=\overline{M}_{n}(x;\hat{\theta}_{n}(x),\hat{g}_{n})-\overline{M}_{n}(x;\theta_{0}(x),\hat{g}_{n})\leavevmode\nobreak\
,\quad\text{and}$ (A.41) $\displaystyle\mathsf{Stoch}^{(2)}(x;\hat{g}_{n})$
$\displaystyle=M(x;\theta_{0}(x),g_{0})-M(x;\theta_{0}(x),\hat{g}_{n})\leavevmode\nobreak\
,$ (A.42)

respectively.

We now give bounds for the terms (A.35), (A.36), (A.37), and (A.38). To handle
(A.35), observe that Parts (i) and (ii) of Assumption A.4 imply that

$\sqrt{\frac{n}{\underline{\lambda}^{2}}}|\mathsf{Bias}(\bm{x}^{(d)};\hat{\theta}_{n}(x),\hat{g}_{n})|\lesssim(1+\|\hat{\theta}_{n}(\bm{x}^{(d)})\|_{\infty})\delta_{n,B}\lesssim\delta_{n,B}\left(1+\frac{\underline{\lambda}^{1/2}}{n^{1/4}}\delta_{n,\theta}\right)$
(A.43)

with probability greater than $1-\rho_{n,\theta}$. Moreover, a Taylor
expansion, second-order smoothness, i.e., Assumption A.1, and Assumption A.4,
Part (ii), give that

$\displaystyle\sqrt{\frac{n}{\underline{\lambda}^{2}}}\mathsf{Nuis}(\bm{x}^{(d)};\theta_{0}(x),\hat{g}_{n})$
$\displaystyle=\sqrt{\frac{n}{\underline{\lambda}^{2}}}\partial_{g}M(\bm{x}^{(d)};\theta_{0}(x),g_{0})[g-g_{0}]$
$\displaystyle+\sqrt{\frac{n}{\underline{\lambda}^{2}}}\partial_{g,g}M(\bm{x}^{(d)};\theta_{0}(x),g_{0})[g-g_{0}]\lesssim\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|g-g_{0}\|_{2,\infty}^{2}\lesssim\delta^{2}_{n,g}\leavevmode\nobreak\
.$ (A.44)

with probability greater than $1-\rho_{n,g}$.

Next, we handle the term (A.36). By Assumption A.2, a Taylor expansion gives

$\displaystyle\mathsf{Stoch}^{(1)}(x;\hat{\theta}_{n}(x),\hat{g}_{n})$
$\displaystyle=(\hat{\theta}_{n}(x)-\theta_{0}(x))\overline{M}_{n}^{(1)}(x;\theta_{0}(x),\hat{g}_{n})$
$\displaystyle+(\hat{\theta}_{n}(x)-\theta_{0}(x))^{2}\overline{H}_{n}(x;\tilde{\theta}_{0}(x),\hat{g}_{n})$
(A.45)

for some, potentially different, $\tilde{\theta}_{0}(x)$ between
$\hat{\theta}_{n}(x)$ and $\theta_{0}(x)$. To bound this term, observe that

$\displaystyle|\overline{M}_{n}^{(1)}(x;\theta_{0}(x),\hat{g}_{n})|\leq|\overline{M}_{n}^{(1)}(x;\theta_{0}(x),g_{0})|+|\overline{M}_{n}^{(1)}(x;\theta_{0}(x),\hat{g}_{n})-\overline{M}_{n}^{(1)}(x;\theta_{0}(x),g_{0})|\leavevmode\nobreak\
.$ (A.46)

Hence, Assumption A.2 and Assumption A.4, Parts (ii) and (iii) imply that

$\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\mathsf{Stoch}^{(1)}(\bm{x}^{(d)};\hat{\theta}_{n}(\bm{x}^{(d)}),\hat{g}_{n})\|_{\infty}\leq\delta_{n,m}\delta_{n,\theta}+\frac{\underline{\lambda}^{1/2}}{n^{1/4}}\delta_{n,\theta}\delta_{n,J}+\delta^{2}_{n,\theta}$
(A.47)

with probability greater than $1-\rho_{n,m}-\rho_{n,\theta}-\rho_{n,J}$.
Moreover, we have that

$\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\mathsf{Stoch}^{(2)}(\bm{x}^{(d)};\hat{g}_{n})\|_{\infty}\leq\delta_{n,S}$
(A.48)

with probability greater than $1-\rho_{n,S}$.

Finally, we handle the terms (A.37) and (A.38). Assumption A.4, Part (ii), and
Assumption A.1 imply that

$\displaystyle\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|(\hat{\theta}_{n}(\bm{x}^{(d)})-\theta_{0}(\bm{x}^{(d)}))^{2}H(\bm{x}^{(d)};\tilde{\theta}_{0}(\bm{x}^{(d)}),\hat{g}_{n})\|_{\infty}$
$\displaystyle\lesssim\delta^{2}_{n,\theta}\quad\text{and}$ (A.49)
$\displaystyle\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|(\hat{\theta}_{n}(\bm{x}^{(d)})-\theta_{0}(\bm{x}^{(d)}))(M^{(1)}(\bm{x}^{(d)};\theta_{0}(x),\hat{g}_{n})-M^{(1)}(\bm{x}^{(d)};\theta_{0}(x),g_{0}))\|_{\infty}$
$\displaystyle\lesssim\delta_{n,\theta}\delta_{n,g}$ (A.50)

with probabilities greater than $1-\rho_{n,\theta}$ and $1-\rho_{n,g}$,
respectively. Putting the pieces together, the decomposition (A.34) and the
lower-boundedness of the Jacobian $M^{(1)}(\cdot;\theta,g)$ imply that

$\displaystyle\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|R_{n}(\bm{x}^{(d)})-U(\bm{x}^{(d)})\|_{\infty}$
(A.51)
$\displaystyle\lesssim\delta^{2}_{n,g}+\delta^{2}_{n,\theta}+\delta_{n,B}+\delta_{n,S}+\delta_{n,\theta}\left(\delta_{n,m}+\delta_{n,g}+\underline{\lambda}^{1/2}n^{-1/4}\left(\delta_{n,B}+\delta_{n,J}\right)\right)$

with probability greater than
$1-\rho_{n,m}+\rho_{n,g}+\rho_{n,\theta}+\rho_{n,S}+\rho_{n,J}$, by the bounds
(A.43), (A.44), (A.47), (A.48), (A.49), and (A.50).

With this in place, consider the decomposition

$\displaystyle R_{n}(x)$
$\displaystyle=\left(\frac{1}{n}\sum_{i=1}^{n}\overline{u}(x,D_{i})-U_{n}(x)\right)-\frac{1}{n}\sum_{i=1}^{n}\overline{u}(x,D_{i})+\Delta_{n}(x),\quad\text{where}$
(A.52) $\displaystyle\Delta_{n}(x)$
$\displaystyle=R_{n}(x)-U(x)\leavevmode\nobreak\ ,$ (A.53)

and we recall that the function $\overline{u}(\cdot,\cdot)$ is defined in
Assumption A.3. Fix a rectangle $\mathsf{R}=[a_{l},a_{u}]$ in $\mathcal{R}$.
Define the normalized functions

$\hat{u}(\bm{x}^{(d)},D_{i})=\Lambda^{-1/2}\overline{u}(\bm{x}^{(d)},D_{i})\quad\text{and}\quad\hat{U}(\bm{x}^{(d)})=\Lambda^{-1/2}U(\bm{x}^{(d)})$

and the analogously normalized rectangle
$\tilde{\mathsf{R}}=[\Lambda^{-1/2}a_{l},\Lambda^{-1/2}a_{u}]$ and the
enlarged rectangle
$\tilde{\mathsf{R}}_{t}=[\Lambda^{-1/2}a_{l}-\mathbf{1}_{d}t,\Lambda^{-1/2}a_{u}+\mathbf{1}_{d}t]$.
Observe that the decomposition (A.52) implies that

$\displaystyle\Big{|}P\left\\{\sqrt{n}R_{n}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}\Big{|}$
$\displaystyle=\Big{|}P\left\\{\sqrt{n}\Lambda^{-1/2}R_{n}(\bm{x}^{(d)})\in\tilde{\mathsf{R}}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}\right\\}\Big{|}$
$\displaystyle\quad\leq\Big{|}P\left\\{\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})\in\tilde{\mathsf{R}}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}_{t}\right\\}\Big{|}$
(A.54)
$\displaystyle\quad\quad+\Big{|}P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}_{t}}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}\right\\}\Big{|}$
(A.55)
$\displaystyle\quad\quad+P\left\\{\sqrt{n}\|\hat{U}_{n}(\bm{x}^{(d)})-\frac{1}{n}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})\|_{\infty}\geq\frac{1}{2}t\right\\}$
(A.56)
$\displaystyle\quad\quad+P\left\\{\sqrt{n}\|\Lambda^{-1/2}\Delta_{n}(\bm{x}^{(d)})\|_{\infty}\geq\frac{1}{2}t\right\\}$
(A.57)

for each $t>0$. The proof is complete by providing appropriate bounds for the
terms (A.54) through (A.57).

We bound the normal approximation term (A.54) through the application of the
following quantitative central limit theorem, stated as Theorem 2.1 of
Chernozhuokov et al., (2022).

######  Lemma A.3 (Chernozhuokov et al.,, 2022, Theorem 2.1, ).

Let $X_{1},\ldots,X_{n}$ be a collection of independent, centered, random
vectors in $\mathbb{R}^{d}$ and let $Z$ be a centered Gaussian random vector
with covariance matrix

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}\left[X_{i}X_{i}^{\top}\right].$ (A.58)

If there exist absolute constants $c$, $C_{1}$, and $\varphi$ such that the
bounds

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}\left[X_{i,j}^{2}\right]\geq
c\leavevmode\nobreak\
,\quad\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}\left[X_{i,j}^{4}\right]\leq
C_{1}\varphi^{2}\leavevmode\nobreak\
,\quad\text{and}\quad\|X_{i,j}\|_{\psi_{1}}\leq\varphi$ (A.59)

hold, then the inequality

$\sup_{\mathsf{R}\in\mathcal{R}}\Big{|}P\left\\{\frac{1}{\sqrt{n}}\sum_{i=1}^{n}X_{i}\in\mathsf{R}\right\\}-P\left\\{\sqrt{n}Z\in\mathsf{R}\right\\}\Big{|}\leq
C_{2}\left(\frac{\varphi^{2}\log^{5}(dn)}{n}\right)^{1/4}$ (A.60)

holds for some constant $C_{2}$ that depends only on $c$ and $C_{1}$.

In particular, observe that

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}\left[\hat{u}^{2}(x^{(j)},D_{i})\right]=1$
(A.61)

by definition and that

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}\left[\hat{u}^{4}(x^{(j)},D_{i})\right]\leq\left(\varphi/\underline{\lambda}\right)^{2}$
(A.62)

by Assumption A.3. Similarly we have that

$\|\hat{u}\left(x^{j}),D_{i}\right)\|_{\psi_{1}}\leq(\varphi/\underline{\lambda})$
(A.63)

by Assumption A.3. Consequently, as

$\operatorname{Var}(\hat{u}(x^{(j)},D_{i}))=\Lambda^{-1/2}\operatorname{Var}\left(Z\right)\Lambda^{-1/2},$
(A.64)

by definition, Lemma A.3 implies that

$\Big{|}P\left\\{\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\hat{u}\left(D_{i}\right)\in\mathsf{\tilde{R}}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{\tilde{R}}_{t}\right\\}\Big{|}\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}\leavevmode\nobreak\
.$ (A.65)

In turn, to bound the term (A.55), we have that

$\Big{|}P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}\right\\}\Big{|}\lesssim
t\sqrt{\log d}\leavevmode\nobreak\ .$ (A.66)

by Lemma A.1. Moreover, to bound the term (A.56), recall that

$P\left\\{\sqrt{n}\|\hat{U}_{n}(\bm{x}^{(d)})-\frac{1}{n}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})\|_{\infty}\geq\delta_{n,u}\right\\}\leq\rho_{n,u}$
(A.67)

by Assumption A.3. Thus, by choosing $t=C\delta_{n}$, the bound (A.51) implies
that the sum of term (A.56) and term (A.57) is upper bounded by $C\rho_{n}$.
Hence, by plugging this choice of $t$ into (A.66) and (A.67), we can conclude
that

$\Big{|}P\left\\{\sqrt{n}R_{n}(\bm{x}^{(d)})\in\mathsf{R}\right\\}-P\left\\{Z_{n}\in\mathsf{R}\right\\}\Big{|}\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+\delta_{n}\sqrt{\log(d)}+\rho_{n}\leavevmode\nobreak\
,$ (A.68)

as required.

####  A.4.2. Proof of Theorem A.2, Part (ii)

The result follows from an argument whose structure is similar to the proof of
Theorem A.2, Part (i). Again, we take $\theta_{0}(x)=0$ for all $x$, without
loss of generality. We are interested in studying the discrepancy

$\displaystyle R_{n}^{*}(x)$
$\displaystyle=\hat{\theta}_{\mathsf{h}}(x)-\hat{\theta}_{n}(x)=(\hat{\theta}_{\mathsf{h}}(x)-\theta_{0}(x))-R_{n}(x).$

We know from (A.52) that

$\displaystyle
R_{n}(x)=\left(\frac{1}{n}\sum_{i=1}^{n}\overline{u}(x,D_{i})-U_{n}(x)\right)-\frac{1}{n}\sum_{i=1}^{n}\overline{u}(x,D_{i})+\Delta_{n}(x)\leavevmode\nobreak\
,$ (A.69)

where $\Delta_{n}(x)$ is defined in (A.53). On the other hand, as
$\hat{\theta}_{\mathsf{h}}(x)$ is constructed with a random half-sample
$\mathsf{h}$ of the data $\mathbf{D}_{n}$, we have that

$\displaystyle\hat{\theta}_{\mathsf{h}}(x)-\theta_{0}(x)=\left(\frac{2}{n}\sum_{i\in\mathsf{h}}^{n}\overline{u}(x,D_{i})-U_{\mathsf{h}}(x)\right)-\frac{2}{n}\sum_{i\in\mathsf{h}}^{n}\overline{u}(x,D_{i})+\Delta_{\mathsf{h}}(x)\leavevmode\nobreak\
,$ (A.70)

where $U_{\mathsf{h}}(x)$ and $\Delta_{\mathsf{h}}(x)$ are constructed with
the half-sample $\mathsf{h}$.

The proof of Theorem A.2, Part (i), worked by giving a high probability bound
for the first and third term in (A.69) and showing that the the second term
satisfies a central limit theorem. Here, as we are interested in giving a
bound conditioned on the data $\mathbf{D}_{n}$, we show that the difference
between the second terms in (A.69) and (A.70) satisfies a central limit
theorem on the event that the first and third terms in (A.69) and (A.70)
satisfy a specified bound, which we show holds with high probability. In
particular, let $\mathcal{F}_{n}(t)$ and $\mathcal{F}_{\mathsf{h}}(t)$ denote
the events that

$\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\Delta_{n}(\bm{x}^{(d)})\|_{\infty}\leq
t/4\quad\text{and}\quad\sqrt{\frac{n}{\underline{\lambda}^{2}}}\|\Delta_{\mathsf{h}}(\bm{x}^{(d)})\|_{\infty}\leq
t/4\leavevmode\nobreak\ ,$ (A.71)

respectively. Similarly, let $\mathcal{H}_{n}(t)$ and
$\mathcal{H}_{\mathsf{h}}(t)$ denote the events that

$\displaystyle\sqrt{n}\|\frac{1}{n}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})-\hat{U}_{n}(\bm{x}^{(d)})\|_{\infty}\leq
t/4\quad\text{and}$ (A.72)
$\displaystyle\sqrt{n}\|\frac{2}{n}\sum_{i\in\mathsf{h}}^{n}\hat{u}(\bm{x}^{(d)},D_{i})-\hat{U}_{\mathsf{h}}(\bm{x}^{(d)})\|_{\infty}\leq
t/4$

respectively, where $\hat{U}_{\mathsf{h}}(\bm{x}^{(d)})$ is defined
analogously to $\hat{U}_{n}(\bm{x}^{(d)})$. Define the event
$\mathcal{E}_{n}(t)=\mathcal{F}_{n}(t)\cap\mathcal{F}_{\mathsf{h}}(t)\cap\mathcal{H}_{n}(t)\cap\mathcal{H}_{\mathsf{h}}(t)$.
Fix a hyper-rectangle $\mathsf{R}$ in $\mathcal{R}$. Recall the definitions of
the transformed rectangle $\tilde{\mathsf{R}}$ and the enlarged transformed
rectangle $\tilde{\mathsf{R}}_{t}$. On the event $\mathcal{E}_{n}(t)$, we have

$\displaystyle|P\left\\{\sqrt{n}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}\mid\mathbf{D}_{n}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}|$
$\displaystyle=|P\left\\{\sqrt{n}\Lambda^{-1/2}R_{n}^{*}(\bm{x}^{(d)})\in\tilde{\mathsf{R}}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}\right\\}|$
$\displaystyle\quad\leq|P\left\\{\frac{2}{\sqrt{n}}\sum_{i\in\mathsf{h}}^{n}\hat{u}(\bm{x}^{(d)},D_{i})-\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})\in\tilde{\mathsf{R}}_{t}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}_{t}\right\\}|$
$\displaystyle\quad+|P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\tilde{\mathsf{R}}_{t}\right\\}|$

for each $t>0$. As the data $D_{\mathsf{h}}$ are drawn independently and
identically with distribution $P$ in $\mathbf{P}$ and we have assumed that
$\delta_{n/2}\lesssim\delta_{n}$ and $\rho_{n/2}\lesssim\rho_{n}$, by setting
$t=C\delta_{n}$, Assumption A.3 and the bound (A.51) imply that the event
$\mathcal{E}_{n}(t)$ occurs with probability greater than $1-C\rho_{n}$. Thus,
Lemma A.1 implies that

$\displaystyle|P\left\\{\sqrt{n}R_{n}^{*}(\bm{x}^{(d)})\in\mathsf{R}\mid\mathbf{D}_{n}\right\\}-P\left\\{Z\in\mathsf{R}\right\\}|$
(A.73)
$\displaystyle\quad\leq\bigg{|}P\left\\{\frac{2}{\sqrt{n}}\sum_{i\in\mathsf{h}}^{n}\hat{u}(\bm{x}^{(d)},D_{i})-\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})\in\tilde{\mathsf{R}}_{t}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{1/2}Z\in\tilde{\mathsf{R}}_{t}\right\\}\bigg{|}$
(A.74) $\displaystyle\quad\quad+\delta_{n}\sqrt{\log(d)}$

with probability greater than $1-C\rho_{n}$. Hence, it suffices to bound the
term (A.74).

To this end, we apply a coupling argument introduced in Yadlowsky et al.,
(2023), which is similar to a Poissonization technique studied in Præstgaard
and Wellner, (1993) (see also Section 3.6.2 of van der Vaart and Wellner,
(2013)). In particular, let $V_{i}$ be a random variable taking the value $1$
when $i$ is an element of the subset $\mathsf{h}$ and taking the value $-1$
otherwise. Observe that

$\displaystyle\frac{2}{n}\sum_{i\in\mathsf{h}}\hat{u}(\bm{x}^{(d)},D_{i})-\frac{1}{n}\sum_{i=1}^{n}\hat{u}(\bm{x}^{(d)},D_{i})=\frac{1}{n}\sum_{i=1}^{n}V_{i}\hat{u}(\bm{x}^{(d)},D_{i})\leavevmode\nobreak\
.$ (A.75)

Let $\tilde{V}_{1},\ldots,\tilde{V}_{n}$ denote a collection of random
variables valued on $\\{-1,1\\}$. We define their joint distribution as
follows. Let $Q_{n}$ denote a random variable with distribution
$\mathsf{Bin}\left(n,1/2\right)$. If $Q_{n}\geq n/2$, then choose $Q_{n}-n/2$
indices $i$ in $[n]$ with $V_{i}=-1$ and set $\tilde{V}_{i}=1$. If
$Q_{n}<n/2$, then choose $n/2-Q_{n}$ indices with $V_{i}=1$ and set
$\tilde{V}_{i}=-1$. Set $\tilde{V}_{i}=V_{i}$ for all other units. Observe
that the collection $\tilde{V}_{i}$ are independent and identically
distributed Rademacher random variables. With this in place, we obtain the
decomposition

$\frac{1}{n}\sum_{i=1}^{n}V_{i}\hat{u}(\bm{x}^{(d)},D_{i})=\frac{1}{n}\sum_{i=1}^{n}\tilde{V}_{i}\hat{u}(\bm{x}^{(d)},D_{i})+\frac{1}{n}\sum_{i=1}^{n}\left(V_{i}-\tilde{V}_{i}\right)\hat{u}(\bm{x}^{(d)},D_{i})\leavevmode\nobreak\
.$ (A.76)

Let $\mathcal{V}(t)$ denote the event that

$\sqrt{n}\Big{\|}\frac{1}{n}\sum_{i=1}^{n}\left(V_{i}-\tilde{V}_{i}\right)\hat{u}(\bm{x}^{(d)},D_{i})\Big{\|}_{\infty}>t\leavevmode\nobreak\
.$ (A.77)

Fix any rectangle $\mathsf{R}^{\prime}=[a^{\prime}_{l},a^{\prime}_{u}]$ and
define the enlarged rectangle
$\mathsf{R}^{\prime}_{t}=[a_{l}-t\mathbf{1}_{d},a_{u}+t\mathbf{1}_{d}]$. On
the event $\mathcal{V}(t)$, the decomposition (A.76) implies that

$\displaystyle\big{|}P\left\\{\frac{1}{\sqrt{n}}\sum_{i=1}^{n}V_{i}\hat{u}(\bm{x}^{(d)},D_{i})\in\mathsf{R}^{\prime}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}^{\prime}\right\\}\big{|}$
$\displaystyle\leq\big{|}P\left\\{\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\tilde{V}_{i}\hat{u}(\bm{x}^{(d)},D_{i})\in\tilde{\mathsf{R}}_{t}\mid\mathbf{D}_{n}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}^{\prime}_{t}\right\\}\big{|}$
(A.78)
$\displaystyle\quad\quad+|P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}^{\prime}_{t}\right\\}-P\left\\{\Lambda^{-1/2}Z\in\mathsf{R}^{\prime}\right\\}|$
(A.79)

for all $t>0$. Lemma A.1 implies that (A.79) is less than $t\sqrt{\log(d)}$.
To handle the term (A.78), we apply the following quantitative central limit
theorem, stated as Lemma 4.6 in Chernozhuokov et al., (2022).

######  Lemma A.4 (Chernozhuokov et al.,, 2022, Lemma 4.6, ).

Consider the setting and assumptions of Lemma A.3. Let
$\bm{X}_{n}=(X_{1},\ldots,X_{n})$ collect the observed data and let
$\tilde{V}_{1},\ldots,\tilde{V}_{n}$ be a collection of independent Rademacher
random variables. We have that

$\sup_{\mathsf{R}\in\mathcal{R}}\Big{|}P\left\\{n^{-1/2}\sum_{i=1}^{n}X_{i}\in\mathsf{R}\right\\}-P\left\\{n^{-1/2}\sum_{i=1}^{n}\tilde{V}_{i}X_{i}\in
R\mid\bm{X}_{n}\right\\}\Big{|}\leq
C_{2}\left(\frac{\varphi^{2}\log^{5}(dn)}{n}\right)^{1/4}.$ (A.80)

with probability greater than

$1-C\frac{\varphi\log^{3/2}(dn)}{n^{1/2}}$ (A.81)

for some constant $C_{2}$ that depends only on the constants $C_{2}$ and $c$
defined in the statement of Lemma A.3.

Thus, on the event $\mathcal{V}(t)$, Lemma A.3 and Lemma A.4 imply that

$\displaystyle\big{|}P\left\\{\frac{1}{n}\sum_{i=1}^{n}V_{i}\hat{u}(\bm{x}^{(d)};D_{i})\in\mathsf{R}^{\prime}\mid\mathbf{D}_{n}\right\\}-P\left\\{\sqrt{n}Z\in\mathsf{R}^{\prime}\right\\}\big{|}\lesssim\left(\frac{\varphi^{2}\log^{5}(dn)}{\underline{\lambda}^{2}n}\right)^{1/4}+t\sqrt{\log(d)}\leavevmode\nobreak\
,$

with probability greater than
$1-Cn^{-1/2}\underline{\lambda}^{-1}\varphi\log^{3/2}(dn)$. Hence, it suffices
to give a high probability bound on $\mathcal{V}(t)$ for a suitable choice of
$t$.

To this end, observe that

$\displaystyle G_{n}$
$\displaystyle=\Big{\|}\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\left(V_{i}-\tilde{V}_{i}\right)\hat{u}(\bm{x}^{(d)};D_{i})\Big{\|}_{\infty}\quad\text{is
equidistributed with}\quad$ (A.82)
$\displaystyle\quad\Big{\|}\frac{2}{\sqrt{n}}\sum_{i=1}^{|Q_{n}-n/2|}\hat{u}(\bm{x}^{(d)};D_{i})\Big{\|}_{\infty}\leavevmode\nobreak\
.$ (A.83)

Consider the decomposition

$\displaystyle P\left\\{G_{n}\geq t\right\\}$ $\displaystyle\leq
P\left\\{G_{n}\geq
t,|Q_{n}-n/2|\leq\delta\frac{n}{2}\right\\}+P\left\\{|Q_{n}-n/2|\geq\delta\frac{n}{2}\right\\}$
$\displaystyle\leq P\left\\{\max_{1\leq
k\leq\delta\frac{n}{2}}\Big{\|}\frac{2}{\sqrt{n}}\sum_{i=1}^{k}\hat{u}\left(D_{i}\right)\Big{\|}_{\infty}\geq
t\right\\}+P\left\\{|Q_{n}-n/2|\geq\delta\frac{n}{2}\right\\},$ (A.84)

for some $\delta>0$ to be chosen. Observe that

$P\left\\{|Q_{n}-n/2|\geq\frac{\delta n}{2}\right\\}\leq
2\exp\left(-\frac{\delta^{2}n}{6}\right)$ (A.85)

by the multiplicative Chernoff bound. To bound the first term in (A.84), we
combine two inequalities. The first inequality is the following standard
Bernstein-type bound, stated in Song et al., (2019).

######  Lemma A.5 (Lemma A.2, Song et al., (2019)).

Let $Z_{1},\ldots,Z_{n}$ be independent, centered, random vectors in
$\mathbb{R}^{d}$. Define the quantity

$\sigma^{2}=\max_{j\in[d]}\sum_{i=1}^{n}\mathbb{E}\left[Z^{2}_{i,j}\right]$
(A.86)

and assume that $\|Z_{ij}\|_{\psi_{1}}\leq\phi$ for all $i\in[n]$ and
$j\in[d]$. The inequality