EconBase
← Back to paper

Doubly Robust Uniform Confidence Band for the Conditional Average Treatment Effect Function

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.

61,935 characters · 16 sections · 69 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.

Doubly Robust Uniform Confidence Band for the Conditional Average Treatment Effect Function

\doublespacing

abstractIn this paper, we propose a doubly robust method to estimate the heterogeneity of the average treatment effect with respect to observed covariates of interest. We consider a situation where a large number of covariates are needed for identifying the average treatment effect but the covariates of interest for analyzing heterogeneity are of much lower dimension. Our proposed estimator is doubly robust and avoids the curse of dimensionality. We propose a uniform confidence band that is easy to compute, and we illustrate its usefulness via Monte Carlo experiments and an application to the effects of smoking on birth weights. Keywords: average treatment effect conditional on covariates, uniform confidence band, double robustness, Gaussian approximation. JEL classification codes: C14, C21

Introduction

In this paper, we propose a doubly robust method to estimate the heterogeneity of the average treatment effect with respect to observed covariates of interest. To describe our methodology, we consider the potential outcome framework. Let $Y_{1}$ and $Y_{0}$ be potential individual outcomes in two states, with treatment and without treatment, respectively. For each individual, the observed outcome $ Y $ is $Y=DY_{1}+(1-D)Y_{0}$, where $D$ denotes an indicator variable for the treatment, with $D=0$ if an individual is not treated and $ D=1$ if an individual is treated. We assume that independent and identically distributed observations $\{(Y_{i},D_{i},\mathbf{Z}_{i}):i=1,\ldots ,n\}$ of $(Y,D,\mathbf{Z})$ are available, where $\mathbf{Z} \in \mathbb{R}^{p}$ denotes a $p$-dimensional vector of covariates.

Suppose that a researcher is interested in evaluating the average treatment effect conditional on only a subset of covariates $\mathbf{X}$, which is of a substantially lower dimension than $\mathbf{Z}$, where $\mathbf{Z}\equiv (\mathbf{X}^{\top },\mathbf{V}^{\top })^{\top }\in \mathbb{R}^{d}\times \mathbb{R}^{m}$, $p \equiv d+m$. That is, we are interested in a case where $d \ll p$.

The main object of interest in this paper is the conditional average treatment effect function (CATEF); namely:

align[align omitted — 96 chars of source]

When $d \geq 3$, it is difficult to plot $g(\mathbf{x})$, not to mention low precision due to the curse of dimensionality. Hence, for practical reasons, we focus on the case that $d=1$ or $d=2$, while $p$ is often of a much higher dimension.

To achieve identification of the CATEF, we assume that $Y_1$ and $Y_0$ are independent of $D$ conditional on $\mathbf{Z}$ (known as the unconfoundedness assumption):

align[align omitted — 63 chars of source]

where $\perp$ denotes the independence. For (ref) to be plausible in applications, applied researchers tend to consider a large number of covariates $\mathbf{Z}$. Note that in our setup, the treatment may be confounded in the sense that the treatment assignment may not be independent of the potential outcome variables given $\mathbf{X}$ only. To satisfy the unconfoundedness condition, a much larger set of conditioning variables $\mathbf{Z}$ needs to be employed.

Different roles of covariates between $\mathbf{X}$ and $\mathbf{V}$ are noted in the recent literature. For example, Ogburn-et-al-JRSSB consider a similar issue in the context of the local average treatment effect (LATE) of Imbens:Angrist:94. Ogburn-et-al-JRSSB emphasize that conditioning on a large number of covariates $\mathbf{Z}$ may be required to make it plausible that the binary instrument is valid. In their empirical example, Ogburn-et-al-JRSSB revisit the analyses of Poterba1995 and Abadie2003 to examine whether participation in 401(k) pension plans increases household savings. In their example, the vector of covariates $\mathbf{Z}$ for the identifying assumption consists of income, age, marital status, and family size, whereas the variable of interest $\mathbf{X}$ is income. Abrevaya-et-al-JBES also consider the case of investigating the effect of smoking during pregnancy on birth weights. They are interested in estimating (ref) with $\mathbf{X}$ being the age of mother; however, as noted in Abrevaya-et-al-JBES, it is unlikely that conditioning only on the age of mother would achieve the unconfoundedness assumption with nonexperimental data. As a result, it is necessary to consider a high-dimensional $\mathbf{Z}$, including the age of the mother.

The fact that a high-dimensional $\mathbf{Z}$ needs to be employed for (ref) to be plausible in an application makes a fully nonparametric estimation approach impractical because of the curse of dimensionality. For example, the propensity score is not nonparametrically estimable in moderately sized samples, if the dimension of $\mathbf{Z}$ is high. One obvious alternative is to use a parametric model for the propensity score; however, it may lead to misleading results if the parametric model is misspecified.

With the aim of providing a practical method and, at the same time, reducing sensitivity to model misspecification, we propose to use a doubly robust method based on parametric regression and propensity score models. Our estimator of the CATEF is doubly robust in the sense that it is consistent when at least one of the regression model and the propensity score model is correctly specified. Specifically, we first estimate CATEF$(\mathbf{Z})$ using a doubly robust procedure: we estimate a parametric regression model of the outcome on $\mathbf{Z}$ for each treatment status and a parametric model for the probability of selecting into the treatment given $\mathbf{Z}$; we then combine the parametric estimation results in a doubly robust fashion to construct an estimate of CATEF$(\mathbf{Z})$. We then obtain an estimate of CATEF$(\mathbf{X})$ by adopting the local linear smoothing of CATEF$(\mathbf{Z})$. As a result, we avoid high-dimensional smoothing with respect to $\mathbf{Z}$ but mitigate the problem of misspecification by both the doubly robust estimation and low-dimensional nonparametric smoothing with respect to $\mathbf{X}$.

We emphasize that we are willing to assume parametric specifications for the propensity score and regression models as functions of $\mathbf{Z}$ to avoid the curse of dimensionality, but not for CATEF$(\mathbf{X})$. One may consider parametric estimation of CATEF$(\mathbf{X})$, as Ogburn-et-al-JRSSB estimate their LATE parameter using least squares approximations. However, note that even if the parametric specification of CATEF$(\mathbf{Z})$ is correct, the resulting specification of CATEF$(\mathbf{X})$ may not be correctly specified since, for example, $\mathbb{E}[ \mathbf{Z}| \mathbf{X}]$ is possibly highly nonlinear. To avoid this misspecification, we estimate CATEF$(\mathbf{X})$ nonparametrically.

Because the CATEF is a functional parameter, as a tool of inference, we propose to use a uniform confidence band for the CATEF. Our construction of the uniform confidence band is based on some analytic approximation of the supremum of a Gaussian process using arguments built on Piterbarg:96, combined with a Gaussian approximation result of CCK:2014 and an empirical process result of Ghosal/Sen/vanderVaart:00. Our method is simple to implement and does not rely on resampling techniques.

