Efficient and Robust Estimation of the Generalized LATE Model

The GLATE model provides a way to conduct causal inference under endogeneity when the treatment is multi-valued and unordered. As mentioned above, the identification result (conditional on the covariates) is first established in Heckman and Pinto (2018a) by using the unordered monotonicity condition. Lee and Salanié (2018) proposes another method of identification in a similar model of multi-valued treatment. Their method is concerned with continuous instruments, while the GLATE is framed in terms of discrete-valued instruments. When the treatment levels are ordered, Angrist and Imbens (1995) derives the identification and estimation results for the causal parameter, which is a weighted average of LATEs across different treatment levels.

The literature on semiparametric efficiency in program evaluation starts with the seminal work of Hahn (1998), which studies the benchmark case of estimating the average treatment effect (ATE) under unconfoundedness. For multi-level treatment, Cattaneo (2010) studies the efficient estimation of causal parameters implicitly defined through over-identified non-smooth moment conditions. In the case where unconfoundedness fails and instruments are present, Frölich (2007) calculates the SPEB for the LATE parameter, and Hong and Nekipelov (2010a) extend to the estimation of parameters implicitly defined by moment restrictions. In a more general framework encompassing missing data, Chen et al. (2004) and Chen et al. (2008) studies semiparametric efficiency bounds and efficient estimation of parameters defined through overidentifying moment restrictions. However, there is currently no theoretical research on semiparametric efficient estimation in models that encompasses endogeneity and unordered multiple treatment levels.

Several ways are available for calculating the EIF for semiparametric estimators, as illustrated by Newey (1990) and Ichimura and Newey (2022). Semiparametric efficiency calculations can be used to construct robust (Neyman orthogonal) moment conditions. This method is illustrated in Newey (1994) and Chernozhukov et al. (2016). Based on the Neyman orthogonality condition, Chernozhukov et al. (2018) introduces the DML method that suits high dimensional settings. This is because Donsker properties and stochastic equicontinuity conditions are no longer required in deriving the asymptotic distribution of the semiparametric estimator.

For testing the GLATE model, Sun (2021) proposes a bootstrap test which is the generalization and improvement of the test studied by Kitagawa (2015) in the binary LATE model.

The GLATE model has received attention in the recent empirical literature due to its ability to model multi-valued treatment. Kline and Walters (2016) evaluate the cost-effectiveness of Head Start, classifying Head Start and other preschool programs as different treatment levels against the control group of no preschool. Galindo (2020) assesses the impact of different childcare choice in Colombia on children’s development. Pinto (2021) studies the neighborhood effects and voucher effects in housing allocations using data from the Moving to Opportunity experiment. Our theoretical analysis of the GLATE model presents important tools for estimation and inference that can be applied to those empirical settings.

We assume a finite collection of instrument values $\mathcal{Z}=\{z_{1},\cdots,z_{N_{Z}}\}$ and a finite collection of treatment values $\mathcal{T}=\{t_{1},\cdots,t_{N_{T}}\}$, where $N_{Z}$ and $N_{T}$ are respectively the total number of instrument and treatment levels. The sets $\mathcal{T}$ and $\mathcal{Z}$ are categorical and unordered. The instrumental variable $Z$ denotes which of the $N_{Z}$ instrument levels is realized. The random variables $T_{z_{1}},\cdots,T_{z_{N_{Z}}}$, each taking values in $\mathcal{T}$, denote the collection of potential treatments under each instrument status. Thus, the observed treatment level is the random variable $T=T_{Z}=\sum_{z\in\mathcal{Z}}\mathbf{1}\{Z=z\}T_{z}$. For each given treatment level $t\in\mathcal{T}$, there is a potential outcome $Y_{t}\in\mathcal{Y}\subset\mathbb{R}$. The observed outcome is denoted by $Y=Y_{T}=\sum_{t\in\mathcal{T}}\mathbf{1}\{T=t\}Y_{t}$. The random vector $X\in\mathcal{X}\subset\mathbb{R}^{d_{X}}$ contains the set of covariates. The observed data is a random sample $(Y_{i},T_{i},Z_{i},X_{i}),1\leq i\leq n$.

