# Double machine learning and design in batch adaptive experiments

Harrison H. Li
Stanford University  Art B. Owen
Stanford University

(September 2023)

###### Abstract

We consider an experiment with at least two stages or batches and $O(N)$
subjects per batch. First, we propose a semiparametric treatment effect
estimator that efficiently pools information across the batches, and show it
asymptotically dominates alternatives that aggregate single batch estimates.
Then, we consider the design problem of learning propensity scores for
assigning treatment in the later batches of the experiment to maximize the
asymptotic precision of this estimator. For two common causal estimands, we
estimate this precision using observations from previous batches, and then
solve a finite-dimensional concave maximization problem to adaptively learn
flexible propensity scores that converge to suitably defined optima in each
batch at rate $O_{p}(N^{-1/4})$. By extending the framework of double machine
learning, we show this rate suffices for our pooled estimator to attain the
targeted precision after each batch, as long as nuisance function estimates
converge at rate $o_{p}(N^{-1/4})$. These relatively weak rate requirements
enable the investigator to avoid the common practice of discretizing the
covariate space for design and estimation in batch adaptive experiments while
maintaining the advantages of pooling. Our numerical study shows that such
discretization often leads to substantial asymptotic and finite sample
precision losses outweighing any gains from design.

##  1 Introduction

In sequential experimentation, we can use earlier observations to adjust our
treatment allocation policy for subsequent observations and thereby gain
improved estimation of causal effects in the overall study. For instance, for
an experiment with one treatment arm and one control arm, Neyman, (1934)
showed that choosing the number of subjects in each arm to be proportional to
the outcome standard deviation of that arm minimizes the variance of the
treatment effect estimate based on the difference in means. While these
standard deviations are unknown, they can be estimated using the initial data.
Then the Neyman allocation can be approximated to improve the sample
efficiency of the remainder of the experiment (Hahn et al.,, 2011; Blackwell
et al.,, 2022; Zhao,, 2023; Dai et al.,, 2023).

We study a version of this design problem for an experiment divided into a
small number of stages or _batches_. The design, or treatment assignment
mechanism, can be updated _adaptively_ for later batches based on the
observations from earlier batches to improve the precision of the causal
estimate computed at the end of the experiment. A salient feature of our
setting is knowledge of pre-treatment covariates that can further improve
precision. Thus, we conceptualize our design problem as choosing a _propensity
score_ for each batch. The propensity score specifies the probability that a
subject receives treatment given their covariates (throughout, we consider the
setting of a binary treatment). The propensity score is well known to be a key
mathematical object to be estimated in the causal analysis of observational
data (e.g. Rosenbaum and Rubin, (1983)). In a randomized experiment it is
known and under the control of the investigator. Hence, it can be exploited
for design.

Our specific design objective is to minimize an appropriate scalarization of
the asymptotic covariance matrix of an estimator that efficiently pools
information across all batches of the experiment. After describing our
mathematical notation and setup in Section 2, we present and study an oracle
version of this pooled estimator in Section 3. That oracle is given knowledge
of some nuisance parameters, typically infinite-dimensional mean or variance
functions. It pertains to a so-called “non-adaptive batch experiment” where
treatment is assigned with possibly varying but nonrandom propensity scores
across batches. In certain cases, our oracle pooled estimator asymptotically
dominates the best possible alternative that aggregates single batch
estimates, regardless of the per-batch propensities. This justifies designing
for the pooled estimator instead of an aggregation-based alternative.

For such a design procedure to be useful, however, we must show that the
targeted asymptotic precision of the oracle pooled estimator is in fact
attainable in a batched experiment where the propensity scores used are
adaptive (data-dependent) and the nuisance parameters need to be estimated. We
address these challenges in Section 4 by extending the framework of double
machine learning formalized by Chernozhukov et al., (2018), hereafter DML, to
what we call a “convergent split batch adaptive experiment,” or CSBAE. In a
CSBAE, observations in each batch are split into $K$ folds, and treatment is
assigned in each fold according to an adaptive propensity score that only
depends on observations from previous batches within the same fold. Within a
given batch, the $K$ adaptive propensity scores are further required to
converge to a common limit at rate $O_{p}(N^{-1/4})$ in root mean square
(RMS). Our DML extension then shows that by plugging in estimated nuisance
functions, we can construct a feasible estimator in a CSBAE that is
asymptotically equivalent to the oracle pooled estimator computed on the
limiting non-adaptive batch experiment. The nuisance function estimates only
need to converge at the rate $o_{p}(N^{-1/4})$.

Section 5 details a finite dimensional concave maximization procedure
(Algorithm 1) that provably constructs a CSBAE for which the limiting
propensities in each batch are sequentially optimal within a function class
satisfying standard complexity conditions. Hence, we can effectively design
for our pooled estimator in a batch adaptive experiment with theoretical
guarantees that our final estimator will indeed attain the targeted optimal
asymptotic precision, even with a fairly flexible propensity score learning
method and nonparametric machine learning estimates for nuisance functions.

To the best of our knowledge, existing work either designs for a less
efficient alternative to our pooled estimator, or discretizes the covariate
space at the design and estimation stages to construct a feasible variant of
the pooled estimator on a batch adaptive experiment. Our simulations in
Section 6 suggest that the latter approach in particular can lead to
substantial precision losses that swamp any gains from design, even with a
moderate number of continuous covariates. Thus, for the practitioner, we
provide an end-to-end design and estimation procedure to efficiently handle
continuous covariates in a batched experiment.

###  1.1 Related work

There has been substantial research interest in adaptive experiment designs in
recent years. In many applications, treatment assignments are updated in an
attempt to maximize the (expected) response values of either those in the
experiment, as in adaptive bandit algorithms (Russo et al.,, 2018; Hao et
al.,, 2020), or those in the superpopulation from which the experimental
subjects are assumed to arrive (Xu et al.,, 2018; Kasy and Sautmann,, 2021).
Inference on data collected from these algorithms can be challenging since the
treatment assignment rules often do not converge (Hadad et al.,, 2021; Zhang
et al.,, 2020, 2021). By contrast, in our setting where the goal is purely
statistical (maximizing asymptotic precision of the treatment estimate), the
design objective is a static propensity score to be learned consistently. An
interesting direction for further study would be to design for a mixture of
both statistical and non-statistical objectives. For example, one might expand
the literature on tie-breaker designs (Owen and Varian,, 2020; Morrison and
Owen,, 2022; Li and Owen,, 2023; Kluger and Owen,, 2023)) to the setting of
batched experiments.

The present work can be viewed as an extension of Hahn et al., (2011) in
several directions. Those authors considered a two batch experiment to
estimate the average treatment effect (ATE) as precisely as possible. Using
data from the first batch to estimate variance functions, they estimate the
asymptotic variance of a pooled version of the semiparametric efficient ATE
estimator of Hirano et al., (2003) for a coarsely discretized covariate. Then,
they learn a propensity score for the second batch that approximately
minimizes this variance. The covariate discretization ensures nuisance
functions and optimal propensities can be estimated at parametric
$O_{p}(N^{-1/2})$ rates without parametric assumptions. Consequently, a
feasible version of the pooled estimator indeed attains the targeted
asymptotic variance on the batch adaptive experiment. We generalize this
pooling construction beyond the setting of ATE estimation and relax these rate
requirements to those described in the previous section. This permits more
efficient handling of continuous covariates through nonparametric nuisance
function estimates and more flexible adaptive propensity scores.

Other approaches to extend the work of Hahn et al., (2011) include Kato et
al., (2020), who consider an online setting where subjects from a stationary
superpopulation enter one at a time, without batches. Similarly, the
literature on covariate-adjusted response-adaptive (CARA) designs has focused
on different but related objectives, both statistical and ethical (Zhang et
al.,, 2007; Zhu and Zhu,, 2023). In the batched setting, Tabord-Meehan, (2022)
proposes a method to learn a variance-minimizing stratification of the
covariate space of fixed size, avoiding the need to discretize the space prior
to observing the data as in Hahn et al., (2011). Cytrynbaum, (2021) showed
that by performing a form of highly stratified treatment assignment called
local randomization, consistent variance function estimates from the first
batch make it possible to attain the semiparametric lower bound for ATE
estimation in the second batch with optimal propensity score without having to
estimate the conditional mean functions. Neither of these approaches, however,
maintains the efficiency advantages of pooling. They also do not immediately
extend beyond ATE estimation.

##  2 Setup and notation

Let $T\geqslant 2$ be the number of batches in the experiment. Each subject
$i=1,\ldots,N_{t}$ in batch $t=1,\ldots,T$ has observed covariates
$X_{ti}\in{}^{d}$ and potential outcomes $Y_{ti}(0),Y_{ti}(1)\in\real$. We
place these in the vectors $S_{ti}=(X_{ti}^{\top},Y_{ti}(0),Y_{ti}(1))^{\top}$
which are exogenous in our model. Let $Z_{ti}\in\\{0,1\\}$ be the binary
treatment indicator for this subject, which is controlled by the investigator.
Under the usual stable unit value treatment assumption (SUTVA), the observed
outcome is

$Y_{ti}=Z_{ti}Y_{ti}(1)+(1-Z_{ti})Y_{ti}(0).$ (1)

Then the available data for the subject is
$W_{ti}=(X_{ti}^{\top},Z_{ti},Y_{ti})^{\top}\in\mathcal{W}$. We assume that
the vectors $S_{ti}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}P^{S}$ for
$t=1,\ldots,T$ and $i=1,\ldots,N_{t}$, for some distribution $P^{S}$. Appendix
A relaxes this assumption to permit certain forms of non-stationarity across
batches, such as covariate shifts. It will be convenient to define the
functions

$m_{0}(z,x)=\mathbb{E}[Y(z)\\!\mid\\!X=x]\quad\text{and}\quad
v_{0}(z,x)=\textnormal{Var}(Y(z)\\!\mid\\!X=x),$ (2)

for $z\in\\{0,1\\}$ and $x\in\mathcal{X}$. These expectations are taken under
$P^{S}$.

Let $N=N_{1}+\ldots+N_{T}$. Then as in Hahn et al., (2011), Che and Namkoong,
(2023) and others, we consider a proportional asymptotic regime

$\lim_{N\rightarrow\infty}\frac{N_{t}}{N}=\kappa_{t}\in(0,1),\quad
t=1,\ldots,T$ (3)

as $N\rightarrow\infty$. In settings where the batch sizes are fully
controlled by the experimenter, it may be theoretically preferable to make
initial batch sample sizes a vanishing fraction of the total sample size
(Zhao,, 2023). However, in many settings the batch sizes are exogenously
constrained to satisfy (3) unless observations are discarded.

For various $q\geqslant 1$ and probability measures $P$ on some space
$\Omega$, it will be useful to consider function norms of the form

$\|f\|_{q,P}=\left(\int|f(w)|^{q}\mathrm{d}P(w)\right)^{1/q}$

for $f\in L^{q}(P)$. We will use propensity scores denoted by $e(\cdot)$ with
various subscripts. A propensity score $e(\cdot)$ specifies
$e(x)=\textnormal{Pr}(Z=1\mid X=x)$, the probability of treatment conditional
on covariates. We will typically require propensity scores to lie in
$\mathcal{F}_{\gamma}$, the set of all measurable functions on $\mathcal{X}$
taking on values in the interval $[\gamma,1-\gamma]$ for some
$\gamma\in[0,1/2)$. We use $\|A\|$ to denote the square root of the sum of the
squared entries of any vector, matrix, or tensor $A$. For any integer
$p\geqslant 1$, $\mathbb{S}_{+}^{p}$ will denote the set of symmetric positive
semidefinite $p\times p$ real matrices, and $\mathbb{S}_{++}^{p}$ will be the
set of symmetric positive definite $p\times p$ real matrices. Finally, for any
real vector $v$, we write $v^{\otimes 2}=vv^{\top}$.

We summarize the preceding requirements for the data generating process in
Assumption 2.1. Assumption 2.1 does not impose any restrictions on the
treatment assignment process, which will be discussed at length in subsequent
sections.

######  Assumption 2.1 (Data generating process).

For some fixed number of batches $T\geqslant 2$, the vectors

$S_{ti}=(X_{ti},Y_{ti}(0),Y_{ti}(1)),\quad 1\leqslant t\leqslant T,\quad
1\leqslant i\leqslant N_{t}$

are independent and identically distributed (i.i.d.) from a distribution
$P^{S}$. Furthermore, the sample sizes $N_{t}$ satisfy (3), and the vector
$W_{ti}=(X_{ti},Z_{ti},Y_{ti})$ is observed where the outcomes $Y_{ti}$
satisfy the SUTVA assumption (1).

###  2.1 Estimands and score equations

Consider the setting where $T=1$ (so we can drop the batch subscript $t$), and
the observations $W_{1},\ldots,W_{N}$ are i.i.d. Suppose additionally that (1)
holds along with the unconfoundedness assumption

$(Y_{i}(0),Y_{i}(1))\perp\\!\\!\\!\perp Z_{i}\mid X_{i},\quad i=1,\ldots,N.$

Then many popular causal estimands $\theta_{0}\in\Theta\subseteq{}^{p}$ are
identified by a score equation

$\mathbb{E}[s(W;\theta_{0},\nu_{0},e_{0})]=0.$

In this score equation, $\nu_{0}$ is a vector of possibly infinite-dimensional
nuisance parameters lying in a nuisance set $\mathcal{N}$, and
$e_{0}=e_{0}(\cdot):\mathcal{X}\rightarrow[0,1]$ is the propensity score.
Following Section 3.1 of Chernozhukov et al., (2018), we will assume for
simplicity that the score $s(\cdot)$ is _linear_ in the sense that

$s(w;\theta,\nu,e)=s_{a}(w;\nu,e)\theta+s_{b}(w;\nu,e),\quad\forall
w\in\mathcal{W},\ \theta\in\Theta,\
(\nu,e)\in\mathcal{N}\times\mathcal{F}_{\gamma}$ (4)

for some $\gamma\in[0,1/2)$,
$s_{a}(\cdot,\nu,e):\mathcal{W}\rightarrow{}^{p\times p}$, and
$s_{b}(\cdot,\nu,e):\mathcal{W}\rightarrow{}^{p}$.

When $T>1$, propensity scores may vary across batches by design or external
constraints. For any propensity $e=e(\cdot)$ and integrable function
$f:\mathcal{W}\rightarrow\real$, we use the subscripted notation
$\mathbb{E}_{e}[f(W)]=\int f(w)\mathrm{d}P_{e}(w)$ where $P_{e}=P_{e}^{W}$ is
the distribution of $W=(X,Z,Y)=(X,Z,ZY(1)+(1-Z)Y(0))$ induced by
$S=(X,Y(0),Y(1))\sim P^{S}$ and $Z\\!\mid\\!X\sim\textnormal{Bern}(e(X))$
under the SUTVA assumption (1). Further let $P^{X}$ be the marginal
distribution of $X$ under $S\sim P^{S}$. Then we will require the following
score equations to hold for some $\gamma\in[0,1/2)$ to identify the our causal
estimand $\theta_{0}$:

$\mathbb{E}_{e}[s(W;\theta_{0},\nu_{0},e^{\prime})]=0,\quad\forall
e,e^{\prime}\in\mathcal{F}_{\gamma}.$ (5)

Note that (5) requires the score $s(\cdot)$ to have mean 0 when _any_
propensity score $e^{\prime}(\cdot)\in\mathcal{F}_{\gamma}$ is plugged in.
This plug-in propensity $e^{\prime}(\cdot)$ may differ from the propensity
$e(\cdot)\in\mathcal{F}_{\gamma}$ used for treatment assignment in the
experiment that generated the observations $W$. Such a robustness property is
satisfied by definition so long as the score $s(\cdot)$ is _doubly robust_. It
is required to ensure the validity of the pooled estimator that we propose in
Section 3. That estimator requires plugging in a mixture propensity score that
averages propensities across all batches $t=1,\ldots,T$. While identification
of $\theta_{0}$ within each batch is possible by only requiring (5) to hold
when $e(\cdot)=e^{\prime}(\cdot)$, this will not be sufficient to ensure
validity of our pooled estimator. We formally restate our requirements on
identification of the estimand $\theta_{0}$ in Assumption 2.2.

######  Assumption 2.2 (Estimand identification).

The estimand $\theta_{0}\in{}^{p}$ of interest satisfies (5) for some
$\gamma\in[0,1/2)$, some nuisance parameters $\nu_{0}$ lying in a known convex
set $\mathcal{N}$, and some score $s(\cdot)$ satisfying (4).

The first estimand we are motivated by is the ATE, given by
$\theta_{0,\textnormal{ATE}}=\mathbb{E}[Y(1)-Y(0)]\in\real$ in our notation.
For an investigator interested in modeling how the treatment effect varies
with $X$, they may instead wish to estimate the regression parameter
$\theta_{0,\textnormal{PL}}\in{}^{p}$ under a linear treatment effect
assumption

$\mathbb{E}[Y(1)-Y(0)\mid X]=\psi(X)^{\top}\theta_{0,\textnormal{PL}}.$ (6)

See Robinson, (1988) for background on semiparametric estimation of
$\theta_{0,\textnormal{PL}}$ under (6), which characterizes the well-known
“partially linear model.” We show next that both $\theta_{0,\textnormal{ATE}}$
and $\theta_{0,\textnormal{PL}}$ are identified by score functions that are
linear in the sense of (4) and robust in the sense of (5), and hence
identifiable according to Assumption 2.2.

######  Example 2.1 (ATE estimation).

Let $\theta_{0}=\theta_{0,\textnormal{ATE}}$ be the estimand of interest. Now
consider the augmented inverse propensity weighting (AIPW) score function

$s_{\textnormal{AIPW}}(W;\theta,\nu,e)=m(1,X)-m(0,X)+\frac{Z(Y-m(1,X))}{e(X)}-\frac{(1-Z)(Y-m(0,X))}{1-e(X)}-\theta$
(7)

for nuisance parameter $\nu=(m(0,\cdot),m(1,\cdot))$. For each $\gamma>0$, it
is well known that
$\mathbb{E}_{e}[s_{\textnormal{AIPW}}(W;\theta_{0},\nu_{0},e^{\prime})]=0$ for
any $e(\cdot),e^{\prime}(\cdot)$ in $\mathcal{F}_{\gamma}$ when
$\nu_{0}=\nu_{0,\textnormal{AIPW}}=(m_{0}(0,\cdot),m_{0}(1,\cdot))$ lies in
the nuisance set
$\mathcal{N}=\mathcal{N}_{\textnormal{AIPW}}=L^{1}(P^{X})\times L^{1}(P^{X})$.
Hence, $s_{\textnormal{AIPW}}(\cdot)$ satisfies the score equation (5). This
score is also linear because
$s_{\textnormal{AIPW}}=s_{\textnormal{AIPW},a}\theta+s_{\textnormal{AIPW},b}$
for

$\displaystyle s_{\textnormal{AIPW},a}(W;\nu,e)$
$\displaystyle=-1,\quad\text{and}$ $\displaystyle
s_{\textnormal{AIPW},b}(W;\nu,e)$
$\displaystyle=m(1,X)-m(0,X)+\frac{Z(Y-m(1,X))}{e(X)}-\frac{(1-Z)(Y-m(0,X))}{1-e(X)}.$

which completes the task of showing that $\theta_{0,\textnormal{ATE}}$
satisfies Assumption 2.2.

######  Example 2.2 (Partially linear model).

Suppose the linear treatment effect assumption (6) holds and
$\theta_{0}=\theta_{0,\textnormal{PL}}$ is the estimand of interest. Now
consider the weighted least squares score

$s_{\textnormal{EPL}}(W;\theta,\nu,e)=w(X;\nu,e)(Z-e(X))(Y-m(0,X)-Z\psi(X)^{\top}\theta)\psi(X)$
(8)

for nonnegative weights

$w(X,\nu,e)=(v(0,x)e(x)+v(1,x)(1-e(x)))^{-1}$

