EconBase
← Back to paper

Structure-agnostic Optimality of Doubly Robust Learning for Treatment Effect Estimation

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

51,241 characters · 12 sections · 50 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Structure-agnostic Optimality of Doubly Robust Learning for Treatment Effect Estimation

abstractAverage treatment effect estimation is the most central problem in causal inference with application to numerous disciplines. While many estimation strategies have been proposed in the literature, the statistical optimality of these methods has still remained an open area of investigation, especially in regimes where these methods do not achieve parametric rates. In this paper, we adopt the recently introduced structure-agnostic framework of statistical lower bounds, which poses no structural properties on the nuisance functions other than access to black-box estimators that achieve some statistical estimation rate. This framework is particularly appealing when one is only willing to consider estimation strategies that use non-parametric regression and classification oracles as black-box sub-processes. Within this framework, we prove the statistical optimality of the celebrated and widely used doubly robust estimators for both the Average Treatment Effect (ATE) and the Average Treatment Effect on the Treated (ATT), as well as weighted variants of the former, which arise in policy evaluation. \footnote{Accepted for presentation at the Conference on Learning Theory (COLT) 2025.}
keywordsCausal inference, semiparametric estimation, minimax lower bounds

Introduction

Estimating the average treatment effect is one of the central problems in causal inference and has found important applications in numerous disciplines such as economics hirano2003efficient,imbens2004nonparametric, education oreopoulos2006estimating, epidemiology little2000causal,wood2008empirical and political science mayer2011does. In view of its practical importance, a large body of work is devoted to developing statistically efficient estimators for the average treatment effect based on regression robins1994estimation,robins1995analysis,imbens2003mean, matching heckman1998matching,rosenbaum1989optimal,abadie2006large and propensity scores rosenbaum1983central,hirano2003efficient as well as their combinations.

Despite the plethora of estimation algorithms for the average treatment effect, little is known about the statistical limits of estimating the average treatment effect within some formal minimax optimality framework. Existing minimax optimality results only apply to statistical quantities that resemble or are natural variants of the average treatment effect (see e.g. robins2009semiparametric,balakrishnan2019hypothesis,kennedy2022minimax,robins2008higher). On the other hand, optimality results for the average effect are only known in the regime where the non-parametric components of the data generating process are estimable at a fast enough rate (typically $n^{-1/4}$). In this regime, the average effect is estimable at root-$n$ and the field of semi-parametric efficiency newey1994asymptotic has provided optimal variance constants that multiply the leading rate. Finally, the prior work of bradic2019minimax, characterizes minimax optimal conditions for root-$n$ estimability, albeit in a model where the effect is assumed to be constant for every unit in the population as well as other linearity assumptions. However, the optimal achievable estimation rate for any estimation quality of the non-parametric components has not been established. We provide the first tight statistical lower bound for the average treatment effect within the structure-agnostic minimax optimality framework balakrishnan2023fundamental, which is an optimality framework that is a natural fit for understanding the limits of estimators that use machine learning algorithms as black-box regression oracles. Within this structure-agnostic paradigm, we show that the celebrated doubly robust estimation algorithm robins1994estimation achieves minimax optimal mean-squared-error rates, up to constant factors.

Given a binary treatment $D\in \{0, 1\}$ and an outcome of interest $Y\in \mathbb{R}$, we let $Y(1), Y(0)$ denote the random potential outcomes that we would have observed from each unit, had we treated them with each possible treatment $d\in \{0,1\}$. Two central problems of causal analysis are the estimation of the average treatment effect (ATE) and the average treatment effect on the treated (ATT) heckman1998matching, defined correspondingly as the causal estimands:

align[align omitted — 155 chars of source]

We consider the case when all potential confounders $X\in \ensuremath{{\cal X}} \subseteq \mathbb{R}^K$ of the treatment and the outcome are observed; a setting that has received substantial attention in the causal inference literature. In particular, we will make the widely used assumption of conditional ignorability:

align[align omitted — 86 chars of source]

We assume that we are given data that consist of samples of the tuple of random variables $(X, D, Y)$, that satisfy the basic consistency property

align[align omitted — 25 chars of source]

Without loss of generality, the data generating process obeys the regression equations:

equation[equation omitted — 206 chars of source]

where $U,V$ are noise variables. The outcome regression function $g_0(d,x)$ and the propensity score $m_0(x)$ are commonly referred to as nuisance functions. Note that when the outcome $Y$ is also binary, then the non-parametric functions $g_0$ and $m_0$, as well as the marginal probability law of the covariates $X$, fully determine the likelihood of the observed data.

Under conditional ignorability, consistency and the overlap assumption that both treatment values are probable conditional on $X$, i.e., $m_0(X)\in [c, 1-c]$ almost surely, for some $c>0$, it is well known that the ATE and ATT are identified by the statistical estimands:

align[align omitted — 137 chars of source]

Our goal is to derive the statistically optimal estimation rates for the ATE and the ATT. We will also be interested in a weighted variant of the average treatment effect (WATE):

align[align omitted — 50 chars of source]

where $w(x)\in \mathbb{R}$ is a given weight function defined on $\ensuremath{{\cal X}}$. Such weighted average effects typically arise in the evaluation of personalized policies, where $w: \ensuremath{{\cal X}}\to \{0,1\}$ corresponds to a personalized treatment policy tao2019doubly,hirano2003efficient. Note that the ATE is a special case of the WATE with $w(X)=1$. Similar to the ATE, the WATE is identified by the statistical estimand:

equation[equation omitted — 105 chars of source]