The description above establishes a random sampling model where the researcher only observes one potential outcome, the one associated with the observed treatment. This implies that the sample of $Y$, observed from an individual with treatment $T=t$, comes from the conditional distribution of $Y_{t}$ given $T=t$ rather than from the marginal distribution of $Y_{t}$. In general, this fact leads to identifications issues and presents challenges for causal inference. To overcome these problems, we impose further structures on the model.

Assumption 1 and 2 provide the multi-valued analog of Assumption 2.1 in Abadie (2003). Assumption 1 restricts that the instrument $Z$ is independent with the potential treatments and outcomes once we condition on $X$. Assumption 2 is the conditional version of the unordered monotonicity condition proposed by Heckman and Pinto (2018a). It means that when we focus on a particular treatment level $t$ and a pair $(z,z^{\prime})$ of instrument values, the binary environment should satisfy the usual monotonicity constraint in the LATE model. Specifically, the unordered monotonicity condition requires that a shift in the instrument moves all agents uniformly toward or against each possible treatment value.²²2As pointed out by Vytlacil (2002), the LATE monotonicity condition is a restriction across individuals on the relationship between different hypothetical treatment choices defined in terms of an instrument.

We define the type $S$ of an individual as the vector of the potential treatments, that is,

By construction, $S$ is not observed. Assumption 2, the unordered monotonicity condition, is essentially a restriction on $\mathcal{S}\equiv\textit{supp}(S)$, the support of $S$. Denote the elements in $\mathcal{S}$ by $s_{1},\cdots,s_{N_{S}}$, where $N_{S}$ is the cardinality of $\mathcal{S}$. A convenient way to characterize $\mathcal{S}$ is by using the $N_{Z}\times N_{S}$ matrix $R\equiv(s_{1},\cdots,s_{N_{S}})$. The matrix $R$ is referred to as the response matrix since it describes how each type of individuals’ treatment choice responds to the instrument.

The role of $S$ is to assist the identification of the counterfactual outcomes by dividing the population into a finite number of groups, where identification can be achieved within specific groups. Those groups are defined as follows. For $k=0,\cdots,N_{Z}$, let $\Sigma_{t,k}$ be the set of types in which the treatment level $t$ appears exactly $k$ times. That is,

where $s[i]$ denotes the $i$th element of the vector $s$. In particular, the collection $\Sigma_{t,k},k=0,\cdots,N_{Z}$ forms a partition of $\mathcal{S}$.

For individuals with type $S$ in the same type set $\Sigma_{t,k}$, their treatment response in terms of $T=t$ is in a way homogeneous. Thus, it is easier intuitively to identify the marginal distribution of the potential outcome $Y_{t}$ within each $\Sigma_{t,k}$. More specifically, we define the local average structural functions (LASF) and the local average structural functions for the treated (LASF-T) as follows.

Before presenting the identification results for the above two classes of parameters, we illustrate the GLATE model in the following two examples.

We introduce some matrix notations related to the type $S$. For each treatment level $t\in\mathcal{T}$, let $B_{t}$ be a binary matrix of the same dimension as the response matrix $R$ with each element of $B_{t}$ signifying whether the corresponding element in the response matrix is $t$. That is, $B_{t}[i,j]$, the $(i,j)$th element of $B_{t}$, is whether $T_{z_{i}}$ equals $t$ for the subpopulation $S=s_{j}$. Define $b_{t,k}\equiv\left(\mathbf{1}\{s_{1}\in\Sigma_{t,k}\},\cdots,\mathbf{1}\{s_{N_{S}}\in\Sigma_{t,k}\}\right)B_{t}^{+},$ where $B_{t}^{+}$ is the Moore-Penrose inverse of $B_{t}$.

For convenience, we also need some notations regarding conditional expectations. Let