with nuisance parameter $\nu=(m(0,\cdot),v(0,\cdot),v(1,\cdot)^{\top}$. Let
$\mathcal{F}(\mathcal{X};I)$ be the set of all measurable functions
$f:\mathcal{X}\rightarrow I$. Then if
$\nu_{0}=(m_{0}(0,\cdot),v_{0}(0,\cdot),v_{0}(1,\cdot))$ lies in the nuisance
set
$\mathcal{N}=\mathcal{N}_{\textnormal{EPL}}=L^{2}(P^{X})\times\mathcal{F}(\mathcal{X};[c,\infty))\times\mathcal{F}(\mathcal{X};[c,\infty))$
for some $c>0$, we have
$\mathbb{E}_{e}[s_{\textnormal{EPL}}(W;\theta_{0},\nu_{0},e^{\prime})]=0$ for
any $e(\cdot),e^{\prime}(\cdot)$ in $\mathcal{F}_{0}$. Furthermore,
$s_{\textnormal{EPL}}(\cdot)$ is linear, because
$s_{\textnormal{EPL}}=s_{\textnormal{EPL},a}\theta+s_{\textnormal{EPL},b}$ for

$\displaystyle s_{\textnormal{EPL},a}(W;\nu,e)$
$\displaystyle=-w(X;\nu,e)Z(Z-e(X))\psi(X)\psi(X)^{\top},\quad\text{and}$
$\displaystyle s_{\textnormal{EPL},b}(W;\nu,e)$
$\displaystyle=w(X;\nu,e)(Z-e(X))(Y-m(0,X))\psi(X).$

Thus, $\theta_{0,\textnormal{EPL}}$ satisfies Assumption 2.2 with $\gamma=0$
and the score $s_{\textnormal{EPL}}(\cdot)$. Note that the score equations (5)
hold for any nonnegative weight functions $w(\cdot,\nu,e)\in L^{1}(P^{X})$,
though the specific choice in $s_{\textnormal{EPL}}(\cdot)$ is
semiparametrically efficient (Chamberlain,, 1992; Ma et al.,, 2006).

Some of our later results pertain to general estimands identified by
Assumption 2.2. Others will be specialized to the settings of Examples 2.1 and
2.2.

##  3 Oracle pooled and aggregate estimation (non-adaptive)

Here we propose and analyze an oracle estimator $\hat{\theta}^{*}$ of a
generic estimand $\theta_{0}$ satisfying Assumption 2.2 with score $s(\cdot)$.
It is an oracle in the sense that it uses the unknown true value of the
nuisance parameter $\nu_{0}$ from the score equations (5). We discuss feasible
estimation of $\theta_{0}$, including estimation of $\nu_{0}$, in Section 4.
The estimator $\hat{\theta}^{*}$ pools observations across all batches
$t=1,\ldots,T$. We then prove a central limit theorem (CLT) for it in the
setting of a _non-adaptive batch experiment_ where treatment in each batch
$t=1,\ldots,T$ is assigned according to a fixed (non-random) propensity score
$e_{t}(\cdot)$:

######  Definition 3.1 (Non-adaptive batch experiment).

A non-adaptive batch experiment has data generating process satisfying
Assumption 2.1 and treatment assignments satisfying

$Z_{ti}=\bm{1}(U_{ti}\leqslant e_{t}(X_{ti})),\quad t=1,\ldots,T,\
i=1,\ldots,N_{t}$

for some _nonrandom_ (i.e., non-adaptive) batch propensity scores
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$ and uniformly distributed random variables
$\\{U_{ti}\mid t=1,\ldots,T,i=1,\ldots,N_{t}\\}$ that are i.i.d. and
independent of the vectors $\\{S_{ti}\mid t=1,\ldots,T,i=1,\ldots,N_{t}\\}$.

Next, we compare the pooled estimator $\hat{\theta}^{*}$ to an alternative
oracle that makes an optimal linear aggregation of per-batch estimates. It
also satisfies a CLT, but we show that our pooling strategy dominates
aggregation in terms of efficiency in the setting of Examples 2.1 and 2.2. We
remark that some authors (e.g. Tabord-Meehan, (2022)) refer to this
aggregation approach as pooling, but we reserve that term for pooling data,
not estimators.

###  3.1 Pooled oracle estimator

The main idea behind the construction of our oracle estimator
$\hat{\theta}^{*}$ is as follows. After collecting the observations from all
batches $t=1,\ldots,T$ in a non-adaptive batch experiment, we ignore the batch
structure and pool together the observations across batches. Now consider a
random draw $W=(X,Z,Y)$ from these pooled observations $\\{W_{ti}\mid
1\leqslant t\leqslant T,1\leqslant i\leqslant N_{t}\\}$. In the notation of
Section 2, it is straightforward to show that the distribution of $W$ is
$P_{e_{0,N}}=\sum_{t=1}^{T}(N_{t}/N)P_{e_{t}}$, where $e_{0,N}(\cdot)$ is the
mixture propensity score

$e_{0,N}(x)=\textnormal{Pr}(Z=1\mid X)=\sum_{t=1}^{T}\frac{N_{t}}{N}e_{t}(x).$
(9)

It will also be helpful to define the limiting mixture propensity score
$e_{0}(\cdot)$ under the proportional asymptotics (3):

$e_{0}(x)=\sum_{t=1}^{T}\kappa_{t}e_{t}(x).$ (10)

When a particular set of nonrandom batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$ is relevant, we omit an additional
subscript by letting $\mathbb{E}_{0,N}[f(W)]=\mathbb{E}_{e_{0,N}}[f(W)]$ and
$\mathbb{E}_{0}[f(W)]=\mathbb{E}_{e_{0}}[f(W)]$. Using this notation, the
oracle pooled estimator $\hat{\theta}^{*}$ is derived by solving the sample
analogue of the mixture score equations
$\mathbb{E}_{0,N}[s(W;\theta_{0},\nu_{0},e_{0,N})]=0$ for $\theta$:

$\hat{\theta}^{*}=-\biggl{(}\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{a}(W_{ti};\nu_{0},e_{0,N})\biggr{)}^{-1}\biggl{(}\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{b}(W_{ti};\nu_{0},e_{0,N})\biggr{)}.$
(11)

This $\hat{\theta}^{*}$ is only defined when the matrix inverse in (11)
exists. That will be the case with probability tending to 1 so long as
$\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]$ is invertible, which we will require
in our theoretical results. One such result is a CLT for $\hat{\theta}^{*}$.
The proofs of all our technical results are provided in Appendix C.

######  Proposition 3.1 (Oracle CLT).

Let $\theta_{0}\in{}^{p}$ be an estimand satisfying Assumption 2.2 for some
$\gamma\in[0,1/2)$, some score $s(\cdot)$, and some nuisance parameters
$\nu_{0}$. Suppose observations $\\{W_{ti}\mid
t=1,\ldots,T,i=1,\ldots,N_{t}\\}$ are collected from a non-adaptive batch
experiment with batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)\in\mathcal{F}_{\gamma}$, and define
$e_{0,N}=e_{0,N}(\cdot)$ and $e_{0}=e_{0}(\cdot)$ as in (9) and (10),
respectively. Further assume the following conditions hold:

  1. 1.

For some sequence $\delta_{N}\downarrow 0$, we have

$\displaystyle\bigl{(}\mathbb{E}_{0}[\|s_{a}(W;\nu_{0},e_{0,N})-s_{a}(W;\nu_{0},e_{0})\|^{2}]\bigr{)}^{1/2}$
$\displaystyle\leqslant\delta_{N},\quad\text{and}$
$\displaystyle\bigl{(}\mathbb{E}_{0}[\|s(W;\theta_{0},\nu_{0},e_{0,N})-s(W;\theta_{0},\nu_{0},e_{0})\|^{2}]\bigr{)}^{1/2}$
$\displaystyle\leqslant\delta_{N}.$

  2. 2.

$\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]$ is invertible and
$\mathbb{E}_{0}[\|s(W;\theta_{0},\nu_{0},e_{0})\|^{2}]<\infty$.

  3. 3.

For some $q>2$ and $C<\infty$ we have
$\mathbb{E}_{0}[\|s(W;\theta_{0},\nu_{0},e_{0,N})\|^{q}]\leqslant C$ for all
sufficiently large $N$.

Then with $\hat{\theta}^{*}$ as defined in (11), we have

$\displaystyle\sqrt{N}(\hat{\theta}^{*}-\theta_{0})$
$\displaystyle\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0})$

where

$\displaystyle V_{0}$
$\displaystyle=\bigl{(}\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]\bigr{)}^{-1}\bigl{(}\mathbb{E}_{0}[s(W;\theta_{0},\nu_{0},e_{0})^{\otimes
2}]\bigr{)}\bigl{(}\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]\bigr{)}^{-1}.$

###### Proof.

See Appendix C.1. ∎

For the estimands $\theta_{0,\textnormal{ATE}}$ and
$\theta_{0,\textnormal{PL}}$, following (11) we use the scores
$s_{\textnormal{AIPW}}(\cdot)$ and $s_{\textnormal{EPL}}(\cdot)$, respectively
to derive the oracle estimates

$\displaystyle\hat{\theta}^{*}_{\textnormal{AIPW}}$
$\displaystyle=\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{\textnormal{AIPW},b}(W_{ti};\nu_{0,\textnormal{AIPW}},e_{0,N})\quad\text{
(recall $s_{\textnormal{AIPW},a}=-1$)},\quad\text{and}$ (12)
$\displaystyle\hat{\theta}^{*}_{\textnormal{EPL}}$
$\displaystyle=-\biggl{(}\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{\textnormal{EPL},a}(W_{ti};\nu_{0,\textnormal{EPL}},e_{0,N})\biggr{)}^{-1}\biggl{(}\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{\textnormal{EPL},b}(W_{ti};\nu_{0,\textnormal{EPL}},e_{0,N})\biggr{)},$
(13)

We now specialize the generic oracle CLT of Proposition 3.1 to these two
estimators under some regularity conditions.

######  Assumption 3.1.

[Regularity for estimating $\theta_{0,\textnormal{ATE}}$] For some $C<\infty$
and $q>2$, we have $(\mathbb{E}[|Y(z)|^{q}])^{1/q}\leqslant C$ and
$\mathbb{E}[Y(z)^{2}\mid X=x]\leqslant C$ for all $z=0,1$ and
$x\in\mathcal{X}$.

######  Assumption 3.2.

[Regularity for estimating $\theta_{0,\textnormal{PL}}$] For some $C<\infty$
and $q>2$, Assumption 3.1 holds. Additionally, $\|\psi(x)\|\leqslant C$ for
all $x\in\mathcal{X}$, and there exists $c>0$ such that $v_{0}(z,x)\geqslant
c$ for all $z=0,1$ and $x\in\mathcal{X}$. Finally, the linear treatment effect
assumption (6) holds.

######  Corollary 3.1 (Oracle CLT for $\hat{\theta}_{\textnormal{AIPW}}^{*}$).

Suppose Assumption 3.1 holds, and let $\\{W_{ti}\mid
t=1,\ldots,T,i=1,\ldots,N_{t}\\}$ be observations from a non-adaptive batch
experiment with batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)\in\mathcal{F}_{\gamma}$ for some $\gamma>0$.
Then
$\sqrt{N}(\hat{\theta}^{*}_{\textnormal{AIPW}}-\theta_{0,\textnormal{ATE}})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0,\textnormal{AIPW}})$
where

$\displaystyle V_{0,\textnormal{AIPW}}$
$\displaystyle=\mathbb{E}\left[\frac{v_{0}(1,X)}{e_{0}(X)}+\frac{v_{0}(0,X)}{1-e_{0}(X)}+(m_{0}(1,X)-m_{0}(0,X)-\theta_{0,\textnormal{ATE}})^{2}\right].$
(14)

###### Proof.

See Appendix C.2. ∎

######  Corollary 3.2 (Oracle CLT for $\hat{\theta}_{\textnormal{EPL}}^{*}$).

Suppose Assumption 3.2 holds, and let $\\{W_{ti}\mid
t=1,\ldots,T,i=1,\ldots,N_{t}\\}$ be observations from a non-adaptive batch
experiment with batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)\in\mathcal{F}_{0}$, where
$\mathbb{E}[e_{0}^{2}(X)(1-e_{0}(X))^{2}\psi(X)\psi(X)^{\top}]\in\mathbb{S}^{p}_{++}$.
Then
$\sqrt{N}(\hat{\theta}^{*}_{\textnormal{EPL}}-\theta_{0,\textnormal{PL}})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0,\textnormal{EPL}})$,
where

$\displaystyle V_{0,\textnormal{EPL}}$
$\displaystyle=\left(\mathbb{E}\left[\frac{e_{0}(X)(1-e_{0}(X))}{v_{0}(0,X)e_{0}(X)+v_{0}(1,X)(1-e_{0}(X))}\psi(X)\psi(X)^{\top}\right]\right)^{-1}.$
(15)

###### Proof.

See Appendix C.3. ∎

###  3.2 The aggregated oracle estimator

When the covariate space $\mathcal{X}$ is finite, the oracle pooled estimator
$\hat{\theta}^{*}_{\textnormal{AIPW}}$ is equivalent to an oracle variant of
the estimator proposed by Hahn et al., (2011). The complexity of constructing
a _feasible_ pooled estimator for an _adaptive_ experiment in more general
settings has led other authors to instead consider single batch estimators
that can lose considerable efficiency. For instance, Cytrynbaum, (2021)
proposes simply discarding the first batch in a two-batch experiment when
computing the final estimate, which is clearly inadmissible in our
proportional asymptotic regime (46). Section 3.2 of Tabord-Meehan, (2022)
suggests instead taking a linear aggregation of estimates computed separately
on each batch, as described above. We will now show that even the best
linearly aggregated oracle estimator is asymptotically dominated by our pooled
estimators $\hat{\theta}^{*}_{\textnormal{AIPW}}$ and
$\hat{\theta}^{*}_{\textnormal{EPL}}$. While Tabord-Meehan, (2022)
hypothesized in their Appendix C.2 that this may be true for ATE estimation,
they do not pursue this further as are unable to construct a feasible pooled
estimator attaining the targeted oracle variance when batch propensities are
chosen adaptively using their stratification trees. By contrast, our design
approach will allow us to construct such an estimator using our extension of
double machine learning in Section 4.

Fix a non-adaptive batch experiment with batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$. For each batch $t=1,\ldots,T$, consider an
(oracle) estimator $\hat{\theta}_{t,\textnormal{AIPW}}^{*}$ for
$\theta_{0,\textnormal{ATE}}$ computed by solving the empirical analogue of
the score equations
$\mathbb{E}_{e_{t}}[s_{\textnormal{AIPW}}(W;\theta_{0,\textnormal{ATE}},\nu_{0,\textnormal{AIPW}},e_{t})]=0$
that averages only those observations in batch $t$:

$\hat{\theta}_{t,\textnormal{AIPW}}^{*}=-\biggl{(}\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}s_{\textnormal{AIPW},a}(W_{ti};\nu_{0,\textnormal{AIPW}},e_{t})\biggr{)}^{-1}\biggl{(}\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}s_{\textnormal{AIPW},b}(W_{ti};\nu_{0,\textnormal{AIPW}},e_{t})\biggr{)}.$

By applying Proposition 3.1 with a single batch, for each $t=1,\ldots,T$ we
obtain the CLT

$\displaystyle\sqrt{N_{t}}(\hat{\theta}_{t,\textnormal{AIPW}}^{*}-\theta_{0,\textnormal{ATE}})$
$\displaystyle\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{t,\textnormal{AIPW}}),$

where

$\displaystyle V_{t,\textnormal{AIPW}}$
$\displaystyle=A_{t,\textnormal{AIPW}}^{-1}B_{t,\textnormal{AIPW}}A_{t,\textnormal{AIPW}}^{-1},\quad\text{for}$
$\displaystyle A_{t,\textnormal{AIPW}}$
$\displaystyle=\mathbb{E}_{e_{t}}[s_{\textnormal{AIPW},a}(W;\nu_{0,\textnormal{AIPW}},e_{t})],\quad\text{and}$
$\displaystyle B_{t,\textnormal{AIPW}}$
$\displaystyle=\mathbb{E}_{e_{t}}[s(W;\theta_{0,\textnormal{ATE}},\nu_{0,\textnormal{AIPW}},e_{t})^{\otimes
2}].$

Now, as stated in Slud et al., (2018), the asymptotically unbiased linear
combination of
$\hat{\theta}_{1,\textnormal{AIPW}}^{*},\ldots,\hat{\theta}_{T,\textnormal{AIPW}}^{*}$
with the smallest asymptotic covariance matrix with respect to the
semidefinite ordering is the inverse covariance weighted estimator

$\hat{\theta}_{\textnormal{AIPW}}^{*,(\mathrm{LA})}=\biggl{(}\,\sum_{t=1}^{T}\kappa_{t}V_{t,\textnormal{AIPW}}^{-1}\biggr{)}^{-1}\sum_{t=1}^{T}\kappa_{t}V_{t,\textnormal{AIPW}}^{-1}\hat{\theta}_{t,\textnormal{AIPW}}^{*}.$

This optimal linearly aggregated estimator
$\hat{\theta}_{\textnormal{AIPW}}^{*,(\mathrm{LA})}$ satisfies the CLT

$\sqrt{N}\bigl{(}\hat{\theta}_{\textnormal{AIPW}}^{*,(\mathrm{LA})}-\theta_{0,\textnormal{ATE}}\bigr{)}\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{\textnormal{AIPW}}^{(\mathrm{LA})}),\quad
V_{\textnormal{AIPW}}^{(\mathrm{LA})}=\biggl{(}\sum_{t=1}^{T}\kappa_{t}V_{t,\textnormal{AIPW}}^{-1}\biggr{)}^{-1}.$

We can similarly define the linearly aggregated oracle estimator
$\hat{\theta}_{\textnormal{EPL}}^{*,(\mathrm{LA})}$ for
$\theta_{0,\textnormal{PL}}$ based on combining per-batch estimates from the
score $s_{\textnormal{EPL}}(\cdot)$, which satisfies a CLT
$\sqrt{N}(\hat{\theta}^{*,(\mathrm{LA})}_{\textnormal{EPL}}-\theta_{0,\textnormal{PL}})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{\textnormal{EPL}}^{(\mathrm{LA})})$.
Here $V_{\textnormal{EPL}}^{(\mathrm{LA})}$ is given by replacing
$s_{\textnormal{AIPW}}(\cdot)$ with $s_{\textnormal{EPL}}(\cdot)$,
$\nu_{0,\textnormal{AIPW}}$ with $\nu_{0,\textnormal{EPL}}$, and
$\theta_{0,\textnormal{ATE}}$ with $\theta_{0,\textnormal{PL}}$ in the
definition of $V_{\textnormal{AIPW}}^{(\mathrm{LA})}$. Our main result is then
that regardless of the batch propensities,
$\hat{\theta}_{\textnormal{AIPW}}^{*,(\mathrm{LA})}$ and
$\hat{\theta}_{\textnormal{EPL}}^{*,(\mathrm{LA})}$ are asymptotically
dominated by our pooled estimators $\hat{\theta}_{0,\textnormal{AIPW}}$ and
$\hat{\theta}_{0,\textnormal{EPL}}$, respectively. This motivates our work in
Section 5 that designs for these estimators.

######  Theorem 3.1 (Pooling dominates linear aggregation).

Under the conditions of Corollary 3.1, $V_{0,\textnormal{AIPW}}\leqslant
V^{(\mathrm{LA})}_{\textnormal{AIPW}}$. Under the conditions of Corollary 3.2,
$V_{0,\textnormal{EPL}}\preccurlyeq V^{(\mathrm{LA})}_{\textnormal{EPL}}$.

###### Proof.

See Appendix C.4. ∎

##  4 Feasible pooled estimation in batch adaptive experiments

The oracle estimator $\hat{\theta}^{*}$ of (11) depends on nuisance parameters
$\nu_{0}$ that are unknown in practice. Additionally, recall that our CLT for
$\hat{\theta}^{*}$ (Proposition 3.1) holds for a non-adaptive batch
experiment. Our goal is to choose propensities adaptively in each batch to
improve precision. Therefore, we would like to develop a feasible estimator
$\hat{\theta}$ that attains the targeted asymptotic variance for experiments
where treatment is assigned _adaptively_ , even when the nuisance parameters
$\nu_{0}$ must be estimated.

As mentioned above, our construction of such a feasible estimator
$\hat{\theta}$ is based on extending the double machine learning (DML)
framework of Chernozhukov et al., (2018). The main requirements for
$\hat{\theta}$ to have the same asymptotic variance as the corresponding
oracle are convergence rate guarantees for both nuisance parameter estimates
and the adaptive propensities. Our DML extension ensures that these rate
requirements can be made sub-parametric, enabling the use of somewhat flexible
machine learning methods.