Since the nuisance functions $g_0$ and $m_0$ in (ref) are unknown and may have complex structures, and since the dimension $K$ of the covariates $X$ can be large relative to the number of data $n$ in many applications, it is extremely suitable to apply modern machine learning (ML) methods for the non-parametric, flexible and adaptive estimation of these nuisance functions, including penalized linear regression methods belloni2014pivotal,van2014asymptotically,chernozhukov2022automatic,zou2005regularization, random forest methods breiman2001random,hastie2009random,biau2008consistency,wager2015adaptive,syrgkanis2020estimation, gradient boosted forests friedman2001greedy,buhlmann2003boosting,zhang2005boosting and neural networks schmidt2020nonparametric,farrell2021deep, as well as ensemble and model selection approaches that combine all the above using out-of-sample cross-validation metrics wolpert1992stacked,zhang1993model,freund1997decision,van2007super,sill2009feature,wegkamp2003model,Arlot2010,chetverikov2021cross.

Motivated by the wide adoption and use of black-box adaptive estimation methods polley2019package,ledell2020h2o,wang2021flaml,karmaker2021automl for these non-parametric components of the data generating process, as well as their superior empirical performance bach2024hyperparameter, even in the context of treatment effect estimation, we will examine the problem of statistical optimality within the structure agnostic minimax framework that was recently introduced in balakrishnan2023fundamental. In particular, the only assumption that we will be making about our data generating process is that we have access to estimates $\hat{g}$ and $\hat{m}$ that achieve some statistical error rate, as measured by the mean-squared error, i.e.

equation[equation omitted — 180 chars of source]

where for any function $v:\ensuremath{{\cal X}} \to \mathbb{R}$, we denote $\|v(X)\|_{P_X,2}:=\sqrt{\mathbb{E}[v(X)^2]}$. Having access to such estimates for these two non-parametric components and imposing the aforementioned estimation error constraints on the data generating process, we resolve the optimal statistical rate achievable by any estimation algorithm for the parameters of interest.

The structure agnostic framework is particularly appealing as it essentially restricts any estimation approach to only use non-parametric regression estimates as a black-box and not tailor the estimation strategy to particular structural assumptions about the regression function or the propensity. These further structural assumptions can many times be brittle and violated in practice, rendering the tailored estimation strategy invalid or low-performing. Hence, the structure agnostic statistical lower bound framework has the benefit that it yields lower bounds that can be matched by estimation procedures that are easy to deploy and robust in their details.

We show that up to constant factors no estimation algorithm for the WATE and ATT can achieve squared error rates for the parameter of interest that are better than:

align[align omitted — 230 chars of source]

respectively. These lower bounds apply even if we impose rate restrictions on stronger nuisance error metrics, e.g. $\|v(X)\|_{P_X,\infty} = \sup_{x\in \mathrm{supp}(P_X)} |v(x)|$. Furthermore, our lower bound constructions apply even when the outcome $Y$ is binary and, in the case of the WATE, they apply even when we know one of the two outcome response functions, i.e. $g_0(0, \cdot)$ or $g_0(1,\cdot)$. Importantly, these lower bounds are well-known to be achievable by the well-established and widely used doubly robust estimators derived from a first-order debiasing scheme, also known as estimators with the mixed bias property rotnitzky2021characterization.

For general non-parametric functional estimation, it has been shown decades ago that if the function possesses certain smoothness properties, then higher-order debiasing schemes can be designed that lead to improved error rates bickel1988estimating,birge1995estimation. Specifically, first-order debiasing methods are suboptimal even when the nuisance function estimators are minimax optimal. Estimators based on higher-order debiasing have also been proposed and analyzed for functionals that arise in causal inference problems robins2008higher,van2014higher,robins2017minimax,liu2017semiparametric,kennedy2022minimax. However, the fast rates of these methods crucially rely on the structure of the underlying function classes. Unlike first-order debiasing methods, higher-order methods are not structure-agnostic, in the sense that their error rates no longer apply to black-box estimators of the nuisance functions and the corresponding estimators are many times cumbersome to deploy in practice.

In contrast, our results show that first-order debiasing is structure-agnostic optimal for estimating both WATE and ATT. Our results extend the recent work of balakrishnan2023fundamental, which proposed the structure agnostic minimax optimality framework and proved the statistical optimality of doubly robust estimators of the expected conditional co-variance functional, defined as $\theta^{\textsc{Cov}} = \mathbb{E}[(D-\mathbb{E}[D\mid X])(Y - \mathbb{E}[Y\mid X])]$. However, the approach in balakrishnan2023fundamental cannot be easily modified to handle the average treatment effect functionals that we study here and which arguably possess a more central role in the causal inference literature.

Technical contributions

Our proof of the lower bounds uses the method of fuzzy hypotheses that reduces our estimation problem to the problem of testing a pair of mixtures of hypotheses. Such methods are widely adopted in establishing lower bounds for non-parametric functional estimation problems tsybakov2008introduction and have been used to address the minimax estimation errors of various causal functionals in different settings, including the expected conditional covariance robins2009semiparametric,balakrishnan2019hypothesis, CATE function at a given point kennedy2022minimax and variance-weighted ATE robins2008higher. Surprisingly however, there is very little understanding of the estimation limit of the vanilla ATE -- arguably a central parameter in causal inference -- even under Holder-smoothness assumptions that are extensively investigated in the literature. Focusing on the structure-agnostic setting introduced in the previous section, this paper takes an important step towards closing this gap.