This paper contributes to the literature on doubly robust estimation by demonstrating that the doubly robust procedures are useful for estimating the CATEF. In this paper, we focus on the so-called augmented inverse probability weighting estimator that was originally proposed by RobinsRotnitzkyZhao94 for the estimation of the mean RobinsRotnitzky95, ScharfsteinRotnitzkyRobins99. Their estimator appears to be the first estimator to be recognized as being doubly robust. Since then, many other alternative doubly robust estimators have been proposed in the literature. For example, the inverse probability weighting regression adjustment estimator KangSchafer07, Wooldridge07, Wooldridge10 is widely known and has been implemented in statistical software packages. See the introduction of Tan10 for a comprehensive summary of other doubly robust estimators. Doubly robust estimators have been advocated for use in many different areas of application: See, for example, LuncefordDavidian04 for medicine, GlynnQuinn10 for political science, Wooldridge10 for economics, and SchaferKang08 for psychology. There are also doubly robust estimators available for different settings including instrumental variables estimation Tan06, OkuiSmallTanRobins12 and estimation under multivalued treatments DeryaUysal:15. It would not be difficult to extend our method to allow these other doubly robust estimators and to consider different settings. However, to keep the analysis simple, in this paper, we focus on the augmented inverse probability weighting estimator of the CATEF.

The CATEF is mathematically equivalent to “$V$-adjusted variable importance” of vanderLaan06, who proposes it as a measure of variable importance in prediction. vanderLaan06 proposes a doubly robust estimator of $V$-adjusted variable importance. Contrary to ours, he considers the projection of the $V$-adjusted variable importance on a parametric working model and does not consider a nonparametric estimation. Moreover, a uniform confidence band is not examined in vanderLaan06.

In a recent paper, Abrevaya-et-al-JBES consider the estimation of the CATEF\footnote{Our paper is independent of Abrevaya-et-al-JBES and it is started without knowing their work.}; however, there are two main differences of this paper relative to Abrevaya-et-al-JBES. First, we propose the doubly robust procedure to estimate the CATEF. Abrevaya-et-al-JBES consider the inverse probability weighting estimator. The inverse probability weighting estimator suffers from model misspecification when the propensity score model is misspecified and from the curse of dimensionality when it is estimated nonparametrically. Second, we present a method to construct a uniform confidence band, whereas Abrevaya-et-al-JBES only provide a pointwise confidence interval.

The remainder of the paper is organized as follows. Section (ref) presents the doubly robust estimation method, Section (ref) gives an informal description of how to construct a two-sided, symmetric uniform confidence band when the dimension of $\mathbf{X}$ is one, and Section (ref) deals with a general case and provides formal theoretical results. In Section (ref), the results of Monte Carlo simulations demonstrate that in finite samples, our doubly robust estimator works well, and the proposed confidence band has desirable coverage properties. Section (ref) gives an empirical application, and Section (ref) concludes. The proofs are contained in Appendix (ref).

Doubly Robust Estimation of the Average Treatment Effect Conditional on Covariates of Interest

In this section, a doubly robust method for estimating the CATEF is proposed. We first estimate the CATEF for all the covariates using a doubly robust method. We then obtain the CATEF for the covariates of interest using a nonparametric approach.

Define:

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

where $\pi (\mathbf{z})$ is the propensity score and $\mu_j (\mathbf{z})$ for $j=0,1$ are called regression functions. Note that $\mu_j (\mathbf{z}) = E( Y_j | \mathbf{Z}= \mathbf{z})$ for $j=0,1$ under unconfoundedness. Let $\pi (\mathbf{z},\beta )\ $and $\mu _{j}(\mathbf{z},\alpha _{j})$ for $j=0,1$ denote parametric models of $\pi (\mathbf{z})$ and $\mu _{j}(\mathbf{z})$, respectively.\footnote{$\mu _{j}(\mathbf{z},\alpha _{j})$ may also be called “marginal structural models” of Robins00.} A doubly robust procedure requires that either $\pi (\mathbf{z})$ or $\mu _{j}(\mathbf{z})$ for $j=0,1$ should be correctly specified, thereby allowing for misspecification in $\pi (\mathbf{z})$ or in $\mu _{j}(\mathbf{z})$. Let $\theta_0 \equiv ( \alpha_{10}^\top, \alpha_{00}^\top,\beta_0^\top)^\top$ denote the vector of true or pseudo-true parameter values that optimize some criterion functions.

We consider the augmented inverse probability weighting approach. Let:

eqnarray*[eqnarray* omitted — 513 chars of source]

where $\mathbf{W} \equiv (Y,\mathbf{Z}^{\top })^{\top }\ $and $\theta \equiv (\alpha _{1}^{\top },\alpha _{0}^{\top },\beta ^{\top })^{\top }.$ The first terms in $\psi_{1}(\mathbf{W},\alpha _{1},\beta ) $ and $\psi_{0}(\mathbf{W},\alpha _{0},\beta ) $ correspond to inverse probability weighting. The second terms are augmented terms that make the procedure doubly robust.

The following lemma gives regularity conditions under which $g(\mathbf{x})$ is identified.

lem[Identification of the CATEF] Assume that (ref) holds and $0 < \pi(\mathbf{Z},\beta_0 ) < 1$ almost surely. Suppose that either $\beta_0$ satisfies $\mathbb{E}\left[D|\mathbf{Z}\right] = \pi (\mathbf{Z},\beta_0 )$ almost surely or $\alpha_{10}$ and $\alpha_{00}$ satisfy $\mathbb{E}\left[Y_1|\mathbf{Z}\right] = \mu _{1}(\mathbf{Z},\alpha _{10})$ and $\mathbb{E}\left[Y_0|\mathbf{Z}\right] = \mu _{0}(\mathbf{Z},\alpha _{00})$ almost surely. Then: \begin{align*} g(\mathbf{x}) = \mathbb{E}\left[ \psi (\mathbf{W},\theta_0 )|\mathbf{X}=\mathbf{x} \right].\end{align*}

Lemma (ref) suggests that one may estimate $g(\mathbf{x})$ by running the nonparametric regression of $\psi (\mathbf{W}, \hat \theta)$ on $\mathbf{X}_i$, where $\hat \theta$ is a consistent parametric estimator of $\theta_0$. Moreover, this lemma implies that the CATEF can be identified through $\psi (\mathbf{W}, \theta_0)$ if either the regression models ($\mu_1 (\mathbf{z}, \alpha_1)$ and $\mu_0 (\mathbf{z}, \alpha_0)$) or the propensity score model ($\pi (\mathbf{z},\beta)$) is correctly specified (or both). That is, even if $\mu_1 (\mathbf{z}, \alpha_1)$ and $\mu_0 (\mathbf{z}, \alpha_0)$ do not represent the true conditional expectation functions, provided that $\pi (\mathbf{z}, \beta)$ is correct, the CATEF is identified. Similarly, even if $\pi (\mathbf{z}, \beta)$ is misspecified, provided that $\mu_1 (\mathbf{z}, \alpha_1)$ and $\mu_0 (\mathbf{z}, \alpha_0)$ are correct, the CATEF is identified.