be the vector of functions that describes the conditional distribution of the instrument $Z$. For each treatment level $t\in\mathcal{T}$, let

be the vector that describes the conditional treatment probabilities given each level of the instrument. Denote

as the vector that contains the conditional outcomes for each treatment level $t$. Notice that the functions $\pi$, $P_{t}$, and $Q_{t}$ are all identified.

Theorem 2.1 identifies $p_{t,k}$, the size of the subpopulation $\Sigma_{t,k}$, and the local structural function for that subpopulation. The only exception when the identification fails is when the type set $\Sigma_{t,0}$, in which case the individual never chooses the treatment $t$. This identification result is a modification of Theorem T-6 in Heckman and Pinto (2018a) that explicitly accounts for the presence of covariates $X$. Bayes rule is applied to convert the conditional result into the unconditional one. The following theorem presents the identification result for the LASF-T.

Let $\mathcal{Z}_{t,k}\subset\mathcal{Z}$ be the set of instrument values that induces the treatment level $t$ in the type set $\Sigma_{t,k}$. That is, $\mathcal{Z}_{t,k}\equiv\left\{z_{i}\in\mathcal{Z}\mathrel{\mathop{\ordinarycolon}}s[i]=t\text{, for all }s\in\Sigma_{t,k}\right\}$, where $s[i]$ denotes the $i$th element of the vector $s$. Then define $\pi_{t,k}\equiv\sum_{z\in\mathcal{Z}_{t,k}}\pi_{z}$ as the total probability of those instrument values.

The identification results are illustrated using the two examples.

For the rest of the paper, we assume that $Y_{t},t\in\mathcal{T}$ have finite second moments. This is necessary since we are studying efficiency. Let $\iota$ denote the column vector of ones and $\zeta(Z,X,\pi)$ the diagonal matrix with the diagonal elements being $\mathbf{1}\{Z=z\}/\pi_{z}(X),z\in\mathcal{Z}$. The following theorem gives the efficient influence function (EIF) and the SPEB for the parameters identified in the preceding section.

The EIF in Theorem 3.1 can be interpreted as the moment condition from the identification results modified by an adjustment term due to the presence of unknown infinite-dimensional parameters. Take $\psi^{\beta_{t,k}}$ as an example, the terms

and

are respectively the adjustment terms due to the presence of $Q_{t}$ and $P_{t}$.

From the expression of $\psi^{\beta_{t,k}}$, we can see that the SPEB would be large when $p_{t,k}$ is small. This is because $p_{t,k}$ measures the size of the subpopulation $S\in\Sigma_{t,k}$ on which the LASF is estimated. When $p_{t,k}$ is small, we run into the weak identification issue. In Section 5, we study inference procedures that are robust against weak identification issues.

One benefit of the EIFs is that we can easily calculate the covariance matrix of different estimators. Consider an example where we are interested in two LASFs $\beta_{1}$ and $\beta_{2}$, whose EIF is given by $\psi_{1}$ and $\psi_{2}$, respectively. If the two estimators $\hat{\beta}_{1}$ and $\hat{\beta}_{2}$ are both semiparametric efficient, then their covariance matrix equals $\mathbb{E}[\psi_{1}\psi_{2}^{\prime}]$.

The derived SPEB helps determine whether an estimation procedure is efficient. In this section, we focus on the condition expectation projection (CEP) estimator.⁴⁴4The terminology “condition expectation projection” is adopted from the papers Chen et al. (2008) and Hong and Nekipelov (2010a), whereas Hahn (1998) refers to these estimators as “nonparametric imputation based estimators.” Define

