EconBase
← Back to paper

Sparsity Double Robust Inference of Average Treatment Effects

The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.

43,369 characters

Sparsity Double Robust Inference of Average Treatment Effects



\maketitle

\begin{abstract}
Many popular methods for building confidence intervals on causal effects under high-dimensional confounding
require strong ``ultra-sparsity'' assumptions that may be difficult to validate in practice. To alleviate this difficulty,
we here study a new method for average treatment effect estimation that yields asymptotically exact
confidence intervals assuming that either the conditional response surface or the conditional probability of
treatment allows for an ultra-sparse representation (but not necessarily both). This guarantee allows us to
provide valid inference for average treatment effect in high dimensions under considerably more generality than
available baselines. In addition, we showcase that our results are semi-parametrically efficient.
\end{abstract}

\section{Introduction}

Average treatment effect estimation is a core problem in causal inference, and has been the topic
of a considerable amount of recent literature \citep{imbens2015causal}. In this paper, we focus on the
task average treatment effect estimation with high-dimensional confounders:  We have access to $n$ \emph{i.i.d.}
samples $(X_i, \, Y_i, \, W_i) \in \xx \times \mathbb{R} \times \cb{0, \, 1}$, where $X_i$ denotes high-dimensional
pre-treatment features ($\xx \subset \mathbb{R}^p$ with $p \gg n$),
$W_i$ is the treatment assignment, and $Y_i$ is our outcome of interest. Causal effects are defined via potential
outcomes $\cb{Y_i(0), \, Y_i(1)}$, such that we observe $Y_i = Y_i(W_i)$ and the average treatment effect is
defined as $\tau = \mathbb{E}{Y_i(1) - Y_i(0)}$ \citep{neyman1923applications,rubin1974estimating}.
Finally, we assume that there are no unmeasured confounders, i.e., the treatment assignment $W_i$ may not
be randomized, but can be treated as such once we control for $X_i$, i.e.,
$\cb{Y_i(0), \, Y_i(1)} \indep W_i \cond X_i$ \citep{rosenbaum1983central}.
Throughout, we also assume overlap, such that $\eta \leq \mathbb{P} \p{W_i \cond X_i = x} \leq 1 - \eta$ for all $x$ and some
$\eta > 0$.


In the low-dimensional case, one of the most prominent approaches to average treatment effect estimation
is via augmented inverse-propensity weighting \citep{robins1994estimation},
\begin{equation}
\label{eq:AIPW}
\hat{\tau} = \frac{1}{n} \sum_{i = 1}^n \p{\hat{\mu}_{(1)}\p{X_i} - \hat{\mu}_{(0)}\p{X_i} + \frac{W_i - \he(X_i)}{\he(X_i)(1 - \he(X_i)}\p{Y_i - \hat{\mu}_{(W_i)}(X_i)}},
\end{equation}
where $e(x) = \mathbb{P} \p{W_i \cond X_i = x}$ is the propensity score, $\mu_{(w)}(x) = \mathbb{E} \p{Y_i(w) \cond X_i = x}$ are
conditional response surfaces, and the quantities above with hats are estimates thereof.
A celebrated property of this estimator is that it is double robust, meaning that it is consistent whenever
either \smash{$\he(x)$} or the \smash{$\hat{\mu}_{(w)}(x)$} are consistent \citep{scharfstein1999adjusting}.
Moreover, \smash{$\hat{\tau}$} is $\sqrt{n}$-consistent and semiparametrically efficient whenever the following
risk bounds hold \citep{farrell2015robust}
\begin{equation}
\label{eq:DR}
\mathbb{E}{\p{\hat{\mu}_{(W)}(X) - \mu_{(W)}(X)}^2} \mathbb{E}{\p{\he(X) - e(X)}^2} = o\p{\frac{1}{n}}.
\end{equation}
This statement is not sensitive to the structure of the estimators \smash{$\he(x)$} or the \smash{$\hat{\mu}_{(w)}(x)$}
provided we use an appropriate type of sample splitting \citep{chernozhukov2016double,zheng2011cross}, and thus
allows for considerable methodological flexibility. For example, \citet{farrell2018deep} establish conditions under
which \eqref{eq:DR} holds when \smash{$\he(x)$} or the \smash{$\hat{\mu}_{(w)}(x)$} are fit using neural networks.
These results on augmented inverse-propensity weighting can also be applied when $X_i$ is
high dimensional; however, in this case, the required risk bound can be difficult to satisfy.
In particular, except in extreme cases, the condition \eqref{eq:DR} effectively requires both
$\mu_{(w)}(x)$ and $e(x)$ to admit very sparse representations.

In this paper, we study a doubly robust construction that is specifically designed for
the high-dimensional case, and can be used for valid inference of $\tau$ under substantially
weaker sparsity assumptions than standard augmented inverse-propensity weighting.
We focus on the case where $\mu_{(w)}(x)$ and $e(x)$ have a high dimensional
linear-logistic specification (we omit intercepts for conciseness of presentation),
\begin{equation}
\label{eq:model}
\mu_{(w)}(x) = x' \beta_{(w)}, \ \ e(x) = 1/\p{1 + \exp(- x' \theta)}, \ \ \beta_{(w)}, \, \theta \in \mathbb{R}^p,
\end{equation}
and consider an estimator that is $\sqrt{n}$-consistent for $\tau$ under the condition that
either $\theta$ or the $\beta_{(w)}$ (but not necessarily both) satisfy the type of sparsity
condition that is usually required for high-dimensional inference
\citep{javanmard2014confidence,van2014asymptotically,zhang2014confidence}.
We refer to this property as sparsity double robustness.


The issue of sparsity doubly robustness has been an open question since the recent development of high-dimensional inference. This literature requires sparsity level $o(\sqrt{n}/\log p)$ for inference, a condition stronger than $o(n/\log p)$ needed for consistent estimation. Such a gap has only been addressed very recently in \cite{javanmard2015biasing}, who found that the sparsity level of only one parameter needs to  satisfy  $o(\sqrt{n}/\log p)$, not both. However, their work only addresses the linear models and heavily relies on the Gaussianity assumption of the design. In this paper, we show that such sparsity doubly robustness result holds true for nonlinear models without Gaussian designs.


