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.
10,140,264 characters · 36 sections · 0 citation commands
Regression Discontinuity Design with Potentially Many Covariates
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.
In causal or treatment effect analysis, discontinuities in regression functions induced by an assignment variable can provide useful information to identify certain causal effects. The regression discontinuity design (RDD) has been widely applied in observational studies to identify the average treatment effect at the discontinuity point. For the RDD, the causal parameters of interest are identified by some contrasts of the left and right limits of the conditional mean functions. See e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), an edited volume by Cattaneo and Escanciano (2017), and references therein.
In the growing literature on the RDD analysis, this paper focuses on the RDDs where covariates are included in the estimation, which are extensively studied by Calonico, Cattaneo, Farrell and Titiunik (2019) (hereafter, CCFT). See also Fr�lich and Huber (2019) for an alternative estimation method based on kernel smoothing after localization around the cutoff. In practice, researchers often augment the regression models for the RDD analysis with various additional predetermined covariates such as demographic or socioeconomic characteristics for data units. For several RDD estimators using covariates based on local polynomial regression methods, CCFT investigated the MSE expansion, asymptotic efficiency, and data-driven bandwidth selection methods. Furthermore, CCFT developed asymptotic distributional approximations for those estimators and proposed valid inference procedures by constructing bias and variance estimators with covariate adjustment. These results may be considered as extensions of the analyses in Calonico, Cattaneo and Titiunik (2014) (hereafter, CCT) combined with robust bias correction methods in Calonico, Cattaneo and Farrell (2018, 2020) to incorporate covariates in the RDD analysis. See also Calonico, Cattaneo, Farrell and Titiunik (2017) for a statistical package on these methods.
In randomized controlled trials, regression adjustment using covariates is a common practice since it is always helpful to improve asymptotic efficiency of the causal effect estimator as far as a full set of treatment-covariate interactions is included (Lin, 2013). Also a recent paper by Lei and Ding (2021) proposed a bias correction method for the regression adjustment estimator with a diverging number of covariates. On the other hand, in the RDD analysis, which is a quasi-experiment setup, the efficiency gain by introducing covariates is not necessarily guaranteed, and CCFT provided a concrete guideline by clarifying the conditions to achieve consistency and efficiency gain for the covariate adjusted RDD estimator. Typically the efficiency improves when the projection coefficients of the covariates on the outcome are equal for both control and treatment groups. Since practitioners also commonly incorporate covariates for the RDD analysis, CCFT's guideline has a large impact in applied research. When we use covariates, it is common to employ their transformations and interactions, and the number of these terms can be pretty large. This paper adds a further guideline for practitioners who also face a large number of covariates. To begin with, the (weighted) OLS estimation in CCFT is not applicable when the number of covariates is larger than the sample size. Also, in the above scenario for efficiency improvement, it is beneficial to augment CCFT's procedure with covariate selection by high-dimensional statistical methods particularly when the regression coefficients for the conditional mean function satisfy certain sparsity.
For point estimation on the causal effect parameter identified by the RDD, we consider the Lasso estimator and its post-selection estimator based on the local linear regression (i.e., eq. (2) of CCFT). The combination of localization using kernel weights and $\ell_{1}$-penalization to deal with high-dimensional covariates is particularly relevant for the RDD analysis, where the effective sample size would be typically small due to the localization so that the effect of dimensionality of covariates becomes severer. Theoretically, we derive the $\ell_{1}$-risk properties of our local Lasso estimator and its post-selection version. Practically, based on our simulation study, we recommend the CCFT estimator after selecting covariates by the $\ell_{1}$-penalization even for a relatively small number of covariates, which exhibits desirable MSE properties and stability across different setups.
For inference, we propose to select covariates with the local Lasso. We show that the inference based on the selected covariates can be implemented in the same manner as in CCFT. We also show that when the effect of the additional covariates on the potential outcomes with or without treatment is invariant, our approach can lead to improved efficiency at the cost of an additional sparsity condition. This sparsity condition is trivially satisfied when the set of active covariates is unknown but fixed. Our simulation results demonstrate that our post-selection confidence interval exhibits robust performances in terms of both coverages and lengths, even for a relatively small number of covariates.
This paper also contributes to the large literature on high-dimensional methods in econometrics and statistics (see, e.g., B�hlmann and van de Geer, 2011, and Belloni et al., 2018, for an overview) by combining the kernel localization with $\ell_{1}$-penalization to handle high-dimensional covariates. Our inference problem can be formulated as the one for low-dimensional parameters in high-dimensional models. In statistics literature, many papers investigated this issue, such as Belloni, Chernozhukov and Hansen (2014), van de Geer, et al. (2014), and Zhang and Zhang (2014). However, these approaches are not directly applicable to the RDD context because the current problem concerns the inference on a jump in a nonparametric regression model.\footnote{A recent paper by Krei\ss\ and Rothe (2023) investigates a similar estimator to ours, and discusses an inference method based on the approach by Armstrong and Koles�r (2018).}
This paper is organized as follows. Section (ref) introduces our basic setup and local Lasso estimator, and presents the $\ell_{1}$-risk properties. In Section (ref), we discuss the validity of CCFT's inference after selecting covariates by our Lasso procedure. Section (ref) provides discussions on some extensions. A step-by-step procedure for implementation of our method is described in Section (ref). To illustrate the proposed method, Section (ref) conducts a simulation study, and Section (ref) presents an empirical example based on the Head Start data.
In this subsection, we present our basic setup and introduce the local Lasso estimator for the RDD with possibly high-dimensional covariates. For each unit $i=1,\ldots,n$, we observe an indicator variable $T_{i}$ for a treatment ($T_{i}=1$ if treated and $T_{i}=0$ otherwise), and outcome $Y_{i}=Y_{i}(0)\cdot(1-T_{i})+Y_{i}(1)\cdot T_{i}$, where $Y_{i}(0)$ and $Y_{i}(1)$ are potential outcomes for $T_{i}=0$ and $T_{i}=1$, respectively. Note that we cannot observe $Y_{i}(0)$ and $Y_{i}(1)$ simultaneously. Our purpose is to make inference on the causal effect of the treatment, or more specifically, some distributional aspects of the difference of the potential outcomes $Y_{i}(1)-Y_{i}(0)$. The RDD analysis focuses on the case where the treatment assignment $T_{i}$ is completely or partly determined by some observable covariate $X_{i}$, called the running variable. For example, to study the effect of class size on pupils' achievements, it is reasonable to consider the following setup: the unit $i$ is school, $Y_{i}$ is an average exam score, $T_{i}$ is an indicator variable for the class size ($T_{i}=0$ for one class and $T_{i}=1$ for two classes), and $X_{i}$ is the number of enrollments. For more examples, see e.g. Imbens and Lemieux (2008), Cattaneo, Titiunik and Vazquez-Bare (2020), Cattaneo and Escanciano (2017), and references therein.
Depending on the assignment rule for $T_{i}$ based on $X_{i}$, we have two cases, called the sharp and fuzzy RDDs. In this section, we focus on the sharp RDD and discuss the fuzzy RDD in Section (ref). In the sharp RDD, the treatment is deterministically assigned as $T_{i}=\mathbb{I}\{X_{i}\geq\bar{x}\}$, where $\mathbb{I}\{\cdot\}$ is the indicator function and $\bar{x}$ is a known discontinuity (cutoff) point. Throughout the paper, we normalize $\bar{x}=0$ to simplify the presentation. A parameter of interest, in this case, is the average causal effect at the discontinuity point:
Since the difference $Y_{i}(1)-Y_{i}(0)$ is unobservable, we need a tractable representation of $\tau$ in terms of quantities that can be estimated by data. If the conditional mean functions $\mathbb{E}[Y_{i}(1)|X_{i}=x]$ and $\mathbb{E}[Y_{i}(0)|X_{i}=x]$ are continuous at the cutoff point $x=0$, then the average causal effect $\tau$ can be identified as a contrast of the left and right limits of the conditional mean $\mathbb{E}[Y_{i}|X_{i}=x]$ at $x=0$, that is
As argued in CCFT, it is usually the case that practitioners have access to additional covariates (denoted by $Z_{i}\in\mathbb{R}^{p}$) and augment their empirical models with $Z_{i}$ to estimate the causal effect $\tau$ of interest. This practically relevant setup is extensively studied in CCFT for the case where $Z_{i}$ is low-dimensional. In this paper, we consider the case of possibly high-dimensional $Z_{i}$, and propose a new point estimation method for $\tau$ and an adjustment of CCFT's inference method.
We examine the case where the additional covariates $Z_{i}$ are predetermined in the sense that $Z_{i}=Z_{i}(0)\cdot(1-T_{i})+Z_{i}(1)\cdot T_{i}$ but $Z_{i}(0)=_{d}Z_{i}(1)$ for the potential covariates $Z_{i}(0)$ and $Z_{i}(1)$ for $T_{i}=0$ and $T_{i}=1$, respectively. Motivated by CCFT's recommended model (in their eq. (2)), we propose the local Lasso estimator $\hat{\theta}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}^{\prime})^{\prime}$ that solves
where $|\gamma|_{1}=\sum_{j=1}^{p}|\gamma_{j}|$ is the $\ell_{1}$-norm of $\gamma$, $\gamma_{j}$ means the $j$-th element of $\gamma$, $K(\cdot)$ is a kernel function, $b_{n}$ is a bandwidth, and $\lambda_{n}$ is a penalty level. Popular choices for $K(\cdot)$ are the uniform and triangular kernels supported on $[-b_{n},b_{n}]$. Based on ((ref)), our point estimator for $\tau$ is given by $\hat{\tau}$.
Our preliminary simulation results suggest that the local Lasso estimator for $\tau$ is somewhat biased in finite samples. Therefore, our recommendation for point estimation is to employ a post-selection method. Let $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ for a non-negative sequence $\{\zeta_{n}\}$, and $Z_{\hat{S},i}$ be a subvector of $Z_{i}$ selected by $\hat{S}$. Then the local post-Lasso estimator $\bar{\theta}=(\bar{\alpha},\bar{\tau},\bar{\beta}_{-},\bar{\beta}_{+},\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$ is defined as a solution of the local least square:
where $h_{n}$ is another bandwidth, and the estimator for $\tau$ is given by $\bar{\tau}$.
Several points are worthy of remark for this estimator. First, without the $\ell_{1}$-penalization, our estimator reduces to the local linear-type estimator recommended by CCFT's eq. (2). Therefore, the proposed estimator is a natural generalization of CCFT's when the dimension of $Z_{i}$ is high. Second, without the kernel weights for localization, our estimator in ((ref)) reduces to the conventional Lasso estimator. However, since our parameter of interest $\tau$ is identified as a local object in ((ref)), it is crucial to introduce such localization to avoid misspecification bias of the conditional mean functions. Third, it is often the case that the kernel function $K(\cdot)$ has bounded support. In this case, the effective sample size would be typically of orders $nb_{n}$ and $nh_{n}$. Thus even if the dimension of $Z_{i}$ is relatively small compared to the original sample size $n$, the $\ell_{1}$-penalization would be useful especially for small values of $b_{n}$ and $h_{n}$. Finally, the trimming term $\zeta_{n}$ to obtain the set $\hat{S}$ is introduced to stabilize numerical results (see ((ref)) below for our recommended choice based on simulation studies), and theoretically we may set as $\zeta_{n}=0$.
We now present risk properties of the local Lasso estimators $\hat{\theta}$ and $\bar{\theta}$. Let $G_{i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{i}^{\prime})^{\prime}$ be the vector of regressors in ((ref)), $G_{i,j}$ be the $j$-th element of $G_{i}$, and $\Theta_{n}=\arg\min_{\theta}\mathbb{E}[K(X_{i}/b_{n})(Y_{i}-G_{i}^{\prime}\theta)^{2}]$ be an argmin set. We impose the following assumptions.
Assumption (ref) defines $\theta_{n}^{*}$ as an approximate linear predictor or the linear projection on the set of included variables in the index set $S^{*}$ since $\mathbb{E}[K(X_{i}/b_{n})G_{i,j}\epsilon_{i}]=0$ for all $j\in S^{*}$. This assumption is general enough to cover the setup in CCFT, which assumes $p$ is fixed. In RDD analyses, it is common to introduce many generated covariates, such as transformations of initial covariates like polynomials, interactions, and various basis functions, without knowing which of them are relevant a priori.\footnote{To motivate the use of generated covariates, it is insightful to note that the asymptotic variance of CCFT's RDD estimator is proportional to $\mathrm{Var}(\{(Y_{i}(1)-Y_{i}(0))-(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma\}^{2}|X_{i}=0)$ for some $\gamma$, which is considered as the (conditional) variance of the linear projection error. Although CCFT considered linear projection due to the constraint on dimensionality, it is clear that the asymptotic variance is minimized by employing the conditional expectation $\mathbb{E}[Y_{i}(1)-Y_{i}(0)|Z_{i}(1)-Z_{i}(0),X_{i}=0]$ instead of the linear projection $(Z_{i}(1)-Z_{i}(0))^{\prime}\gamma$. Therefore, it is natural to extend CCFT's approach to high-dimensional settings by employing generated covariates or series approximation for the conditional mean.}
Although it is beyond the scope of this paper, the definition of $\theta_{n}^{*}$ could be modified to be an approximate minimizer which does not belong to $\Theta_{n}$ but gets closer to it at some suitable rate. Since it complicates the exposition and derivation as in Krei\ss\ and Rothe (2023) or Belloni, Chernozhukov and Hansen (2014), we maintain this exact sparsity assumption. For example, such an extension for approximate sparsity will be useful to allow the situation where the conditional mean satisfies $E[Y|X,Z]=E[Y|X,Z_{\mathcal{S}}]$ for some sparse set $\mathcal{S}$ but the conditional mean function $E[Y|X,Z_{\mathcal{S}}]$ is nonlinear in $Z_{\mathcal{S}}$ so that the exact sparsity assumption typically fails.
Assumption (ref) (1) contains a set of moment conditions, which are introduced to verify a local-type Bernstein inequality in Lemma (ref). The extra factor $b_{n}$ is due to the presence of the kernel weight $K(X_{i}/b_{n})$. Assumption (ref) (2) is a localized version of the compatibility condition. A sufficient condition for this is the so-called restricted eigenvalue condition. More specifically, $\min_{\beta:|\beta|_{0}\leq s^{*}}\frac{1}{nb_{n}}\sum_{i=1}^{n}K\left(\frac{X_{i}}{b_{n}}\right)\frac{\beta^{\prime}G_{i}G_{i}^{\prime}\beta}{\beta^{\prime}\beta}$, where $|\beta|_{0}$ denotes the cardinality of $\beta$, provides a lower bound for the compatibility constant $(\phi^{*})^{2}$ (see, e.g., Section 6.13 of B�hlmann and van der Geer, 2011). If CCFT's method is feasible with each subset of covariates of dimension $2s^{*}$, then the restricted eigenvalue condition is indeed satisfied. Assumption (ref) contains assumptions on the kernel $K$ and bandwidth $b_{n}$, which are standard in the literature of nonparametric methods. Note that since $\theta_{n}^{*}$ and $S^{*}$ depend on $b_{n}$, Assumption (ref) should be satisfied along each sequence $\{b_{n}\}$. Also, the deviation bounds on the prediction and estimation errors of $\hat{\theta}$ will be given as functions of $s^{*}$. While a precise condition on $s^{*}$ is hard to specify and depends on the sampling distribution, it would be typically smaller order than $\sqrt{nb_{n}}$ to satisfy the compatibility condition in Assumption (ref) (2).
Let $\hat{\gamma}_{\hat{S}}$ be the subvector of $\hat{\gamma}$ selected by $\hat{S}$, $\hat{\theta}_{\hat{S}}=(\hat{\alpha},\hat{\tau},\hat{\beta}_{-},\hat{\beta}_{+},\hat{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $Z_{\hat{S},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}$, $G_{\hat{S},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S},i}^{\prime})^{\prime}$, and $m_{n}=\lambda_{\min}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}\right)^{-1}$, where $\lambda_{\min}(A)$ means the minimum eigenvalue of a matrix $A$. Also let $\hat{S}_{1}=\{i:0<|\hat{\gamma}_{i}|<\zeta_{n}\}$, $S_{n}=\hat{S}\cup\hat{S}_{1}=\{i:\hat{\gamma}_{i}\neq0\}$, $Z_{\hat{S}_{1},i}$ be the subvector of $Z_{i}$ selected by $\hat{S}_{1}$, $G_{\hat{S}_{1},i}=(1,T_{i},X_{i},T_{i}X_{i},Z_{\hat{S}_{1},i}^{\prime})^{\prime}$, and $m_{1n}=\lambda_{\max}\left(\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/h_{n})G_{\hat{S}_{1},i}G_{\hat{S}_{1},i}^{\prime}\right)$, where $\lambda_{\max}(A)$ means the maximum eigenvalue of a matrix $A$. The $\ell_{1}$-risk properties of the local Lasso and post-Lasso estimators (for the case of $h_{n}=b_{n}$) are obtained as follows.
The proof of this theorem is presented in Appendix (ref). This theorem characterizes the risk properties of the estimators $\hat{\theta}$ and $\bar{\theta}$ around $\theta_{n}^{*}$. The risk bound of $\hat{\theta}$ depends on the tuning parameter $\lambda_{n}$, the number of non-zero coefficients $s^{*}$, and the compatibility constant $\phi^{*}$. Note that the decay rate of $\lambda_{n}$ is bounded from below by $\sqrt{\log p/(nb_{n})}$. Thus, the risk bound of $\hat{\theta}$ gets worse as the number of covariates $p$ increases or the effective sample size $nb_{n}$ due to the kernel localization decreases. The result ((ref)) for the post-selection estimator $\bar{\theta}$ shows that the deviation from the original Lasso estimator $\hat{\theta}$ is small when tuning parameter $\lambda_{n}$ or the number of selected covariates $|S_{n}|$ is small, or the minimum eigenvalue of $\frac{1}{nh_{n}}\sum_{i=1}^{n}K(X_{i}/b_{n})G_{\hat{S},i}G_{\hat{S},i}^{\prime}$ is large. If the trimming parameter $\zeta_{n}$ is of similar magnitude of $\lambda_{n}$, then the three terms in the bounds are of similar magnitude. If $\zeta_{n}$ is smaller order of magnitude than $\lambda_{n}$, then the first term will dominate. We suggest some practical choice of the trimming term $\zeta_{n}$ in Section (ref) based on our simulation studies.
The above theorem is on estimation of the coefficients of the best linear predictor $\theta_{n}^{*}$ defined in Assumption (ref) (1). Additionally suppose that the assumptions of Lemma 1 of CCFT hold true, and the covariates $Z_{i}$ are predetermined. Then we can guarantee that the second element of $\theta_{n}^{*}$ coincides with the average causal effect $\tau$ in ((ref)) so that Theorem (ref) provides the conditions for the consistency and convergence rate of $\hat{\tau}$ to $\tau$. If $Z_{i}$ are not predetermined (i.e., $Z_{i}(0)\neq_{d}Z_{i}(1)$), then $\hat{\tau}$ typically converges to $\tau$ minus some bias component, which is obtained as a limit of CCFT's bias term in their Lemma 1.
Our estimators and above theorem can be extended to other regression models that contain the covariates $\{T_{i}Z_{i},(1-T_{i})Z_{i}\}$, $(Z_{i}-\bar{Z})$, or $\{T_{i}(Z_{i}-\bar{Z}),(1-T_{i})(Z_{i}-\bar{Z})\}$ as in CCFT. However, as shown in Lemma 1 of CCFT, such estimators require more stringent conditions to guarantee the consistency for $\tau$. Furthermore, the local Lasso regression ((ref)) can be extended to incorporate polynomials of $X_{i}$ and $T_{i}X_{i}$ even though this paper focuses on the local linear model.
Finally, we discuss the choices of the localization bandwidths $b_{n}$ and $h_{n}$ and regularization parameter $\lambda_{n}$. We can use the MSE-optimal bandwidth based on the suggestion by CCFT and the regularization parameter $\lambda_{n}$ using cross-validation by Friedman, Hastie and Tibshirani (2010) or a data-driven choice by Belloni, Chernozhukov and Hansen (2014) among others. See Theorem (ref) in the next subsection for their justification, and Section (ref) for a detail on our practical recommendation.
We next consider interval estimation and hypothesis testing on the average causal effect $\tau$. For finite or low-dimensional $Z_{i}$, we recommend to use CCFT's bias corrected inference method. This subsection argues that we can still apply CCFT's inference procedure for high-dimensional $Z_{i}$, provided that CCFT's conditions remain valid for $S^{*}$ and the subvector $\theta_{n,S^{*}}^{*}$ of $\theta_{n}^{*}$ selected by $S^{*}$.
In this subsection, we specify the tuning constant $\zeta_{n}$ to obtain $\hat{S}=\{j:|\hat{\gamma}_{j}|\ge\zeta_{n}\}$ as
where we set $\varrho_{n}=\log\log\log n$. This choice of $\varrho_{n}$ is based on the simulation experiments in Section (ref), and it is not shown to be optimal but works reasonably well. Based on the selected covariates by $\hat{S}$ with $\zeta_{n}$ in ((ref)), we apply CCFT's bias corrected t-ratio to conduct statistical inference on the causal effect parameter $\tau$.
Consider the local post-Lasso estimator $\bar{\tau}$ defined by ((ref)). As shown in Appendix (ref), the dominant term of $\bar{\tau}$ can be characterized as
where $e_{2}=(0,1,0,0)^{\prime}$, $G_{1i}=(1,T_{i},X_{i},T_{i}X_{i})^{\prime}$, $\xi_{i}=T_{i}\xi_{i}(1)+(1-T_{i})\xi_{i}(0)$ with $\xi_{i}(t)=Y_{i}(t)-Z_{i}(t)^{\prime}\gamma_{Y}$, and
with $\tilde{Y}=Y(1)-Y(0)-\mathbb{E}[Y(1)-Y(0)|X=0]$. Indeed, the asymptotic linear form in ((ref)) is analogous to the one derived for the case of fixed dimensional $Z_{i}$ in CCFT (except that $\gamma_{Y,j}=0$ for $j\notin S^{*}$). Therefore, the pre-asymptotic bias and variance of $\bar{\tau}$ can be analogously written as
respectively, where $e_{1}=(1,0)^{\prime}$, $R=\left[
\right]^{\prime}$, $q=(1,-\gamma^{*\prime})^{\prime}$, and
with $\mathbf{Y}(0)=(Y_{1}(0),\ldots,Y_{n}(0))^{\prime}$, $\mathbf{Y}(1)=(Y_{1}(1),\ldots,Y_{n}(1))^{\prime}$, $\mathbf{X}=(X_{1},\ldots,X_{n})^{\prime}$, $\mathbf{Z}^{*}(0)=(Z_{S^{*},1}(0),\ldots,Z_{S^{*},n}(0))^{\prime}$, and $\mathbf{Z}^{*}(1)=(Z_{S^{*},1}(1),\ldots,Z_{S^{*},n}(1))^{\prime}$.
By estimating the unknown components, the pre-asymptotic bias and variance can be estimated as
respectively, where $\bar{q}^{\prime}=(1,-\bar{\gamma}_{\hat{S}}^{\prime})^{\prime}$, $\bar{\mu}_{-}^{(2)}$ and $\bar{\mu}_{+}^{(2)}$ are local polynomial estimators of $\mu_{-}^{(2)}$ and $\mu_{+}^{(2)}$ for the elements corresponding to $Z_{\hat{S},i}$, respectively, and $\bar{\Sigma}_{+}$ and $\bar{\Sigma}_{-}$ are conditional variance estimators of $\Sigma_{-}$ and $\Sigma_{+}$, respectively, such as the nearest neighborhood or plug-in estimator in Section 7.9 of CCFT's supplement. Based on these estimators, the t-ratio for $\tau$ is obtained as
which is exactly the same as the t-ratio in Theorem 2 of CCFT but using the selected covariates $Z_{\hat{S},i}$. By extending the theoretical developments in CCFT, we obtain the following result.
The results in ((ref)) and ((ref)) are analogous to CCFT's Theorems 1 and 2, respectively. This theorem theoretically supports to employ the bias correction and bandwidth selection methods by CCFT based on the selected covariates $Z_{\hat{S},i}$. See Section (ref) below for our practical recommendation. The assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$ is natural for predetermined covariates but may be relaxed by introducing additional regularity conditions (see, Krei\ss\ and Rothe, 2023). Other assumptions except for the last one are also imposed in CCFT. The assumption $(\sqrt{\log p}+\sqrt{nh_{n}}h_{n}^{2})(h_{n}^{2}+b_{n}^{2}+\lambda_{n}+\zeta_{n}+\zeta_{n}^{2}/\lambda_{n})s^{*}\to0$ is used to control the remainder term in ((ref)), and can be considered as a sparsity assumption to restrict the growth rate of $s^{*}$.\footnote{In the standard Lasso literature, we typically impose $\frac{s^{*}\log p}{\sqrt{n}}\to0$ and the minimal penalty level requirement on $\lambda_{n}$. For comparison, consider the following standard setting for tuning parameters, where $b_{n}\sim h_{n}\sim n^{-1/5}$, $\lambda_{n}=a_{n}\sqrt{\log p/(nb_{n})}$ with a slowly diverging $a_{n}$ and $\zeta_{n}=O(\lambda_{n})$. Then the condition on $s^{*}$ reduces to $\frac{s^{*}a_{n}\log p}{\sqrt{nh_{n}}}\to0$, which is analogous to the standard case.} Although the conditions $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ and $\frac{\bar{\mathcal{V}}}{\mathcal{V}}\overset{p}{\to}1$ are high level, these are typically satisfied for the bias and variance estimators discussed in CCFT.\footnote{Under the assumption $\partial^{2}\mathbb{E}[Z_{i}(0)|X_{i}=x]/\partial x^{2}=\partial^{2}\mathbb{E}[Z_{i}(1)|X_{i}=x]/\partial x^{2}$, the components $q^{\prime}\mu_{-}^{(2)}$ and $q^{\prime}\mu_{+}^{(2)}$ in $\mathcal{B}$ become $\mu_{Y-}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(0)|X_{i}=x]/\partial x^{2}\right|_{x=0}$ and $\mu_{Y+}^{(2)}=\left.\partial^{2}\mathbb{E}[Y_{i}(1)|X_{i}=x]/\partial x^{2}\right|_{x=0}$, respectively. Thus, in this case, the convergence rates of the conventional local polynomial estimators for $\mu_{Y-}^{(2)}$ and $\mu_{Y+}^{(2)}$ guarantee $\sqrt{\frac{nh_{n}^{5}}{\mathcal{V}}}(\bar{\mathcal{B}}-\mathcal{B})\overset{p}{\to}0$ (see, e.g., Fan and Gijbels, 1992, and Ruppert and Wand, 1994).} See Remark (ref) below for a specific example of the variance estimator $\bar{\mathcal{V}}$.