The CEP procedure first estimates $\pi_{z}$, $h_{Y,t,z}$, and $h_{t,z}$ by using nonparametric estimators $\hat{\pi}_{z}$, $\hat{h}_{Y,t,z}$, and $\hat{h}_{t,z}$ respectively. These estimators can be constructed based on series or local polynomial estimation. Then $Q_{t,z}$ and $P_{t,z}$ are estimated using $\hat{Q}_{t,z}=\hat{h}_{Y,t,z}/\hat{\pi}_{z}$ and $\hat{P}_{t,z}=\hat{h}_{t,z}/\hat{\pi}_{z}$. The vectors of estimators $\hat{Q}_{t}$ and $\hat{P}_{t}$, $\hat{\pi}$ are stacked in an obvious way. Let $\hat{\pi}_{{t,k}}=\sum_{z\in\mathcal{Z}_{t,k}}\hat{\pi}_{z}$. The CEP estimators for the structural parameters are defined by

The next proposition shows that the CEP estimators are semiparametrically efficient. The result is similar in style to Hahn’s (1998) Proposition 4 that the low-level regularity conditions are omitted. Instead, the proposition assumes the high-level condition that the CEP estimators are asymptotically linear, which means they are asymptotically equivalent to sample averages. More formally, an estimator $\hat{\beta}$ of $\beta$ is asymptotically linear if it admits an influence function. That is, there exists an iid sequence $\psi_{i}$ with zero mean and finite variance such that

Since each element of the conditional expectations $h_{Y,t,z}$, $h_{t,z}$, and $\pi_{z}$ can be considered as coming from a binary LATE model, the regularity conditions in Hong and Nekipelov (2010b) should work with little modification.

The reason that this type of estimator is efficient is well explained in Ackerberg et al. (2014). The estimation problem here falls into their general semiparametric model, where the finite-dimensional parameter of interest is defined by unconditional moment restrictions. They show that the semiparametric two-step optimally weighted GMM estimators, the CEP estimators in this case, achieve the efficiency bound since the parameters of interest are exactly identified. Discussions related to this phenomenon can also be found in Chen and Santos (2018).

We next examine the efficient estimation of other policy-relevant parameters that can be derived from the parameters $\left(\beta_{t,k},\gamma_{t,k},p_{t,k},q_{t,k}\right)$. As an example, consider the type set $\Sigma_{t}\equiv\cup_{k=1}^{N_{Z-1}}\Sigma_{t,k}$, which is referred to as $t$-switchers. This subpopulation contains individuals who switch between $t$ and other treatments when given different levels of instruments. It is a generalization of the concept of compliers in the binary LATE framework.⁵⁵5Recall that switchers are also illustrated in Example 2. The LASF for the subpopulation $\Sigma_{t}$ is given by

Similarly, one can also define

which represents the LASF-T for the subpopulation of $t$-treated $t$-switchers.

For some subpopulations, a treatment effect can be identified. This point is already illustrated with Example 2 in the discussion of the identification of the usual LATE parameter. We further illustrate this point with Example 2.

To summarize the above examples using a general expression, let $\phi=\phi(\underline{p},\underline{q},\underline{\beta},\underline{\gamma})$ be a finite-dimensional parameter, where $\phi(\cdot)$ is a known continuously differentiable function, and $\underline{p}$ is the vector containing all identifiable $p_{t,k}$’s, that is, $\underline{p}\equiv\{p_{t,k}\mathrel{\mathop{\ordinarycolon}}t\in\mathcal{T},1\leq k\leq N_{Z}\}$. Let $\underline{q},\underline{\beta}$, and $\underline{\gamma}$ be defined analogously. A natural estimator can be defined through the CEP estimates, $\phi(\hat{\underline{p}},\hat{\underline{q}},\hat{\underline{\beta}},\hat{\underline{\gamma}})$. The delta method can help calculate the efficiency bound of $\phi$ and show the efficiency of $\phi(\hat{\underline{p}},\hat{\underline{q}},\hat{\underline{\beta}},\hat{\underline{\gamma}})$. In fact, by Theorem 25.47 of van der Vaart (1998), we immediately have the following corollary, which shows that plug-in estimators are efficient.