The typical DML setting assumes access to a single sample $W_{1},\ldots,W_{N}$
of i.i.d. observations. An example of this setting is a non-adaptive batch
experiment with $T=1$ and propensity $e_{0}(\cdot)$. Then a standard DML
estimator is based on two ingredients: a Neyman orthogonal score and cross-
fitting. Neyman orthogonality of the score $s(\cdot)$ at $(\nu_{0},e_{0})$
means a local insensitivity of the score equations to perturbations in
$(\nu_{0},e_{0})$ in any direction:

$\frac{\partial}{\partial\lambda}\mathbb{E}_{e_{0}}\bigl{[}s(W_{i};\theta_{0},\nu_{0}+\lambda(\nu-\nu_{0}),e_{0}+\lambda(e-e_{0}))\bigr{]}=0,\quad\forall(\nu,e)\in\mathcal{N}\times\mathcal{F}_{\gamma}.$

It is well known that the scores $s_{\textnormal{AIPW}}(\cdot)$ and
$s_{\textnormal{EPL}}(\cdot)$ are Neyman orthogonal (Chernozhukov et al.,,
2018). Given a Neyman orthogonal score $s(\cdot)$, DML proceeds by
constructing an estimator $\hat{\theta}$ by cross-fitting. In cross-fitting,
the indices $1,\ldots,N$ are partitioned into $K$ (roughly) equally sized
folds $\mathcal{I}_{1},\ldots,\mathcal{I}_{K}$. Then $\hat{\theta}$ is
computed as the solution to the empirical score equations
$N^{-1}\sum_{k=1}^{K}\sum_{i\in\mathcal{I}_{k}}s(W_{i};\hat{\theta},\hat{\nu}^{(-k)},\hat{e}^{(-k)})=0$,
where for each $k=1,\ldots,K$, $\hat{\nu}^{(-k)}$ and $\hat{e}^{(-k)}(\cdot)$
are estimates of $\nu_{0}$ and $e_{0}(\cdot)$, respectively. Each pair
$(\hat{\nu}^{(-k)},\hat{e}^{(-k)}(\cdot))$ depends only on the observations
$\\{W_{i}\mid i\notin\mathcal{I}_{k}\\}$ outside fold $k$. The sample
splitting ensures that for each $k=1,\ldots,K$, the estimates
$(\hat{\nu}^{(-k)},\hat{e}^{(-k)})$ are independent of the observations in
fold $k$. By the arguments of Chernozhukov et al., (2018), such independence
is key to guarantee that the feasible estimator $\hat{\theta}$ is equivalent
(up to first order asymptotics) to the oracle $\hat{\theta}^{*}$ solving
$N^{-1}\sum_{k=1}^{K}\sum_{i\in\mathcal{I}_{k}}s(W_{i};\hat{\theta},\nu_{0},e_{0})=0$
(cf. (11)), even when the estimates $(\hat{\nu}^{(-k)},\hat{e}^{(-k)})$
converge at sub-parametric rates.

To maintain this independence in an _adaptive_ batched experiment, we require
sample splitting at the _design_ stage, as illustrated in Figure 1, along with
convergence of the adaptive propensity scores. Our notion of a convergent
split batch adaptive experiment (CSBAE) in Definition 4.1 formalizes this. The
main idea is to split the observations in _every_ batch $t=1,\ldots,T$ into
$K$ folds. Re-using the notation above from the standard DML setting with
$T=1$, we let $\mathcal{I}_{k}$ denote the set of batch and observation
indices $(t,i)$ assigned to fold $k=1,\ldots,K$. Then the adaptive propensity
used to assign treatment to a subject in batches $t=2,\ldots,T$ is allowed to
only depend on observations in previous batches from the same fold as this
subject. To ensure that the adaptivity does not introduce any additional
variability (up to first-order asymptotics) into the final estimator, a CSBAE
requires these adaptive propensities to converge to nonrandom limits
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$ at RMS rate $O_{p}(N^{-1/4})$. While this
convergence requirement may appear restrictive, in Section 5 we show how it
can be ensured by design by solving an appropriate finite-dimensional concave
maximization procedure. Moreover, the limiting propensity scores from this
procedure will be provably optimal, in a sense we make more precise in Section
5.

Figure 1: A graphical representation of the dependencies among the
observations in two batches of a CSBAE (Definition 4.1). Note the lack of
vertical arrows, indicating independence between the observations in fold $j$
and fold $k$ for any distinct $j$, $k$ in $\\{1,\ldots,K\\}$.

######  Definition 4.1 (Convergent split batch adaptive experiment).

A convergent split batch adaptive experiment (CSBAE) is an experiment with
data generating process satisfying Assumption 2.1 where each observation index
$(t,i)$ is assigned to one of $K$ folds
$\mathcal{I}_{1},\ldots,\mathcal{I}_{K}$. The fold assignments are such that
$n_{t,k}=|\\{(t,i)\in\mathcal{I}_{k}\mid i=1,\ldots,N_{t}\\}|$, the number of
observations in batch $t$ assigned to fold $k$, satisfies
$|n_{t,k}-N_{t}/K|\leqslant 1$ for all $t=1,\ldots,T$, $k=1,\ldots,K$. Now let
$P_{N,t}^{X,(k)}$ be the empirical distribution on
$\\{X_{ti}\mid(t,i)\in\mathcal{I}_{k},1\leqslant i\leqslant N_{t}\\}$, the
covariates in batch $t$ and fold $k$, and define $\mathcal{S}_{t}^{X,(k)}$ to
be the $\sigma$-algebra generated by the covariates
$\\{X_{ti}\mid(t,i)\in\mathcal{I}_{k},1\leqslant i\leqslant N_{t}\\}$ in batch
$t$ and fold $k$ along with the observations
$\\{W_{ui}\mid(u,i)\in\mathcal{I}_{k},u=1,\ldots,t-1,i=1,\ldots,N_{u}\\}$ in
fold $k$ and any of the previous batches $1,\ldots,t-1$. We further require
the following for each batch $t=1,\ldots,T$ and fold $k=1,\ldots,K$:

  1. 1.

Treatment is assigned according to an adaptive propensity
$\hat{e}_{t}^{(k)}(\cdot)$ that is measurable with respect to
$\mathcal{S}_{t}^{X,(k)}$. That is, the treatment indicators can be
represented as

$Z_{ti}=\bm{1}(U_{ti}\leqslant\hat{e}_{t}^{(k)}(X_{it})),\quad(t,i)\in\mathcal{I}_{k}$
(16)

where $\\{U_{ti}:1\leqslant t\leqslant T,1\leqslant i\leqslant N_{t}\\}$ is a
collection of i.i.d. uniformly distributed random variables independent of the
vectors $\\{S_{ti}\mid t=1,\ldots,T,i=1,\ldots,N_{t}\\}$.

  2. 2.

For some nonrandom propensity $e_{t}(\cdot)$, the adaptive propensity
$\hat{e}_{t}^{(k)}(\cdot)$ satisfies

$\|\hat{e}_{t}^{(k)}-e_{t}\|_{2,P_{N,t}^{X,(k)}}=O_{p}(N^{-1/4}),\quad
t=1,\ldots,T,\ k=1,\ldots,K.$ (17)

######  Remark 4.1.

The left-hand side of equation (17) uses an $L^{2}$ norm on the empirical
distribution $P_{N,t}^{X,(k)}$ of the covariates of the subjects that will be
assigned treatment according to the learned propensity
$\hat{e}_{t}^{(k)}(\cdot)$. These covariates will also be used to learn
$\hat{e}_{t}^{(k)}(\cdot)$ itself in our propensity learning procedure of
Section 5. Thus, we can interpret (17) as a rate requirement on the “in-
sample” convergence of $\hat{e}_{t}^{(k)}(\cdot)$.

Given a CSBAE, our feasible estimator is

$\hat{\theta}=-\biggl{(}\frac{1}{N}\sum_{k=1}^{K}\sum_{(t,i)\in\mathcal{I}_{k}}s_{a}(W_{ti};\hat{\nu}^{(-k)},\hat{e}^{(-k)})\biggr{)}^{-1}\biggl{(}\frac{1}{N}\sum_{k=1}^{K}\sum_{(t,i)\in\mathcal{I}_{k}}s_{b}(W_{ti};\hat{\nu}^{(-k)},\hat{e}^{(-k)})\biggr{)}.$
(18)

As in the standard (single batch) DML setting, for each $k=1,\ldots,K$,
$\hat{\nu}^{(-k)}$ and $\hat{e}^{(-k)}(\cdot)$ are estimates of the nuisance
parameters $\nu_{0}$ and the mixture propensity $e_{0,N}(\cdot)$ defined in
(9), respectively, that depend only on the observations
$\\{W_{ti}\mid(t,i)\notin\mathcal{I}_{k}\\}$ outside fold $k$. These
observations are fully independent of the observations in fold $k$ (across
_all_ batches $t=1,\ldots,T$) by the construction of a CSBAE. As in the single
batch case, given $o_{p}(N^{-1/4})$ convergence of the estimators
$\hat{\nu}^{(-k)}$ to $\nu_{0}$, this independence along with (17) enable a
DML-style argument that $\hat{\theta}$ is asymptotically equivalent to the
oracle $\hat{\theta}^{*}$ under (17) computed on a counterfactual non-adaptive
batch experiment with propensities $e_{1}(\cdot),\ldots,e_{T}(\cdot)$. This
argument proceeds by coupling the treatment indicators $Z_{ti}$ in the CSBAE
with counterfactual treatment indicators
$\tilde{Z}_{ti}=\bm{1}(U_{ti}\leqslant e_{t}(X_{ti}))$.

######  Assumption 4.1 (Score properties and convergence rates for estimating
nuisance parameters and the mixture propensity in a CSBAE).

Observations $\\{W_{ti}:1\leqslant t\leqslant T,1\leqslant i\leqslant
N_{t}\\}$ are collected from a CSBAE with limiting batch propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$. Additionally, the estimand $\theta_{0}$ of
interest is identified as in Assumption 2.2 by some score $s(\cdot)$, nuisance
parameters $\nu_{0}\in\mathcal{N}$, and $\gamma\in[0,1/2)$, such that the
propensity collection $\mathcal{F}_{\gamma}$ contains
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$. Defining $W_{ti}(z)=(Y_{ti}(z),X_{ti},z)$
for $z=0,1$, the score $s(\cdot)$ has the following properties:

  1. (a)

The matrix $\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]$ is invertible and
$\mathbb{E}_{0}\bigl{[}\|s(W;\theta_{0},\nu_{0},e_{0})\|^{2}\bigr{]}<\infty$.

  2. (b)

The mapping
$\lambda\mapsto\mathbb{E}_{0,N}[s(W;\theta_{0},\nu_{0}+\lambda(\nu-\nu_{0}),e_{0,N}+\lambda(e-e_{0,N})]$
is twice continuously differentiable on $[0,1]$ for each
$(\nu,e)\in\mathcal{T}=\mathcal{N}\times\mathcal{F}_{\gamma}$.

  3. (c)

All propensities $e(\cdot)\in\mathcal{F}_{\gamma}$ satisfy

$\mathbb{E}\bigl{[}s(W_{ti}(1);\theta_{0},\nu_{0},e)-s(W_{ti}(0);\theta_{0},\nu_{0},e)\\!\mid\\!X_{ti}\bigr{]}=0,$
(19)

for all $t=1,\dots,T$ and $i=1,\dots,N_{t}$.

Also, there exist estimators $\hat{\nu}^{(-k)}$ and $\hat{e}^{(-k)}(\cdot)$ of
$\nu_{0}$ and the mixture propensity $e_{0,N}(\cdot)$ defined in (9),
respectively, that depend only on the observations outside fold $k$ of the
CSBAE. Next, there are nonrandom subsets $\mathcal{T}_{N}\subseteq\mathcal{T}$
containing $(\nu_{0},e_{0,N}(\cdot))$, such that for all $k=1,\ldots,K$,
$\textnormal{Pr}((\hat{\nu}^{(-k)},\hat{e}^{(-k)}(\cdot))\in\mathcal{T}_{N})\rightarrow
1$ as $N\rightarrow\infty$. The sets $\mathcal{T}_{N}$ shrink quickly enough
for the following to hold for all $(\nu,e)\in\mathcal{T}_{N}$, all
$\lambda\in(0,1)$ and all $z\in\\{0,1\\}$ when $N$ is sufficiently large:

$\displaystyle\bigg{\|}\frac{\partial}{\partial\lambda}\mathbb{E}_{0,N}[s(W;\theta_{0},\nu_{0}+\lambda(\nu-\nu_{0}),e_{0,N}+\lambda(e-e_{0,N}))]\Big{|}_{\lambda=0}\bigg{\|}$
$\displaystyle\leqslant N^{-1/2}\delta_{N}$ (20)
$\displaystyle\bigg{\|}\frac{\partial^{2}}{\partial\lambda^{2}}\mathbb{E}_{0,N}[s(W;\theta_{0},\nu_{0}+\lambda(\nu-\nu_{0}),e_{0,N}+\lambda(e-e_{0,N}))]\bigg{\|}$
$\displaystyle\leqslant N^{-1/2}\delta_{N}$ (21)
$\displaystyle\left(\mathbb{E}_{0}\bigl{[}\|s_{a}(W;\nu,e)-s_{a}(W;\nu_{0},e_{0})\|^{2}\bigr{]}\right)^{1/2}$
$\displaystyle\leqslant\delta_{N}$ (22)
$\displaystyle\left(\mathbb{E}_{0}\bigl{[}\|s(W;\theta_{0},\nu,e)-s(W;\theta_{0},\nu_{0},e_{0})\|^{2}\bigr{]}\right)^{1/2}$
$\displaystyle\leqslant\delta_{N}$ (23)
$\displaystyle\bigl{(}\mathbb{E}_{0}\bigl{[}\|s_{a}(W(z);\nu,e)\|^{q}\bigr{]}\bigr{)}^{1/q}$
$\displaystyle\leqslant C$ (24)
$\displaystyle\bigl{(}\mathbb{E}_{0}\bigr{[}\|s(W(z);\theta_{0},\nu,e\|^{q}\bigr{]}\bigr{)}^{1/q}$
$\displaystyle\leqslant C.$ (25)

Finally, letting $\mathcal{S}^{(-k)}$ be the $\sigma$-algebra generated by the
observations $\\{W_{ti}:(t,i)\notin\mathcal{I}_{k}\\}$ outside fold $k$ across
all batches 1 through $T$, we require

$S_{t}^{(-k)}(z)=o_{p}(N^{-1/4}),\quad z\in\\{0,1\\},\quad k=1,\ldots,K$ (26)

where

$S_{t}^{(-k)}(z)=\sqrt{\frac{1}{n_{t,k}}\sum_{i:(t,i)\in\mathcal{I}_{k}}\Bigl{\|}\mathbb{E}\bigl{[}s(W_{ti}(z);\hat{\nu}^{(-k)},\hat{e}^{(-k)})-s(W_{ti}(z);\nu_{0},e_{0,N})\bigm{|}\mathcal{S}^{(-k)},X_{ti}\bigr{]}\Bigr{\|}^{2}}.$

######  Theorem 4.1 (Feasible CLT for a CSBAE).

Suppose Assumption 4.1 holds. Then for $\hat{\theta}$ defined in (18) there
exists a non-adaptive batch experiment with propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$ for which

$\hat{\theta}=\hat{\theta}^{*}+o_{p}(N^{-1/2}).$

Here $\hat{\theta}^{*}$ is the oracle (11) computed on this non-adaptive batch
experiment. Then
$\sqrt{N}(\hat{\theta}-\theta_{0})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0})$
where $V_{0}$ is the limiting covariance derived in Proposition 3.1.

###### Proof.

See Appendix C.5. ∎

In Assumption 4.1, the conditions (a) and (b) along with the inequalities (20)
through (25) are direct extensions of Assumptions 3.1 and 3.2 in Chernozhukov
et al., (2018) for ordinary DML ($T=1$). The equations (19) and (26) are
additional requirements that enable the dependence across batches in a CSBAE
to be sufficiently weak so that $\hat{\theta}$ computed on the CSBAE is
asymptotically equivalent to a version of it computed on the limiting non-
adaptive batch experiment.

We can show that Assumption 4.1 is satisfied for estimating
$\theta_{0,\textnormal{ATE}}$ and $\theta_{0,\textnormal{PL}}$ with
$s_{\textnormal{AIPW}}(\cdot)$ and $s_{\textnormal{EPL}}(\cdot)$, respectively
under simple rate conditions on nuisance parameter and propensity estimation
rates that mirror those in Section 5 of Chernozhukov et al., (2018) for
single-batch DML estimators. Then we apply Theorem 4.1 to construct feasible
pooled estimators $\hat{\theta}_{\textnormal{AIPW}}$ and
$\hat{\theta}_{\textnormal{EPL}}$ as special cases of equation (18), which
attain the oracle asymptotic variances $V_{0,\textnormal{AIPW}}$ and
$V_{0,\textnormal{EPL}}$ defined in (14) and (15), respectively.

######  Corollary 4.1 (Feasible pooled estimation of
$\theta_{0,\textnormal{ATE}}$ in a CSBAE).

Suppose observations are collected from a CSBAE for which the regularity
conditions of Assumption 3.1 hold for some $q>2$ and $C<\infty$ and the
limiting batch propensities $e_{1}(\cdot),\ldots,e_{T}(\cdot)$ are in
$\mathcal{F}_{\gamma}$ for some $\gamma>0$. Additionally, for each
$k=1,\ldots,K$, suppose we have estimates $\hat{m}^{(-k)}(\cdot)$ and
$\hat{e}^{(-k)}(\cdot)$ of the mean function $m_{0}(\cdot)$ and mixture
propensity $e_{0,N}(\cdot)$, respectively, both depending only on the
observations outside fold $k$, such that the following are true:

  1. 1.

$\|\hat{m}^{(-k)}(z,\cdot)-m_{0}(z,\cdot)\|_{2,P^{X}}=o_{p}(N^{-1/4}),\quad
z=0,1,$

  2. 2.

$\|\hat{e}^{(-k)}-e_{0,N}\|_{2,P^{X}}=O_{p}(N^{-1/4}),$

  3. 3.

$\|\hat{m}^{(-k)}(z,\cdot)-m_{0}(z,\cdot)\|_{q,P^{X}}\leqslant C,\quad z=0,1$,

  4. 4.

$\hat{e}^{(-k)}(\cdot)\in\mathcal{F}_{\gamma}$ with probability tending to 1.

Then
$N^{1/2}(\hat{\theta}_{\textnormal{AIPW}}-\theta_{0,\textnormal{ATE}})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0,\textnormal{AIPW}})$.

###### Proof.

See Appendix C.6. ∎

######  Corollary 4.2 (Feasible estimation of $\theta_{0,\textnormal{PL}}$ in
a CSBAE).

Suppose observations are collected from a CSBAE for which the regularity
conditions of Assumption 3.2 hold for some $q>2$ and $0<c<C<\infty$, and for
which the limiting batch propensities $e_{1}(\cdot),\ldots,e_{T}(\cdot)$ are
in $\mathcal{F}_{0}$ with
$\mathbb{E}[e_{0}^{2}(X)(1-e_{0}(X))^{2}\psi(X)\psi(X)^{\top}]\in\mathbb{S}^{p}_{++}$.
For each $k=1,\ldots,K$, assume we have estimates $\hat{m}^{(-k)}(0,\cdot)$,
$\hat{v}^{(-k)}(\cdot,\cdot)$, and $\hat{e}^{(-k)}(\cdot)$ of the mean
function $m_{0}(0,\cdot)$, the variance function $v_{0}(\cdot,\cdot)$, and the
mixture propensity $e_{0,N}(\cdot)$, respectively, all depending only on the
observations outside fold $k$, such that the following are true:

  1. 1.

$\|\hat{m}^{(-k)}(0,\cdot)-m_{0}(0,\cdot)\|_{2,P^{X}}=o_{p}(N^{-1/4})$,

  2. 2.

$\|\hat{v}^{(-k)}(z,\cdot)-v_{0}(z,\cdot)\|_{2,P^{X}}=o_{p}(1),\quad z=0,1$,

  3. 3.

$\|\hat{e}^{(-k)}-e_{0,N}\|_{2,P^{X}}=O_{p}(N^{-1/4})$,

  4. 4.