remIn this paper, we focus on cases in which $\mathbf{X}$ is continuous. When $\mathbf{X}$ is discrete, the CATEF can be estimated by the sample average of $\psi (\mathbf{W}, \hat \theta)$ using the sub-sample for each possible value of $\mathbf{X}$ and an estimator $\hat \theta$ of $\theta_0$. Moreover, constructing a confidence band is standard when $\mathbf{X}$ takes a finite number of values.

Parametric Estimation of $\theta$

For concreteness, we consider the following estimation procedure for $\theta_0$. However, how $\theta_0$ is estimated does not alter our results provided that the rate of convergence is sufficiently fast so that Assumption (ref)(7) given below is satisfied. For each $j = 0,1$, we estimate $\alpha_j$ by least squares:

align[align omitted — 161 chars of source]

We estimate $\beta$ by maximum likelihood (e.g., probit or logit):

align[align omitted — 174 chars of source]
remWhen the dimension of $\mathbf{Z}$ is not too high, an alternative to parametric estimation of $\psi (\mathbf{W},\theta_0 )$ is to estimate its nonparametric counterpart via local polynomial estimators as in Rothe/Firpo:16. However, this would not work when the dimension of $\mathbf{Z}$ is sufficiently high (see related remarks in Rothe/Firpo:16). The latter is the case we focus on in the paper.

Local Linear Estimation of $g$

We consider a local linear estimator of $g(\mathbf{x})$. Assume that $g(\mathbf{x})$ is twice continuously differentiable. For each $\mathbf{x}=(x_{1},\ldots ,x_{d})$, the local linear estimator of $g(\mathbf{x})$ can be obtained by minimizing:

equation*[equation* omitted — 236 chars of source]

with respect to $\gamma \equiv (\gamma_0, \gamma_1^\top) ^\top \in \mathbb{R}^{d+1}$, where $\mathbf{K}(\cdot)$ is a kernel function on $\mathbb{R}^d$ and $h_n$ is a sequence of bandwidths. More specifically, let $ \hat{g}( \mathbf{x})=\mathbf{e}_{1}^{\top }\hat{\gamma}(\mathbf{x})$, where $ \hat{\gamma}(\mathbf{x})\equiv \arg \min_{\gamma \in \mathbb{R} ^{d+1}}S_{n}(\gamma )$ and $\mathbf{e}_{1}$ is a column vector whose first entry is one, and the rest are zero.

Effect of First Stage Estimation

In our setting, we can carry out inference as if $\theta_0$ were known. This result would not be a surprise given that our first-stage estimation is parametric and our second-stage estimation is nonparametric: the rate of the convergence in the first-stage estimation is faster than that of the second stage. This feature of no first-order effect of the first-stage estimation in the second stage turns out to be more general than our setup. It is indeed closely related to doubly robustness.

If we model $g(\mathbf{x})$ parametrically or more generally approximate $g(\mathbf{x})$ by linear projection, it can be estimated by running an OLS of $\psi (\mathbf{W},\hat \theta )$ on $\mathbf{X}$. Because of the built-in feature of double robustness, it can be shown that the limiting distribution of the OLS estimator of $\psi (\mathbf{W},\hat \theta )$ on $\mathbf{X}$ is equivalent to that of the infeasible OLS estimator of $\psi (\mathbf{W}, \theta_0 )$ on $\mathbf{X}$. Furthermore, even if we estimate $\pi(\cdot)$ and $\mu_j(\cdot)$ $(j=0,1)$ nonparametrically when the dimension of $\mathbf{Z}$ is moderate, there will be no estimation effect from the first stage as well. For example, see chen2008, Rothe/Firpo:16 and CEIN2016 among others for related results.

An Informal Description of a Uniform Confidence Band

In this section, we provide an informal description of how to construct a two-sided, symmetric uniform confidence band. For simplicity, we focus on the leading case where $d=1$. Let $\mathcal{I} \equiv [a,b]$ denote an interval of interest for which we build a uniform confidence band. Assume that $\mathcal{I}$ is a subset of the support of $X$. We use nonbold $x$ to mean that $x$ is one-dimensional.

algorithm[algorithm omitted — 2,386 chars of source]

We make some remarks on the proposed algorithm. In step (1), the factor $n^{1/5} \times n^{-2/7}$ is multiplied in the definition of $h_n$ to ensure that the bias is asymptotically negligible by undersmoothing. In step (2), one can estimate $f_X$ and $\sigma^2(x) \equiv \text{Var}\left[ \psi (\mathbf{W},\theta_0 )|X=x \right]$ using the standard kernel density and regression estimators with the same kernel function $K(\cdot)$ and the same bandwidth $h$ and also with an estimator of $\theta_0$. In step (3), we may restrict the bandwidth such that $h_{n} \ll (b-a)$ (which is satisfied asymptotically), thereby imposing the condition that $\log (h_{n}^{-1}(b-a))$ is positive. The critical value proposed in step (3) is strictly positive if $\alpha$ is not too close to one or if $n$ is large enough.