Our method starts with a functional form that closely resembles \eqref{eq:AIPW}. However, we choose
our estimators of $\mu_{(w)}(x)$ and $e(x)$ in  ways that carefully exploit the geometry of sparseness
in \eqref{eq:model} and are thus able to improve on its performance. A closely related estimator has been independently studied by \cite{tan2018model}, who considered potentially misspecified models but did not provide results on sparsity doubly robustness. Our main construction is as follows,
modulo some algorithmic tweaks (including a type of sample splitting):
\begin{align}
\label{eq:cbps}
&\hat{\theta}_{(w)} = \argmin_{\theta}\cb{\frac{1}{n} \sum_{i = 1}^n \p{ \mathds{1}\{W_i \neq w\} X_i'\theta + \mathds{1} \{W_i = w\} \exp\left(- X_i' \theta}\right) + \lambda_\theta \Norm{\theta}_1} \\
\label{eq:lasso}
&\hat{\beta}_{(w)} = \argmin_{\beta} \cb{  \frac{1}{n} \sum_{W_i = w} \exp(- X_i' \hat{\theta}_{(w)})\p{Y_i - X_i' \beta}^2  + \lambda_\beta \Norm{\beta}_1} \\
\label{eq:tauhat}
&\hat{\tau} = \frac{1}{n} \sum_{i = 1}^n \Big[\p{X_i'\hat{\beta}_{(1)} + W_i  [1+\exp(-X_i'\hat{\theta}_{(1)})] (Y_i-X_i'\hat{\beta}_{(1)})  }  \\
& \qquad\qquad\qquad -  \p{X_i'\hat{\beta}_{(0)} + (1-W_i)  [1+\exp(-X_i'\hat{\theta}_{(0)})] (Y_i-X_i'\hat{\beta}_{(0)})  } \Big].\nonumber
\end{align}
As discussed in Section \ref{sec:sdr}, we can study this estimator from two different perspectives.
If $\beta_{(w)}$ is very sparse, then the solution to \eqref{eq:lasso} converges at a fast rate, while the
solution to the propensity model \eqref{eq:cbps} effectively debiases \smash{$\hat{\beta}_{(w)}$} even if
\smash{$\hat{\theta}_{(w)}$} is not particularly accurate. Meanwhile, if $\theta_{(w)}$ is very sparse, then the
converse holds. Our proof exploits this idea to establish sparsity double robustness.

The idea of  fitting a propensity model that can also leverage the shape of the conditional response surface
has generated considerable interest in recent years. The key observation here is that, in addition to being a
consistent estimator when $\theta$ is very sparse, \eqref{eq:cbps} also ``balances'' the inverse-propensity weighted features among
the treated and control samples in finite samples \citep{chan2015globally,hainmueller,imai2014covariate,tan2017regularized,zhao2016covariate}
\begin{equation}
\label{eq:balance}
\frac{1}{n} \sum_{i = 1}^n X_i \approx \frac{1}{n} \sum_{W_i = w} \frac{X_i}{1 + \exp(- X_i' \hat{\theta}_{(w)})}.
\end{equation}
The advantage of balancing is that, if the linear model for $Y$ is well specified, then balancing as in
\eqref{eq:balance} is sufficient for eliminating confounding, even when \smash{$\hat{\theta}_{(w)}$} itself
may be inconsistent or misspecified \citep{athey2016approximate,hirshberg2017balancing,kallus2018balanced,zhao2017entropy,zubizarreta2015stable}.
Note that, here, we estimate separate models for \smash{$\mathbb{P} \p{W_i = 0 \cond X_i = x}$}
and \smash{$\mathbb{P} \p{W_i = 1 \cond X_i = x}$}, parametrized by \smash{$\theta_{(0)}$} and \smash{$\theta_{(1)}$} respectively. This parametrization is based on (\ref{eq:model}) and reads
$$\mathbb{P} \p{W_i=w\cond X_i=x}=1/(1+\exp(-x'\theta_{(w)})) \qquad \text{for}\qquad w\in\{0,1 \}. $$
Notice that by (\ref{eq:model}), we have that $\theta_{(1)}=\theta $ and $ \theta_{(0)}=-\theta $.
Asymptotically, we expect both parameter vectors to be consistent, \smash{$-\hat{\theta}_{(0)}, \, \hat{\theta}_{(1)} \approx \theta$},
but finite-sample differences between \smash{$\hat{\theta}_{(0)}$} and \smash{$\hat{\theta}_{(1)}$} play a key role in enabling the
balance \citep{imai2014covariate}.



Our main finding is that an estimator  constructed via the above ``balancing'' principle  achieves sparsity double robustness, meaning that it attains
$\sqrt{n}$-consistency given strong enough sparsity assumptions $o(\sqrt{n}/\log p)$ on either $\theta$ or the $\beta_{(w)}$, but not necessarily both. As discussed further below, this property is considerably stronger than the standard double robustness property \eqref{eq:DR} in the high-dimensional setup \eqref{eq:model}.







\subsection{Related Work}

Double robust and/or semiparametrically efficient estimation has a long tradition in the literature on causal inference
\citep{chernozhukov2016double,farrell2015robust,hahn1998role,hirano2003efficient,newey2018cross,
robins1,robins1994estimation,scharfstein1999adjusting,tan2010bounded,van2006targeted}.
More recently, it has been shown that with high dimensional confounders, we can improve the behavior
of double-robust-type estimators by having them directly exploit the geometry of sparsity.

As one of the first result in this direction, \citet{athey2016approximate} showed,
given sufficient sparsity on the outcome function in \eqref{eq:model}, $\lVert \beta_{(w)} \rVert_0 \ll \sqrt{n} / \log(p)$,
we can achieve $\sqrt{n}$-consistency without any assumptions on the propensity score  beyond
overlap by simply using weights that balance moments as follows
(the \smash{$\hat{\beta}_{(w)}$} are estimated via the lasso):
\begin{equation}
\label{eq:arb}
\begin{split}
&\hat{\tau} = \frac{1}{n} \sum_{i = 1}^n X_i' \p{\hat{\beta}_{(1)} - \hat{\beta}_{(0)}} + \hat{\gamma}_i(W_i) (2W_i - 1) \p{Y_i - X_i'\hat{\beta}_{(W_i)}}, \\
& \hat{\gamma}(w) = \argmin_\gamma \cb{\frac{1}{n^2} \sum_{W_i = w} \gamma_i^2 + \Norm{\frac{1}{n} \sum_{i = 1}^n \p{1 - \gamma_i \, \mathds{1}\cb{W_i = w} } X_i}_\infty^2}.
\end{split}
\end{equation}
Conceptually, this approach is related to several papers that stress the important of covariate balance for accurate estimation
of treatment effects \citep{chan2015globally,imai2014covariate,kallus2018balanced,zhao2016covariate,zubizarreta2015stable}.
\citet{hirshberg2017balancing} establish conditions under which this estimator is efficient.

The main downside of the approximate residual balancing estimator \eqref{eq:arb} is that it always
requires sparsity of the outcome model, and cannot use a well specified and sparse propensity model
to compensate for a complex outcome model. Our sparsity double robustness result, which only requires
strong sparsity of either $\theta$ or the $\beta_{(w)}$ in \eqref{eq:model} directly addresses this limitation; and,
as shown in our experiments, yields substantial gains in accuracy when $\theta$ is in fact sparse.

Our result is most closely related to a recent proposal by \citet{chernozhukov2018double}, who studied any linear functional whose Riesz representer admits an (approximate) linear representation. In another paper, \citet{chernozhukov2016double} considers theoretical results for estimators based on
learning conditional mean function and the propensity score. In both papers, the key condition is that the product of $\ell_2$-loss for learning the two nuisance parameters is $o(n^{-1/2})$, a condition referred to as rate double robustness; see Definition 2 in \cite{Smucler2019unifying}. Sufficient conditions for rate double robustness have been provided in these works in terms of sparsity levels. For example, Remark 5.2 of \citet{chernozhukov2016double} shows that rate double robustness is guaranteed when the product of two sparsity levels is $o(n)$, while Remark 7 of \citet{chernozhukov2018double} points out that under the assumption of bounded $\ell_1$-norm of both parameters, rate double robustness holds whenever one of the sparsity levels is $o(\sqrt{n}/\log p)$.

The sparsity doubly robustness in this paper contributes to the literature by providing a different perspective. We show that efficient estimation is also possible in certain cases in which rate double robustness might not hold. One such example is  when the logistic parameter has bounded $\ell_1$-norm and has sparsity level $o(\sqrt{n}/\log p)$ and the conditional parameter has sparsity level $o(n^{3/4}/\log p)$ with potentially large $\ell_1$-norm. In this example,  we can still derive $1/\sqrt{n}$-consistency although we are not aware of any results that can guarantee rate double robustness.

In addition, our work is also different from \citet{chernozhukov2018double} in terms of specification. In the context of average treatment effect estimation, the formulation in \citet{chernozhukov2018double} means
that we need there to exist (potentially sparse) vectors $\xi_{(0)}$ and $\xi_{(1)}$  whose $\ell_1$-norms are bounded (see Definition 3 or 4 therein) as well as such that $|1/(1 - e(x))  - x' \xi_{(0)}| \approx 0$ and $| 1/e(x) -x' \xi_{(1)}|\approx 0$ uniformly across $x$. This may be a reasonable assumption
if $x$ was in fact constructed as a basis expansion of some simpler measured features; however, it appears to be
difficult to justify more generally. One contribution of this paper relative to \citet{chernozhukov2018double} is that we achieve sparsity double robustness using the natural linear-logistic specification \eqref{eq:model}.




We also note two recent papers that consider estimators that resemble ours. \citet{ning2017high} consider an
estimator that, in the spirit of \citet{belloni2014inference}, first fit a penalized covariate-balancing propensity model,
and then re-fit without penalty those coefficients that correspond to features that are relevant to outcome modeling.
Meanwhile, \citet{tan2018model}  augments a penalized covariate-balancing propensity model in an outcome regression; it turns out that his covariate-balancing mechanism designed to address the issue of misspecification is also helpful for relaxing sparsity requirements.  Neither paper, however, achieves sparsity double robustness as discussed here; rather, they require both the outcome
parameter vector $\beta$ and the propensity parameter vector $\theta$ to be ultra-sparse---or, if there is misspecification
they require the population minimizers of both the outcome and propensity loss functions to be ultra-sparse.  Under the framework of \cite{Smucler2019unifying,Rotnizky2019mixed}, \citet{ning2017high,tan2018model} are classified as examples of model double robustness, which means that one of the models (either conditional mean or propensity score) is misspecified. Rate double robustness requires that the product of the $\ell_2$-norms of the estimation errors in two models is of the order $o(n^{-1/2})$.

\section{Sparsity Double Robust Estimation}
\label{sec:sdr}


Whenever a parameter is identified through a moment condition, like \eqref{eq:AIPW}, a direct loss minimization that does not take into account this moment condition may not guarantee desirable properties.  Controlling inferential features of high-dimensional estimates is extremely difficult; most, if not all, require strict sparsity conditions.
We aim to control optimality at estimation by directly embedding the leading term of the bias into a constraint of  newly designed estimators.

The main idea behind our construction is that we use estimators \smash{$\hat{\beta}_{(0)}$}, etc.,
of $\beta_{(0)}$, etc., that have two complementary properties. When the underlying
parameter $\beta_{(0)}$ is ultra-sparse, then \smash{$\hat{\beta}_{(0)}$} converges to $\beta_{(0)}$
in $\ell_1$-norm. Furthermore, even when $\beta_{(0)}$ is not ultra-sparse, \smash{$\hat{\beta}_{(0)}$}
still has a useful covariate-balancing property implied by its Karush-Kuhn-Tucker (KKT) conditions
that can be put to good use (and similar guarantees hold for \smash{$\hat{\beta}_{(1)}$}, \smash{$\hat{\theta}_{(0)}$}
and \smash{$\hat{\theta}_{(1)}$}).
We then provide two separate consistency and asymptotic normality proofs for our estimator: One that
assumes that $\beta_{(w)}$ is ultra-sparse and relies on  KKT conditions for the \smash{$\hat{\theta}_{(w)}$} estimator to
debias a very accurate \smash{$\hat{\beta}_{(w)}$} estimator, and a second that assumes that $\theta$
is ultra-sparse and relies KKT conditions for the \smash{$\hat{\beta}_{(w)}$} estimator to
debias a very accurate \smash{$\hat{\theta}_{(w)}$} estimator. Of course, only one of these
arguments needs to hold for us to achieve asymptotic normality, and thus our estimator is sparsity double robust.
This argument was inspired by the one used by \citet{chernozhukov2018double}; however, as discussed in the
related works section, \citet{chernozhukov2018double} make the somewhat unusual assumption that $1/e(x)$ can
 be approximated by a sparse linear model (rather than the assumption we make here, i.e., a sparse logistic model
 for $e(x)$).



 \subsection{Doubly Robust Balancing via Moment Targeting} \label{sec:kkt}

In this section, we briefly sketch the argument behind our main formal result, and use it to motivate
the form of our estimator. One of the main ingredients is moment targeting: We design estimators  such that they satisfy certain moment conditions that help reduce the bias at estimation. The  construction for estimators for $\beta_{(w)}$ and $\theta_{(w)} $ is based on the structure of the bias in the final estimator for $\mathbb{E} \mu_{(w)}(X_i) $.  We emphasize that the argument here is only heuristic; formal arguments
are given in the appendix.


Given these preliminaries, observe that the treatment effect estimator under consideration can be written in a familiar form
\[
\hat \tau = \hat \mu_{(1)} - \hat \mu_{(0)},
\]
where $\hat\mu_{(w)}$ is an estimate of $\mu_{(w)} = \mathbb{E}{Y_i(w)}$. We use
\[
\hat{\mu}_{(w)}=\frac{1}{n}\sum_{i=1}^{n}  X_{i}' \hat{\beta}_{(1)}  +  \hat \gamma_ i (w) \mathds{1}\{W_i=w\} (Y_i - X_{i}' \hat{\beta}_{(w)}), \ \ \ \ \hat \gamma_i (w) =   1 + \exp(- X_i'\hat{\theta}_{(w)}),
\]
and construct \smash{$\hat{\mu}_{(0)}$} analogously. Our goal is to choose  $\hat{\theta}_{(w)}$ (and hence $\hat{\gamma}_{i}(w)$)
such as to control the errors \smash{$\hat{\mu}_{(w)} - \mu_{(w)}$} under flexible sparsity conditions.

To motivate our choice of $\hat{\theta}_{(1)}$, let us first consider the case where
$\beta_{(1)}$ is very sparse, i.e., $\|\beta_{(1)}\|_0\ll \sqrt{n}/\log p $. Notice that
\begin{align*}
\hat{\mu}_{(1)}-\mu_{(1)}  & =n^{-1}\sum_{i=1}^{n}(X_{i}'\beta_{(1)}-\mu_{(1)})+n^{-1}\sum_{i=1}^{n}W_{i}\varepsilon_{i,(1)}\hat{\gamma}_{i}(1)\\
&\qquad +n^{-1}\sum_{i=1}^{n}\left[1-W_{i}\hat{\gamma}_{i}(1)\right]X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right),
\end{align*}
where $\varepsilon_{i,(w)}=Y_i(w)-X_i'\beta_{(w)} $. The first two terms on the right hand side are asymptotically normal with mean zero under weak
consistency conditions on $\hat{\theta}_{(1)}$ that only require a moderate amount of sparsity on $\theta$.
Meanwhile, the last term can be bounded using Holder's inequality,
\begin{align*}
& \left| n^{-1}\sum_{i=1}^{n}\left[1-W_{i}\hat{\gamma}_{i}(1)\right]X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right) \right|\\
& =\left| n^{-1}\sum_{i=1}^{n}\left[1-W_{i}(1+\exp(-X_{i}'\hat{\theta}_{(1)}))\right]X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right) \right| \\
&\leq \left\Vert n^{-1}\sum_{i=1}^{n}\left[1-W_{i}(1+\exp(-X_{i}'\hat{\theta}_{(1)}))\right]X_{i} \right\Vert_{\infty} \left\Vert \hat{\beta}_{(1)}-\beta_{(1)}\right\Vert_{1}.
\end{align*}
Under sparsity assumption $\|\beta_{(1)}\|_{0}=o(\sqrt{n}/\log p)$, we can typically obtain $\|\hat{\beta}_{(1)}-\beta_{(1)}\|_{1}=o_P(1/\sqrt{\log p})$ via sparse methods \citep[e.g.,][]{negahban2012unified}.
Meanwhile, the KKT conditions for the estimator in \eqref{eq:cbps} with $w=1$ automatically yields \citep{tan2017regularized}
$$
\left\Vert n^{-1}\sum_{i=1}^{n}\left[1-W_{i}(1+\exp(-X_{i}'\hat{\theta}_{(1)}))\right]X_{i} \right\Vert_{\infty}=O_P(\sqrt{n^{-1}\log p}),
$$
thus bounding the bias to the order of $o_P(n^{-1/2})$. This is the first example of moment targeting. The KKT condition of the estimator provides a convenient moment condition for the purpose of bias reduction.

The above argument closely mirrors the argument used by \citet{athey2016approximate} to obtain
$\sqrt{n}$-consistent estimates of $\tau$ when $\|\beta_{(w)}\|_0\ll \sqrt{n}/\log p$. The main difference
  with our approach is that \citet{athey2016approximate} do not fit a model for $\theta$, but instead
directly optimize the weights \smash{$\hat{\gamma}$} via quadratic programming as in
\citet{javanmard2014confidence} and \citet{zubizarreta2015stable}. That in turn, leads to somewhat loss of flexibility whenever the outcome model is not sparse.

Here, the fact that we also model $\theta$
enables us to alternatively exploit sparsity in $\theta$ and correspondingly relax assumptions on $\beta_{(1)}$.
To do so, note that
\begin{align*}
\hat{\mu}_{(1)}-\mu_{(1)}  & =n^{-1}\sum_{i=1}^{n}(X_{i}'\beta_{(1)}-\mu_{(1)})+n^{-1}\sum_{i=1}^{n}W_{i}\varepsilon_{i,(1)}\hat{\gamma}_{i}(1)\\
&\qquad +n^{-1}\sum_{i=1}^{n}\left[1-W_{i}(1+\exp(-X_{i}'\theta))\right]X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right)\\
&\qquad +n^{-1}\sum_{i=1}^{n}W_{i}\left[\exp(-X_{i}'\theta)-\exp(-X_{i}'\hat{\theta}_{(1)})\right]X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right).
\end{align*}
Again, the sum of the first three terms above is asymptotically Gaussian on the $\sqrt{n}$-scale
under only weak assumptions on $\hat{\beta}_{(1)}$.
To handle the last term, we can use Taylor expansion to argue that (we will make this rigorous in the proof of our main result)
\begin{multline*}
\left| n^{-1}\sum_{i=1}^{n}W_{i}\exp(-X_{i}'\hat{\theta}_{(1)})X_{i}'(\hat{\theta}_{(1)}-\theta)X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right)\right|\\
\lesssim \left\Vert n^{-1}\sum_{i=1}^{n}W_{i}\exp(-X_{i}'\hat{\theta}_{(1)})X_{i}X_{i}'\left(\hat{\beta}_{(1)}-\beta_{(1)}\right)\right\Vert_{\infty} \left\Vert \hat{\theta}_{(1)}-\theta \right\Vert_{1}.
\end{multline*}
Now, given sufficient sparsity on $\theta$, i.e., $\|\theta\|_0\ll \sqrt{n}/\log p $ we can verify that
$\|\hat{\theta}_{(1)}-\theta \|_1=o_P(1/\sqrt{\log p})$. Meanwhile, the
the KKT condition for the estimator $\hat{\beta}_{(1)}$ in \eqref{eq:lasso} automatically yields that
the first component above is $O_P(\sqrt{n^{-1}\log p})$.
Thus, we also expect $\hat{\mu}_{(0)}$ to be accurate when $\theta$ is very sparse, even if $\beta$ is not. This is another example of moment targeting in that the KKT condition for $\hat{\beta}_{(1)} $ again provides a convenient bound for bounding the bias.




\subsection{Sample splitting for optimality}


The above discussion provides some helpful conceptual guidance on how to pick good estimators of the unknown $\beta_{(w)}$ and $\theta_{(w)}$.
To achieve optimality in most general terms, we invoke a special scheme of sample-splitting similar to cross-fitting. Under the usual cross-fitting scheme, the influence function is evaluated on observations that are not used to estimate the nuisance parameters (in our case $\beta_{(w)}, \theta_{(w)} $). Cross-fitting has been used to reduce bias terms in many semiparametric and high-dimensional models \citep[see, e.g.,][]{chernozhukov2016double,newey2018cross,schick1986asymptotically,zheng2011cross}. Here, our approach requires us to only cross-fit $\hat \beta_{(w)}$, but not $\hat \theta_{(w)}$.





The entire sample is divided into two parts $\mathcal{J}$ and $\mathcal{J}^c$.   For $F\in\{\mathcal{J},\mathcal{J}^c,\}$, estimators trained using the sample $F$ are denoted with $\hat \beta_{(w),F}$ and $\hat \theta_{(w),F}$, respectively. For expositional simplicity, we assume that $|\mathcal{J}|=|\mathcal{J}^c|=n/2 $.
Then, for $ (w,F)\in \{0,1\}\times \{\mathcal{J},\mathcal{J}^c \}$, we define the estimator of the mean
\begin{equation} \label{eq: mu_w_J}
\hat{\mu}_{(w),F }=\frac{1}{ |F|} \sum_{i \in F}  \p{ X_{i}' \hat{\beta}_{(w),F^c}    +  \hat \gamma_ i (w, F)  \mathds{1} \{W_i=w\} (Y_i - X_{i}' \hat{\beta}_{(w),F^c}) }
\end{equation}
where
 the weight function is defined in-sample
\[
\hat \gamma_i (w, F  ) =  1 +  \exp(- X_i'\hat{\theta}_{(w),F})  ,
\]
  $\hat \theta_{(w), F}$ is defined in Algorithm \ref{alg1}, and  $\hat \beta_{(w),F}$ is given by
\begin{equation}\label{eq: est beta}
\hat{\beta}_{(w),F} = \arg\min_{\beta} \cb{  \frac{1}{|F| } \sum_{ i \in F} \mathds{1}\{W_i = w\} \exp(-  X_i' \hat{\theta}_{(w),F})\p{Y_i - X_i'\beta}^2  + \lambda_\beta \Norm{\beta}_1}.
\end{equation}
Algorithm \ref{alg1} presents details of the propensity estimation.
The   loss functions in \eqref{eq: est beta} and \eqref{eq:loss2} were recently utilized in \citet{tan2018model} but the  proposed average treatment effects estimator therein does not achieve sparsity double robustness.





  \begin{algorithm} [H]
\caption{Optimistic penalized covariate-balancing propensity estimation}
\label{alg1}
\begin{algorithmic}
    \REQUIRE - a training sample $F\in\{\mathcal{J},\mathcal{J}^c \} $, a treatment status indicator $w \in \{0,1\}$ a tuning parameter  $\lambda_\theta \asymp \sqrt{\log(p)/n}$ and a pre-defined constant $\kappa$
    \STATE  Compute
    \begin{equation}\label{eq:loss2}
\check \theta_{(w), F}  \leftarrow  \arg\min_{\theta}\cb{\frac{1}{|F|} \sum_{i \in F}\left[ \mathds{1}\{W_i\neq w\}X_i' \theta +\mathds{1}\{W_i=w\}\exp\p{  -X_i'\theta  } \right] + \lambda_\theta \Norm{\theta}_1 }
    \end{equation}
    \IF{$ \| \check \theta_{(w), F} \|_1 > \kappa$,}
        \STATE
         \[
\hat  \theta_{(w), F}  \leftarrow  \arg\min_{\theta}\cb{ \| \theta\|_1, \mbox{ s.t. }  \left \|  \frac{1}{|F|} \sum_{i \in F} \Big[  1- \mathds{1}\{W_i=w\} \bigl(  1+\exp(-X_i'  \theta)\bigl)\Big] X_{i}  \right \|_{\infty}  \leq \lambda_\theta }
    \]
    \ELSE
        \STATE $\hat  \theta_{(w), F} \leftarrow \check \theta_{(w), F}$
          \ENDIF
           \RETURN $\hat  \theta_{(w), F} $
\end{algorithmic}
\end{algorithm}


The method presented here splits the sample into two subsamples. Although one can easily follow the same principle and split the sample  into multiple subsamples, we  do  not pursue this option here  for notational simplicity. We now define
\[
\hat \mu _{(1)} =(\hat{\mu}_{(1),\mathcal{J} }+\hat{\mu}_{(1),\mathcal{J}^c })/2.
\]
Similarly, we can define $\hat \mu _{(0)} $.
Then, the average treatment effect estimator is defined as
\begin{equation}\label{eq:tau}
\hat \tau = \hat \mu _{(1)}-\hat \mu _{(0)}.
\end{equation}


\section{Formal Results}

We now turn to a formal characterization of the average treatment effect estimator in (\ref{eq:tau}), with the aim of providing asymptotic Gaussianity whenever one, but not both, of $\theta$, $\beta$ is estimated consistently. We begin by listing some theoretical assumptions necessary for the development of the theoretical guarantees.

First, we assume that the covariate space and the parameter space are both
subsets of Euclidean space; specifically, we assume that $X \in [a,b]^p$ and $\theta \in \mathcal{B}_1 (r) \subset \mathbb{R}^p$ for some bounded $r>0$,
where $\mathcal{B}_1(r)$ is an $\ell_1$ ball with radius $r$.

 The results discussed below hold whenever,
 the tuning parameters (in Algorithm \ref{alg1}) are chosen to be proportional to $\sqrt{\log(p)/n}$.
  Moreover, we also assume that $\kappa$ is chosen to be larger than $\| \theta_{(w)}\|_1$. Our procedure is not particularly sensitive to the choice of $\kappa$;   in practice it suffices to choose  a large enough number.



 \begin{assumption}[Eigenvalue]\label{assu: RE condition}
 The minimum and maximum eigenvalues of  $\mathbb{E}[ X_iX_i']$ are contained in a bounded interval that does not contain zero.
 \end{assumption}


 Our next assumption controls the regularity properties of the errors within both models \eqref{eq:model}.
 Let $\varepsilon_{i,(w)} = Y_i - X_i'\beta_{(w)}$ and $v_{i,(w)} = \mathds{1}\{W_i=w\} - e_{(w)}(X_i)$. Note that in the context of models \eqref{eq:model} the unconfoundedness assumption implies  $\varepsilon_{i,(w)} \perp v_{i,(w)} | X_i$ and from now on we will work with this slightly weaker assumption.

  \begin{assumption}[Model]\label{assu: unconfoundedness}
  $X_i $ has a bounded sub-Gaussian norm. Moreover, for $w\in\{0,1\}$, $ \varepsilon_{i,(w)}$ is sub-Gaussian.
 \end{assumption}


 Now, observe that Assumption \ref{assu: unconfoundedness} is  very weak and in particular it is not implying consistent estimation in  the outcome  model. The boundedness of $\|X_i\|_\infty$ and $\|\theta\|_1$ guarantees the overlap condition.


 Finally, in the context of the average treatment effects, in order to provide confidence intervals an estimate of the   asymptotic variance of $\hat \tau$  is needed. We show   an asymptotic variance  of $\hat \tau$ takes the form of
 \begin{align*}
  & \mathbb{E} \left[	X_{i}'  \left(\beta_{(1)}  - \beta_{(0)} \right)    - \tau \right]^2  +\mathbb{E}\left[ W_{i}\varepsilon_{i,(1)}  \gamma_i (1) \right]^2
+\mathbb{E}\left[ (1-W_{i})\varepsilon_{i,(0)}  \gamma_i (0) \right]^{2}
\\
& \qquad : = \Omega + V_{(1)} + V_{(0)}.
\end{align*}
Observe that $\Omega$ is the variability  induced primarily from the variability of the design $X$. The other two terms  can be viewed as properly normalized unexplained variance of the models \eqref{eq:model}.

To define variance estimates, we define estimates of $\Omega$, $V_{(1)}$ and  $V_{(0)}$ separately.
We set
\begin{equation}\label{eq:Omega}
\hat \Omega = n^{-1} \sum_{i=1}^n \left( X_{i}'   (\hat \beta_{(1)}  - \hat \beta_{(0)}  ) - \hat \tau\right)^2,
\end{equation}
as well as
\begin{equation}\label{eq:V_w}
\begin{split}
&\hat V_{(w)} =n^{-1} \sum_{i=1}^n   (2W_i -1)^2 \hat \varepsilon_{i,(w)}^2  \hat \gamma_i^2 (w)  \ind\{W_i =w\}, \\
&\hat \varepsilon_{i,(w)} = Y_i - X_i' \hat \beta_{(w)}, \qquad \hat \gamma_i(w) = 1+e^{- (2W_i -1) X_i' \hat \theta_{(w)}}.
\end{split}
\end{equation}
In the above display, $\hat \beta_{(w)}$ and $\hat \theta_{(w)}$ could be the ones computed on one sample, $\mathcal{J}$ or could be different ones; for example to borrow strength across samples we consider
\[
 \hat \beta_{(w)}   = (\hat \beta_{(w),\mathcal{J}}  +\hat \beta_{(1),\mathcal{J}^c} )/2, \qquad   \hat \theta_{(w)} = (\hat \theta_{(w),\mathcal{J}} +\hat \theta_{(w),\mathcal{J}^c})/2 .
\]
Now, we define the variance estimate as
\begin{equation}\label{eq:hatV}
\hat V = \hat \Omega+\hat V_{(0)}+\hat V_{(1)}
\end{equation}
for $\hat \Omega$ and $\hat V_{(w)}$ defined in \eqref{eq:Omega} and \eqref{eq:V_w}, respectively.
We show that the construction above is appropriate for such circumstances and leads to asymptotically optimal confidence sets.


 \begin{thm}\label{thm: main}
 Let Assumptions \ref{assu: RE condition} and \ref{assu: unconfoundedness} hold. Then,
 as long as one of the following two conditions holds,
\begin{itemize}
\item[\mbox{\rm (i)}] {\mbox{\rm (Ultra-sparse outcome model)}} $\|\beta_{(w)}\|_{0}=o(\sqrt{n}/\log p)$
and $\|\theta\|_{0}=o(n/\log p)$,
\item[\mbox{\rm (ii)}] {\mbox{\rm (Ultra-sparse propensity model)}} $\|\theta\|_{0}=o(\sqrt{n}/\log p)$
and $\|\beta_{(w)}\|_{0}=O(n^{3/4}/\log p)$,
\end{itemize}
 we have the following representation of $\hat \tau$ as defined in \eqref{eq:tau}
 \[
\sqrt{n}(\hat{\tau}-\tau)=n^{-1/2}\sum_{i=1}^{n} \psi(Y_{i},W_{i},X_{i},\tau,e_{(1)}(X_{i})) +o_{P}(1),
\]
where
\begin{equation} \label{eq: influence fun}
\begin{split}
&\psi(Y_{i},W_{i},X_{i},\tau,e_{(1)}(X_{i}))= \\
&\ \ \ \ \ \ X_{i}'\p{\beta_{(1)} - \beta_{(0)}} + W_i \p{\frac{Y_{i} - X_{i}'\beta_{(1)}}{e_{(1)}(X_{i})}} - (1-W_{i})\p{\frac{Y_{i} - X_{i}'\beta_{(0)}}{1-e_{(1)}(X_{i})}}-\tau.
\end{split}
\end{equation}
In particular, under these assumptions we have
 \[
\sqrt{n}(\hat{\tau}-\tau)  \overset{\text{d}}{\rightarrow}  \mathcal{N}(0,V_*), \ \ \ \ \
V_* = \mathbb{E} \left[ \psi^2(Y_{i},W_{i},X_{i},\tau,e_{(1)}(X_{i}))\right].
 \]
Moreover,  under the same assumptions, for $\hat V$ defined in \eqref{eq:hatV} we have
\[
\hat{V}=V_*+o_{P}(1),
\]
and in turn $\sqrt{n} (\hat \tau - \tau) /\sqrt{\hat V}   \overset{\text{d}}{\rightarrow}  \mathcal{N}(0,1 )$.
 \end{thm}

By well known results \citep[e.g.,][]{hahn1998role,newey1994asymptotic,robins1}, $\psi$ in (\ref{eq: influence fun}) is the efficient influence function and $V_*$  is the semiparametric efficiency lower bound. Therefore, our estimator $\hat{\tau}$ in (\ref{eq:tau}) is a semiparametrically efficient estimator. As discussed before, the sparsity requirement in Theorem \ref{thm: main} is considerably weaker that needed by existing estimators in the high-dimensional linear-logistic model \eqref{eq:model}, including the methods discussed in \citet{athey2016approximate},
\citet{belloni2014inference}, \citet{farrell2015robust}, \citet{ning2017high} and \citet{tan2018model},
in that we only need either the outcome model or the propensity model to be ultra sparse (but not both).




\section{Numerical Experiments}

In this section we present numerical work where we contrast the behavior of the introduced method with the existing approaches.
We consider the following design setting,   $X_{i}\sim N(0,\Sigma)$ with $\Sigma_{i,j}=\rho^{|i-j|}$
and $\rho=0.6$. We set the sample size  and the number of covariates to be $n=500$ and   $p=600$, respectively. The following structure for the parameters is used.
We consider the propensity model where
$$\theta=a_{\theta}(1,0,1,0,1,0,...,0,1,0,0...,0)'$$
with $\|\theta\|_{0}=s_{\theta}$.  We vary the value of $s_{\theta}$ and  set $a_{\theta}$ such that $\sqrt{\theta'\Sigma_{X}\theta}=1$.

Similarly,  for  the outcome model we consider
$$\beta_{(1)}=a_{\beta}(1,0,1,0,1,0,...,0,1,0,0...,0)'$$
and set $\beta_{(0)}=-\beta_{(1)}$. In other words, non-zero entries
appear only on indices with odd numbers.  For $a_\beta$, we consider two cases.
In the first, we have homoskedastic Erros: The error term is generated from a centered $\chi^{2}(1)$ distribution (since it is light tailed and asymmetric). We set  $\sqrt{\beta_{(1)}'\Sigma_{X}\beta_{(1)}}=\sqrt{2}$. We use $\sqrt{2}$ to get an $R^{2}$ of	50\%  (because the error has variance 2).
The second has heteroskedastic errors: $\varepsilon_{i,(1)} $ is generated according to
	$$\p{4\times\mathbf{1}\{e(X_{i})\leq0.5\}+\mathbf{1}\{e(X_{i})>0.5\}}\xi_i$$ with $\xi_i$ being a centered $\chi^2(1)$ variable independent of $X_i$. Observe that $\varepsilon_{i,(0)}$ is still a centered  $\chi^2(1)$ variable.
 We also consider other values of $a_\beta$ such that the $R^2$ in the homoskedastic case is 10\%. We report the mean squared error (MSE) and coverage probability of 95\% confidence
interval (CP). We compare our methods with two popular alternatives:
\begin{itemize}
	\item \textbf{AIPW.} Augmented inverse propensity weighting \citep{robins1994estimation} is a popular method for estimating the average treatment effect. We implement this using the hdm package in R.
	\item \textbf{ARB.} Approximate residual balancing was proposed by \cite{athey2016approximate}. This method can handle cases in which the propensity score is hard to estimate.
\end{itemize}

\begin{table}[H]
	\caption{Mean squared error (MSE) and coverage probability (CP) across two models both with R-squared= $0.5$. Comparison includes augmented inverse propensity weighting (AIPW), approximate residual balancing (ARB), and  sparsity double robust estimation (SDR) in (\ref{eq:tau}). Parameters $s_\theta$ and $s_\beta$ denote sparsity of the propensity and outcome models, respectively.}
	\label{tab: 1}
\begin{center}
		\begin{tabular}{cr@{\extracolsep{0pt}.}lr@{\extracolsep{0pt}.}lcr@{\extracolsep{0pt}.}lr@{\extracolsep{0pt}.}l}
		& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
		& \multicolumn{9}{c}{Homoscedastic errors}\tabularnewline
		  \cline{4-9} \\
		& \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=30$}\tabularnewline
		& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
		AIPW & 0&041 & 0&954 &  & 0&063 & 0&950\tabularnewline
		ARB & 0&039 & 0&928 &  & 0&043 & 0&916\tabularnewline
		SDR & 0&042 & 0&952 &  & 0&037 & 0&960\tabularnewline
		& \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=30$}\tabularnewline
		& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
		AIPW & 0&058 & 0&932 &  & 0&097 & 0&954\tabularnewline
		ARB & 0&039 & 0&882 &  & 0&043 & 0&924\tabularnewline
		SDR & 0&035 & 0&972 &  & 0&038 & 0&968\tabularnewline
				& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
		& \multicolumn{9}{c}{Heteroscedastic errors}\tabularnewline
			  \cline{4-9} \\
		& \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=30$}\tabularnewline
		& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
		AIPW & 0&223 & 0&920 &  & 0&191 & 0&958\tabularnewline
		ARB & 0&145 & 0&922 &  & 0&129 & 0&950\tabularnewline
		SDR & 0&132 & 0&972 &  & 0&081 & 0&990\tabularnewline
		& \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=30$}\tabularnewline
		& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
		AIPW & 0&245 & 0&934 &  & 0&258 & 0&962\tabularnewline
		ARB & 0&113 & 0&948 &  & 0&102 & 0&946\tabularnewline
		SDR & 0&080 & 0&992 &  & 0&074 & 0&994\tabularnewline
		& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
	\end{tabular}
\end{center}


\end{table}


The results are reported in  Tables \ref{tab: 1} and \ref{tab: 2}. Table \ref{tab: 1} indicated that in the baseline case with extremely sparse $\beta_{(w)}$ and $\theta$, AIPW, approximate residual balancing method and SDR perform very similarly. Note that this is as expected;  all of the methods should be achieving the same asymptotic variance. However, when either the propensity score model or the conditional mean function or both are not extremely sparse, the SDR method delivers smaller MSE. This confirms our theoretical results, which state that our method is guaranteed to provide efficient estimation even if there is lack of extreme sparsity.

Perhaps a more direct way of formalizing this intuition is through analyzing the rate of the remainder similar to the discussion in \citep{newey2018cross}; one way is to express the rate of the remainder for the asymptotic expansion in Theorem \ref{thm: main} in terms of $\|\beta_{(w)}\|_0$ and $\|\theta\|_0$. One can use our technical arguments to show that, compared to AIPW, the remainder of the SDR estimator has the same order or smaller order of magnitude. In Table \ref{tab: 2}, we also report the results by setting $a_\beta$ such that in the homoscedasticity case we have an R-squared of 10\%. The pattern is quite similar; we observe comparable performance when we have extreme sparsity in both models and the proposed method has lower MSE in the absence of such sparsity in either model.


\begin{table}[H]
	\caption{Mean squared error (MSE) and coverage probability (CP) across two models both with R-squared= $0.1$. Comparison includes augmented inverse propensity weighting (AIPW), approximate residual balancing (ARB), and  sparsity double robust estimation (SDR) in (\ref{eq:tau}). Parameters $s_\theta$ and $s_\beta$ denote sparsity of the propensity and outcome models, respectively.}
	\label{tab: 2}
	\begin{center}
		\begin{tabular}{cr@{\extracolsep{-0.2pt}.}lr@{\extracolsep{-0.2pt}.}lcr@{\extracolsep{-0.20pt}.}cr@{\extracolsep{-0.20pt}.}l}
			& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
			& \multicolumn{9}{c}{Homoscedastic errors}\tabularnewline
				  \cline{4-9} \\
			& \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=30$}\\
			& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\\
			AIPW & 0&043 & 0&954 &  & 0&055 & 0&952\tabularnewline
			ARB & 0&038 & 0&915 &  & 0&039 & 0&926\tabularnewline
			SDR & 0&042 & 0&948 &  & 0&035 & 0&965\tabularnewline
			& \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=30$}
			\tabularnewline
			& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
			AIPW & 0&048 & 0&947 &  & 0&086 & 0&962\tabularnewline
			ARB & 0&039 & 0&900 &  & 0&042 & 0&920\tabularnewline
			SDR & 0&035 & 0&977 &  & 0&039 & 0&974
			\\
					& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
	& \multicolumn{9}{c}{Heteroscedastic errors}\tabularnewline
		  \cline{4-9} \\
	& \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=2$, $s_{\beta}=30$}\tabularnewline
	& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
	AIPW & 0&206 & 0&936 &  & 0&211 & 0&932\tabularnewline
	ARB & 0&116 & 0&950 &  & 0&132 & 0&946\tabularnewline
	SDR & 0&073 & 0&980 &  & 0&079 & 0&990\tabularnewline
	& \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=2$} &  & \multicolumn{4}{c}{$s_{\theta}=30$, $s_{\beta}=30$}\tabularnewline
	& \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP} &  & \multicolumn{2}{c}{MSE} & \multicolumn{2}{c}{CP}\tabularnewline
	AIPW & 0&207 & 0&934 &  & 0&204 & 0&936\tabularnewline
	ARB & 0&105 & 0&942 &  & 0&122 & 0&928\tabularnewline
	SDR & 0&057 & 0&992 &  & 0&060 & 0&996\tabularnewline
	& \multicolumn{2}{c}{} & \multicolumn{2}{c}{} &  & \multicolumn{2}{c}{} & \multicolumn{2}{c}{}\tabularnewline
\end{tabular}
\end{center}
\end{table}



\newpage