$\|\hat{m}^{(-k)}(0,\cdot)-m_{0}(0,\cdot)\|_{q,P^{X}}\leqslant C,$

  5. 5.

$\inf_{x\in\mathcal{X}}\hat{v}^{(-k)}(z,x)\geqslant c,\quad z=0,1$.

Then
$N^{1/2}(\hat{\theta}_{\textnormal{EPL}}-\theta_{0,\textnormal{PL}})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{0,\textnormal{EPL}})$.

###### Proof.

See Appendix C.7. ∎

We now compare the rate requirements in Corollaries 4.1 and 4.2 with those
needed to prove a feasible CLT for the linearly aggregated estimators
discussed in Section 3. Consider a batch adaptive experiment _without_ sample
splitting, so that we assign treatment using the propensities
$\hat{e}_{1}(\cdot),\ldots,\hat{e}_{T}(\cdot)$ where each $\hat{e}_{t}(\cdot)$
is possibly random but can only depend on the observations in batches
$1,\ldots,t-1$, and converges to some nonrandom $e_{t}(\cdot)$ at _any_ rate;
in particular, this rate may be slower than the $O_{p}(N^{-1/4})$ rate of
(17). Then the standard DML results in Section 5 of Chernozhukov et al.,
(2018) show that we can construct feasible estimators
$\hat{\theta}_{t,\textnormal{AIPW}}$ and $\hat{\theta}_{t,\textnormal{EPL}}$
that are asymptotically equivalent to the oracle single-batch estimators
$\hat{\theta}_{t,\textnormal{AIPW}}^{*}$ and
$\hat{\theta}_{t,\textnormal{EPL}}^{*}$ in Section 3, so long as we plug in
the true propensity score $\hat{e}_{t}(\cdot)$ used to assign treatment and
use cross-fitting at the _estimation_ stage, even if all rates in Corollaries
4.1 and 4.2 are weakened to $o_{p}(1)$. Unfortunately, we cannot extend this
construction to our pooled estimator $\hat{\theta}$ when $T\geqslant 2$, since
$\hat{\theta}$ plugs in estimates of the mixture propensity $e_{0,N}(\cdot)$,
which does not correspond (in general) to a propensity actually used for
treatment in any batch. Our numerical studies in Section 6 suggest, however,
that the weaker rate requirements for a feasible CLT for the linearly
aggregated estimators (compared to our pooled estimators) make little
difference in practice. Indeed, there we also see finite sample advantages to
pooling, beyond those predicted by the asymptotics.

##  5 Batch adaptive learning of the optimal propensity score

We now discuss how to learn adaptive propensity scores
$\hat{e}_{t}^{(k)}(\cdot)$ that satisfy (17) with limiting propensity scores
$e_{t}(\cdot)$ that maximize asymptotic precision of the final estimators
$\hat{\theta}_{\textnormal{AIPW}}$ and $\hat{\theta}_{\textnormal{EPL}}$
constructed in the previous section, as measured by their asymptotic
covariance matrices $V_{0,\textnormal{AIPW}}$ and $V_{0,\textnormal{EPL}}$.
This generates a CSBAE on which, by Theorem 4.1, feasible estimators
$\hat{\theta}_{0,\textnormal{AIPW}}$ and $\hat{\theta}_{0,\textnormal{EPL}}$
achieve the targeted asymptotic variances $V_{0,\textnormal{AIPW}}$ and
$V_{0,\textnormal{EPL}}$.

While $V_{0,\textnormal{AIPW}}$ is scalar and so there is no ambiguity in what
it means to maximize asymptotic precision of
$\hat{\theta}_{\textnormal{AIPW}}$, when $\theta_{0,\textnormal{PL}}$ is
multivariate ($p>1$), $V_{0,\textnormal{EPL}}$ is a matrix. To handle the
multivariate setting, we follow classical literature on experiment design in
regression models by scalarizing an _information matrix_
$\mathcal{I}\in\mathbb{S}_{+}^{p}$ (typically an inverse covariance matrix)
using an information function
$\Psi:\mathbb{S}_{+}^{p}\rightarrow\real\cup\\{-\infty\\}$. See textbooks such
as Atkinson et al., (2007); Pukelsheim, (2006) for more background. We
generically write $\mathcal{I}=\mathcal{I}(e,\eta)$ to emphasize that in our
setting, $\mathcal{I}$ will be indexed by a propensity score $e=e(\cdot)$ in
some function class $\mathcal{F}_{*}$ along with some unknown nuisance
functions $\eta=\eta(\cdot)$ to be estimated. Then the design objective is to
learn an optimal propensity $e^{*}(\cdot)$, in the sense of maximizing
$\Psi(\mathcal{I})$:

$e^{*}(\cdot)\in\operatorname*{arg\,max}_{e\in\mathcal{F}_{*}}\Psi(\mathcal{I}(e;\eta)).$
(27)

The nuisance functions $\eta$ in the information matrix are distinct from the
nuisance functions $\nu$ in the score function. Examples of $\eta$ for
information matrices based on $V_{0,\textnormal{AIPW}}$ and
$V_{0,\textnormal{EPL}}$ are given in Section 5.2.

###  5.1 Generic convergence rates with concave maximization

We start by considering a feasible procedure to learn the propensity
$e^{*}(\cdot)$ of (27) when the information matrix $\mathcal{I}$ and the
function class $\mathcal{F}_{*}$ being optimized over generically satisfy
Assumption 5.1 below. An explanation of how to apply this generic procedure to
construct an appropriate CSBAE that designs for the estimators
$\hat{\theta}_{\textnormal{AIPW}}$ and $\hat{\theta}_{\textnormal{EPL}}$ is
deferred to Section 5.2. Here we instead focus on a precise exposition of our
propensity learning procedure and the technical assumptions needed to ensure
convergence, without explicitly invoking any notation or setup from previous
sections. We begin with some generic structure on the information matrix
$\mathcal{I}$ and the function class $\mathcal{F}_{*}$. Roughly, we require
the information matrix to be strongly concave in $e(\cdot)$ and the function
class $\mathcal{F}_{*}$ to be not too complex. We also allow for interval
constraints on the per-batch expected proportion of subjects treated.

######  Assumption 5.1 (Generic optimization setup).

The information matrix $\mathcal{I}=\mathcal{I}(e,\eta)$ in (27) takes the
form

$\mathcal{I}(e;\eta)=\mathbb{E}_{X\sim P}\left[f(e(X),\eta(X))\right]$ (28)

where $\eta:\mathcal{X}\rightarrow\mathcal{W}$ is a vector of possibly unknown
functions taking values in some compact set $\mathcal{W}\subseteq{}^{r}$, and
$P$ is a covariate distribution from which i.i.d. observations
$X_{1},\ldots,X_{n}\in\mathcal{X}$ are drawn. Additionally, the function
$f:[0,1]\times\mathcal{W}\rightarrow\mathbb{S}_{+}^{p}$ appearing in (28)
satisfies the following properties:

  * •

There exists an extension of $f(\cdot,\cdot)$ with continuous second partial
derivatives to an open neighborhood of r+1 containing
$(\delta,1-\delta)\times\mathcal{W}$ for some $\delta<0$.

  * •

For some $c>0$ we have

$-f^{\prime\prime}(e,w)\in\mathbb{S}_{+}^{p},\quad\text{and}\quad\textnormal{tr}(f^{\prime\prime}(e,w))\leqslant-c,\quad\forall(e,w)\in[0,1]\times\mathcal{W}$
(29)

where $f^{\prime\prime}(e,w)$ denotes the second partial derivative of
$f(\cdot,\cdot)$ with respect to the first argument, evaluated at $(e,w)$.

Next, the collection of propensity scores $\mathcal{F}_{*}$ to be optimized
over takes the form
$\mathcal{F}_{*}=\mathcal{F}_{*}(m_{L},m_{H};P)=\\{e\in\mathcal{E}\mid
m_{L}\leqslant\mathbb{E}_{P}[e(X)]\leqslant m_{H}\\}$, for some base
propensity class $\mathcal{E}\subseteq\mathcal{F}_{0}$ and some known budget
constraints $(m_{L},m_{H})$ with $0\leqslant m_{L}\leqslant m_{H}\leqslant 1$.
The base propensity class $\mathcal{E}$ is convex and closed in $L^{2}(P)$ and
additionally satisfies the following properties for all $n\geqslant 1$:

  1. 1.

There exists $C<\infty$ for which

$\int_{0}^{1}\sqrt{\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))}\,\mathrm{d}\epsilon\leqslant
C\quad\forall n\quad\text{w.p.1}$ (30)

where $\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))$ is the metric
entropy of the base propensity class $\mathcal{E}$ in $L^{2}(P_{n})$, as
defined in Appendix B.1, and $P_{n}$ is the empirical distribution on the
observations $X_{1},\ldots,X_{n}$.

  2. 2.

Given any $x_{1},\ldots,x_{n}\in\mathcal{X}$, there exists a convex set
$E_{n}\subseteq[0,1]^{n}$, possibly depending on $x_{1},\ldots,x_{n}$, such
that

     * •

For every $e\in\mathcal{E}$, we must have
$(e(x_{1}),e(x_{2}),\ldots,e(x_{n}))\in E_{n}$, and

     * •

For every $(e_{1},e_{2},\ldots,e_{n})\in E_{n}$, there exists
$e\in\mathcal{E}$ with $e(x_{i})=e_{i}$ for all $i=1,\ldots,n$.

  3. 3.

There exists $e_{*}(\cdot)\in\mathcal{F}_{*}$ such that
$\mathbb{E}_{P}[f(e_{*}(X),\eta(X))]\succcurlyeq cI$ for some $c>0$.

  4. 4.

There exist $e_{L}(\cdot),e_{H}(\cdot)\in\mathcal{E}$ such that
$\mathbb{E}_{P}[e_{L}(X)]>m_{L}$ and $\mathbb{E}_{P}[e_{H}(X)]<m_{H}$.

######  Remark 5.1.

If (29) only holds for $(e,w)$ in $[\gamma,1-\gamma]\times\mathcal{W}$ for
some $\gamma\in(0,1/2)$ (instead of on all of $[0,1]\times\mathcal{W}$) and
the base propensity class $\mathcal{E}$ is chosen to lie within
$\mathcal{F}_{\gamma}$, then one can rewrite (28) as

$\mathcal{I}(e;\eta)=\mathbb{E}_{X\sim
P}[\tilde{f}(e(X),\eta(X))],\quad\tilde{f}(e,w)=f(L_{\gamma}(e),w)$

where $L_{\gamma}(e)=\gamma+(1-2\gamma)e$ is an invertible linear mapping.
Then the remainder of Assumption 5.1 holds as stated with $f$ replaced by
$\tilde{f}$ and the base propensity class $\mathcal{E}$ replaced by
$\tilde{\mathcal{E}}=\\{L_{\gamma}^{-1}(e(\cdot))\mid
e(\cdot)\in\mathcal{E}\\}$, and so all results below depending on Assumption
5.1 hold by considering optimization over $\tilde{\mathcal{E}}$ instead of
$\mathcal{E}$.

Similar to the idea of empirical risk minimization in supervised learning
(ERM, e.g. Vapnik, (1991) and Montanari and Saeed, (2022)), when Assumption
5.1 is satisfied, our learning procedure replaces the unknown population
expectation appearing in the objective (27) by a sample average over
observations
$X_{1},\ldots,X_{n}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}P$, where the
generic sample size $n$ diverges. In our designed CSBAE, these observations
will be the covariates of those subjects in a given batch $t=1,\ldots,T$ of a
CSBAE within a given fold $k=1,\ldots,K$, so $n$ will be identified with the
quantity $n_{t,k}$ of Definition 4.1. The main computational procedure is then
a finite dimensional optimization over the values of the propensity score at
the $n$ points $X_{1},\ldots,X_{n}$:

$(\hat{e}_{1},\ldots,\hat{e}_{n})\in\operatorname*{arg\,max}_{(e_{1},\ldots,e_{n})\in
F_{n}}\Psi\biggl{(}n^{-1}\sum_{i=1}^{n}f(e_{i},\hat{\eta}(X_{i}))\biggr{)}.$
(31)

Here $\hat{\eta}(\cdot)$ is an estimate of the nuisance function $\eta(\cdot)$
in (28), and the optimization set $F_{n}=F_{n}(m_{L},m_{H})\subseteq{}^{n}$ is
defined by

$F_{n}(m_{L},m_{H})=\biggl{\\{}(e_{1},\ldots,e_{n})\in E_{n}\mid
m_{L}\leqslant\frac{1}{n}\sum_{i=1}^{n}e_{i}\leqslant
m_{H}\biggr{\\}}\subseteq{}^{n}$ (32)

where $(m_{L},m_{H})$ are budget constraints as in Assumption 5.1 and $E_{n}$
is as in the numbered condition 2 of that assumption.

We convert the vector $(\hat{e}_{1},\ldots,\hat{e}_{n})$ in (31) to a
propensity score $\hat{e}(\cdot)$ by taking $\hat{e}(\cdot)$ to be any member
of the base propensity class $\mathcal{E}$ of Assumption 5.1 with
$\hat{e}(X_{i})=\hat{e}_{i}$ for each $i=1,\ldots,n$. The existence of such a
propensity $\hat{e}(\cdot)$ is guaranteed by the numbered condition 2 in
Assumption 5.1. Then as in the ERM literature, we use empirical process
arguments to guarantee the learned propensity $\hat{e}(\cdot)$ converges to
$e^{*}(\cdot)$ at rate $O_{p}(n^{-1/4})$ under appropriate restrictions on the
complexity of the base propensity class $\mathcal{E}$. Our main such
restriction is the finite entropy integral requirement in (30). Examples of
function classes satisfying this condition can be found in the literature on
empirical processes. A sufficient condition is that
$\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))\leqslant
K\epsilon^{-2+\delta}$ for some $\delta>0$, some $K<\infty$, and all
$\epsilon>0$. Examples 5.1 through 5.4 below show that this requirement is
loose enough to admit some relatively rich function classes.

######  Example 5.1 (Monotone).

Suppose $d=1$ and let $\mathcal{E}$ be the set of nondecreasing functions in
$\mathcal{F}_{0}$. By Lemma 9.11 of Kosorok, (2008), we know that
$\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))\leqslant K/\epsilon$ for
each $\epsilon>0$ and some positive universal constant $K<\infty$.

######  Example 5.2 (Lipschitz).

Again suppose $d=1$ and let $\mathcal{E}$ be the set of $L$-Lipschitz
functions in $\mathcal{F}_{0}$ for some fixed $L>0$. If $\mathcal{X}$ is a
bounded closed interval, then the discussion preceding Example 5.11 of
Wainwright, (2019) shows that
$\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))\leqslant K/\epsilon$ for
each $\epsilon>0$ and some positive universal constant $K<\infty$ (which may
depend on $L$).

######  Example 5.3.

(VC-subgraph class) Let $\mathcal{E}$ be any subset of $\mathcal{F}_{0}$ that
is closed and convex in $L^{2}(P)$, and whose subgraphs are a Vapnik-
Chervonenkis (VC) class, meaning they have a finite VC dimension $V$. A
special case is a fully parametric class like
$\\{x\mapsto\theta^{\top}\xi(x)\mid\theta^{\top}\bm{1}_{p}\leqslant
1,\theta\succcurlyeq 0\\}$ where $\xi(x)\in[0,1]^{p}$ is a known set of basis
functions and $\bm{1}_{p}:=(1,\ldots,1)\in{}^{p}$. Then by Theorem 2.6.7 of
van der Vaart and Wellner, (1996),
$\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))\leqslant K(1/\epsilon)^{2V-2}$
for some universal $K>0$ depending on the VC dimension $V$ of the subgraphs.
Note that $K$ may depend on $p$.

######  Example 5.4.

(Symmetric convex hull of VC-subgraph class) Let $\mathcal{E}_{0}$ be a VC-
subgraph class of functions. The symmetric convex hull of $\mathcal{E}_{0}$ is
defined as

$\textnormal{sconv}(\mathcal{E}_{0})=\biggl{\\{}\,\sum_{i=1}^{m}\omega_{i}e_{i}\Bigm{|}e_{i}\in\mathcal{E}_{0},\sum_{i=1}^{m}|\omega_{i}|\leqslant
1\biggr{\\}}.$

Now suppose $\mathcal{E}$ is contained within
$\overline{\textnormal{sconv}(\mathcal{E}_{0})}$, the pointwise closure of
$\textnormal{sconv}(\mathcal{E}_{0})$. Let $V<\infty$ be the VC dimension of
the collection of subgraphs of functions in $\mathcal{E}_{0}$. Then by Theorem
2.6.9 of van der Vaart and Wellner, (1996) we have
$\log\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P_{n}))\leqslant
K(1/\epsilon)^{2(1-1/V)}$ for all $\epsilon>0$. For example, with
$\textnormal{expit}(x)=\exp(x)/(1+\exp(x))$ we can take

$\mathcal{E}=\biggl{\\{}\,\sum_{i=1}^{m}\omega_{i}\textnormal{expit}(\theta_{i}^{\top}\varphi(x))\Bigm{|}\omega_{i}\geqslant
0,\sum_{i=1}^{m}\omega_{i}\leqslant 1\biggr{\\}}$ (33)

where $\varphi(\cdot)$ is any vector of $p$ real-valued basis functions, $m$
can be made arbitrarily large, and $\theta_{1},\ldots,\theta_{m}\in{}^{p}$ are
arbitrary. This choice of $\mathcal{E}$ is evidently a closed and convex
subset of $\textnormal{sconv}(\mathcal{E}_{0})$ with
$\mathcal{E}_{0}=\\{\textnormal{expit}(\theta^{\top}\varphi(x))\mid\theta\in{}^{p}\\}$.
Note the collection $\mathcal{E}_{0}$ is indeed a VC-subgraph class by Lemmas
2.6.15 and 2.6.17 of van der Vaart and Wellner, (1996), as each function in
$\mathcal{E}_{0}$ is the composition of the monotone function
$\textnormal{expit}(\cdot)$ with the $p$-dimensional vector space of functions
$\\{x\mapsto\theta^{\top}\varphi(x)\mid\theta\in{}^{p}\\}$.

Next, we construct sets $E_{n}$ that satisfy condition 2 of Assumption 5.1 for
the base propensity classes in Examples 5.1 to 5.4. For the set of monotone
functions in one dimension (Example 5.1) we can take

$E_{n}=\\{(e_{1},\ldots,e_{n})\mid 0\leqslant e_{\pi(1)}\leqslant
e_{\pi(2)}\leqslant\cdots\leqslant e_{\pi(n)}\leqslant 1\\}$

where $\pi(\cdot)$ is the inverse of the function that maps each
$i\in\\{1,\ldots,n\\}$ to the rank of $x_{i}$ among $x_{1},\ldots,x_{n}$ (with
any ties broken in some deterministic way). For the set of $L$-Lipschitz
functions (Example 5.2) we can take

$E_{n}=\bigl{\\{}(e_{1},\ldots,e_{n})\bigm{|}|e_{\pi(i)}-e_{\pi(i-1)}|\leqslant
L(x_{\pi(i)}-x_{\pi(i-1)}),\,i=2,\dots,n\bigr{\\}}.$

For the parametric class in Example 5.3 we can take

$E_{n}=\bigl{\\{}(\theta^{\top}\xi(x_{1}),\ldots,\theta^{\top}\xi(x_{n}))\mid\theta^{\top}\bm{1}_{p}\leqslant
1,\theta\succcurlyeq 0\bigr{\\}}.$

Finally, for the class (33) in Example 5.4 we take

$E_{n}=\biggl{\\{}\,\sum_{i=1}^{m}\omega_{i}\bigl{(}\textnormal{expit}(\theta_{i}^{\top}\varphi(x_{1})),\ldots,\textnormal{expit}(\theta_{i}^{\top}\varphi(x_{n}))\bigr{)}\Bigm{|}\omega_{i}\geqslant
0,\sum_{i=1}^{m}\omega_{i}\leqslant 1\biggr{\\}}.$

The convergence rates of $\hat{e}(\cdot)$ to $e^{*}(\cdot)$ will be proven
using strong concavity of the design objective (27) on the space
$\mathcal{F}^{*}$. This is ensured by Assumption 5.1 along with the following
conditions on the information function $\Psi(\cdot)$:

######  Assumption 5.2 (Information function regularity).