remIt is straightforward to modify the algorithm above to construct one-sided symmetric confidence bands. Define a new critical value by \begin{equation*} c_{one-sided}(1-\alpha) \equiv \left( a_n^2 - 2 \log \{ \log [(1-\alpha )^{-1}] \} \right)^{1/2}. \end{equation*} Then, for each $x \in \mathcal{I}$, we set the one-sided symmetric confidence bands: \begin{align*} \hat{g}(x) - c_{one-sided}(1-\alpha) \frac{\hat{s}(x)}{\sqrt{nh_n}} \leq g(x), \end{align*} or \begin{align*} g(x) \leq \hat{g}(x) + c_{one-sided}(1-\alpha) \frac{\hat{s}(x)}{\sqrt{nh_n}}. \end{align*}
remWhen $x$ is more than one dimension, the algorithm may be revised as follows. Obviously, we need to use multivariate kernels and pointwise standard errors should be adjusted because the rate of convergence becomes $nh_n^d$. The value of $\lambda$ stays the same when we use a product kernel. For example, if $\mathbf{K}$ is the product Gaussian kernel, then $\lambda =0.5$. The formulas of $a_n$ and $c(1-\alpha)$ need to be changed. $a_n$ is the largest solution to the following equation: \begin{align*} mes(\mathcal{I}) {h_n}^{-d} \lambda^{d/2} (2\pi)^{-(d+1)/2 } a_n^{d-1} \exp(-a_n^2/2) = 1, \end{align*} where $\text{mes}(\mathcal{I})$ is the Lebesgue measure of $\mathcal{I}$. When $d =2$, the critical value has the form $c(1-\alpha) \equiv a_n + c/a_n$, where $c$ is the smallest value that satisfies \begin{align*} \exp \left( -2 e^{-c -c^2 / 2 a_n^2}\right)\left( 1+ \frac{c}{a_n^2}\right) \ge 1-\alpha. \end{align*} When $d=3$, we have that $c(1-\alpha) \equiv a_n + c/a_n$, where $c$ is the smallest value that satisfies \begin{align*} \exp \left( -2 e^{-c -c^2 / 2 a_n^2}\right)\left( \left( 1+ \frac{c}{a_n^2}\right)^2 - 2 \frac{1}{a_n^2} \right) \ge 1- \alpha. \end{align*} We note that in this paper, we assume that $d <4$ (see Assumption (ref)).
remWe may compare our proposal with the critical value based on the $(1-\alpha)$ quantile of the Gumbel distribution, which is given by: \begin{equation*} c_{\infty}(1-\alpha) \equiv a_{n} + \frac{-\log \{ \log [(1-\alpha)^{-1/2}] \}}{a_{n}}. \end{equation*} Note that: $$ c_{\infty}(1-\alpha) - c(1-\alpha) = \left[\frac{-\log \{ \log [(1-\alpha)^{-1/2}] \}}{a_{n}}\right]^2, $$ which is strictly positive for small $\alpha$ but converges to zero as $a_n$ diverges. Hence, we expect that in finite samples, the confidence band based on $c_{\infty}(1-\alpha)$ is too wide and has a higher coverage probability than the nominal level. It is shown in the next section that the critical value based on the Gumbel distribution is accurate only up to the logarithmic rate, where our proposed critical value is precise in a polynomial rate. This is because our proposal uses a higher-order expansion of Piterbarg:96, whose approximation error is of a polynomial rate. See Theorem (ref) in Section (ref) for details.
remOur construction of critical values is based on a simple analytic method that is easy to compute. Alternatively, one may rely on bootstrap methods to compute critical values for the uniform confidence band. For example, see claeskens2003 for smoothed bootstrap confidence bands and CLR for multiplier bootstrap confidence bands. chernozhukov2013 show that in general settings including high dimensional models, Gaussian multiplier bootstrap methods yield critical values for which the approximation error decreases polynomially in the sample size. Roughly speaking, both our simple analytic correction and multiplier bootstrap methods yield critical values that are accurate at polynomial rates. A refined theoretical analysis is necessary to determine which type of the critical value is better asymptotically.
remThe proposed confidence band can be used to test whether the CATEF is constant. Suppose that our null hypothesis is that $g(\mathbf{x})$ is constant in $ \mathcal{I}$. This null hypothesis can be written as $g(\mathbf{x}) = g_{\mathcal{I}}$, where $g_{\mathcal{I}} = \mathbb{E}[ g(\mathbf{x}) | \mathbf{x} \in \mathcal{I}]$. Since $g_{\mathcal{I}}$ can be estimated at the parametric ($\sqrt{n}\,$) rate and the estimator thus converges faster than $\hat g(\mathbf{x})$, we can ignore the estimation error for $g_{\mathcal{I}}$. We reject the constancy of $g(\mathbf{x})$, if the confidence band does not include the estimate of $g_{\mathcal{I}}$ for some $\mathbf{x} \in \mathcal{I}$.

Asymptotic Theory

In this section, we establish asymptotic theory. Let $U \equiv \psi (\mathbf{W},\theta_0 ) - g(\mathbf{X})$ and let $U_i \equiv \psi (\mathbf{W}_i,\theta_0 ) - g(\mathbf{X}_i)$ for $i=1,\ldots,n$. Let $\hat s^2 (\mathbf{x})$ be the estimator of the asymptotic variance of $\hat g(\mathbf{x})$. Let $s_n^2(\mathbf{x})$ denote the population version of the asymptotic variance of the estimator:

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

Assume that the $d$-dimensional kernel function is the product of $d$ univariate kernel functions.\ That is, $\mathbf{K}(\mathbf{s}) = \prod_{j=1}^d K(s_j)$, where $\mathbf{s} \equiv (s_1,\ldots,s_d)$ is a $d$-dimensional vector and $K(\cdot)$ is a kernel function on $\mathbb{R}$. Let $\rho_d (\mathbf{s}) = \prod_{j=1}^d \rho(s_j)$, where:

align[align omitted — 91 chars of source]

for each $j$. We make the following assumptions.

assumLet $d < 4$. \begin{enumerate} • $\mathcal{I} \equiv \prod_{j=1}^d [a_j,b_j]$, where $a_j < b_j, j=1,\ldots,d$, and $\mathcal{I}$ is a strict subset of the support of $\bm{X}$. • The distribution of $\mathbf{X}$ has a bounded Lebesgue density $f_{\mathbf{X}}(\cdot)$ on $\mathbb{R}^d$. Furthermore, $f_{\mathbf{X}}(\cdot)$ is bounded below from zero with continuous derivatives on $\mathcal{I}$. • The density of $U$ is bounded, $\mathbb{E}[U^2|\mathbf{X}=\mathbf{x}]$ is continuous on $\mathcal{I}$, and $\sup_{\mathbf{x} \in \mathbb{R}^d} \mathbb{E}[U^4|\mathbf{X}=\mathbf{x}] < \infty$. • $g(\cdot)$ is twice continuously differentiable on $\mathcal{I}$. • $\mathbf{K}(\mathbf{s}) = \prod_{j=1}^d K(s_j)$, where $K(\cdot)$ is a kernel function on $\mathbb{R}$ that has finite support on $[-1,1]$, $\int_{-1}^1 u K(u) du=0$, $\int_{-1}^1 K(u) du=1$, symmetric around zero, and six times differentiable. • $h_n = C n^{-\eta}$, where $C$ and $\eta$ are positive constants such that $\eta < 1/(2d)$ and $\eta > 1/(d+4)$. • $\inf_{n \geq 1} \inf_{\mathbf{x} \in \mathcal{I}} s_n(\mathbf{x}) > 0$ and $s_n(\mathbf{x})$ is continuous for each $n \geq 1$. Furthermore, $\mathbf{x} \mapsto \mathbb{E}\left[ U^2 | \mathbf{X} =\mathbf{x} \right]f_{\mathbf{X}}(\mathbf{x})$ is Lipschitz continuous. • There exists an estimator $\hat s^2 (\mathbf{x})$ such that \begin{align*} \sup_{\mathbf{x} \in \mathcal{I}} \left| \hat s^2 (\mathbf{x}) - s_n^2(\mathbf{x}) \right| = O_p( n^{-c}) \end{align*} for some constant $c > 0$. • $\max\left\{ (n h_n^d)^{1/2}|\psi (\mathbf{W}_{i},\hat \theta) - \psi (\mathbf{W}_{i},\theta_0) |: i=1,\dots,n \right\} = O_p( n^{-c})$ for some constant $c > 0$. \end{enumerate}