Our main technical contribution is a collection of carefully-designed constructions of the hypotheses that are built on asymmetric perturbations in the space of nuisance functions. We construct the perturbations in a sequential manner, with the perturbation of one nuisance estimate being dependent on the other nuisance. We note that due to the more complicated relationships between the estimand and the data distribution, existing ways to construct composite hypotheses robins2009semiparametric,kennedy2022minimax,balakrishnan2023fundamental are no longer applicable to our setting, as we expand next.

In balakrishnan2023fundamental, the authors investigate the estimation problem of three functionals: quadratic functionals in Gaussian sequence models, quadratic integral functionals and the expected conditional covariance. They establish their lower bound by reducing it to lower-bounding the error of a related hypothesis testing problem. The error is then lower-bounded by constructing priors (mixtures) of the composite null and alternate distribution. The priors they construct are based on adding or subtracting bump functions on top of a fixed hypothesis in a symmetric manner, which is a standard proof strategy for functional estimation problems ingster1994minimax,robins2009semiparametric,arias2018remember,balakrishnan2019hypothesis. The reason why the proof strategy of balakrishnan2023fundamental fails for WATE and ATT is that the functional relationships between the nuisance parameters and these target parameters take significantly different forms. Specifically, the target parameters that balakrishnan2023fundamental investigates are all in the form of

equation[equation omitted — 93 chars of source]

where $f,g$ are unknown nuisance parameters that lie in some Hilbert space ${\mathcal{H}}$. To be concrete, consider the example of the expected conditional covariance $\theta^{\textsc{Cov}}$. Let $\mu_0(x) = \mathbb{E}\left[ Y \mid X=x\right]$, then we have that $\theta^{\textsc{Cov}} = \mathbb{E}[DY] - \int m_0(x)\mu_0(x)\text{d} p_X(x)$ where $p_X$ is the marginal density of $X$. The first term, $\mathbb{E}[DY]$, can be estimated with a standard ${\mathcal{O}}(n^{-{1}/{2}})$ rate, so what remains to be estimated is the second term which is exactly in the form of (ref). However, the ATE and ATT functionals do not take this inner product form. Instead they are, respectively, of the form:

align*[align* omitted — 309 chars of source]

Stepping outside of the realm of inner product functionals is the major challenge in extending existing approaches of establishing lower bounds to the problem of estimating WATE and ATT, and very different constructions are required, which is our main technical innovation.

Notation

We use $P_X$ to denote the marginal distribution of the confounding factors $X$ in the model ((ref)). For any function $f:\mathbb{R}^n\mapsto\mathbb{R}^k$ and distribution $P$ over $\mathbb{R}^n$, we define its $L^r$-norm as $\left\Vert f\right\Vert _{P,r} = \left(\int \left\Vert f\right\Vert^r \text{d}P \right)^{{1}/{r}},\quad r\in(0,+\infty)$ and $\left\Vert f\right\Vert_{P,\infty} = \text{ess sup} \left\{ f(X): X\sim P\right\}$. We also slightly abuse notation and use $\left\Vert f\right\Vert _{r}$ instead, when the distribution is clear from context. For two sequences $(a_n)_{n\geq 1}$ and $(b_n)_{n\geq 1}$, we write $a_n={\mathcal{O}}(b_n)$ if there exists a constant $C>0$ such that $|a_n|\leq C|b_n|,\forall n\geq 1$, and we write $a_n=\Omega(b_n)$ if there exists a constant $c>0$ such that $|a_n|\geq c|b_n|,\forall n\geq 1$.

Structure-agnostic estimation of average treatment effect

To analyze the statistical limit of estimating treatment effect without making assumptions on regularity properties of nuisance functions, we adopt the structure-agnostic framework introduced by balakrishnan2023fundamental. Specifically, we assume the existence of black-box estimates $\hat{m}(x)$ and $\hat{g}(d,x)$ of $m(x)$ and $g(d,x)$ that are accurate in the sense of $L^2$ distance:

equation[equation omitted — 293 chars of source]

where $e_n, e_n'$ and $f_n$ are arbitrary positive numbers that depend on the sample size $n$ used to estimate the nuisance functions. Note that here we assume that the estimators $\hat{m},\hat{g}$ are already known to the statistician rather than a part of the estimation process. The reason for considering this setup is that we do not want to open the black box of how these estimators are obtained. In practice, these estimators can be obtained by leveraging estimation methods such as Lasso bickel2009simultaneous, random forest syrgkanis2020estimation, deep neural networks chen1999improved,schmidt2020nonparametric,farrell2021deep among others.

Apart from the nuisance estimators, we also assume access to i.i.d. data $\{(X_i,D_i,Y_i)\}_{i=1}^n$ that are also independent of the data used to obtain the nuisance estimators. In this way, we fully disentangle the complete estimation procedure into a learning phase where any estimation method can be used to obtain nuisance estimators from a portion of the data, and an second phase that leads to the final parameter estimate. While estimation of WATE/ATT does not necessarily need to follow this procedure, this is a typical pipeline implemented in practice, since it allows a flexible combination of black-box machine learning methods and estimators that cleverly leverage the structural properties of the model. As discussed before, the doubly robust estimators achieve error rates in the form of (ref). Surprisingly, however, it has long been unknown whether one can actually do better than doubly robust estimators, which we address next.