The information function
$\Psi:\mathbb{S}_{+}^{p}\rightarrow\real\cup\\{-\infty\\}$ is concave,
continuous, and nondecreasing with respect to the semidefinite ordering
$\succcurlyeq$ on p×p and satisfies the following conditions:

  * (a)

For every $k>0$, $\inf_{B\succcurlyeq
kI}\Psi(B)>\sup_{A\in\mathbb{S}_{+}^{p}\setminus\mathbb{S}_{++}^{p}}\Psi(A)=:\Psi_{0}$.

  * (b)

$\Psi(\cdot)$ is twice continuously differentiable on $\mathbb{S}_{++}^{p}$,
such that for all $0<k<K$, there exists $C>0$ such that
$\|\nabla\Psi(A)-\nabla\Psi(B)\|\leqslant C\|A-B\|$ whenever $KI\succcurlyeq
A\succcurlyeq kI$ and $KI\succcurlyeq B\succcurlyeq kI$.

  * (c)

For every $0<k<K$, $kI\preccurlyeq A\preccurlyeq KI$ implies
$\tilde{k}I\preccurlyeq\nabla\Psi(A)\preccurlyeq\tilde{K}I$ for some
$0<\tilde{k}<\tilde{K}$.

  * (d)

For every $K<\infty$ and $\tilde{\Psi}_{0}>\Psi_{0}$, there exists $k>0$ such
that for all $0\preccurlyeq A\preccurlyeq KI$ with
$\Psi(A)\geqslant\tilde{\Psi}_{0}$, we have $A\succcurlyeq kI$.

We can show that Assumption 5.2 is satisfied by two common information
functions: the “$A$-optimality” function
$\Psi_{a}(\cdot)=-\textnormal{tr}((\cdot)^{-1})$ with $\Psi_{a}(M):=-\infty$
whenever $M$ is singular, and the “$D$-optimality” function
$\Psi_{d}(\cdot)=\log(\det(\cdot))$. The $A$-optimality criterion corresponds
to minimizing the average (asymptotic) variance of the components of the
estimand, while $D$-optimality corresponds to minimizing the volume of the
ellipsoid spanned by the columns of the (asymptotic) covariance matrix.

######  Lemma 5.1.

The information functions $\Psi_{d}(\cdot)$ and $\Psi_{a}(\cdot)$ satisfy
Assumption 5.2.

###### Proof.

See Appendix C.8. ∎

There are some common information functions that do not satisfy Assumption
5.2. For example, the “$E$-optimality” function
$\Psi_{e}(\cdot)=\lambda_{\min}(\cdot)$, where $\lambda_{\min}(M)$ refers to
the smallest eigenvalue of $M\in{}^{p\times p}$, is not differentiable.
Similarly the function $\Psi_{c}(\cdot)=-c^{\top}(\cdot)^{-1}c$ (for some
fixed $c\in{}^{p}$), corresponding to “$c$-optimality,” does not satisfy
condition (c). We leave open the question of whether the $O_{p}(n^{-1/4})$
convergence rate of $\hat{e}(\cdot)$ to $e^{*}(\cdot)$ in Lemma 5.2 can be
extended to these (and other) information functions using different
techniques, and now prove this rate under Assumptions 5.1 and 5.2.

######  Lemma 5.2 (Convergence of generic concave maximization routine).

Suppose Assumption 5.1 holds for some information matrix $\mathcal{I}$ of the
form (28), covariate distribution $P$, and budget constraints $(m_{L},m_{H})$.
Further assume that for some sequence $\alpha_{n}\downarrow 0$, we have an
estimate $\hat{\eta}(\cdot)$ of $\eta(\cdot)$ satisfying

$\hat{\eta}(x)\in\mathcal{W},\ \forall
x\in\mathcal{X},\quad\text{and}\quad\|\hat{\eta}-\eta\|_{2,P_{n}}=O_{p}(\alpha_{n})$
(34)

where $\eta(\cdot)$ is defined by $\mathcal{I}$ by (28). Then for any
information function $\Psi(\cdot)$ satisfying Assumption 5.2, the following
statements are true:

  1. 1.

There exists an optimal propensity function $e^{*}(\cdot)$ satisfying (27)
which is unique $P$-almost everywhere.

  2. 2.

There exist optimal finite sample treatment probabilities
$(\hat{e}_{1},\ldots,\hat{e}_{n})\in[0,1]^{n}$ satisfying (31), where
$F_{n}=F_{n}(m_{L},m_{H})$ is defined as in (32).

  3. 3.

For any such optimal probabilities $(\hat{e}_{1},\ldots,\hat{e}_{n})$, there
exists a propensity score $\hat{e}(\cdot)\in\mathcal{E}$ for which
$\hat{e}(X_{i})=\hat{e}_{i}$ for each $i=1,\ldots,n$. Any such function
$\hat{e}(\cdot)$ satisfies both
$\|\hat{e}-e^{*}\|_{2,P}=O_{p}(n^{-1/4}+\alpha_{n})$ and
$\|\hat{e}-e^{*}\|_{2,P_{n}}=O_{p}(n^{-1/4}+\alpha_{n})$.

###### Proof.

See Appendix C.9. ∎

###  5.2 Convergence of batch adaptive designs

We now leverage Lemma 5.2 to develop a procedure (Algorithm 1) that can learn
adaptive propensities $\hat{e}_{t}^{(k)}(\cdot)$ with the convergence
guarantees (17) so that when used for treatment assignment, they lead to a
CSBAE with limiting propensities that optimize objectives of the form (27)
with information matrices based on $V_{0,\textnormal{AIPW}}^{-1}$ and
$V_{0,\textnormal{EPL}}^{-1}$. This shows we can effectively design for the
estimators $\hat{\theta}_{\textnormal{AIPW}}$ and
$\hat{\theta}_{\textnormal{EPL}}$.

Algorithm 1 CSBAE for estimating $\theta_{0,\textnormal{ATE}}$ resp.
$\theta_{0,\textnormal{PL}}$

0: Base propensity class $\mathcal{E}$ satisfying conditions in Assumption
5.1, initial propensity $e_{1}(\cdot)\in\mathcal{F}_{\epsilon_{1}}$ for some
$\epsilon_{1}>0$, information function $\Psi(\cdot)$ satisfying Assumption
5.2, number of folds $K\geqslant 2$

1: for batch $t=1,\ldots,T$ do

2: Observe subject covariates $X_{t1},\ldots,X_{tN_{t}}$

3: Split subject indices $i=1,\ldots,N_{t}$ into $K$ folds
$\mathcal{I}_{1},\ldots,\mathcal{I}_{K}$, of size as equal as possible

4: for fold $k=1,\ldots,K$ do

5: Label the covariates in batch $t$, fold $k$ by
$X_{t1}^{(k)},\ldots,X_{tn_{t,k}}^{(k)}$

6: if $t=1$ then

7: Set
$(\hat{e}_{t1}^{(k)},\ldots,\hat{e}_{tn_{t,k}}^{(k)})=\big{(}e_{1}(X_{t1}^{(k)}),\ldots,e_{1}(X_{tn_{t,k}}^{(k)})\big{)}$

8: else

9: Compute $(\hat{e}_{t1}^{(k)},\ldots,\hat{e}_{tn_{t,k}}^{(k)})$ using right-
hand side of (40) resp. (43) given batch $t$ budget constraints
$m_{L,t}\leqslant m_{H,t}$

10: end if

11: Draw
$(U_{t1}^{(k)},\ldots,U_{tn_{t,k}}^{(k)})\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}\textnormal{Unif}(0,1)$

12: Assign treatments according to
$Z_{ti}^{(k)}=\bm{1}(U_{ti}^{(k)}\leqslant\hat{e}_{ti}^{(k)})$,
$i=1,\ldots,n_{t,k}$

13: Choose and store any $\hat{e}_{t}^{(k)}(\cdot)\in\mathcal{E}$ with
$\hat{e}_{t}^{(k)}(X_{ti}^{(k)})=\hat{e}_{ti}^{(k)},i=1,\ldots,n_{t,k}$

14: end for

15: Observe outcomes $\\{W_{ti}\mid 1\leqslant i\leqslant N_{t}\\}$ in batch
$t$

16: Compute and store $N^{1/4}$-consistent estimates
$\big{(}\hat{v}_{1:t}^{(1)}(\cdot,\cdot),\ldots,\hat{v}_{1:t}^{(K)}(\cdot,\cdot)\big{)}$
of $v_{0}(\cdot,\cdot)$ where for each fold $k=1,\ldots,K$,
$\hat{v}_{1:t}^{(k)}(\cdot,\cdot)$ depends only on the observations
$\\{W_{ui}:(u,i)\in\mathcal{I}_{k},1\leqslant u\leqslant t\\}$ in fold $k$ and
batches $1,\ldots,t$

17: end for

18: Compute final estimator $\hat{\theta}$ via (18), using
$(s_{a,\textnormal{AIPW}},s_{b,\textnormal{AIPW}})$, resp.
$(s_{a,\textnormal{PL}},s_{b,\textnormal{PL}})$.

For simplicity, we assume treatment in the first batch is assigned according
to a non-random propensity $e_{1}(\cdot)\in\mathcal{F}_{\epsilon_{1}}$ for
some $\epsilon_{1}>0$. We let
$\hat{e}_{1}^{(k)}(\cdot)=e_{1}(\cdot),k=1,\ldots,K$. Then for later batches
$t=2,\ldots,T$, the target propensities are taken to be one of the following,
for an information function $\Psi(\cdot)$ satisfying Assumption 5.2:

$e_{t,\textnormal{AIPW}}^{*}(\cdot)\in\operatorname*{arg\,max}_{e_{t}(\cdot)\in\mathcal{F}_{*,t}}\Psi(V_{0:t,\textnormal{AIPW}}^{-1})\quad\text{or}\quad
e_{t,\textnormal{EPL}}^{*}(\cdot)\in\operatorname*{arg\,max}_{e_{t}(\cdot)\in\mathcal{F}_{*,t}}\Psi(V_{0:t,\textnormal{EPL}}^{-1}).$
(35)

Above, $V_{0:t,\textnormal{AIPW}}$ and $V_{0:t,\textnormal{EPL}}$ are the
asymptotic variances of the oracle pooled estimators
$\hat{\theta}_{\textnormal{AIPW}}$ and $\hat{\theta}_{\textnormal{EPL}}$ of
(12) and (13), respectively, _when computed using observations in a non-
adaptive batch experiment with only $t$ batches and propensities
$e_{1}(\cdot),\ldots,e_{t}(\cdot)$_. By (14) we can compute

$V_{0:t,\textnormal{AIPW}}=V_{0:t,\textnormal{AIPW}}(e_{t};\eta_{0,\textnormal{AIPW}})=\mathbb{E}_{P^{X}}\bigg{[}\frac{v_{0}(1,X)}{e_{0:t}(X)}+\frac{v_{0}(0,X)}{1-e_{0:t}(X)}+(\tau_{0}(X)-\theta_{0})^{2}\bigg{]}$
(36)

where $\eta_{0,\textnormal{AIPW}}(x)$ includes the components
$(v_{0}(0,x),v_{0}(1,x),\tau_{0}(x),\theta_{0})$. Similarly, by (15) we have

$V_{0:t,\textnormal{EPL}}(e_{t};\eta_{0,\textnormal{EPL}})=\Bigg{(}\mathbb{E}_{P^{X}}\bigg{[}\frac{e_{0:t}(X)(1-e_{0:t}(X))}{v_{0}(0,X)e_{0:t}(X)+v_{0}(1,X)(1-e_{0:t}(X))}\psi(X)\psi(X)^{\top}\bigg{]}\Bigg{)}^{-1}$
(37)

where $\eta_{0,\textnormal{EPL}}(x)$ includes the components
$(v_{0}(0,x),v_{0}(1,x))$. In both of the preceding equations, the dependence
on the batch $t$ propensity score $e_{t}(\cdot)$ is through the mixture
$e_{0:t}(\cdot)$ given by

$e_{0:t}(x):=\left(\sum_{u=1}^{t}\kappa_{u}\right)^{-1}\left(\sum_{u=1}^{t-1}\kappa_{u}e_{u}(x)+\kappa_{t}e_{t}(x)\right),\quad
x\in\mathcal{X}.$ (38)

Finally, the optimization set $\mathcal{F}_{*,t}$ in (35) is

$\mathcal{F}_{*,t}=\mathcal{F}_{*}(m_{L,t},m_{H,t};P^{X})=\\{e(\cdot)\in\mathcal{E}\mid
m_{L,t}\leqslant\mathbb{E}_{P^{X}}[e(X)]\leqslant m_{H,t}\\}$ (39)

which satisfies all the conditions in Assumption 5.1 with covariate
distribution $P^{X}$, budget constraints $(m_{L,t},m_{H,t})$, and base
propensity class $\mathcal{E}$.

We do not target the final covariances
$V_{0,\textnormal{AIPW}}=V_{0:T,\textnormal{AIPW}}$ and
$V_{0,\textnormal{EPL}}=V_{0:T,\textnormal{EPL}}$ in our discussion here since
the sample sizes and budget constraints in future batches may not be known. If
they are known in advance, then we can learn propensities for all future
batches simultaneously at the time batch 2 covariates are observed, and indeed
target $V_{0:T,\textnormal{AIPW}}$ or $V_{0:T,\textnormal{EPL}}$ at that
stage.

Now suppose we split our observations into $K$ folds as in a CSBAE, and for
notational simplicity we re-index the covariates in each batch $t=1,\ldots,T$,
fold $k=1,\ldots,K$ as $X_{t1}^{(k)},\ldots,X_{tn_{t,k}}^{(k)}$. Then
following (31), we can estimate $e_{t,\textnormal{AIPW}}^{*}(\cdot)$ for each
batch $t\geqslant 2$ within each fold $k=1,\ldots,K$ by computing

$\Big{(}\hat{e}_{t1,\textnormal{AIPW}}^{(k)},\ldots,\hat{e}_{tn_{t,k},\textnormal{AIPW}}^{(k)}\Big{)}\in\operatorname*{arg\,max}_{(e_{1},\ldots,e_{n_{t,k}})\in
F_{n_{t,k}}}\Psi\biggl{(}\Big{(}\hat{V}_{0:t,\textnormal{AIPW}}^{(k)}\Big{)}^{-1}\bigg{)},\quad
k=1,\ldots,K.$ (40)

Here $F_{n_{t,k}}=F_{n_{t,k}}(m_{L,t},m_{H,t})$ is defined as in (32), and the
estimate $\hat{V}_{0:t,\textnormal{AIPW}}^{(k)}$ of
$V_{0:t,\textnormal{AIPW}}$ is given by

$\displaystyle\hat{V}_{0:t,\textnormal{AIPW}}^{(k)}$
$\displaystyle=\hat{V}_{0:t,\textnormal{AIPW}}^{(k)}\Bigl{(}e_{1},\ldots,e_{n_{t,k}};\hat{e}_{1}^{(k)}(\cdot),\ldots,\hat{e}_{t-1}^{(k)}(\cdot),\hat{v}_{1:(t-1)}^{(k)}(\cdot,\cdot)\Bigr{)}$
$\displaystyle=-\frac{1}{n_{t,k}}\sum_{i=1}^{n_{t,k}}\frac{\hat{v}_{1:(t-1)}^{(k)}\big{(}1,X_{ti}^{(k)}\big{)}}{\hat{e}_{(0:t)i}^{(k)}}+\frac{\hat{v}_{1:(t-1)}^{(k)}\big{(}0,X_{ti}^{(k)}\big{)}}{1-\hat{e}_{(0:t)i}^{(t)}}$
(41)

where each estimate $\hat{v}_{1:(t-1)}^{(k)}(z,\cdot)$ of the variance
function $v_{0}(z,\cdot)$ is computed using only the observations from batches
$u=1,\ldots,t-1$ within fold $k$. Note that any plug-in estimate of
$\tau_{0}(\cdot)$ and $\theta_{0}$ does not affect the optimization (40) and
so can be omitted. The dependence of $\hat{V}_{0:t,\textnormal{AIPW}}^{(k)}$
on the optimization variables $e_{1},\ldots,e_{n_{t,k}}$ is through the
mixture quantities

$\hat{e}_{(0:t)i}^{(k)}:=\frac{1}{N_{1:t}}\Bigg{(}\,\sum_{u=1}^{t-1}N_{u}\hat{e}_{u}^{(k)}(X_{ti}^{(k)})+N_{t}e_{i}\Bigg{)},\quad
i=1,\ldots,n_{t,k}$ (42)

where for $u=1,\ldots,t-1$, $\hat{e}_{u}^{(k)}(\cdot)$ is the (possibly
adaptive) propensity used to assign treatment in batch $u$, fold $k$. By
comparing (42) and (38), we see that for each batch $u=1,\ldots,t$, $N_{u}/N$
is being used as a plug-in estimate of $\kappa_{u}$ and the adaptive
propensity score $\hat{e}_{u}^{(k)}(\cdot)$ is used as a plug-in estimate of
its limit $e_{u}(\cdot)$. Finally, as in the conclusion of Lemma 5.2, the
learned adaptive propensity $\hat{e}_{t,\textnormal{AIPW}}^{(k)}(\cdot)$ to be
used for treatment assignment in batch $t$, fold $k$ is taken to be any choice
in the predetermined base collection $\mathcal{E}$ with
$\hat{e}_{t,\textnormal{AIPW}}^{(k)}\left(X_{ti}^{(k)}\right)=\hat{e}_{ti,\textnormal{AIPW}}^{(k)}$
for all $i=1,\ldots,n_{t,k}$.

Learning $e_{t,\textnormal{EPL}}^{*}(\cdot)$ is exactly analogous; first we
compute

$\Big{(}\hat{e}_{t1,\textnormal{EPL}}^{(k)},\ldots,\hat{e}_{tn_{t,k},\textnormal{EPL}}^{(k)}\Big{)}\in\operatorname*{arg\,max}_{(e_{1},\ldots,e_{n_{t,k}})\in
F_{n_{t,k}}}\Psi\bigg{(}\Big{(}\hat{V}_{0:t,\textnormal{EPL}}^{(k)}\Big{)}^{-1}\bigg{)}$
(43)

where

$\displaystyle\hat{V}_{0:t,\textnormal{EPL}}^{(k)}$
$\displaystyle=\hat{V}_{0:t,\textnormal{EPL}}^{(k)}\Bigl{(}e_{1},\ldots,e_{n_{t,k}};\hat{e}_{1}^{(k)}(\cdot),\ldots,\hat{e}_{t-1}^{(k)}(\cdot),\hat{v}_{1:(t-1)}^{(k)}(\cdot,\cdot)\Bigr{)}$
(44)
$\displaystyle=\frac{1}{n_{t,k}}\sum_{i=1}^{n_{t,k}}\frac{\hat{e}_{(0:t)i}^{(k)}\big{(}1-\hat{e}_{(0:t)i}^{(k)}\big{)}}{\hat{v}_{1:(t-1)}^{(k)}\big{(}0,X_{ti}^{(k)}\big{)}\hat{e}_{(0:t)i}^{(k)}+\hat{v}_{1:(t-1)}^{(k)}\big{(}1,X_{ti}^{(k)}\big{)}\big{(}1-\hat{e}_{(0:t)i}^{(k)}\big{)}}\psi\big{(}X_{ti}^{(k)}\big{)}\psi\big{(}X_{ti}^{(k)}\big{)}^{\top}.$
(45)

Then we assign treatment with any propensity
$\hat{e}_{t,\textnormal{EPL}}^{(k)}(\cdot)$ in the base propensity class
$\mathcal{E}$ satisfying
$\hat{e}_{t,\textnormal{EPL}}\big{(}X_{ti}^{(k)}\big{)}=\hat{e}_{ti,\textnormal{EPL}}^{(k)}$
for all $i=1,\ldots,n_{t,k}$.

The main additional regularity condition required to ensure the adaptive
propensities $\hat{e}_{t,\textnormal{AIPW}}^{(k)}(\cdot)$ and
$\hat{e}_{t,\textnormal{EPL}}^{(k)}$ above converge at the desired
$O_{p}(N^{-1/4})$ RMS rate to $e_{t,\textnormal{AIPW}}^{*}$ and
$e_{t,\textnormal{EPL}}^{*}$ of (35) is the same rate of convergence in the
estimates $\hat{\eta}^{(k)}_{\textnormal{AIPW}}(\cdot)$ and
$\hat{\eta}^{(k)}_{\textnormal{EPL}}(\cdot)$ of the nuisance parameters
$\eta_{0,\textnormal{AIPW}}(\cdot)$ and $\eta_{0,\textnormal{EPL}}(\cdot)$. We
also strengthen the sample size asymptotics (3) by requiring