Most of the assumptions are standard. Condition (2) of Assumption (ref) rules out discrete covariates. If all regressors are discrete, then the estimation problem reduces to a parametric estimation problem. In this case, one may consider a multiple testing approach as in Lee:Shaikh:2014 by defining subpopulations with observed cells of discrete covariates. If some covariates are discrete and others are continuous, then one may use a smoothing approach proposed in Li:Racine:2004.

Condition (5) assumes that the kernel function has finite support. This assumption is for the simplicity of the paper and can be dropped at the expense of complicated proofs. It also assumes that the kernel function is differentiable. This assumption is crucial and excludes, for example, the uniform kernel. One of the bandwidth conditions in $h_n$ (that is, $\eta > 1/(d+4)$ in condition (6)) imposes undersmoothing, so that we can ignore the bias asymptotically. The rule-of-thumb bandwidth proposed in Section (ref) satisfies the required rate conditions.

remNote that $d < 4$ is necessary to ensure that $\eta < 1/(2d)$ and $\eta > 1/(d+4)$ can hold jointly. It is possible to extend our asymptotic theory to the case that $d \geq 4$ using a higher-order local polynomial estimator under stronger smoothness conditions. In this paper, we limit our attention to the local linear estimator since we are mainly interested in low dimensional $\bm{x}$.
remAn estimator of $\hat s^2 (\mathbf{x})$ is readily available. For example, we may consider \begin{align} \hat s^2(\mathbf{x}) = \frac{\hat{\sigma}^2(\mathbf{x})}{ \hat f_\mathbf{X}(\mathbf{x})} \int \mathbf{K}^2 \left( \mathbf{u} \right) d \mathbf{u}, \end{align} where $\hat{f}_\mathbf{X}(\cdot)$ is the kernel density estimator and $\hat{\sigma}^2(\mathbf{x})$ is a nonparametric estimator of ${\sigma}^2(\mathbf{x})$ using $\{ (\hat{U}_i^2, \mathbf{X}_i): i=1,\ldots,n \}$ with $\hat{U}_i \equiv \psi (\mathbf{W}_i,\hat{\theta} ) - \hat{g}(\mathbf{X}_i)$. Recall (ref) for its concrete form for the one-dimensional case. Alternatively, we may set \begin{align*} \hat s^2(\mathbf{x}) = \frac{1}{nh_n^{d}} \sum_{i=1}^n \frac{\hat{U}_i^2}{\hat{f}^2_\mathbf{X}(\mathbf{x})} \mathbf{K}^2 \left( \frac{\mathbf{X} - \mathbf{x}}{h_n} \right). \end{align*} For either estimator, it is straightforward to verify condition (8) of Assumption (ref) using the standard results in kernel estimation.
remNote that Condition (9) of Assumption (ref) is merely a sufficient (but not necessary) condition. This condition is satisfied, for example, if $\| \hat \theta - \theta_0 \| = O_p( n^{-1/2})$, functions $\beta \mapsto \pi (\mathbf{Z},\beta)$ and $\alpha_j \mapsto \mu _{j}(\mathbf{Z},\alpha _{j}), j =0,1,$ are Lipschitz continuous, $\pi (\mathbf{Z},\beta_0 )$ is bounded between $\epsilon$ and $1-\epsilon$ for some constant $\epsilon > 0$, provided that some weak moment conditions on $(Y, \mathbf{Z})$ hold.

Let $a_n \equiv a_n(\mathcal{I})$ be the largest solution to the following equation:

align[align omitted — 122 chars of source]

where $\text{mes}(\mathcal{I})$ is the Lebesgue measure of $\mathcal{I}$; that is, $\text{mes}(\mathcal{I}) = \prod_{j=1}^d (b_j-a_j)$ and:

equation[equation omitted — 107 chars of source]

The following is the main theoretical result of our paper.

thmLet Assumption (ref) hold. Then there exists $\kappa > 0$ such that, uniformly in $t$, on any finite interval: \begin{align} \begin{split} &\mathbb{P}\left( a_n \left[ \max_{\mathbf{x} \in \mathcal{I}} \left| \frac{\hat{g}(\mathbf{x}) - g(\mathbf{x})}{\hat{s}(\mathbf{x})} \right| - a_n \right] < t \right) = \\ &\exp \left( - 2e^{-t - t^2/2a_n^2} \right) \sum_{m=0}^{\lfloor (d-1)/2 \rfloor} h_{m,d-1} a_n^{-2m} \left( 1 + \frac{t}{a_n^2} \right)^{d-2m-1} + O( n^{-\kappa}), \end{split} \end{align} as $n \rightarrow \infty$, where $h_{m,d-1} \equiv \frac{(-1)^m (d-1)!}{ m! 2^m (d-2m-1)!}$ and $\lfloor \cdot \rfloor$ is the integer part of a number.

Notice that the approximation error is of a polynomial rate. As a result, a critical value based on the leading term of the right-hand side of (ref) provides a better approximation than one based on the Gumbel approximation. The result in Theorem (ref) may be of independent interest for constructing the uniform confidence band in nonparametric regression beyond the scope of estimating the CATEF in our context.

remIn a setting different from here, Lee/Linton/Whang:06 propose analytic critical values based on Piterbarg:96 in order to test for stochastic monotonicity, compare its performance with the bootstrap critical values in their Monte Carlo experiments, and find that both perform well in finite samples. However, the discussions in Lee/Linton/Whang:06 are informal and rely on the results of Monte Carlo experiments without the formal proof of establishing the polynomial approximation error.

The conclusion of the theorem can be simplified for special cases. In particular, if $d=1$, then:

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

where $a_n$ is the largest solution to $\text{mes}(\mathcal{I}) {h_n}^{-1} \lambda^{1/2} (2\pi)^{-1 }\exp(-a_n^2/2) = 1$. Also, if $d=2$, then:

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

where $a_n$ is the largest solution to $\text{mes}(\mathcal{I}) {h_n}^{-2} \lambda^{2} (2\pi)^{-3/2 } a_n \exp(-a_n^2/2) = 1$.