comment\section{Structure-agnostic high-probability upper bounds} We first introduce achievable structure-agnostic upper bounds that are quite well-understood in the existing literature. In particular, revisiting the asymptotic normality and root-$n$ consistency results in chernozhukov2017double,chernozhukov2018double,foster2023orthogonal, we argue that the sample-splitting variants of the doubly robust estimators presented in (ref) and (ref), achieve high-probability error bounds, even in the regime where the parameter error rate is not root-$n$. These high probability squared-error bounds depend on the product of the errors in the relevant outcome regressions and the propensity (i.e. $e_n\cdot f_n$), plus the typical parametric rate of $1/n$. When the product of nuisance errors decays faster then $1/n$, then the nuisance estimation error can be asymptotically ignored and the estimate can also be shown to be asymptotically normally distributed. However, the estimation error bounds are valid even when the product of the errors does not lie in this “fast” parametric regime. This result follows from standard arguments and is deferred to the online appendix. \begin{theorem} Suppose that there exists a constant $c\in(0,1)$ such that $c\leq \hat{m}(x)\leq 1-c, \forall x\in\mathrm{supp}(X)$ and $|Y|\leq G$ a.s., for some constant $G$. Then for any $\delta>0$, there exists a constant $C_{\delta}$ such that the doubly robust estimator of the WATE (defined in (ref)) achieves estimation error \begin{equation} \notag \left|\hat{\theta}^{WATE}-\theta^{WATE}\right|^2 \leq C_{\delta}\left( \max\{e_n,e_n'\}\cdot f_n\cdot \left\Vert w\right\Vert_{P_X,\infty}^2 + \frac{1}{n}\left\Vert w\right\Vert_{P_X,2}^2\right). \end{equation} with probability $\geq 1-\delta$. Moreover, the doubly robust estimator of the ATT (defined in (ref)) achieves estimation error \begin{equation} \left|\hat{\theta}^{ATT}-\theta^{\textsc{ATT}}\right|^2 \leq C_{\delta}\left(e_n\cdot f_n+\frac{1}{n}\right) \end{equation} with probability $\geq 1-\delta$. \end{theorem} Theorem (ref) implies that with high probability, the estimation error of the debiased estimator ((ref)) is upper-bounded by the sum of the oracle error which equals ${1}/{n}$ multiplied by the $L^2$ norm of weight function $w$, and the product of the error in estimating nuisance functions $m_0$ and $g_0$, multiplied by the $L^{\infty}$ norm of $w$. Similarly, for estimating the \textsc{ATT}, Theorem (ref) implies that with high probability, the error of the doubly robust estimator ((ref)) is upper-bounded by the sum of the oracle error ${1}/{n}$ and the product of the error in estimating nuisance functions $m_0$ and $g_0(0,\cdot)$. Given the high-probability upper bounds in Theorem (ref), it is natural to ask whether these structure-agnostic guarantees achieved by the doubly robust estimators can be further improved. Our results in the next section provide a negative answer (up to constants).

Main results

In this section, we present our main results that lower-bound the estimation errors in the structural-agnostic setting.

We restrict ourselves to the case of binary outcomes:

assumptionThe outcome variable $Y$ is binary, i.e., $Y\in\{0,1\}$.

Given that the black-box nuisance function estimators satisfy (ref), we define the following constraint set

equation[equation omitted — 445 chars of source]

where

equation[equation omitted — 73 chars of source]

Note that introducing Assumption (ref) and constraints on $P_X$ in (ref) only strengthens the lower bound that we are going to prove, since they provide additional information on the ground-truth model. Moreover, the constraints $0\leq m(x), g(d,x)\leq 1$ naturally holds due to the fact that both the treatment and outcome variables are binary. We then define the minimax $(1-\gamma)$-quantile risk of estimating $\theta^{\textsc{WATE}}$ over a function space ${\mathcal{F}}$ as

equation[equation omitted — 356 chars of source]

where ${\mathcal{Q}}_{P,\gamma}(X)=\inf\left\{x\in\mathbb{R}: P[X\leq x]\geq\gamma\right\}$ denotes the quantile function of a random variable $X$, and $P_{m^*,g^*}$ is the joint distribution of $\{(X_i,D_i,Y_i)\}_{i=1}^n$ which is uniquely determined by the functions $m^*$ and $g^*$. Specifically, let $\mu$ be the uniform distribution on ${\mathcal{X}}\times{\mathcal{D}}\times{\mathcal{Y}}=[0,1]^K\times\{0,1\}\times\{0,1\}$, then the density $p_{m^*,g^*}={\text{d}P_{m^*,g^*}}/{\text{d}\mu}$ can be expressed as

equation[equation omitted — 105 chars of source]

According to (ref), $\mathfrak{M}_{n,\gamma}^{\textsc{WATE}}\left({\mathcal{F}}\right) \geq \rho$ would imply that for any estimator $\hat{\theta}$ of WATE, there must exist some $(m^*,g^*)\in{\mathcal{F}}$, such that under the induced data distribution, the probability of $\hat{\theta}$ having estimation error $\geq\rho$ is at least $1-\gamma$. This provides a stronger form of lower bound compared with the minimax expected risk defined in balakrishnan2023fundamental, in the sense that the lower bound $\mathfrak{M}_{n,\gamma}^{\textsc{WATE}}\left({\mathcal{F}}\right) \geq \rho$ implies a lower bound $(1-\gamma)\rho$ of the minimax expected risk, but the converse does not necessarily hold.

Similarly, one can define the minimax quantile risk for estimating ATT as

equation[equation omitted — 404 chars of source]

The main objective of this section is to derive lower bounds for $\mathfrak{M}_{n,\gamma}^{\textsc{WATE}}\left({\mathcal{F}}_{e_n,e_n',f_n}\right)$ and $\mathfrak{M}_{n,\gamma}^{\textsc{ATT}}\left({\mathcal{F}}_{e_n,e_n',f_n}\right)$ in terms of $e_n,e_n',f_n$ and $n$. We also need to assume that the estimators $\hat{m}(x): [0,1]^K\mapsto[0,1]$ and $\hat{g}(d,x):\{0,1\}\times[0,1]^K\mapsto[0,1]$ are bounded away from $0$ and $1$.

assumptionThere is a constant $c\in(0,1/2)$ such that $c \leq \hat{m}(x), \hat{g}(d,x) \leq 1-c, \forall d\in\{0,1\}$, $x\in[0,1]^K$.

The assumption that $c \leq \hat{m}(x)\leq 1-c$ is common in deriving upper bounds for doubly robust estimators. On the other hand, the assumption that $c\leq \hat{g}(d,x)\leq 1-c$ is typically not needed for deriving upper bounds, but it is also made in prior works for proving lower bounds for estimating the expected conditional covariance $\mathbb{E}\left[\mathrm{Cov}(D,Y\mid X) \right]$ robins2009semiparametric,balakrishnan2023fundamental. Now we are ready to state our main results.

theoremFor any constant $\gamma\in\left({1}/{2},1\right)$ and estimators $\hat{m}(x)$ and $\hat{g}(d,x)$ that satisfy Assumption (ref), for any given weight function $w$, the minimax risk of estimating the WATE is \begin{equation} \notag \mathfrak{M}_{n,\gamma}^{WATE}\left({\mathcal{F}}_{e_n,e_n',f_n}\right) = \Omega\left( \max\{e_n,e_n'\} f_n\cdot\Vert w\Vert _{P_X,\infty}^2 + \Vert w\Vert_{P_X,2}^2/n\right). \end{equation}
remarkIf we only assume that $c\leq\hat{m}(x), \hat{g}(1,x) \leq 1-c$ in Assumption (ref), then we would have the lower bound \begin{equation} \notag \mathfrak{M}_{n,\gamma}^{WATE}\left({\mathcal{F}}_{e_n,e_n',f_n}\right) = \Omega\left( e_n' f_n\cdot\Vert w\Vert _{P_X,\infty}^2 + \Vert w\Vert_{P_X,2}^2/n\right). \end{equation} Furthermore, this lower bound still holds in the case where we know the baseline response, i.e., $\hat{g}(0,x)=g_0(0,x)=0$.
theoremFor any constant $\gamma\in\left({1}/{2},1\right)$ and estimators $\hat{m}(x)$ and $\hat{g}(d,x)$ that satisfy Assumption (ref), the minimax risk of estimating the ATT is given by \begin{equation} \notag \mathfrak{M}_{n,\gamma}^{ATT}\left({\mathcal{F}}_{e_n,e_n',f_n}\right) = \Omega\left( e_n f_n + 1/n\right). \end{equation}
remarkAs discussed in Section (ref), generic machine learning estimators are typically guaranteed to have small $L^2$ errors. However, the lower bounds presented in Theorem (ref) and (ref) still hold even if we replace the $L^2$ norm constraints in (ref) are replaced with stronger $L^r (2\leq r\leq +\infty)$ constraints.

Theorems (ref) and (ref) provide lower bounds of the minimax estimation errors for the WATE and ATT, in terms of the sample size and the estimation error of the black-box nuisance function estimators. Our lower bounds exactly matches the well-known upper bounds attained by the doubly robust estimators (see Section (ref)), indicating that doubly robust estimators are minimax optimal in the structural-agnostic setup.

Proof of Theorem (ref)

In this section, we give the proof outline of our main result, Theorem (ref), for the lower bound of estimating WATE. Omitted details in the proof can be found in The proof of Theorem (ref) can be found in Section (ref) in the appendix. We first introduce some preliminary results that our proof will rely on.

Preliminaries

Let $H$ be the Hellinger distance defined as $H(P,Q) = \frac{1}{2}\int \big(\sqrt{P(\text{d}x)} - \sqrt{Q(\text{d}x)}\big)^2$ for any distributions $P,Q$. The first result that we will introduce is due to robins2009semiparametric and upper-bounds the Hellinger distance between two mixtures of product measures.

Formally, let $\mathcal{X}=\cup_{j=1}^m \mathcal{X}_j$ be a measurable partition of the sample space. Given a vector $\lambda=\left(\lambda_1, \ldots, \lambda_m\right)$ in some product measurable space $\Lambda=\Lambda_1 \times \cdots \times \Lambda_m$, let $P$ and $Q_{\lambda}$ be probability measures on $\mathcal{X}$ such that the following statements hold:

enumerate[1.] • $P\left(\mathcal{X}_j\right)=Q_\lambda\left(\mathcal{X}_j\right)= p_j$ for every $\lambda \in \Lambda$, and • The probability measures $P$ and $Q_{\lambda}$ restricted to $\mathcal{X}_j$ depend on the $j$-th coordinate $\lambda_j$ of $\lambda$ only.

Let $p$ and $q_{\lambda}$ be the densities of the measures $P$ and $Q_\lambda$ that are jointly measurable in the parameter $\lambda$ and the observation $x$, and $\pi$ be a probability measure on $\Lambda$. Define $b= m\max _j \sup _\lambda \int_{\mathcal{X}_j} \left(q_\lambda-p\right)^2/p d \mu$ and the mixed density $q=\int q_\lambda d \pi(\lambda)$, then we have the following result.

lemma(robins2009semiparametric, Theorem 2.1, simplified version) Suppose that the mixed densities are equal, i.e. that $q=p$, and that $n\max\{1,b\}\max_j p_j \leq A$ for all $j$ for some positive constant $A$, then there exists a constant $C$ that depends only on $A$ such that, for any product probability measure $\pi=\pi_1 \otimes \cdots \otimes \pi_m$, $ H\left(P^{\otimes n}, \int Q_\lambda^{\otimes n} d \pi(\lambda)\right) \leq \max_j p_j\cdot Cn^2b^2. $
remarkTheorem (ref) considers a special case of robins2009semiparametric, Theorem 2.1. The original variant of the theorem considers a more general setting where the measures $p$ are also indexed by $\lambda$, i.e. $p_{\lambda}$ and where $p$ is the mixture density. Here, we only need the special cases where all $P_{\lambda}$'s are equal to $P$. The original version of the theorem also required that all $p_{\lambda}$ satisfy that $\underline{B} \leq p \leq \bar{B}$ for some constants $\underline{B}, \bar{B}$. In our special case, we no longer need to assume that. The only step in the proof of robins2009semiparametric that makes use of this assumption is that $\max_j \sup_{\lambda} \int_{{\mathcal{X}}_j} \frac{p^2}{p_{\lambda}}\frac{\text{d}\mu}{p_j} \leq \frac{\bar{B}}{\underline{B}}$ (see the arguments following their proof of Lemma 5.2). However, in our setting this term is simply \begin{equation} \notag \max_j \sup_{\lambda} \int_{{\mathcal{X}}_j} \frac{p^2}{p_{\lambda}}\frac{d\mu}{p_j} = \max_j p_j^{-1} \int_{{\mathcal{X}}_j} p d\mu = \max_j p_j^{-1} P({\mathcal{X}}_j) = 1. \end{equation}
lemma(tsybakov2008introduction, Theorem 2.15) Let $\pi$ be a probability distribution on a set (measure space) of distributions $\mathcal{P}$ with common support ${\mathcal{X}}$, which induce the distribution $ Q_1(A)=\int Q^{\otimes n}(A) d \pi(Q), \quad \forall A \subset {\mathcal{P}}. $ Suppose that there exists $P\in{\mathcal{P}}$ and a functional $T: {\mathcal{P}}\mapsto\mathbb{R}$ which satisfies \begin{equation} T(P)\leq c, \quad \pi(\{Q: T(Q) \geq c+2 s\})=1 \end{equation} for some $s>0$. If $H^2\left(P^{\otimes n}, Q_1\right) \leq \delta<2$, then $\inf_{\hat{T}: {\mathcal{X}}^n\mapsto\mathbb{R}} \sup_{P \in \mathcal{P}} P\left[\left|\hat{T}-T(P)\right|\geq s\right] \geq {\big(1-\sqrt{\delta(1-\delta / 4)}\big)}/{2}.$ Then it holds that $\inf_{\hat{T}: {\mathcal{X}}^n\mapsto\mathbb{R}} \sup_{P \in \mathcal{P}} {\mathcal{Q}}_{P,1-\gamma}\left(\left|\hat{T}-T(P)\right|^2\right) \geq s^2$, where $\gamma = {\big(1+\sqrt{\delta(1-\delta / 4)}\big)}/{2}$,

Partitioning the covariate space

The following lemma states that for an arbitrary weight function $w(x)$, one can always partition the domain into two subsets that have the same amount of weights.

definitionWe say that a set ${\mathcal{S}}\subseteq\mathbb{R}^K$ is a hyperrectangle collection if it can be partitioned into a finite number of disjoint hyperrectangles in $\mathbb{R}^K$.
lemmaLet ${\mathcal{S}}\subseteq\mathbb{R}^K$ be a hyperrectangle collection and $w(x): [0,1]^K\mapsto\mathbb{R}$ be a non-negative Lebesgue-integrable function such that $\int_{[0,1]^K} w(x) \text{d} \mu_L(x) > 0$, then ${\mathcal{S}}$ can be partitioned into two hyperrectangle collections ${\mathcal{S}}_1, {\mathcal{S}}_2$ such that $\mu_L({\mathcal{S}}_1) = {\mu_L({\mathcal{S}})}/2$ and \begin{equation} \notag \int_{{\mathcal{S}}_1} w(x) d \mu_L(x) = \frac{1}{2} \int_{{\mathcal{S}}} w(x) d \mu_L(x), \end{equation} where $\mu_L$ is the Lebesgue measure on $\mathbb{R}^K$.

Let $P_X$ be the uniform distribution on $\mathrm{supp}(X)=[0,1]^K$ and $p_X$ be its density. It is easy to see that $P_X\big[|w(X)|>\Vert w\Vert _{P_X, \infty}/2\big] > 0$. Assuming without loss of generality that $P_X\big[w(X)>\Vert w\Vert _{P_X, \infty}/2\big] > 0$ (otherwise we consider $-w$ instead of $w$), we can define the "truncated" weight function $\hat{w}(x)=w(x)\mathbbm{1}\left\{w(x)>\Vert w\Vert_{P_X, \infty}/2\right\}$. Applying Lemma (ref) to $\hat{w}\cdot w$, recursively, for $m$ times, with $m\in\mathbb{Z}_{+}$, we can partition $[0,1]^K$ into $M=2^m$ hyperrectangle collections $B_1, B_2,\cdots, B_M$, such that $\mu_L(B_j) = 1/M$ and

equation[equation omitted — 143 chars of source]

Since $P_X$ is the uniform distribution on $[0,1]^K$, the above implies that

equation[equation omitted — 173 chars of source]

Let $\lambda_i, i=1,2,\cdots, M$ be i.i.d. Rademacher random variables taking values $+1$ and $-1$ both with probability $0.5$. We define

equation[equation omitted — 182 chars of source]
propositionWe have \begin{subequations} \begin{align} &\mathbb{E}_{\lambda} \Delta(\lambda,x) = \sum_{j=1}^{M/2} \mathbb{E}\lambda_j\left(\mathbbm{1}\left\{x\in B_{2j}\right\}-\mathbbm{1}\left\{x\in B_{2j-1}\right\}\right) = 0,\quad\forall x\in[0,1]^K \\ &\mathbb{E}_X w(X)\hat{w}(X) \Delta(\lambda,X) = 0,\quad\forall\lambda\in\{0,1\}^{M/2} \\ &\Delta(\lambda,x)^2 = \sum_{j=1}^{M/2} \left(\mathbbm{1}\left\{x\in B_{2j}\right\}-\mathbbm{1}\left\{x\in B_{2j-1}\right\}\right)^2 = 1,\quad\forall x\in[0,1]^K,\lambda\in\{0,1\}^{M/2}. \end{align} \end{subequations}
remarkThe construction of bump functions $\Delta(\lambda,x)$ in the form of (ref) has also been used in a line of prior works for proving minimax lower bounds balakrishnan2023fundamental. However, here we need to carefully construct the partition $B_j$ of the whole domain to handle non-uniform weights. We note that if we only wanted to deal with an ATE and not a WATE, then we would have simply chosen $B_j$ to be an equi-partition of the $[0,1]^K$ space and the above constructions of the regions $B_j$, related to balancing the given weights, would not be needed.

Core part of lower bound construction

Having completed all preparation steps, we are now ready to present our proof for Theorem (ref). The remaining part of Section (ref) is organized as follows. In Section (ref) and (ref), we first establish our lower bound $\Omega\left(e_n'f_n\Vert w\Vert _{P_X,\infty}^2\right)$ under the following weaker version of Assumption (ref), as previously mentioned in Remark (ref):

assumptionThere exists a constant $c>0$ such that $c \leq \hat{m}(x), \hat{g}(1,x) \leq 1-c, \forall x\in[0,1]^K$.

We separately present our proof of this lower bound for the two cases $e_n'\geq f_n$ and $e_n'<f_n$. Interestingly, these two cases need to be handled separately using different constructions of the composite hypotheses. In Section S.2 in the supplementary material, we show how the lower bound $\Omega\left(e_nf_n\Vert w\Vert _{P_X,\infty}^2\right)$ can be derived in a similar fashion. To conclude our proof, it remains to prove the lower bound ${\mathcal{O}}\left(n^{-1}\Vert w\Vert _{P_X,2}^2\right)$, which is the standard oracle error and can be found in the supplementary material.

Case 1: $e_n\geq f_n$

In this case, we define

equation[equation omitted — 366 chars of source]

where $\alpha,\beta>0$ are constants that will be specified later in Lemma (ref), where we will verify that $(m_{\lambda},g_{\lambda})$ belongs to the constrained set ${\mathcal{F}}_{e_n,e_n',f_n}$ and thus are valid probabilities in particular. Compared with standard approaches for constructing the composite hypotheses ingster1994minimax,robins2009semiparametric,arias2018remember,balakrishnan2019hypothesis, we employ an asymmetric construction which means that the nuisance functions are non-linear in the Rademacher variables $\lambda$ (in particular the function $g_{\lambda}$ depends non-linearly in $\lambda$ due to the dependence on $m_{\lambda}$ in the denominator). As discussed in Section (ref), such type of non-standard constructions are necessary since the functional that we need to estimate has a different structure than those handled in previous works. We first prove some basic properties of our construction.

propositionFor all $x\in[0,1]^K$, we have \begin{subequations} \begin{align} &\mathbb{E}_{\lambda} m_{\lambda}(x) = \hat{m}(x) - \hat{m}(x)\frac{\beta}{\hat{g}(1,x)}\hat{w}(x)\mathbb{E}_{\lambda}\Delta(\lambda,x) = \hat{m}(x) \\ &\mathbb{E}_{\lambda} \left[g_{\lambda}(1,x)m_{\lambda}(x)\right] = \hat{m}(x)\left(\hat{g}(1,x)+\alpha\hat{w}(x)\mathbb{E}_{\lambda}\Delta(\lambda,x)\right) = \hat{g}(1,x)\hat{m}(x). \end{align} \end{subequations}

We start by bounding the $L^2$ distance from $g_{\lambda}, m_{\lambda}$ to $\hat{g},\hat{m}$.

lemmaAssuming that $\beta\leq{c\left\Vertw\right\Vert_{P_X,\infty}^{-1}}/2$ where $c$ is the constant introduced in Assumption (ref), then the following holds for all $0<r\leq +\infty$: \begin{equation} \notag \begin{aligned} \left\Vertg_{\lambda}(1,X)-\hat{g}(1,X)\right\Vert_{P_X,r} &\leq 2(\alpha+c^{-1}\beta)\Vert \hat{w}(X)\Vert _{P_X,r},\\ \left\Vertm_{\lambda}(X)-\hat{m}(X)\right\Vert_{P_X,r} &\leq c^{-1}\beta\Vert \hat{w}(X)\Vert _{P_X,r}. \end{aligned} \end{equation}

Let $Q_{\lambda}$ be the joint distribution of $(X,D,Y)$ induced by $g_{\lambda}$ and $m_{\lambda}$ and $\mu$ be the uniform distribution on $[0,1]^K\times\{0,1\}\times\{0,1\}$. Define $q_{\lambda} = {\text{d} Q_{\lambda}}/{\text{d} \mu}$. Similarly, let $\hat{P}$ be the joint distribution of $(X,D,Y)$ induced by $\hat{g}$ and $\hat{m}$, and $\hat{p}={\text{d} \hat{P}}/{\text{d}\mu}$. The next lemma states that the mixture of $Q_{\lambda}$ with prior $\pi(\lambda)$ is exactly equal to $\hat{P}$.

lemmaLet $Q = \int Q_{\lambda} \text{d} \pi(\lambda)$ and $q = {\text{d} Q}/{\text{d} \mu} = \int q_{\lambda} \text{d} \pi(\lambda)$, then $\hat{p}=q$.

The following lemma implies that the Hellinger distance between the empirical distribution under $\hat{P}$ and $Q_{\lambda}$ with prior $\pi(\lambda)$ can be made arbitrarily small, as long as the domain $\mathrm{supp}(X)$ is partitioned into sufficiently many pieces.

lemmaFor any $\delta>0$, as long as $M \geq \max\{n,32Cn^2/(c^4\delta)\}$ where $c$ is the constant introduced in Assumption (ref) and $C$ is the constant implied by Lemma (ref) for $A=4c^{-2}$, we have $H^2\left(\hat{P}^{\otimes n},\int Q_{\lambda}^{\otimes n}\text{d} \pi(\lambda)\right) \leq \delta.$

As the final building block for establishing our lower bound, we prove the following lemma, which implies that with proper choices of $\alpha$ and $\beta$, $m_{\lambda}, g_{\lambda}$ are close (in the sense of $L^2$-distance) to $\hat{m}$ and $\hat{g}$ respectively, and that the separation condition ((ref)) holds with distance $s=\Omega\left(\sqrt{e_n f_n}\left\Vertw\right\Vert_{P_X,\infty}\right)$.

lemmaLet $\alpha = \sqrt{e_n'}/\big(4\Vert \hat{w}(X)\Vert _{P_X,2}\big),\quad \beta = c\sqrt{f_n}/\big(4\Vert \hat{w}(X)\Vert _{P_X,2}\big)$, then for sufficiently large $n$, we have $(m_{\lambda}, g_{\lambda})\in{\mathcal{F}}_{e_n,e_n',f_n}$ and \begin{equation} \mathbb{E}_X\left[w(X) g_{\lambda}(1,X)\right] \geq \mathbb{E}\left[w(X) \hat{g}(1,X) \right] + \frac{1}{2}\alpha\beta\mathbb{E}\left[\frac{w(X)\hat{w}(X)^2}{\hat{g}(1,X)}\right],\forall\lambda\in\{0,1\}^{M/2}. \end{equation}

We are now ready to prove Theorem (ref) in the case when $e_n'\geq f_n$. For any $\gamma>{1}/{2}$, there exists some $\delta\in(0,2)$ such that ${(1+\sqrt{\delta(1-\delta/4)})}/{2}=\gamma$. We choose $M \geq \max\{n,32Cn^2/{c^4\delta}\}$ and ${\mathcal{P}}=\{\hat{P}\}\cup\left\{Q_{\lambda}:\lambda\in\{0,1\}^{M/2}\right\}$, $P=\hat{P}$, $\pi$ be the discrete uniform distribution on $\left\{Q_{\lambda}:\lambda\in\{0,1\}^{M/2}\right\}$, $s=\frac{1}{4}\alpha\beta \mathbb{E}\big[{w(X)\hat{w}(X)^2}/{\hat{g}(1,X)}\big]$ in the context of Lemma (ref). Then Lemma (ref) and (ref) imply that all the listed conditions are satisfied for the WATE functional $T(P) = \theta^{\textsc{WATE}}(P) = \mathbb{E}_P\left[ w(X)\left( g(1,X)-g(0,X)\right)\right]$. Therefore, by Lemma (ref), we have

equation[equation omitted — 965 chars of source]

Case 2: $f_n > e_n'$

In this case, we consider a different construction as follows:

equation[equation omitted — 406 chars of source]

where $\Delta(\lambda,x)$ is defined in (ref) and $\alpha,\beta>0$ are constants that need to be specified later. The remaining steps follow a similar reasoning as the previous subsection. Due to space limit, we present the proof in Section (ref).

Conclusion

We investigated the statistical limit of treatment effect estimation in the structural-agnostic regime, which is an appropriate lower bound technique when one wants to only consider estimation strategies that use generic black-box estimators for the various nuisance functions involved in the estimation of treatment effects. We establish the minimax optimality of the celebrated and widely used in practice doubly robust learning strategies via reducing the estimation problem to a hypothesis testing problem, and lower bound the error of any estimation algorithm via non-standard constructions of the fuzzy hypotheses. Our results show that these estimators are optimal, in the structure agnostic sense, even in the slow rate regimes, where the implied rate for the functional of interest is slower than root-$n$. Hence, any improvements upon these estimation strategies need to incorporate elements of the structure of the nuisance functions and cannot simply invoke generic adaptive regression approaches as black-box sub-processes. While the focus of this paper is on treatment effect estimation, we believe that our techniques can be extended to address structure agnostic minimax lower bounds of more general functional estimation problems.

\acks{VS is supported by NSF Award IIS-2337916. JJ is partially supported by NSF Award IIS-2337916.}