$\frac{N_{t}}{N}=\kappa_{t}+O(N^{-1/4}),\quad t=1,\ldots,T.$ (46)

######  Theorem 5.1 (Convergence of Algorithm 1).

Suppose $T\geqslant 2$, fix a batch $t\in\\{2,\ldots,T\\}$ and suppose
Assumption 2.1 holds along with (46). Further assume treatment in batches
$1,\ldots,t-1$ is assigned according to a CSBAE where the batch 1 propensities
are
$\hat{e}_{1}^{(1)}(\cdot)=\ldots=\hat{e}_{1}^{(K)}(\cdot)=e_{1}(\cdot)\in\mathcal{F}_{\epsilon_{1}}$
for some $\epsilon_{1}>0$, and that for each fold $k=1,\ldots,K$:

  * •

There exists an estimator $\hat{v}^{(k)}(\cdot)$ of the variance function
$v_{0}(\cdot)$ depending only on the observations $\\{W_{ui}^{(k)}\mid
1\leqslant u\leqslant t-1\\}$ in batches $1,\ldots,t-1$ assigned to fold $k$,
such that
$\|\hat{v}^{(k)}(z,\cdot)-v_{0}(z,\cdot)\|_{2,P^{X}}=O_{p}(N^{-1/4})$ for
$z=0,1$.

  * •

There are universal constants $0<c<C<\infty$ for which

$c\leqslant\inf_{(z,x)\in\\{0,1\\}\times\mathcal{X}}\min\big{(}\hat{v}^{(k)}(z,x),v_{0}(z,x)\big{)}\leqslant\sup_{(z,x)\in\\{0,1\\}\times\mathcal{X}}\max\big{(}\hat{v}^{(k)}(z,x),v_{0}(z,x)\big{)}\leqslant
C.$

  * •

The information function $\Psi(\cdot)$ satisfies Assumption 5.2.

Let $\mathcal{E}\subseteq\mathcal{F}_{0}$ be any base propensity class
satisfying the conditions of Assumption 5.1, and define $P_{N,t}^{(k),X}$ to
be the empirical distribution on the covariates
$X_{t1}^{(k)},\ldots,X_{tn_{t,k}}^{(k)}$ in batch $t$, fold $k$. Then for any
budget constraints $0\leqslant m_{L,t}\leqslant m_{H,t}\leqslant 1$ and each
fold $k=1,\ldots,K$, the following holds:

  1. 1.

(Design for $\hat{\theta}_{\textnormal{AIPW}}$) There exists a target
propensity $e_{t,\textnormal{AIPW}}^{*}(\cdot)\in\mathcal{E}$ satisfying (35)
that is unique $P^{X}$-almost surely. Additionally, there exists a solution
$(\hat{e}_{t1,\textnormal{AIPW}}^{(k)},\ldots,\hat{e}_{tn_{t,k},\textnormal{AIPW}}^{(k)})$
to (40); any such solution has the property that any propensity
$\hat{e}_{t}^{(k)}(\cdot)\in\mathcal{E}$ with
$\hat{e}_{t}^{(k)}(X_{ti}^{(k)})=\hat{e}_{ti,\textnormal{AIPW}}^{(k)}$ for
$i=1,\ldots,n_{t,k}$ satisfies

$\|\hat{e}_{t}^{(k)}-e_{t,\textnormal{AIPW}}^{*}\|_{2,P^{X}}+\|\hat{e}_{t}^{(k)}-e_{t,\textnormal{AIPW}}^{*}\|_{2,P_{N,t}^{X}}=O_{p}(N^{-1/4}).$

If additionally, the linear treatment effect assumption (6) holds for some
basis function $\psi(X)$ containing an intercept, then:

  1. 2.

(Design for $\hat{\theta}_{\textnormal{EPL}}$) There exists a target
propensity $e_{t,\textnormal{EPL}}^{*}(\cdot)\in\mathcal{E}$ satisfying (35)
that is unique $P^{X}$-almost surely. Furthermore, there exists a solution
$(\hat{e}_{t1,\textnormal{EPL}}^{(k)},\ldots,\hat{e}_{tn_{t,k},\textnormal{EPL}}^{(k)})$
to (43); any such solution has the property that any propensity
$\hat{e}_{t}^{(k)}(\cdot)\in\mathcal{E}$ with
$\hat{e}_{t}^{(k)}(X_{ti}^{(k)})=\hat{e}_{ti,\textnormal{EPL}}^{(k)}$ for
$i=1,\ldots,n_{t,k}$ satisfies

$\|\hat{e}_{t}^{(k)}-e_{t,\textnormal{EPL}}^{*}\|_{2,P^{X}}+\|\hat{e}_{t}^{(k)}-e_{t,\textnormal{EPL}}^{*}\|_{2,P_{N,t}^{X}}=O_{p}(N^{-1/4}).$

###### Proof.

See Appendix C.10. ∎

##  6 Numerical simulations

We implement Algorithm 1 to construct some synthetic CSBAE’s that illustrate
the finite sample performance of our proposed methods. For simplicity we
consider $T=2$ batches throughout. Our evaluation metric is the average mean
squared error (AMSE) of the estimators $\hat{\theta}_{\textnormal{AIPW}}$ and
$\hat{\theta}_{\textnormal{EPL}}$ computed at the end of each CSBAE. As a
baseline, we also compute feasible variants of the linearly aggregated
estimators $\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ and
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ computed on a “simple RCT”:
that is, a non-adaptive batch experiment with a constant propensity score in
each batch. We additionally consider the approach of Hahn et al., (2011) for
design and estimation of $\theta_{0,\textnormal{ATE}}$. This is equivalent to
using Algorithm 1 for design (without sample splitting, i.e. $K=1$) and using
the pooled $\hat{\theta}_{\textnormal{AIPW}}$ as the final estimator, but with
the covariates $X$ replaced everywhere by a coarse discretization $S=S(X)$. As
a hybrid we also consider using the discretized covariates $S$ for design but
computing the final estimates $\hat{\theta}_{\textnormal{AIPW}}$ and
$\hat{\theta}_{\textnormal{EPL}}$ using the full original covariate $X$. To
separately attribute efficiency gains to design and pooling, we also consider
the linearly aggregated estimators
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ and
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ when a modification of
Algorithm 1 that targets $V_{2,\textnormal{AIPW}}$ and
$V_{2,\textnormal{EPL}}$ (the asymptotic variances of the estimators
$\hat{\theta}_{2,\textnormal{AIPW}}$ and $\hat{\theta}_{2,\textnormal{EPL}}$
given in Section 3, depending only on observations in batch 2) is used for
design.

We consider four data generating processes (DGPs), distinguished by whether
the covariate dimension $d$ is 1 or 10 and whether the conditional variance
functions $v_{0}(\cdot,\cdot)$ are homoskedastic (with
$v_{0}(0,x)=v_{0}(1,x)=1$ for all $x\in\mathcal{X}$) or heteroskedastic (with
$v_{0}(0,x)=v_{0}(1,x)/2=\exp((\bm{1}_{d}^{\top}x)/(2\sqrt{d}))$. The scaling
by the covariate dimension $d$ in the heteroskedastic variance functions
ensures the variance in $v_{0}(z,X)$ is independent of the covariate dimension
$d$. In all of the DGPs the covariates are i.i.d. spherical Gaussian, i.e.
$P^{X}=\mathcal{N}(0,I_{d})$. The outcome mean functions are taken to be
$m_{0}(0,x)=m_{0}(1,x)=\bm{1}_{d}^{\top}x$ for
$\bm{1}_{d}=(1,\ldots,1)^{\prime}\in{}^{d}$. For estimating
$\theta_{0,\textnormal{PL}}$ we use the basis functions
$\psi(x)=(1,x^{\top})^{\top}\in{}^{p}$ where $p=d+1$. Note
$\theta_{0,\textnormal{ATE}}=\theta_{0,\textnormal{PL}}=0$. The potential
outcomes $Y(0),Y(1)$ are generated as follows:

$\displaystyle(\epsilon(0),\epsilon(1))\mid X$
$\displaystyle\sim\mathcal{N}((0,0)^{\top},\textnormal{diag}(v_{0}(0,X),v_{0}(1,X)))$
$\displaystyle Y(z)$ $\displaystyle=m_{0}(z,X)+\epsilon(z),\quad z=0,1.$

For each DGP we run Algorithm 1 with $K=2$ folds, batch sample sizes
$N_{1}=N_{2}=1000$, information function $\Psi=\Psi_{a}(\cdot)$ (corresponding
to $A$-optimality), and treatment fraction constraints $m_{L,t}=m_{H,t}=0.2$
for $t=1,2$. In Appendix D, we present additional simulation results where
$m_{L,1}=m_{H,1}$ remain at 0.2 but $m_{L,2}=m_{H,2}=0.4$, so that the
treatment budget for the second batch has increased. The initial propensity
score $e_{1}(\cdot)$ is taken to be constant (i.e. $e_{1}(x)=0.2$ for all
$x\in\mathcal{X}$), and the base propensity class $\mathcal{E}$ is the set of
all $1$-Lipschitz functions taking values in $[0,1]$ when $d=1$ (cf. Example
5.2). When $d=10$, we take $\mathcal{E}$ as in (33), with
$\\{\theta_{1},\ldots,\theta_{m}\\}$ the collection of vectors
$(a_{1},\ldots,a_{11})^{\prime}\in{}^{11}$ with each coordinate
$a_{i}\in\\{-2,-1,0,1,2\\}$ and no more than two of the $a_{i}$’s nonzero.

The discretization used to implement the approach of Hahn et al., (2011)
partitions d into four bins based on the quartiles of $\bm{1}_{d}^{\top}X$.
Note that this partition is along the single dimension along which the
variance functions $v_{0}(\cdot,\cdot)$ vary in the heteroskedastic DGPs, so
we would expect this to perform better than in practice, where the structure
of the variance functions is not known (it could possibly be learned, as in
Tabord-Meehan, (2022)). The choice of four bins is based on the experiments of
Hahn et al., (2011), which find minimal performance difference between two and
six bins for the DGP’s they consider. We denote their “binned” AIPW estimator
by $\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$. Recall, as indicated
above, that this is equivalent to $\hat{\theta}_{\textnormal{AIPW}}$ when the
only available covariate is the binned $S(X)$ and no cross-fitting is used.
For $\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ the conditional means
$\tilde{m}_{0}(z,s)=\mathbb{E}(Y(z)\mid S=s),z=0,1$ are estimated
nonparametrically by sample outcome means among all units with $Z=z$ and $S=s$
in the appropriate batch(es) and fold(s). Similarly, the conditional variances
$\tilde{v}_{0}(z,s)=\textnormal{Var}(Y(z)\mid S=s)$ are estimated by sample
outcome variances.

For all simulations, the concave maximizations are performed using the CVXR
software (Fu et al.,, 2017) and the MOSEK solver (MOSEK ApS,, 2022). All
estimates $\hat{m}(z,\cdot)$ of the mean function $m_{0}(z,\cdot),z=0,1$ are
computed by fitting a generalized additive model (GAM) to the outcomes $Y$ and
covariates $X$ from the observations with treatment indicators equal to $z$ in
the appropriate fold(s) and batch(es). The GAMs use a thin-plate regression
spline basis (Wood,, 2003) and the degrees of freedom are chosen using the
generalized cross-validation procedure implemented in the mgcv package in R
(Wood,, 2004). Variance function estimates $\hat{v}(z,\cdot)$ are computed by
first computing $\hat{m}(z,\cdot)$ as above on the appropriate observations,
then fitting a GAM on these same observations to predict the squared residuals
$(Y-\hat{m}(z,X))^{2}$ from $X$.

###  6.1 Average treatment effect

Table 1: Simulated and asymptotic relative efficiencies (as defined in Section 6 of the various design and estimation approaches for $\theta_{0,\textnormal{ATE}}$ that we study numerically under each of the four DGP’s described in the text.  DGP | Estimator | Design | Sim. rel. eff. (90% CI) | Asymp. rel. eff.
---|---|---|---|---
$d=1$, Homoskedastic | $\hat{\theta}_{\textnormal{AIPW}}$ | Flexible | 0.989 (0.965, 1.013) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}$ | Binned | 1.011 (0.977, 1.046) | 0.999
$\hat{\theta}_{\textnormal{AIPW}}$ | Simple RCT | 1.016 (1.003, 1.029) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ | Flexible | 0.984 (0.971, 0.998) | 0.999
$\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ | Binned | 0.919 (0.877, 0.961) | 0.885
$d=1$, Heteroskedastic | $\hat{\theta}_{\textnormal{AIPW}}$ | Flexible | 1.016 (0.968, 1.064) | 1.048
$\hat{\theta}_{\textnormal{AIPW}}$ | Binned | 1.046 (0.997, 1.096) | 1.041
$\hat{\theta}_{\textnormal{AIPW}}$ | Simple RCT | 0.997 (0.979, 1.016) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ | Flexible | 1.008 (0.971, 1.045) | 1.024
$\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ | Binned | 1.001 (0.946, 1.057) | 0.963
$d=10$, Homoskedastic | $\hat{\theta}_{\textnormal{AIPW}}$ | Flexible | 1.039 (0.990, 1.089) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}$ | Binned | 1.011 (0.959, 1.065) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}$ | Simple RCT | 1.081 (1.052, 1.110) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ | Flexible | 1.047 (1.018, 1.077) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ | Binned | 0.477 (0.439, 0.519) | 0.417
$d=10$, Heteroskedastic | $\hat{\theta}_{\textnormal{AIPW}}$ | Flexible | 1.081 (1.023, 1.143) | 1.036
$\hat{\theta}_{\textnormal{AIPW}}$ | Binned | 1.050 (0.986, 1.117) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}$ | Simple RCT | 1.102 (1.070, 1.135) | 1.000
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ | Flexible | 0.991 (0.952, 1.030) | 1.024
$\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ | Binned | 0.651 (0.597, 0.708) | 0.595

Table 1 shows the performance of the various design and estimation procedures
for $\theta_{0,\textnormal{ATE}}$ that we consider. Each entry is a “relative
efficiency”: that is, a ratio of the MSE of the baseline approach (which
computes the linearly aggregated
$\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$ on a simple RCT) to the MSE
of the relevant approach. The simulated relative efficiencies in the table
estimate these MSE’s by averaging the squared error of each estimator over
1,000 simulations. The 90% confidence intervals for the true finite sample
relative efficiency are computed using 10,000 bootstrap replications of these
1,000 simulations. Finally, the asymptotic relative efficiencies in Table 1
are computed by estimating the asymptotic variance of each estimator using the
appropriate formula, i.e. $V_{0,\textnormal{AIPW}}$ for
$\hat{\theta}_{\textnormal{AIPW}}$ and $V_{\textnormal{AIPW}}^{(\mathrm{LA})}$
for $\hat{\theta}_{\textnormal{AIPW}}^{(\mathrm{LA})}$. The computation is
based on a non-adaptive batch experiment with second batch propensity
$e_{2}(\cdot)$ equal to the average of the learned propensities from the
relevant approach across the 1,000 simulations, as a closed form solution for
the limiting $e_{2}(\cdot)$ is not easily obtained in general. All
expectations over the covariate distribution are computed using Monte Carlo
integration.

For both of the homoskedastic DGPs, it is straightforward to show using
Jensen’s inequality that $e_{2}^{*}(x)=0.2$ for all $x$. It can be further
shown that the unpooled and pooled estimators are asymptotically equivalent.
Thus, for the homoskedastic DGPs, there is no asymptotic efficiency gain to be
had. With unequal budget constraints there is some asymptotic benefit to
pooling with homoskedastic variance functions (Appendix D.1); however the
optimal design will still be the simple RCT. Nonetheless, we do observe some
finite sample benefits to pooling in Table 1. For example, the simulated
relative efficiency of $\hat{\theta}_{\textnormal{AIPW}}$ on the simple RCT is
significantly larger than 1 for both $d=1$ and $d=10$. We attribute this to
improved nuisance estimates for the pooled estimator, as discussed further in
Appendix D.2. As one might expect, this finite sample improvement from pooling
is apparently offset by variance in both the flexible and binned design
procedures. Still, in the homoskedastic DGPs, our adaptive approaches do not
show significant finite sample performance decline relative to the baseline,
which here is an oracle. With unequal treatment constraints, we further obtain
asymptotic efficiency gains from pooling (Appendix D.1).

We notice that using the discretized covariate $S(X)$ in place of $X$ for both
design and estimation, as in Hahn et al., (2011), leads to a substantial loss
of efficiency. Indeed, when $d=10$ the (asymptotic and simulated) variance of
the estimator $\hat{\theta}_{\textnormal{AIPW}}^{(\text{bin})}$ is more than
double that of our baseline under the homoskedastic DGP, both asymptotically
and in our finite sample simulations. This efficiency loss occurs because the
discretized $S(X)$ explains much less of the variation in the potential
outcomes $Y(z)$ than the original $X$. We expect greater precision losses from
this discretization at the estimation stage when the variance functions
$v_{0}(z,\cdot)$ vary substantially within the strata defined by $S(X)$.

For the heteroskedastic DGPs, we see modest asymptotic efficiency gains from
both pooling and design. Design using the flexible base propensity class leads
to about a 2.4% asymptotic efficiency gain for both $d=1$ and $d=10$, while
pooling provides an additional 1–2% gain. These small asymptotic gains appear
to be largely canceled out at our sample sizes by the finite sample
variability in learning propensity scores, limiting the net finite sample
gains from design. Of course, with greater heteroskedasticity and/or
differences between $v_{0}(0,\cdot)$ and $v_{0}(1,\cdot)$, we would expect
greater efficiency gains from design; in our simulations we have chosen to
keep these differences within common ranges in social science studies as per
Blackwell et al., (2022).

###  6.2 Partially linear model

Unlike for estimating $\theta_{0,\textnormal{ATE}}$, for estimating
$\theta_{0,\textnormal{PL}}$, we see clear efficiency gains over the baseline
from _design_ when $d=1$ (Table 2). For instance, the linearly aggregated
estimator $\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ exhibits a 5.6%
asymptotic efficiency gain as a result of the flexible design; replacing this
with the pooled estimator $\hat{\theta}_{\textnormal{EPL}}$ then yields a
total asymptotic gain of 10.0% over the baseline, even in the homoskedastic
DGP. The analogous gains for the heteroskedastic DGP are slightly larger. We
once again observe a substantial finite sample benefit to pooling, with the
simulated relative efficiency of the approaches using the pooled
$\hat{\theta}_{\textnormal{EPL}}$ tending to be larger than the asymptotic
relative efficiency. We attribute this to both improved use of nuisance
estimates by the pooled estimator (as in the ATE case) as well as a more
fundamental finite sample efficiency boost due to the fact that the asymptotic
variance $V_{0,\textnormal{EPL}}$ is not exact for the oracle pooled
$\hat{\theta}_{\textnormal{EPL}}^{*}$ in finite samples, whereas
$V_{0,\textnormal{AIPW}}$ is exact for $\hat{\theta}_{\textnormal{AIPW}}^{*}$;
see Appendix D.2. When $d=10$, design introduces some more salient finite
sample variance from the errors in the concave maximization procedure. This is
offset in both the homoskedastic and heteroskedastic DGPs by the asymptotic
gains from the flexible design, so that when $\hat{\theta}_{\textnormal{EPL}}$
is used, the flexible design ultimately performs similarly to the simple RCT
in finite samples.

Even if we use $\hat{\theta}_{\textnormal{EPL}}$ and hence the original
covariate(s) $X$ in the final estimation step, we see that the “binned” design
which uses only $S(X)$ for choosing the second batch propensity score
struggles to learn a substantially better propensity score than the simple RCT
in all DGPs. Indeed, in the heteroskedastic DGP with $d=1$, we see a 2.1%
asymptotic efficiency _loss_ relative to the baseline from using the binned
design. By comparison, there is an 11.2% asymptotic efficiency gain from using
the flexible design. The efficiency loss can occur with the binned design
because the objective in (27) changes when working in terms of the discretized
covariate $S(X)$ instead of $X$. In other words, even though the simple RCT
propensity $e_{2}(x)=0.2$ is within the class of propensities that can be
chosen by the binned design, it is worse according to the binned objective
based on $S(X)$, but not according to the objective based on the original $X$.