remIt is standard to obtain pointwise confidence intervals based on normal approximations. Recall that our two-sided symmetric uniform confidence band has the form: \begin{align} \hat{g}(\mathbf{x}) - c(1-\alpha) \frac{\hat{s}(\mathbf{x})}{\sqrt{nh_n^d}} \leq g(x) \leq \hat{g}(\mathbf{x}) + c(1-\alpha) \frac{\hat{s}(\mathbf{x})}{\sqrt{nh_n^d}}, \end{align} where $c(1-\alpha)$ is obtained from Theorem (ref). To obtain two-sided symmetric pointwise confidence intervals, we just need to replace $c(1-\alpha)$ with the usual normal critical value $\Phi^{-1}(1-\alpha/2)$, where $\Phi(\cdot)$ is the standard normal cumulative distribution function. The pointwise confidence interval given in (ref) is different from the one resulting from Abrevaya-et-al-JBES in terms of the formula for $\hat{s}(\mathbf{x})$ in (ref). In their case, they need to estimate $\hat{\sigma}^2(\mathbf{x})$ using $\{ (\tilde{U}_i^2, \mathbf{X}_i): i=1,\ldots,n \}$ with \begin{align*} \tilde{U}_i \equiv \frac{D_i Y_i}{\pi (\mathbf{Z}_i,\hat{\beta} )}- \frac{\left( 1-D_i\right) Y_i}{1-\pi (\mathbf{Z}_i,\hat{\beta})} \hat{g}(\mathbf{X}_i). \end{align*}
remA one-sided version of the uniform confidence band is readily available. Combining Theorems 14.1 and 14.2 of Piterbarg:96 with arguments identical to those used in the proof of Theorem (ref) yields the following proposition. Under Assumption (ref), there exists $\kappa > 0$ such that, uniformly in $t$, on any finite interval: \begin{align} \begin{split} &\mathbb{P}\left( a_n \left[ \max_{\mathbf{x} \in \mathcal{I}} \frac{\hat{g}(\mathbf{x}) - g(\mathbf{x})}{\hat{s}(\mathbf{x})} - a_n \right] < t \right) = \\ &\exp \left( - e^{-t - t^2/2a_n^2} \right) \sum_{m=0}^{\lfloor (d-1)/2 \rfloor} h_{m,d-1} a_n^{-2m} \left( 1 + \frac{t}{a_n^2} \right)^{d-2m-1} + O( n^{-\kappa}), \end{split} \end{align} as $n \rightarrow \infty$. Note that the only differences between (ref) and (ref) are that (i) there is no absolute value on the left side of the equation in (ref) and (ii) there is no factor 2 inside the exponential function in (ref). Hence, for example, if $d=1$, then: \begin{align*} \mathbb{P}\left( a_n \left[ \max_{\mathbf{x} \in \mathcal{I}} \frac{\hat{g}(\mathbf{x}) - g(\mathbf{x})}{\hat{s}(\mathbf{x})} - a_n \right] < t \right) = \exp \left( - e^{-t - t^2/2a_n^2} \right) + O( n^{-\kappa}). \end{align*}

Construction of critical values

We use the leading term on the right-hand side of (ref) as a distribution-like function to construct a uniform confidence band. For example, if $d=1$, we may construct a critical value $c(1-\alpha)$ that satisfies:

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

where $F_{n,1}(t) \equiv \exp \left( - 2e^{-t - t^2/2a_n^2} \right)$. This yields the critical value presented in the Algorithm of Section (ref). Similarly, if $d=2$, we can use:

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

where $F_{n,2}(t) \equiv \exp \left( - 2e^{-t - t^2/2a_n^2} \right)\left( 1 + \frac{t}{a_n^2} \right)$. In finite samples, it might be useful to impose monotonicity of $F_{n,j}(\cdot)$ by rearrangement (see, e.g., CF-VG:09).

remTheorem (ref) implies that: \begin{align*} &\lim_{n \rightarrow \infty} \mathbb{P}\left( a_n \left[ \max_{\mathbf{x} \in \mathcal{I}} \left| \frac{\hat{g}(\mathbf{x}) - g(\mathbf{x})}{\hat{s}(\mathbf{x})} \right| - a_n \right] < t \right) = \exp \left( - 2e^{-t } \right). \end{align*} Thus, one may construct analytical critical values based on the Gumbel distribution. However, this approximation is accurate only up to the logarithmic rate in view of Theorem (ref).

Monte Carlo Experiments

In this section, we present the results of Monte Carlo experiments. These experiments are conducted to see the finite sample performances of the proposed doubly robust estimator and the proposed uniform confidence band. The simulations are conducted by R 3.3.1 with Windows 10. The number of replications is $5000$.

Data generating process

The data generating process follows the potential outcome framework. The notations for the variables are the same as those used in the theoretical part of the paper. We consider cases with $p=10, 30$ and $N= 500, 2000$.

The data generating process is the following. The vector of covariates $\mathbf{Z}= (X_1, \dots, X_p)^{\top} $ is generated by: $\mathbf{Z} \sim N(0,I_p)$, where $I_p$ is the $p$-dimensional identity matrix. The potential outcomes are generated by:

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

where $v \sim N(0, 1)$ and $v$ is independent of $\mathbf{Z}$. The treatment status $D$ is generated by:

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

where $U \sim U[0,1]$, $U$ is independent of $(\mathbf{Z}^{\top}, v)$ and $\Lambda$ is the logistic function. Thus, the propensity score is $\pi(\mathbf{Z}) = \Lambda \left(\sum_{k= p/2}^p X_k/ \sqrt{p/2}\right)$. The observed outcome is $Y=DY_1$.

The parameter of interest is the CATEF for $\mathbf{X} =X_1$. In our specification, the CATEF can be written as:

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

We examine the performance of various statistical procedures regarding this CATEF.

Model specification

To estimate and conduct statistical inferences on $CATEF(x_1)$ using our doubly robust procedure, we need to specify a model for the regression $\mu_j (\mathbf{z})$ for $j=0,1$ and a model for the propensity score $\pi (\mathbf{z})$. We consider two regression models and two propensity score models. One of two models is correctly specified, but the other model is misspecified. We note that our doubly robust procedure is predicted to work well provided that at least one of the regression model and the propensity score model is correctly specified.

We first discuss the model specifications for the regression part. The first regression model is:

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

This model is correctly specified. The coefficients are estimated by OLS using $(1, X_1, \dots, X_p)$ as the explanatory variable. The second regression model, which is misspecified, is:

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

The model is estimated by OLS using $(1, X_1, \dots, X_{p/2})$ as explanatory variables. This model is misspecified because it suffers from sample selection bias introduced by omitting the second half of the regressors which affects the treatment status.

We also consider two models for propensity score. The model for propensity score is:

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

The misspecified model is:

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

Similarly to the case of the regression part, misspecification is introduced by omitting the second half of the regressors. The models for propensity score are estimated by maximum likelihood.