We focus on the LASF $\beta_{t,k}$. The same analysis can be applied to the other parameters. To avoid notational burden in the main text, we drop the subscript $(t,k)$ in $\beta_{t,k}$, $p_{t,k}$, and $b_{t,k}$, and the subscript $t$ in $P_{t}$ and $Q_{t}$.⁶⁶6The full subscripts are kept in the Appendices. It is straightforward to verify that the EIF $\psi^{\beta}$ has zero mean. However, we do not want to use $\psi^{\beta}$ itself as the estimating equation since it contains $1/p$ as a factor. To deal with this problem, we simply multiply $\psi^{\beta}$ by $p$ and define

The corresponding moment condition is

This moment condition is doubly robust, as demonstrated in the following proposition.

The above proposition divides the nonparametric nuisance parameters into two groups, $\pi$ and $(Q,P)$. The doubly robust moment condition is valid if either of these two groups of nuisance parameters is true. On the other hand, if the researcher uses parametric models for these nuisance parameters, then the structural parameter $\beta$ can be recovered provided that at least one of the working nuisance models is correctly specified. Therefore, the doubly robust moment condition is “less demanding” on the researcher’s ability to devise a correctly specified model for the nuisance parameters. The double robustness result in Proposition 4.1 can be seen as the GLATE extension of the existing results in the binary LATE literature (e.g., Tan, 2006; Okui et al., 2012).

The second robustness property is Neyman orthogonality. Moment conditions with this property have reduced sensitivity with respect to the nuisance parameters. Formally, Neyman orthogonality means that the moment condition has zero Gateaux derivative with respect to the nuisance parameters. The result is presented in the following proposition.

In many econometrics models, double robustness and Neyman orthogonality come in pairs. Discussions about their general relationships can be found in Chernozhukov et al. (2016). In practice, double robustness is often used for parametric estimation, as previously explained, whereas Neyman orthogonality is used in estimation with the presence of possibly high-dimensional nuisance parameters.

Next, we apply the double/debiased machine learning (DML) method developed by Chernozhukov et al. (2018) to the moment condition (5). This estimation method works even when the nuisance parameter space is complex enough that the traditional assumptions, e.g., Donsker properties, are no longer valid.⁷⁷7In two-step semiparametric estimations, Donsker properties are usually required so that a suitable stochastic equicontinuity condition is satisfied. See, for example, Assumption 2.5 in Chen et al. (2003). The implementation details are explained below.

The nuisance parameters $Q$, $P$, and $\pi$ are estimated using a cross-fitting method: Take an $L$-fold random partition of the data such that the size of each fold is $n/L$. For $l=1,\cdots,L$, let $I_{l}$ denote the set of observation indices in the $l$th fold and $I^{c}_{l}=\bigcup_{l^{\prime}\neq l}I_{l^{\prime}}$ the set of observation indices not in the $l$th fold. Define $\check{Q}^{l}$, $\check{P}^{l}$, and $\check{\pi}^{l}$ to be the estimates constructed by using data from $I_{l}^{c}$. The DML estimator of $\beta$ is constructed following the moment condition (5):⁸⁸8This is the DML2 estimator defined in Chernozhukov et al. (2018). Another estimator, the DML1 estimator, is proposed in the same paper. We do not study the DML1 estimator since it is asymptotically equivalent to DML2, and the authors generally recommend DML2.

To conduct inference, we also need an estimate for the asymptotic variance of $\check{\beta}$, which we denote by $\sigma^{2}$. The asymptotic variance equals to the expectation of the squared efficient influence function: $\sigma^{2}=\mathbb{E}\left[\psi^{\beta}\right]^{2}=\mathbb{E}[\psi^{2}]/p^{2}$. We first estimate $p$ by using the cross-fitting method, which is essentially given by the denominator of (6):

Then the asymptotic variance can be estimated by

We want to establish the convergence results for the DML estimator uniformly over a class of data generating processes (DGPs) defined as follows. For any two constants $c_{1}>c_{0}>0$, let $\mathcal{P}(c_{1},c_{0})$ be the set of joint distributions of $(Y,T,Z,X)$ such that