Table 2: Same as Table 1, but for estimating $\theta_{0,\textnormal{EPL}}$.
DGP | Estimator | Design | Sim. rel. eff. (90% CI) | Asymp. rel. eff.
---|---|---|---|---
$d=1$, Homoskedastic | $\hat{\theta}_{\textnormal{EPL}}$ | Flexible | 1.150 (1.103, 1.198) | 1.100
$\hat{\theta}_{\textnormal{EPL}}$ | Binned | 1.057 (1.009, 1.107) | 1.026
$\hat{\theta}_{\textnormal{EPL}}$ | Simple RCT | 1.038 (1.018, 1.058) | 1.000
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ | Flexible | 1.049 (1.016, 1.083) | 1.056
$d=1$, Heteroskedastic | $\hat{\theta}_{\textnormal{EPL}}$ | Flexible | 1.310 (1.214, 1.412) | 1.112
$\hat{\theta}_{\textnormal{EPL}}$ | Binned | 1.049 (0.972, 1.133) | 0.979
$\hat{\theta}_{\textnormal{EPL}}$ | Simple RCT | 1.084 (1.014, 1.159) | 1.000
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ | Flexible | 0.983 (0.880, 1.076) | 1.073
$d=10$, Homoskedastic | $\hat{\theta}_{\textnormal{EPL}}$ | Flexible | 1.229 (1.203, 1.256) | 1.023
$\hat{\theta}_{\textnormal{EPL}}$ | Binned | 1.177 (1.147, 1.206) | 1.000
$\hat{\theta}_{\textnormal{EPL}}$ | Simple RCT | 1.228 (1.204, 1.251) | 1.000
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ | Flexible | 1.023 (1.011, 1.035) | 1.019
$d=10$, Heteroskedastic | $\hat{\theta}_{\textnormal{EPL}}$ | Flexible | 0.997 (0.963, 1.032) | 1.071
$\hat{\theta}_{\textnormal{EPL}}$ | Binned | 0.932 (0.899, 0.965) | 1.002
$\hat{\theta}_{\textnormal{EPL}}$ | Simple RCT | 0.982 (0.952, 1.014) | 1.000
$\hat{\theta}_{\textnormal{EPL}}^{(\mathrm{LA})}$ | Flexible | 0.969 (0.942, 0.996) | 1.060

##  7 Discussion

We view our primary technical contribution in this paper to be a careful
extension of the double machine learning framework that enables both
estimation and design in batched experiments based on pooled treatment effect
estimators. This allows the investigator to take advantage of the efficiency
gains from pooling and design without needing to make strong parametric
assumptions or to discretize their covariates. As our numerical study in
Section 6 shows, the latter can more than wipe out any efficiency gains from
design.

Related to our work is the extensive literature on combining observational
data with a (single batch) randomized experiment. In that setting a primary
concern is mitigating bias from unobserved confounders in the observational
data (Rosenman and Owen,, 2021; Gagnon-Bartsch et al.,, 2023). By contrast, in
our setting unconfoundedness holds by design in each batch of the experiment.
It would be useful to examine if the ideas from the present work can be
extended to the observational setting where confounding bias is a concern.

Acknowledgments: The authors thank Stefan Wager, Lihua Lei, and Kevin Guo for
comments that improved the content of this paper. H.L. was partially supported
by the Stanford Interdisciplinary Graduate Fellowship (SIGF). This work was
also supported by the NSF under grant DMS-2152780.

##  Appendix A Nonstationary batches

Assumption 2.1 in the main text supposes the distribution $P^{S}$ of the
covariates and potential outcomes $S_{ti}=(X_{ti},Y_{ti}(0),Y_{ti}(1))$ is
stationary across batches $t=1,\ldots,T$. Here we relax that assumption. Let
$P^{S}_{t}$ be the distribution of the vector $S_{ti}$ in batch
$t=1,\ldots,T$, now allowed to vary across batches (in the main text it is
assumed that $P^{S}_{1}=\dots=P^{S}_{T}=P^{S}$). Then we have the following
relaxation of Assumption 2.1:

######  Assumption A.1 (Relaxation of Assumption 2.1).

For some fixed number of batches $T\geqslant 2$, the vectors

$S_{ti}=(X_{ti},Y_{ti}(0),Y_{ti}(1)),\quad 1\leqslant t\leqslant T,\quad
1\leqslant i\leqslant N_{t}$

are mutually independent such that for each batch $t=1,\ldots,T$, we have
$S_{ti}\sim P_{t}^{S}$. Furthermore, the sample sizes $N_{t}$ satisfy (3), and
the vector $W_{ti}=(X_{ti},Z_{ti},Y_{ti})$ is observed where the outcomes
$Y_{ti}$ satisfy the SUTVA assumption (1).

Now letting $P_{0}^{S}$ be the mixture distribution
$\sum_{t=1}^{T}\kappa_{t}P_{t}^{S}$, we introduce the notation
$\mathbb{E}_{t,e}[f(W)]$, which denotes an expectation under the distribution
$P_{t,e}=P_{t,e}^{W}$ on $W=(X,Z,Y)$ induced by $S=(X,Y(0),Y(1))\sim
P_{t}^{S}$ and $Z\mid X\sim\textnormal{Bern}(e(X))$ for any propensity
$e(\cdot)$ and $t=0,1,\ldots,T$. The notation $P_{t}^{X}$ refers to the
corresponding marginal distribution of the covariates $X$. Then the score
equation (5) will be generalized to

$\mathbb{E}_{t,e}[s(W;\theta_{0},\nu_{0},e^{\prime})]=0,\quad\forall
e,e^{\prime}\in\mathcal{F}_{\gamma},\quad t=1,\ldots,T,$ (47)

which we will require to identify $\theta_{0}$ in each batch:

######  Assumption A.2 (Relaxation of Assumption 2.2).

The estimand $\theta_{0}\in{}^{p}$ of interest satisfies (47) for some
$\gamma\in[0,1/2)$, some nuisance parameters $\nu_{0}$ lying in a known convex
set $\mathcal{N}$, and some score $s(\cdot)$ satisfying (4).

Equation (47) encodes a requirement that the same parameters $\theta_{0}$ and
$\nu_{0}$ satisfy the score equations for all batches $t=1,\ldots,T$; to
ensure this, we require those parameters to be stationary across batches. For
instance, for ATE estimation with the score $s_{\textnormal{AIPW}}(\cdot)$, we
require the conditional mean functions $\mathbb{E}[Y(z)\mid X=x]=m_{0}(z,x)$
to be stationary across batches. For estimation under the partially linear
model with the score $s_{\textnormal{EPL}}$, we also require the outcome
variance functions $\textnormal{Var}_{t}(Y(z)\mid X=x)=v_{0}(z,x)$ to remain
stationary. However, in both cases the covariate distribution can otherwise
vary arbitrarily across batches, as can higher moments of the conditional
distributions of the potential outcomes $Y(z)$ given the covariates $X$.

Due to the possibility of covariate shift, the relevant mixture propensity
scores are now

$e_{0,N}(x)=\sum_{t=1}^{T}\frac{N_{t}}{N}e_{t}(x)\frac{\mathrm{d}P_{t}^{X}}{\mathrm{d}P_{0}^{X}}(x)\quad\text{and}\quad
e_{0}(x)=\sum_{t=1}^{T}\kappa_{t}e_{t}(x)\frac{\mathrm{d}P_{t}^{X}}{\mathrm{d}P_{0}^{X}}(x).$
(48)

These definitions generalize (9) and (10). Here $P_{t}^{X}$ denotes the
marginal distribution of $X$ when $S\sim P_{t}^{S}$. When there is no
covariate shift, we have $\mathrm{d}P_{t}^{X}/\mathrm{d}P_{0}^{X}(x)=1$ for
all $x$ and we recover (9) and (10) in the main text. The expressions in (48)
are derived using Bayes’ rule as the conditional probability that $Z=1$ given
$X=x$ when $(X,Z,Y)$ is drawn uniformly at random from the pooled collection
of observations $\\{W_{ti}\mid 1\leqslant t\leqslant T,1\leqslant i\leqslant
N_{t}\\}$ in a non-adaptive batch experiment with propensities
$e_{1}(\cdot),\ldots,e_{T}(\cdot)$.

Where indicated, we are able to generalize various results in Sections 3 and
4. The generalized results are as stated in the main text, if we make the
following changes to the notation and assumptions:

  1. 1.

Assumptions 2.1 and 2.2 are replaced by Assumptions A.1. and A.2, respectively

  2. 2.

Any references to the mixture propensities in $e_{0,N}(\cdot)$ and
$e_{0}(\cdot)$ correspond to the more general definitions in (48), rather than
(9) and (10).

  3. 3.

Any expectations of the form $\mathbb{E}[f(W)]$ without subscripts are
interpreted as being taken under the distribution $P_{0,e_{0}}$ on $W$.

  4. 4.

Any expectations of the form $\mathbb{E}_{t}[f(W)],t=0,1,\ldots,T$ are
interpreted as being taken under the distribution $P_{t,e_{t}}$ on $W$, and
any expectations of the form $\mathbb{E}_{0,N}[f(W)]$ are interpreted as being
taken under the distribution $P_{0,e_{0,N}}$ on $W$.

##  Appendix B Technical lemmas

Here we give some technical lemmas used in our proofs.

###  B.1 Asymptotics

######  Definition B.1.

For any sequence of random vectors $\\{X_{n}:n\geqslant 1\\}$ and constants
$a_{n}\downarrow 0$, we write $X_{n}=O_{p}(a_{n})$ if
$\lim_{M\rightarrow\infty}\limsup_{n\rightarrow\infty}\textnormal{Pr}(\|X_{n}\|>Ma_{n})=0$.
We write $X_{n}=o_{p}(a_{n})$ if for every $M>0$,
$\limsup_{n\rightarrow\infty}\textnormal{Pr}(\|X_{n}\|>Ma_{n})=0$.

######  Lemma B.1.

Let $X_{n}$ be a sequence of random vectors and $\\{\mathcal{F}_{n},n\geqslant
1\\}$ be a sequence of $\sigma$-algebras such that
$\mathbb{E}\bigl{[}\|X_{n}\|\mid\mathcal{F}_{n}\bigr{]}=o_{p}(1)$. Then
$X_{n}=o_{p}(1)$.

###### Proof of Lemma B.1.

Fixing $M>0$, we have $M\bm{1}(\|X_{n}\|>M)\leqslant\|X_{n}\|$ for all $n$.
Taking conditional expectations given $\mathcal{F}_{n}$ on both sides we have

$P\bigl{(}\|X_{n}\|>M\mid\mathcal{F}_{n}\bigr{)}\leqslant
M^{-1}\mathbb{E}\bigl{[}\|X_{n}\|\mid\mathcal{F}_{n}\bigr{]}.$ (49)

Thus if $\mathbb{E}\bigl{[}\|X_{n}\|\mid\mathcal{F}_{n}\bigr{]}=o_{p}(1)$ we
have $\textnormal{Pr}(\|X_{n}\|>M\mid\mathcal{F}_{n})=o_{p}(1)$ as well. But
$\textnormal{Pr}(\|X_{n}\|>M\mid\mathcal{F}_{n})$ is uniformly bounded so its
expectation converges to zero, i.e., $\textnormal{Pr}(\|X_{n}\|>M)=o(1)$.
Since $M>0$ was arbitrary we conclude that $X_{n}=o_{p}(1)$. ∎

######  Lemma B.2.

$X_{n}=O_{p}(a_{n})$ if and only if for every sequence $b_{n}\uparrow\infty$
we have $\textnormal{Pr}(\|X_{n}\|>b_{n}a_{n})\rightarrow 0$ as
$n\rightarrow\infty$.

###### Proof of Lemma B.2.

Fix $b_{n}\uparrow\infty$ and $\epsilon>0$. If $X_{n}=O_{p}(a_{n})$ then there
exists $M<\infty$ such that
$\limsup_{n\rightarrow\infty}\textnormal{Pr}(\|X_{n}\|>Ma_{n})<\epsilon$.
Since $b_{n}>M$ eventually we conclude
$\limsup_{n\rightarrow\infty}\textnormal{Pr}(\|X_{n}\|>b_{n}a_{n})<\epsilon$
as well. With $\epsilon>0$ arbitrary, the result follows. Conversely now
suppose we do not have $X_{n}=O_{p}(a_{n})$. Then there exists $\epsilon>0$
such that
$\limsup_{n\rightarrow\infty}\textnormal{Pr}(\|X_{n}\|>Ma_{n})\geqslant\epsilon$
for all $M<\infty$. Defining $n_{0}=1$, this ensures that for each
$k=1,2,\ldots$, there exists $n_{k}>n_{k-1}$ so that
$\textnormal{Pr}(\|X_{n_{k}}\|>ka_{n})\geqslant\epsilon$. But then for
$b_{n}=\max\\{k\geqslant 0:n_{k}\leqslant n\\}$, we have $b_{n}\uparrow\infty$
yet $\textnormal{Pr}(\|X_{n}\|>b_{n}a_{n})\geqslant\epsilon$ for all $n\in
n_{1},n_{2},\ldots$ so $\textnormal{Pr}(\|X_{n}\|>b_{n}a_{n})$ does not
converge to 0 as $n\rightarrow\infty$. ∎

We now state an important result in empirical process theory in our proof of
Lemma 5.2 above. For a metric space $(\mathcal{M},d)$, let
$B_{\epsilon}(m_{0})=\\{m\in\mathcal{M}\mid d(m,m_{0})\leqslant\epsilon\\}$ be
the $\epsilon$-ball around $m_{0}\in\mathcal{M}$. For any set
$\mathcal{S}\subseteq\mathcal{M}$, the $\epsilon$ covering number
$\mathcal{N}(\epsilon,\mathcal{S},d)$ of $\mathcal{S}$ is then defined as the
smallest number of $\epsilon$-balls in $\mathcal{M}$ whose union contains
$\mathcal{S}$. For a subset $\mathcal{F}$ of the space $L^{2}(\mathcal{X};P)$
of $P$-square integrable real-valued functions on $\mathcal{X}$ with
$\sup_{f\in\mathcal{F}}\|f\|_{\infty}\leqslant\bar{F}<\infty$, control over
the logarithm of the covering numbers of $\mathcal{F}$ (the _metric entropy_)
over a variety of radii $\epsilon$ under the random metric $L^{2}(P_{n})$
given by
$L^{2}(P_{n})(f_{1},f_{2})=\|f_{1}-f_{2}\|_{2,P_{n}}=(\int(f_{1}(x)-f_{2}(x))^{2}\mathrm{d}P_{n}(x))^{1/2}$
implies control of the empirical process $\sup_{f\in\mathcal{F}}|(P_{n}-P)f|$
where $Qf:=\int f(x)\mathrm{d}Q(x)$ for any measure $Q$ on $\mathcal{X}$. Here
$P_{n}$ is the empirical probability measure on observations
$X_{1},\ldots,X_{n}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}P$. This
result is due to a “chaining” argument of Dudley, (1967). We restate a more
direct version of this result below, which is Lemma A.4 of Kitagawa and
Tetenov, (2018).

######  Lemma B.3.

For a class $\mathcal{F}$ of $P$-measurable functions
$f:\mathcal{X}\rightarrow\real$ with
$\sup_{f\in\mathcal{F}}\|f\|_{\infty}\leqslant\bar{F}$, there exists a
universal constant $K<\infty$ such that for $P_{n}$ the empirical distribution
of $X_{1},\ldots,X_{n}\stackrel{{\scriptstyle\mathrm{iid}}}{{\sim}}P$,

$\sup_{f\in\mathcal{F}}\big{|}(P_{n}-P)f\big{|}\leqslant
K\bar{F}n^{-1/2}\int_{0}^{1}\sqrt{\log\mathcal{N}(\epsilon,\mathcal{F},L^{2}(P_{n}))}\,\mathrm{d}\epsilon.$

######  Remark B.1.

Often, control of the right-hand side in the previous display is shown by
controlling

$\sup_{Q}\int_{0}^{1}\sqrt{\log\mathcal{N}(\epsilon,\mathcal{F},L^{2}(Q)}\,\mathrm{d}\epsilon$

where the supremum is taken over all finitely supported probability measures
$Q$. See Sections 2.5 and 2.6 of van der Vaart and Wellner, (1996) for further
discussion.

We also make use of the following elementary results on covering numbers.

######  Lemma B.4.

Let $\mathcal{E}$ be a collection of functions contained in the class
$\mathcal{F}_{0}$. Let $\mathcal{G}=\\{x\mapsto g(e(x),\eta(x))\mid
e\in\mathcal{E}\\}$ for some $\eta:\mathcal{X}\rightarrow\mathcal{W}$ and
$g:[0,1]\times\mathcal{W}\rightarrow\real$ with $g(\cdot,w)$ continuous on
$[0,1]$ and $\sup_{k\in[0,1],w\in\mathcal{W}}|g^{\prime}(k,w)|\leqslant C$ for
some $C<\infty$, where $g^{\prime}(\cdot,\cdot)$ denotes the partial
derivative of $g(\cdot,\cdot)$ with respect to the first argument. Then for
all probability distributions $P$ on $\mathcal{X}$ and $\epsilon>0$, we have
$\mathcal{N}(C\epsilon,\mathcal{G},L^{2}(P))\leqslant\mathcal{N}(\epsilon,\mathcal{E},L^{2}(P))$.

###### Proof of Lemma B.4.

Fix $\epsilon>0$ and a probability distribution $P$ on $\mathcal{X}$. For each
$e\in\mathcal{E}$ let $h^{(e)}(x)=g(e(x),\eta(x))$ for each $x\in\mathcal{X}$,
so that $h^{(e)}\in\mathcal{G}$. Suppose $\\{e_{1},\ldots,e_{N}\\}$ is an
$\epsilon$-cover of $\mathcal{E}$ in the $L^{2}(P)$ norm. WLOG we can assume
that each $e_{k}$ is a member of $\mathcal{F}_{\gamma}$. Then for each
$k=1,\ldots,N$, by the uniform bound on $g^{\prime}$ we have

$\|h^{(e)}-h^{(e_{k})}\|_{2,P}^{2}=\mathbb{E}_{P}\bigl{[}|g(e(X),\eta(X))-g(e_{k}(X),\eta(X))|^{2}\bigr{]}\leqslant
C^{2}\|e-e_{k}\|_{2,P}^{2}$

so that $\\{g^{(e_{1})},\ldots,g^{(e_{N})}\\}$ is a $C\epsilon$ cover of
$\mathcal{G}$ in the $L^{2}(P)$ norm. ∎

######  Lemma B.5.

Let $\mathcal{E}$ be a collection of functions contained in the class
$\mathcal{F}_{0}$. Define
$\mathcal{E}_{2}^{-}=\\{(f-g)^{2}:f\in\mathcal{E},g\in\mathcal{E}\\}$. Then
for every $\epsilon>0$ and probability measure $P$ on $\mathcal{X}$ we have

$\displaystyle\mathcal{N}(\epsilon;\mathcal{E}_{2}^{-},L^{2}(P))\leqslant\mathcal{N}\Bigl{(}\frac{\epsilon}{4};\mathcal{E},L^{2}(P)\Bigr{)}^{2}.$
(50)

###### Proof.

Fix $\epsilon>0$ and a probability distribution $P$ on $\mathcal{X}$. Define
the collection $\mathcal{E}^{-}=\\{f-g:f\in\mathcal{E},g\in\mathcal{E}\\}$.
Suppose $\\{f_{1},\ldots,f_{N}\\}$ is a $\epsilon/4$ cover of $\mathcal{E}$ in
the $L^{2}(P)$ norm. Then the collection $D=\\{d_{ij}=f_{i}-f_{j}:1\leqslant
i,j\leqslant N\\}$ is a $\epsilon/2$ cover of $\mathcal{E}^{-}$ in the
$L^{2}(P)$ norm, since for any $f-g\in\mathcal{E}^{-}$ there exist $i,j$ such
that $\|f-f_{i}\|_{2,P}\vee\|g-f_{j}\|_{2,P}\leqslant\epsilon/4$ and so

$\|(f-g)-d_{ij}\|_{2,P}\leqslant\|f-f_{i}\|_{2,P}+\|g-f_{j}\|_{2,P}\leqslant\frac{\epsilon}{2}$

showing that
$\mathcal{N}(\epsilon/2,\mathcal{E}^{-},L^{2}(P))\leqslant(\mathcal{N}(\epsilon/4,\mathcal{E},L^{2}(P)))^{2}$.
But by applying Lemma B.4 with $\mathcal{E}=\mathcal{E}^{-}$ and
$g(e,w)=e^{2}$ (hence we can take $C=2$) we have
$\mathcal{N}(\epsilon,\mathcal{E}_{2}^{-},L^{2}(P))\leqslant\mathcal{N}(\epsilon/2,\mathcal{E}^{-},L^{2}(P))$
for all $\epsilon>0$. Chaining together the inequalities preceding two
sentences establishes (50). ∎

###  B.2 Miscellaneous

Here we have some standalone technical lemmas. Their proofs do not depend on
any of our other results.

######  Lemma B.6.

Suppose $X$ and $Y$ are mean zero random vectors in p with finite second
moments where $\textnormal{Cov}(Y)$ has full rank. Then

$\textnormal{Cov}(X)\succcurlyeq\textnormal{Cov}(X,Y)(\textnormal{Cov}(Y))^{-1}\textnormal{Cov}(X,Y)^{\color[rgb]{0,0,1}\top}.$

###### Proof of Lemma B.6.

For any matrix $A\in{}^{p\times p}$ we have $(X+AY)(X+AY)^{\top}\succcurlyeq
0$, hence

$\mathbb{E}[(X+AY)(X+AY)^{\top}]=\textnormal{Cov}(X)+\textnormal{Cov}(X,Y)A^{\top}+A\textnormal{Cov}(X,Y)^{\top}+A\textnormal{Cov}(Y)A^{\top}\succcurlyeq
0$

Taking $A=-\textnormal{Cov}(X,Y)(\textnormal{Cov}(Y))^{-1}$ yields the desired
result. ∎

######  Lemma B.7.

If $Y=\sum_{i=1}^{n}Y_{i}$ is a random vector where $Y_{1},\ldots,Y_{n}$ are
independent with mean 0 and finite second moments, then
$\mathbb{E}[\|Y\|^{2}]=\sum_{i=1}^{n}\mathbb{E}[\|Y_{i}\|^{2}]$.

###### Proof of Lemma B.7.

By assumption we have
$\mathbb{E}[Y_{i}^{\top}Y_{j}]=\mathbb{E}[Y_{i}]^{\top}\mathbb{E}[Y_{j}]=0$
for $i\neq j$, so

$\mathbb{E}[\|Y\|^{2}]=\mathbb{E}\Biggl{[}\bigg{\|}\sum_{i=1}^{n}Y_{i}\bigg{\|}^{2}\Biggr{]}=\mathbb{E}\Biggl{[}\biggl{(}\sum_{i=1}^{n}Y_{i}\biggr{)}^{\top}\biggl{(}\sum_{i=1}^{n}Y_{i}\biggr{)}\Biggr{]}=\sum_{i=1}^{n}\mathbb{E}[Y_{i}^{\top}Y_{i}]=\sum_{i=1}^{n}\mathbb{E}[\|Y_{i}\|^{2}].$

∎

##  Appendix C Proofs

Here we collect proofs of the formal results stated in the main text; some are
generalized to allow for some nonstationarities across batches, as described
in Appendix A.

###  C.1 Proof of Proposition 3.1: CLT for oracle $\hat{\theta}^{*}$

The CLT of Proposition 3.1 for our oracle pooled estimator $\hat{\theta}^{*}$
holds under the numbered generalizations in Appendix A without further
restrictions; here we prove this more general proposition. It is helpful to
begin by noting that $P_{0,e_{0}}=\sum_{t=1}^{t}\kappa_{t}P_{t,e_{t}}$, which
shows, for instance, that for any $P_{0}$-integrable function $f$

$\mathbb{E}_{t}[|f(W)|]=\mathbb{E}_{0}\left[|f(W)|\frac{\mathrm{d}P_{t,e_{t}}}{\mathrm{d}P_{0,e_{0}}}(W)\right]\leqslant\kappa_{t}^{-1}\mathbb{E}_{0}[|f(W)|]$

By score linearity (4) we can write

$\sqrt{N}(\hat{\theta}^{*}-\theta_{0})=-\biggl{(}\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{a}(W_{ti};\nu_{0},e_{0,N})\biggr{)}^{-1}\biggl{(}\frac{1}{\sqrt{N}}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})\biggr{)}$
(51)