We estimate $CATEF(x_1)$ for $x_1 \in \{ -1, -0.5, 0, 0.5, 1\}$ and compute the mean bias (“MEAN”), standard deviation (“SD”), the average of standard error for $\widehat{CATEF} (x_1)$ (“ASE”), and the root mean squared error (“RMSE”). The local linear regression is conducted with the Gaussian kernel, and the preliminary bandwidth ($\hat h$ in Algorithm (1)) is chosen by the method of Ruppert:Sheather:Wand:1995. We also compute the “BIAS”, “SE” and “RMSE” of the corresponding inverse probability weighting estimators and the regression adjustment estimators. Note that the difference between the proposed method and those alternative methods arises only in the estimation of $\psi (\mathbf{W}, \theta_0)$ and the other steps are the same.

We examine the coverage probability of the uniform confidence band for $CATEF (x_1)$ for the range $-1 \le x_1 \le 1$. The nominal coverage probabilities that we consider are 99%, 95% and 90%. We compute the empirical coverage (“CP”), the mean critical value (“Mcri”), and the standard deviation of critical value (“Sdcri”). We also compute the coverage probabilities of the confidence band based on the critical values computed by the Gumbel approximation (“GCP”).

Results

Tables (ref) and (ref) summarize the results on the properties of the estimators. In both tables, DR refers to our doubly robust method, whereas IPW and RA correspond to the inverse probability weighting and regression adjustment methods, respectively. The proposed doubly robust estimator of the CATEF works well in finite samples. As the theory indicates, the proposed estimator exhibits small bias provided that at least one of the regression model and the propensity score model is correctly specified. We find that the regression adjustment estimator is very precise when the regression model is correctly specified. However, it suffers from substantial bias when the regression model is misspecified. The inverse probability weighting estimator also suffers from model misspecification. Moreover, its standard deviation is much larger than those of the doubly robust and regression adjustment estimators. When both models are misspecified, all three estimators suffer from heavy bias. The inverse probability weighting estimator has the largest RMSE because its distribution is more diverse than those of the other two estimators. All the estimators have larger standard deviations when $x=1$ and $x=-1$ compared to those in other points. This is because the number of observations around $x=1$ or $x=-1$ is expected to be smaller than that around, for example, $x=0$ which is the center of the distribution. On the other hand, the magnitude of the bias does not vary much across data points. The standard error for the proposed doubly robust estimator is slightly smaller than the standard deviation, but the difference is not large.

Tables (ref) and (ref) summarize the finite sample properties of the proposed uniform confidence band. The results show that our uniform confidence band has a reasonably good coverage property provided that one of the models is correctly specified. When both models are misspecified, the size distortion is heavy. We observe that the size distortion is heavier when the regression model is misspecified than that in the case of propensity score misspecification. This result indicates that we should carefully model the regression part in order to obtain reliable confidence bands. The average values of the 95% critical values are around 2.75. Because the pointwise critical value is 1.96 and is much smaller than the uniformly valid critical value, it demonstrates the importance of the uniform property of confidence band. The standard deviations of the critical values are small because they change only if the bandwidth changes. The confidence band based on the Gumbel approximation is very conservative.

The results of the Monte Carlo simulation confirm that the proposed doubly robust estimator indeed works well in finite samples provided that one of the regression and propensity score models is correctly specified. The proposed uniform confidence band also has good coverage properties.

An Empirical Application

We apply our uniformly valid confidence band for the CATEF for the effect of maternal smoking on birth weight where the argument of the CATEF is the mother's age. Our aim here is to illustrate our confidence band in comparison with alternative confidence bands. We first discuss the background of this application and the datasets used. We use two different data sets: the dataset from Pennsylvania and that from North Carolina. We then compute various confidence bands for the CATEF and discuss the results.

While the purpose of this application is to illustrate our uniformly valid confidence band and not to present new insights on the effect of smoking, it is still informative to discuss the background of this application. Many studies document that low birth weight is associated with prolonged negative effects on health and educational or labor market outcomes throughout life, although there has been a debate over its magnitude. See, e.g., AlmondCurrie11 for a review. Maternal smoking is considered to be the most important preventable negative cause of low birth weight Kramer87. There are many studies that evaluate the effect of maternal smoking on low birth weight AlmondCurrie11. The program evaluation approach is employed by, for example, AlmondChayLee05, daVeigaWilder08 and WalkerTekinWallace09, and panel data analysis is carried out by Abrevaya06 and AbrevayaDahl08. Here, we are interested in how the effect of smoking changes across different age groups of mothers. WalkerTekinWallace09 examine whether the effect of smoking is larger for teen mothers than for adult mothers and find mixed evidence. Abrevaya-et-al-JBES also consider this problem in their application.

Pennsylvania data

The first dataset consists of observations from white mothers in Pennsylvania in the USA. The dataset is an excerpt from Cattaneo10 and is obtained from the STATA website (“\url{http://www.stata-press.com/data/r13/cattaneo2.dta}”). Note that the dataset was originally used in AlmondChayLee05. We restrict our sample to white and non-Hispanic mothers, and the sample size is 3754. The outcome of interest ($Y$) is infant birth weight measured in grams. The treatment variable ($D$) is a binary variable that is equal to 1 if the mother smokes and 0 otherwise. The set of covariates $\mathbf{Z}$ includes the mother's age, an indicator variable for alcohol consumption during pregnancy, an indicator for the first baby, the mother's educational attainment, an indicator for the first prenatal visit in the first trimester, the number of prenatal care visits, and an indicator for whether there was a previous birth where the newborn died. We are interested in how the effect of smoking varies across different values of the mother's age. Therefore, $\mathbf{X}$ is mother's age in this application.

To estimate the CATEF, we use linear regression models for the regression part and a logit model for propensity score. The explanatory variables used in the regression models and the logit model consist of all the elements of $\mathbf{Z}$, the square of the mother's age, and the interaction terms between the mother's age and all other elements of $\mathbf{Z}$. We estimate the CATEF in the interval between ages 15 and 35.

We compute the following three 95% confidence bands for the CATEF. “Our CB” is the confidence band proposed in this paper. Because $\mathbf{X}$ is univariate in this application, we follow the algorithm in Section 3. We use the Gaussian kernel. The preliminary bandwidth ($\hat h$) is chosen by the method of Ruppert:Sheather:Wand:1995. “Gumbel CB” is the confidence band in which $c(1-\alpha)$ in the algorithm is replaced by that based on the Gumbel approximation (see Remark (ref)). “PW CB” is a pointwise valid confidence band where we replace $c(1-\alpha)$ in the algorithm by the corresponding value from the standard normal distribution (i.e., $1.96$). This provides a valid confidence interval for each point of the CATEF. However, its uniform coverage rate would be smaller than 95%.

Figure 1 plots the estimated CATEF and the three 95% confidence bands for the range between 15 and 35 years of age. The figure also contains the average treatment effect estimate (AIPW estimate) for a reference.

The widths of the three confidence bands are substantially different. The confidence band based on the Gumbel approximation provides the widest band and may not be very informative. The confidence band that is valid only in a pointwise sense gives the narrowest band. This band is not uniformly valid and so may provide misleading information about the CATEF. On the other hand, this provides valuable information if we are interested at a particular point of the CATEF. The confidence band we propose lay between “Gumbel CB” and “PW CB”. While this band is wider than “PW CB”, it is much narrower than “Gumbel CB” and is uniformly asymptotically valid. We see from this figure that our confidence band is informative while being uniformly valid.

The estimated CATEF is decreasing from 15 to around 25 years of age. It is rather stable for the range above 25 years of age. All confidence bands indicate that the CATEF is estimated imprecisely near the ends of the range. Nonetheless, the estimated CATEF indicates that smoking may not have a strong impact when the mother is young. The CATEF is estimated relatively precisely in the middle of the range. For the range between 20 and 30 years of age, even the band based on the Gumbel approximation, which is the widest, does not contain 0. This result provides robust evidence that smoking has a negative impact on birth weight at least for mothers who are 20 to 30 years old. In this particular dataset, the statistical evidence against a constant smoking effect is somewhat weak. The confidence band that is valid only in a pointwise sense may provide an impression that the smoking effect depends on the mother's age. However, the uniformly valid confidence band that we propose marginally contains the straight line that is equal to the ATE estimate. This result illustrates that there is a caveat when we use pointwise confidence intervals, as well as the importance of using uniformly valid confidence bands.

North Carolina data

The second dataset is based on the records between 1988 and 2002 by the North Carolina State Center Health Services. This dataset is used in Abrevaya-et-al-JBES and obtained from Robert Lieli's website (“\url{http://www.personal.ceu.hu/staff/Robert_Lieli/cate-birthdata.zip}). We restrict our sample to white and first--time mothers, and the sample size is 433,558. As in the case of the Pennsylvania data, the outcome is infant birth weight measured in grams and the treatment variable is an indicator for smoking status. The set of covariates $\mathbf{Z}$ includes those used in the analysis of the Pennsylvania data, except an indicator for the first baby because we focus on first--time mothers, and in addition, it includes indicators for gestational diabetes, hypertension, amniocentesis and ultra sound exams. Again, $\mathbf{X}$ is mother's age in this application. The specification for the estimation of the CATEF is the same as before.

The purpose of using this much larger dataset is to examine the effect of the sample size. Our method involves nonparametric kernel regression and it might require a large sample size to yield a reliable result. For example, the result from the Pennsylvania data indicates that the effect of smoking is very small for very young mothers. One might argue that such a result is an artifact of small sample size. The other issue is that the confidence bands obtained using the Pennsylvania data are somewhat wide. We hope that using this larger dataset provides us with narrower confidence bands and more informative statistical results.

Figure 2 plots the estimated CATEF and the three 95% confidence bands for the range between 15 and 35 years of age. Note that the scale of the vertical axis is different from Figure 1. We now obtain much narrower confidence bands. The widths of the three (uniform, point-wise and Gumbel) confidence bands are still different. The estimated CATEF for young mothers is negative and statistically different from 0. The previous result that it is close to 0 may be considered as an artifact of small sample size. The estimated CATEF is decreasing from around 17 to around 29 years of age. For the range above 30 years of age, we obtain relatively wide confidence bands. We reject the null hypothesis of no effect of smoking on birth weights uniformly over 15-35 years of age. These confidence bands do not support the hypothesis that the CATEF is constant because the ATE line exceeds the confidence bands.

One might argue that the difference in the results may stem from the fact that the North Carolina data contains richer information and we use a larger set of covariates. We reexamine the North Carolina data based on the same set of covariates as that for the Pennsylvania data, except an indicator for the first baby. Figure 3 plots the estimated CATEF and confidence bands obtained using this set of covariates. The results in Figure 3 are qualitatively very similar to those in Figure 2. We thus believe that the difference between the results from the Pennsylvania data and the North Carolina data are not from the difference in the covariates but from the difference in the sizes of these two samples.

We thus interpret our findings to indicate that the different results come from the difference in sample size yet our confidence bands reasonably quantify the uncertainty from small sample size. While two data-sets yield different estimates of CATEF, the confidence bands from the Pennsylvania data include the estimated CATEF and the confidence bands from the North Carolina data.

While we use the same data set as that used in Abrevaya-et-al-JBES, it is somewhat difficult to compare their results with ours because of differences in the implementations. In particular, the bandwidths are very different. Our choice of bandwidth is around 0.2, while theirs are between 1.4--11.2. Nonetheless, we make several remarks. Using small bandwidths (1.4 and 2.8), Abrevaya-et-al-JBES observe almost no effect for young mothers and a large negative effect for 25--30 years old mothers. We do not observe such a large difference in the effect across different age groups. Our confidence band is as tight as their confidence band obtained with bandwidth 11.2 even though we use a much smaller bandwidth and our confidence band is uniform. This is possibly because we use an AIPW method which yields a more efficient estimate than an IPW method does.

Conclusion

In this paper, we propose a doubly robust method for estimating the CATEF. We consider the situation where a high-dimensional vector of covariates is needed for identifying the average treatment effect but the covariates of interest are of much lower dimension. Our proposed estimator is doubly robust and does not suffer from the curse of dimensionality. We propose a uniform confidence band that is easy to compute, and we illustrate its usefulness via Monte Carlo experiments and an application to the effects of smoking on birth weights.

There are a few topics to be explored in the future. First, it would be useful to consider the issue of asymptotic biases of the proposed estimator without relying on undersmoothing. For example, it might be possible to extend the approach of hall2013 that avoids undersmoothing for our purposes. Second, it would be an interesting exercise to develop a method for estimating the quantile treatment effects conditional on covariates. Third, it is possible to extend our approach to the local average treatment effect. As mentioned in the Introduction, Ogburn-et-al-JRSSB consider conditioning on $\mathbf{Z}$ to achieve identification, but they estimate the local average treatment effect, say LATE($\mathbf{X}$), as a function of $\mathbf{X}$. However, their specification of LATE($\mathbf{X}$) is parametric. Our approach can be adapted to specify LATE($\mathbf{X}$) nonparametrically and to develop a corresponding uniform confidence band. Fourth, this paper does not cover marginal treatment effects that can be identified using the method of local instrumental variables developed by Heckman:Vytlacil:99,Heckman:Vytlacil:05. It would be interesting to develop a uniform confidence band for the marginal treatment effects.