whenever $\hat{\theta}^{*}$ exists. Defining

$r_{N}=\sum_{t=1}^{T}\frac{N_{t}}{N}\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}(s_{a}(W_{ti};\nu_{0},e_{0,N})-s_{a}(W_{ti};\nu_{0},e_{0})).$

we have from (3) and the law of large numbers that

$\displaystyle\frac{1}{N}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s_{a}(W_{ti};\nu_{0},e_{0,N})$
$\displaystyle=\sum_{t=1}^{T}\frac{N_{t}}{N}\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}s_{a}(W_{ti};\nu_{0},e_{0,N})$
$\displaystyle=r_{N}+\sum_{t=1}^{T}\frac{N_{t}}{N}(\mathbb{E}_{t}[s_{a}(W;\nu_{0},e_{0})]+o_{p}(1))$
$\displaystyle=\sum_{t=1}^{T}\kappa_{t}\mathbb{E}_{t}[s_{a}(W;\nu_{0},e_{0})]+o_{p}(1)$
$\displaystyle=\mathbb{E}_{0}\left[\sum_{t=1}^{T}\kappa_{t}\frac{\mathrm{d}P_{t}}{\mathrm{d}P_{0}}(W)s_{a}(W;\nu_{0},e_{0})\right]+o_{p}(1)$
$\displaystyle=\mathbb{E}_{0}[s_{a}(W;\nu_{0},e_{0})]+o_{p}(1).$ (52)

Then the third equality follows because

$\|r_{N}\|\leqslant\sum_{t=1}^{T}\frac{N_{t}}{N}\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}\|s_{a}(W_{ti};\nu_{0},e_{0,N})-s_{a}(W_{ti};\nu_{0},e_{0})\|=o_{p}(1)$

by Lemma B.1 and condition 1 of the Proposition, noting that

$\displaystyle\mathbb{E}\left[\frac{1}{N_{t}}\sum_{i=1}^{N_{t}}\|s_{a}(W_{ti};\nu_{0},e_{0,N})-s_{a}(W_{ti};\nu_{0},e_{0})\|\right]$
$\displaystyle\leqslant\kappa_{t}^{-1}\mathbb{E}_{0}[\|s_{a}(W;\nu_{0},e_{0,N})-s_{a}(W;\nu_{0},e_{0})\|]$
$\displaystyle\leqslant\kappa_{t}^{-1}(\mathbb{E}_{0}[\|s_{a}(W;\nu_{0},e_{0,N})-s_{a}(W;\nu_{0},e_{0})\|^{2}])^{1/2}$
$\displaystyle\leqslant\kappa_{t}^{-1}\delta_{N}$

for $t=1,\ldots,T$. Invertibility of $\mathbb{E}_{0}[s_{a}]$ from condition 2
ensures that $\hat{\theta}^{*}$ is well-defined with probability tending to 1
by (52).

Next, we fix $c\in{}^{p}$ with $\|c\|=1$ and a batch $t\in\\{1,\ldots,T\\}$.
Define

$U_{N,t,i}=\frac{c^{\top}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})}{(N_{t}\cdot
c^{\top}V_{t,N}c)^{1/2}},\quad i=1,\ldots,N_{t}$

where $V_{t,N}=\mathbb{E}_{t}[s(W;\theta_{0},\nu_{0},e_{0,N})^{\otimes 2}]$.
Evidently the random variables $U_{N,1},\ldots,U_{N,N_{t}}$ are independent,
with $\mathbb{E}[U_{N,t,i}]=0$ for all $i$ by (47), and
$\sum_{i=1}^{N_{t}}\mathbb{E}[U_{N_{t,i}}^{2}]=1$. Furthermore we have

$\lim_{N_{t}\rightarrow\infty}\sum_{i=1}^{N_{t}}\mathbb{E}[|U_{N,t,i}|^{q}]\leqslant
N_{t}\cdot\frac{\mathbb{E}_{t}[\|s(W;\theta_{0},\nu_{0},e_{0,N})\|^{q}]}{(N_{t}\cdot
c^{\top}V_{t,N}c)^{q/2}}=O(N_{t}^{1-q/2})=o(1)$

for all sufficiently large $N$ and some $q>2$ by condition 3. Then by the
Lyapunov CLT we have

$\sum_{i=1}^{N_{t}}U_{N,t,i}=c^{\top}\biggl{(}(c^{\top}V_{t,N}c)^{-1/2}\frac{1}{\sqrt{N_{t}}}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})\biggr{)}\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,1).$
$V_{t,N}=\mathbb{E}_{t}[s(W;\theta_{0},\nu_{0},e_{0,N})^{\otimes
2}]\rightarrow\mathbb{E}_{t}[s(W;\theta_{0},\nu_{0},e_{0})^{\otimes 2}]\equiv
V_{t}$

as $N\rightarrow\infty$. Therefore

$c^{\top}\biggl{(}\frac{1}{\sqrt{N_{t}}}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})\biggr{)}\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,c^{\top}V_{t}c)$

and since $c$ was arbitrary,

$\frac{1}{\sqrt{N_{t}}}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}(0,V_{t}),\quad
t=1,\ldots,T.$

With the left-hand side of the preceding display independent across batches
$t=1,\ldots,T$, we have

$\displaystyle\frac{1}{\sqrt{N}}\sum_{t=1}^{T}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})$
$\displaystyle=\sum_{t=1}^{T}\sqrt{\frac{N_{t}}{N}}\frac{1}{\sqrt{N_{t}}}\sum_{i=1}^{N_{t}}s(W_{ti};\theta_{0},\nu_{0},e_{0,N})$
$\displaystyle\stackrel{{\scriptstyle\mathrm{d\,\,}}}{{\rightarrow}}\mathcal{N}\biggl{(}0,\sum_{t=1}^{T}\kappa_{t}V_{t}\biggr{)}.$
(53)

With

$\displaystyle\sum_{t=1}^{T}\kappa_{t}V_{t}=\sum_{t=1}^{T}\kappa_{t}\mathbb{E}_{t}[s(W;\theta_{0},\nu_{0},e_{0})^{\otimes
2}]$
$\displaystyle=\mathbb{E}_{0}\left[\sum_{t=1}^{T}\kappa_{t}\frac{\mathrm{d}P_{t}}{\mathrm{d}P_{0}}(W)s(W;\theta_{0},\nu_{0},e_{0})^{\otimes
2}\right]$
$\displaystyle=\mathbb{E}_{0}[s(W;\theta_{0},\nu_{0},e_{0})^{\otimes 2}]$

the result of the Proposition follows by (51) and (52).

###  C.2 Proof of Corollary 3.1: CLT for
$\hat{\theta}^{*}_{\textnormal{AIPW}}$

Corollary 3.1, which applies Proposition 3.1 to prove a CLT for the oracle
estimator $\hat{\theta}_{\textnormal{AIPW}}^{*}$ of
$\theta_{0,\textnormal{ATE}}$, holds under the numbered generalizations of
Appendix A under one additional condition: that the mean functions are
stationary, meaning $\mathbb{E}_{t}[Y(z)\mid X=x]=m_{0}(z,x)$ for all
$t=1,\ldots,T$, $z=0,1$, and $x\in\mathcal{X}$. This condition is needed to
ensure Assumption A.2 holds, as discussed in Appendix A.

Our proof proceeds by showing that the conditions of the Corollary imply the
conditions of generalized Proposition 3.1 proven in the previous section with
$\theta_{0}=\theta_{0,\textnormal{ATE}}$ and $s=s_{\textnormal{AIPW}}(\cdot)$.
That is, first we show Assumption A.2 is satisfied with
$\theta_{0}=\theta_{0,\textnormal{ATE}}$ and $s=s_{\textnormal{AIPW}}(\cdot)$.
Then we show the three numbered conditions in Proposition 3.1.

First, for brevity let
$\nu_{0}=\nu_{0,\textnormal{AIPW}}(\cdot)=(m_{0}(0,\cdot),m_{0}(1,\cdot))$,
and note that for any $e$, $e^{\prime}\in\mathcal{F}_{\gamma}$ we have

$\displaystyle\mathbb{E}_{t,e}[s_{\textnormal{AIPW}}(W;\theta_{0,\textnormal{ATE}},\nu_{0},e^{\prime})]$
$\displaystyle=\mathbb{E}_{t,e}[m_{0}(1,X)-m_{0}(0,X)-\theta_{0,\textnormal{ATE}}]$
$\displaystyle\phantom{=}+\mathbb{E}_{t,e}\left[\frac{Z(Y(1)-m_{0}(1,X))}{e^{\prime}(X)}-\frac{(1-Z)(Y(0)-m_{0}(0,X))}{1-e^{\prime}(X)}\right]$
$\displaystyle=\mathbb{E}_{t,e}\left[\frac{e(X)}{e^{\prime}(X)}\mathbb{E}_{t,e}[Y(1)-m_{0}(1,X)\mid
X]\right]$
$\displaystyle\phantom{=}-\mathbb{E}_{t,e}\left[\frac{1-e(X)}{1-e^{\prime}(X)}\mathbb{E}_{t,e}[Y(0)-m_{0}(0,X)\mid
X]\right]$ $\displaystyle=0$

using unconfoundedness and stationarity of the mean function. All the
necessary expectations exist by our assumption that $\gamma>0$. Hence
Assumption A.2 is satisfied. Next we show the conditions of Proposition 3.1:

  1. 1.

Trivially we have

$\mathbb{E}_{0}[|s_{\textnormal{AIPW},a}(W;\nu_{0},e_{0,N})-s_{\textnormal{AIPW},a}(W;\nu_{0},e_{0})|^{2}]=\mathbb{E}_{0}[|-1-(-1)|^{2}]=0.$

Next we compute the following for each
$(\nu,e)\in\mathcal{N}\times\mathcal{F}_{\gamma}$:

$\displaystyle\phantom{=}\,s_{\textnormal{AIPW}}(W;\theta_{0},\nu,e)-s_{\textnormal{AIPW}}(W;\theta_{0},\nu_{0},e_{0})$
$\displaystyle=\left(1-\frac{Z}{e_{0}(X)}\right)(m(1,X)-m_{0}(1,X))$
$\displaystyle\phantom{=}+\left(1-\frac{1-Z}{1-e_{0}(X)}\right)(m_{0}(0,X)-m(0,X))$
$\displaystyle\phantom{=}+Z(Y-m(1,X))(e(X)^{-1}-e_{0}(X)^{-1})$
$\displaystyle\phantom{=}-(1-Z)(Y-m(0,X))((1-e(X))^{-1}-(1-e_{0}(X))^{-1}).$

Plugging in $(\nu,e)=(\nu_{0},e_{0,N})$ gives, by Minkowski’s inequality, that

$\displaystyle(\mathbb{E}_{0}[|s_{\textnormal{AIPW}}(W;\theta_{0},\nu_{0},e_{0,N})-s_{\textnormal{AIPW}}(W;\theta_{0},\nu_{0},e_{0})|^{2}])^{1/2}\leqslant
A_{0}+B_{0}$

where

$\displaystyle A_{0}$
$\displaystyle=(\mathbb{E}_{0}[Z^{2}(Y(1)-m_{0}(1,X))^{2}(e_{0,N}(X)^{-1}-e_{0}(X)^{-1})^{2}])^{1/2}$
$\displaystyle\leqslant\gamma^{-2}(\mathbb{E}_{0}[(Y(1)-m_{0}(1,X))^{2}(e_{0,N}(X)-e_{0}(X))^{2}])^{1/2}$
$\displaystyle=\gamma^{-2}(\mathbb{E}_{0}[(e_{0,N}(X)-e_{0}(X))^{2}v_{0}(1,X)])^{1/2}$
$\displaystyle\leqslant C\gamma^{-2}\|e_{0,N}-e_{0}\|_{2,P_{0}^{X}}.$

By an analogous computation

$\displaystyle B_{0}$
$\displaystyle=(\mathbb{E}_{0}[(1-Z)^{2}(Y(0)-m_{0}(0,X))^{2}((1-e_{0,N}(X))^{-1}-(1-e_{0}(X))^{-1})^{2}])^{1/2}$
$\displaystyle\leqslant C\gamma^{-2}\|e_{0,N}-e_{0}\|_{2,P_{0}^{X}}.$

The result now follows because

$\displaystyle|e_{0,N}(x)-e_{0}(x)|$
$\displaystyle=\Bigg{|}\sum_{t=1}^{T}\left(\frac{N_{t}}{N}-\kappa_{t}\right)e_{t}(x)\frac{\mathrm{d}P_{t}^{X}}{\mathrm{d}P_{0}^{X}}(x)\Bigg{|}$
$\displaystyle\leqslant\left(\sup_{1\leqslant t\leqslant
T}\kappa_{t}^{-1}\Big{|}\frac{N_{t}}{N}-\kappa_{t}\Big{|}\right)\sum_{1\leqslant
t\leqslant T}e_{t}(x)$ $\displaystyle\leqslant T\left(\sup_{1\leqslant
t\leqslant
T}\kappa_{t}^{-1}\Big{|}\frac{N_{t}}{N}-\kappa_{t}\Big{|}\right)=o(1)$

and hence

$\|e_{0,N}-e_{0}\|_{2,P_{0}^{X}}^{2}\leqslant\sup_{x\in\mathcal{X}}|e_{0,N}(x)-e_{0}(x)|^{2}=o(1).$

This bound only uses the fact that propensities are bounded between 0 and 1,
so

$\sup_{x\in\mathcal{X},e_{1}(\cdot),\ldots,e_{T}(\cdot)\in\mathcal{F}_{0}}|e_{0,N}(x)-e_{0}(x)|=o(1).$
(54)

  2. 2.

Evidently $\mathbb{E}_{0}[s_{\textnormal{AIPW},a}(W;\nu_{0},e_{0})]=-1$ is
invertible. Additionally for $z=0,1$ we have $\mathbb{E}[Y(z)^{2}]\leqslant C$
by the moment conditions in Assumption 3.1 so

$\mathbb{E}_{0}[m_{0}(z,X)^{2}]=\mathbb{E}_{0}[(\mathbb{E}_{0}(Y(z)\mid
X))^{2}]\leqslant\mathbb{E}_{0}[Y(z)^{2}]<\infty.$

Now
$s_{\textnormal{AIPW}}(W;\theta_{0,\textnormal{ATE}},\nu_{0,\textnormal{AIPW}},e_{0})$
is the sum of the following terms:

$\displaystyle m_{0}(1,X)\left(1-\frac{Z}{e_{0}(X)}\right),$ $\displaystyle
m_{0}(0,X)\left(\frac{1-Z}{1-e_{0}(X)}-1\right)-\theta_{0,\textnormal{ATE}},\quad\text{and}$
$\displaystyle Y\left(\frac{Z}{e_{0}(X)}-\frac{1-Z}{1-e(X)}\right).$

These are all square integrable, because $Z$, $(e_{0}(X))^{-1}$, and
$(1-e_{0}(X))^{-1}$ are all uniformly bounded.

  3. 3.

From $(\mathbb{E}_{0}[|Y(z)|^{q}])^{1/q}\leqslant C$ for $z=0,1$ by Assumption
3.1, we have

$\mathbb{E}_{0}[|m_{0}(z,X)|^{q}]=\mathbb{E}[|\mathbb{E}_{0}[Y(z)\mid
X]|^{q}]\leqslant\mathbb{E}_{0}[|Y(z)|^{q}]\leqslant C^{q}.$

Then each of

$\displaystyle\mathbb{E}_{0}\left[|m_{0}(1,X)|^{q}\Big{|}1-\frac{Z}{e_{0}(X)}\Big{|}^{q}\right],$
$\displaystyle\mathbb{E}_{0}\left[|m_{0}(0,X)|^{q}\Big{|}\frac{1-Z}{1-e_{0}(X)}-1\Big{|}^{q}\right],\quad\text{and}$
$\displaystyle\mathbb{E}_{0}\left[|Y|^{q}\Big{|}\frac{Z}{e_{0}(X)}-\frac{1-Z}{1-e(X)}\Big{|}^{q}\right]$

is at most $[C(1+\gamma^{-1})]^{q}$ and the desired condition holds by
Minkowski’s inequality.