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.
74,561 characters · 22 sections · 31 citation commands
Weak instrumental variables due to nonlinearities in panel data: A Super Learner Control Function estimator
A triangular structural panel data model with additive separable individual-specific effects is used to model the causal effect of a covariate on an outcome variable when there are unobservable confounders, and some of them are time-invariant. In this setup, it is common to assume that the structural is linear and to use a linear projection in the reduced form equation. However, a linear reduced-form equation might be problematic when the conditional mean of the endogenous covariate and the instrumental variables is nonlinear. The reason is that ignoring the nonlinearity could lead to problems due to weak instruments since the strength of the linear relationship might not be high. This paper proposes to exploit the nonlinearties in the relationship between the endogenous covariate and the instrumental variables as a way to avoid problems of weak instruments due to nonlinearities. \footnote{Instruments are weak when they are weakly correlated to the endogenous covariate. In a cross-sectional setup, weak IVs causes that the median of the 2SLS estimator is biased towards the OLS estimator, 2SLS t-test presents size inflation KeaneARE2024. In panel data, Within-2SLS estimator is consistent when the time dimension grows CaiER2012.}
An example is the relationship between air pollution and children's educational outcomes, where both variables are influenced by unobserved confounders that can be either time-varying or time-invariant. For instance, neighborhood characteristics and family income are often unobserved in publicly available data. While neighborhood characteristics typically remain constant over time, family income fluctuates. Both factors are relevant to a child's academic achievement and exposure to air pollution, leading to endogeneity in pollution exposure. To address this, fixed effects can be used to control for unobserved neighborhood characteristics, while an instrumental variable can help account for time-varying confounders. A commonly used instrumental variable is wind direction, which has a well-established nonlinear relationship with air pollution levels ZabrockiN2022.
This paper studies a triangular structural panel data model in which the reduced form equation is nonlinear and individual specific effects are present. It proposes the use of Super Learning \footnote{It is feasible to use other machine learning methods instead of Super Learning. However, Super Learning automatically combines different methods in an optimal way achieving a convergence rate that is lower or equal to the best learner included in the library of base learners.} for the estimation of the nonlinearities in the model. For simplicity consider the following triangular model without covariates for $i \in \{1,...,N\}$, $t \in \{1,...,T\}$:
where $y_{it}$, and $x_{it}$ are continuous and bounded, $u_i' = [u_{i1}, u_{i2}, ..., u_{iT}]$, $z_i' = [z_{i1}, z_{i2}, ..., z_{iT}]$ are vectors stacking up all the observations available per individual. $\alpha_{i,y}$, and $\alpha_{i,1x}$ represent individual specific effects. \footnote{$\alpha_{i,y}$, and $\alpha_{i,1x}$ can be equal as in the cited example but for generality we allow them to be different.}
The model comprises a structural equation ((ref)) with an endogenous regressor $x_{1it}$ that presents a non-linear relationship with the available instrumental variable $z_{it}$ ((ref)). Both equations present individual effects that represent time-invariant unobserved confounders. In this setup, conditioning on the individual effects allows us to control for time-invariant unobserved confounders. But, estimating all the individual effects leads to incidental parameter bias as the number of parameters to estimate grows with the sample size. Then, in order to eliminate the individual effects, we propose to transform the model (e.g. first-differencing, or a within transformation which is equivalent to obtain the deviations of the original variables from their individual mean). On the other hand, the presence of time-varying unobserved confounders is modeled through a linear relationship between the error terms in both equations ((ref)). We deal with this issue by exploiting available exclusion restrictions (instrumental variables) along with a Control Function approach on the transformed model.
To estimate the causal parameter of interest $\beta_{1_o}$, we provide a two-step estimator. In the first step, we use a super learner VanSAMB2007 to learn the nonlinear relationship of the transformed endogenous variable with the instrumental variables and the observable confounders. After this, we obtain the residuals of the transformed endogenous variable. In the second step, we use the estimated residuals in the first stage as a control function for the endogenous variable in the structural equation. We call this estimator Super Learner Control Function (SLCF) estimator. The main results of the paper are that it the Super Learner Control Function estimator is consistent and asymptotically normal with a convergence rate equal to $\sqrt{N_T}$ \footnote{$N_T = \sum_i^N T_i$.}. An extensive Monte Carlo experiment is performed to test the small sample properties of the proposed estimator. We conclude that the estimator performs well provided that we can accurately learn the nuisance parameter in the first stage. An R package is available, if you are interested please contact the author.
Our paper contributes to the panel data literature in several ways. To the best of our knowledge, it is the first one to study and establish the following results for panel data. First, this paper studies the use of machine learning to estimate panel data models with an endogenous covariate (arising from unobserved confounders or measurement error) and additive individual fixed effects using a control function approach and exploiting nonlinearities. Second, it provides general theoretical results that apply to different panel data transformations \footnote{Panel data transformations are required to eliminate the individual specific effects} while using machine learning methods and control function. Third, it shows that, under suitable conditions, one can identify the transformed nuisance parameter and as a result the transformed control function. Fourth, it shows that a control function estimator with a parametric assumption in the structural equation is a double machine learning estimator. Finally, it shows that the control function estimator is not equivalent to a plug-in IV estimation.
In contrast, the available papers using machine learning to estimate panel data models focus on exogenous covariates (time-varying confounders are observed) and the presence of individual fixed effects and focus on specific panel data transformations. For instance, klosinWP2022 proposes a debiased LASSO estimator to estimate heterogeneous treatment effects of a continuous treatment, which is modeled using a structural nonparametric panel data model with additive fixed effects and strictly exogenous covariates. They account for the individual specific effects by taking the first-difference of the original model. SemenovaQE2023 focuses on a structural equation model for panel data allowing for heterogenous treatment effects with unconfoundedness on observables (sequentially exogenous covariates) and individual specific fixed effects. In order to deal with the individual fixed effects, the authors use a Mundlak approach. They provide an algorithm composed of three steps. In the first step, using cross-fitting they estimate the nuisance functions of the structural and treatment equations to obtain the residuals. In the second step, they estimate the conditional average treatment effect using lasso estimators. In the third step, they provide debiased inference.
This paper also contributes to the literature on weak instruments in panel data by being the first to focus on a setup with short time dimension in the presence of nonlinearities. While available papers on weak instrumental variables in panel data focus on the asymptotic properties of the Within 2SLS estimator showing that it is consistent even with weak instruments when the time dimension grows to infinity. CaiER2012 show that the Within 2SLS is consistent when the individual and time dimension grow to infinity. They also present asymptotical pivotal tests (Anderson-Rubin test, Kleibergen test) for longitudinal data.
In Section (ref), we present motivating examples. In Section (ref), we present the baseline model. In Section (ref), we describe the identification strategy. In Section (ref), we present the SLCF estimator. In Section (ref), we present the large sample properties of the SLCF estimator. In Section (ref), we present a comparison of the proposed estimator with IV and 2SLS plugin estimators. In Section (ref), we present the simulation experiment and its results, in Section (ref) we present the results of the empirical application, and in Section (ref) we present the conclusions of the paper.
This paper is related to the literature on semiparametric panel data models. RodriguezPoo2017JES explains that one can either use a first-differencing approach or a profiling technique to estimate a semi-parametric model with an additive separable individual specific effect and a disturbance term. Li1996EL explains that kernel estimation of the conditional expectation of the first-differenced endogenous regressor is not feasible when the dimension of the conditioning set is larger than 5. As a solution, BaltagiASF2002 presents an additive estimation of the first-differenced non-linear function.
In addition, this paper is related to the literature on orthogonal/debiased estimation, which focuses on providing $\sqrt{N}$ consistent estimators of low-dimensional parameters in the presence of high-dimensional nuisance parameters. Chernozhukov2018EJ defines an orthogonal score as the one that presents a vanishing Gateaux derivative with respect to the nuisance parameters when evaluated at the true finite-dimensional parameter values. The orthogonal score is closely related to the score proposed by Robinson1988E.
Finally, our work relates to the literature of mixed effects models and machine learning. Our model is related to a mixed effects model assuming that covariates are exogenous conditioning on the cluster random effects. EmmeneggerSJS2023 uses double machine learning to estimate a partially linear mixed-effects model. Young2024JRSSB proposes a sandwich boosting for accurately estimating partially linear models for grouped data.
The Super Learner was proposed by VanSAMB2007, it is a heterogeneous ensemble method within Supervised Learning. It consists of different base learners that are trained on the training data set, then their predictions are combined using a meta-learner by minimizing the cross-validated risk. The meta-learner can be a linear regression of the predictions of each learner on the outcome variable. VanSAMB2007 shows that the Super Learner performs asymptotically as well as or better than learner in the library of base learners.
As mentioned in the Introduction, measuring the causal effect of air pollution on health and educational outcomes is not straightforward as it is not randomly assigned. More specifically, time-varying unobservables such as family income might affect air pollution and health or educational outcomes. In addition, unobserved time-invariant variables, such as neighborhood characteristics or genetics, might also be present. While time-invariant unobservables could be dealt with by including individual fixed effects in the model, the presence of time-varying unobservable confounders requires instrumental variables. In this setup, wind direction is a commonly used as an instrumental variable. While it is well documented that the relationship between air pollution and wind direction is nonlinear ZabrockiN2022, it is common to see that the nonlinearity is ignored or modeled with a parametric assumption.
GilraineNBER2022 studies the effect of air pollution on educational outcomes. They specify the following linear regression for subject $s$, cohort $c$, district $d$, year $t$:
where $y_{sdct}$ represents the average grade of subject $s$, cohort $c$, district $d$ at year $t$. $PM2.5_{dt}$ represents the average level of $PM2.5$ at district $d$ in year $t$, $X_{sdct}$ is a vector of control variables, $W_{dt}$ is a vector of weather variables, $\omega_s$ represents subject specific effects, $\alpha_{d}$ represents district fixed effects, $\theta_{c}$ represents cohort fixed effects, $\lambda_t$ are time fixed effects, and $\epsilon_{scdt}$ is the disturbance term.
As mentioned before, $PM2.5_{dt}$ is not randomly assigned and subject to measurement error. Thus, the authors use year-to-year changes in coal production and a shift-share instrument equal to the interaction of fuel shares used for nearby production with national growth rates along with a linear first-stage regression.
As an alternative, following DeryuginaAER2019 we could use wind direction as an instrument of air pollution, but with a nonlinear first-stage equation as follows:
where $WD_{dt}$ is the average wind direction in district $d$ in year $t$, and $\gamma_{d}$ represent district specific effects.
DeryuginaNBERWP2023 studies the effect of air pollution on mortality rates. More precisely, they specify the following linear regression for county $c$, day $d$, month $m$ and year $y$:
where $y_{cdmy}^k$ represents the cumulative mortality rate $k$ days after day $d$, $SO2_{cdmy}$ represents the level of $SO2$ at day $d$ in county $c$, $X_{cdmy}$ is a vector of control variables that include precipitation, wind speed, and temperature, $\alpha_{cm}$ represents county-month fixed effects, $\alpha_{my}$ month-by-year fixed effects, and $\epsilon_{cdmy}$ is the disturbance term.
For identification, the authors use wind direction changes as an instrument for air pollution. The authors assume that the first-stage regression presents a parametric nonlinear specification.
where:
In this equation, $\mathbf{1}\{G_c=g\}$ is an indicator function that is equal to 1 if county c is member of group $g$ and 0 otherwise. $\theta_{cdmy}$ is local wind direction measured in radians. There are 100 instrumental variables.
The authors show that their results are robust to different parametric specifications of the first-stage regression but do not allow for a nonparametric specification. However, the concern is that the parametric specifications do not fully capture the nonlinearity in the first-stage equation. Thus, we propose to use a super learner to estimate the nonlinear first-stage regression as follows:
where $h(\cdot)$ is an unknown nonlinear function.
We study a model composed of a linear structural equation, the outcome variable ($y_{it}$) explained by an endogenous variable ($x_{1it}$), $K-1$ exogenous regressors ($\tilde{x}_{it}$), and an additive individual-specific effect ($\alpha_{i,y}$). The model presents an external instrumental variable ($z_{it}$) for the endogenous regressor. The endogenous regressor is mean dependent on the exogenous regressors, the instrumental variable, and the individual specific effects $\alpha_{i,x}$. We assume that: \sloppy
Assumption (ref) states that each individual $i$ presents $T_i$ observations. The total number of observations is denoted by $N_T = \sum_i T_i$.
Under assumption (ref), $\tilde{x}_{it}$ are strictly exogenous in the structural equation conditional on $\alpha_{i,y}$. In addition, $\{z_{is}, \quad s=1,\dots,T \}$ are uncorrelated with the structural error term $\varepsilon_{it}$ conditional on $\alpha_{i,y}$.
Under assumption (ref), $\tilde{x}_{it}$, $z_{it}$ are strictly exogenous in the reduced form equation conditional on the individual specific effect $\alpha_{i,1x}$. In addition, this assumption states that the endogenous regressor has a nonlinear relationship with the exogenous regressors $\tilde{x}_{it}$, and the variable $z_{it}$.
As a consequence of Assumptions (ref) and (ref), $z_{it}$ is an instrumental variable.
The decomposition of $x_{1it}$ into a conditional expectation and an additive disturbance term is appropriate under the assumption that $x_{1it}$ is continuous (Assumption (ref)). Assumption (ref) needs to be modified when $x_{1it}$ is a non-continuous endogenous explanatory variable.
As a consequence of Assumption (ref), $\omega_{it}$ is uncorrelated with $u_{it}$.
As consequence of Assumptions (ref) and (ref), the structural error term $\varepsilon_{it}$ is correlated with $x_{1it}$. The reason is that $\mathbb{E}[x_{1it}\varepsilon_{it}]=\mathbb{E}[(g_o(\tilde{x}_{it},z_{it})+\alpha_{i,x}+u_{it})\varepsilon_{it}]$ by Assumption (ref). The last expression is equal to $\mathbb{E}[(g_o(\tilde{x}_{it},z_{it})+\alpha_{i,x}+u_{it})(\rho_o u_{it}+\omega_{it})]=\rho_o\mathbb{E}[u_{it}^2]$ by Assumptions (ref) and (ref).
The error term $\omega_{it}$ has zero mean conditional on $x_{1it}$, $\tilde{x}_{it}$, $\alpha_{i,y}$, and $u_{it}$ ($\mathbb{E}[\omega_{it}|x_{1it},\tilde{x}_{it},u_{it},\alpha_{i,y}]=0$) because $\mathbb{E}[ \omega_{it} |x_{1i}, \tilde{x}_{i}, z_i, \alpha_{i,y}] =\mathbb{E}[\omega_{it}|u_{1i}, \tilde{x}_{i},z_i,\alpha_{i,y}]=\mathbb{E}[\omega_{it}\mid u_i]=0$. The first equality is a consequence of the fact that $x_{1it}$ is a one-to-one function with $u_{it}$ conditional on $\tilde{x}_{it}$, $z_{it}$, $\alpha_{i,1x} \subseteq \alpha_{i,y}$ Wooldridge2010Book. The second equality is a consequence of the mean independence of $\varepsilon_{it}$ on the exogeneous covariates, the instrumental variables, and the individual specific effects $\alpha_{i,y}$ (Assumption (ref)).
The identification of the parameter of interest $\beta_{1_o}$ relies on the presence of an instrumental variable $z_{it}$ conditioning on the individual-specific effects (Assumptions (ref) and (ref)). Since we run into an incidental parameter problem, we get rid of the individual-specific effects by using an appropriate transformation of the data $\tau: \mathbb{R} \rightarrow \mathbb{R}$. The typical transformations of the data are first-differencing and within transformation. The first-difference transformation is equal to $y_{it}-y_{it-1}$, and the Within transformation is equal to $y_{it}-\sum_{t=1}^{T_i} y_{it}/T_i$.
where $t_a$ is the first time period available after transforming the data. In the case of first-differencing $t_a=2$, and for a Within-transformation $t_a=1$.
Now, it is tempting to replace equation (ref) into (ref) to obtain:
In equation (ref), one can believe that the estimation of the unknown parameters is straightforward since the error term $\tau u_{it}\beta_{1_o}+\tau \varepsilon_{it}$ in (ref) has zero mean conditional, on an appropriate information set, as a consequence of assumptions (ref), and (ref). Moreover, one can think that the best estimation procedure is to estimate $\tau g_o(\tilde{x}_{it},z_{it})$ in the first stage and in the second stage, plug it in equation (ref) and perform simple OLS. But this is not possible because the regularization bias of the first stage estimator of $\tau g_o(\tilde{x}_{it},z_{it})$ contaminates the estimation of $\beta_{1_o}$ and $\beta_{2_o}$ Chernozhukov2018EJ, Robinson1988E, GuoWP2022 and produces estimates that are not $\sqrt{N_T}$-consistent (See Section (ref) for more details).
In order to deal with this problem, we can use a control function approach and augment the structural equation by controlling for the unobservable $\tau u_{it}$ as follows:
In model (ref), under assumptions (ref) - (ref) we can guarantee zero correlation between $\tau x_{1it}$, and $\tau \omega_{it}$ thanks to the presence of an instrumental variable in the reduced form equation that is mean independent of the error term $\omega_{it}$.
Stacking up the time observations and letting $M_{\tau_i}$ denote the transformation operator mapping $\mathbb{R}^{T_i} \rightarrow \mathbb{R}^{T_{i,a}}$, with $T_{i,a}$ equal to $T_i$ under the within transformation and $T_i-1$ under the first-difference transformation, we obtain the following moment conditions:
where $\phi(O_i;\theta_o,\tau g_o)$ is a vector of score functions equal to $(M_{\tau_i} H_{i})' V_i^{-1} M_{\tau_i} \omega_{i} $, $\theta_o=[\beta_{1_o} \quad \beta_{2_o}' \quad \rho_o]'$, $M_{\tau_i} \omega_{i}=M_{\tau_i} y_{i}-M_{\tau_i} x_{1i}\beta_{1_o}-\sum_{k=2}^K M_{\tau_i} x_{ki}\beta_{2k_o}-\rho_o M_{\tau_i} u_{i}$, $\textbf{0}=[0\quad \cdots \quad 0\quad 0]'$, and $M_{\tau_i} H_{i}=[M_{\tau_i} x_{1i}\quad M_{\tau_i}\tilde{x}_{i}\quad M_{\tau_i} u_{i}]$, and $V_i = \mathbb{E}[(M_{\tau_i}\omega_i)(M_{\tau_i}\omega_i)'| H_i, z_i]=\sigma^2_{\omega}M_{\tau_i}M_{\tau_i}'=\sigma^2_{\omega} \tilde{V}_i $.
Finally, we identify $\tau u_{it}$ since it has zero mean conditional on an appropriate information set $I_{t}=\{\tilde x_{it-j}, z_{it-j}, L_i \leq j \leq U_i \}$ with $L_i \geq t_i - T_i$ and $U_i \leq t_i-1$ as a consequence of Assumption (ref). This identification strategy was first proposed by AvilaMarquez2023, and it was presented at the 28th International Panel Data Conference at the University of Amsterdam in 2023.
\begin{example_un}{ First-difference transformation:}
We could follow BaltagiASF2002 by first-differencing the structural and reduced-form equations such that we obtain:
Now, we can use the proposed control function approach and augment the structural equation by controlling for the unobservable $\Delta u_{it}$ as follows:
Next, we set up the following moment conditions:
where $D_i \omega_{i}=D_i y_{i}-D_i x_{1i}\beta_{1_o}-D_i \tilde{x}_{i}'\beta_{2_o}-\rho_o D_i u_{i}$, $\textbf{0}=[0\quad \cdots \quad 0\quad 0]'$, $D_i H_{i}=[D_i x_{1i}\quad D_i\tilde{x}_{i}'\quad D_i u_{i}]'$, and $V_i = \sigma^2_{\omega} D_iD_i'$.
Finally, we identify $\Delta u_{it}$ since it has zero mean conditional on $\tilde{x}_{it}$, $\tilde{x}_{it-1}$, $z_{it}$, $z_{it-1}$ as a consequence of Assumption (ref).
\end{example_un}
In Section (ref), we propose an estimation method based on the Control Function approach.
As described in the previous section, we can write the transformed structural equation for observation $i$ at period $t$ as follows:
In this transformed equation, the parameters of interest are identified because $\tau x_{1it}$, $\tau \tilde{x}_{it}$ and $\tau u_{it}$ are uncorrelated with error term $\tau\omega_{it}$ under assumptions (ref) - (ref).
In addition, it is clear that if $\tau u_{it}$ is observed, we could estimate the three parameters $\beta_{1_o}$, $\beta_{2_o}$, $\rho_o$ using the sample counterparts of the population moment conditions (ref). But since $\tau u_{it}$ is not observed, we propose the following estimation procedure:
First Step:
We partition the set $\{1, 2, ..., N\}$ in B subsets $S_1$, $S_2$, ..., $S_B$, and denote $n_{T,b} = \sum_{i \in S_b} T_i$. This sample splitting is equal to the one proposed by EmmeneggerSJS2023 and it was developed without knowledge of it. Alternatively, if the cluster structure in the data is known (individuals exhibit correlation within clusters but not across them, and both cluster membership and the number of clusters are identified) the data should be partitioned accordingly. Then, for each $b \in{B}$ we estimate $\tau g_o(\tilde{x}_{it},z_{it}) = \mathbb{E}[\tau x_{1it}|I_t]$ using a super learner with partition $S_b^c= \{i \in S_g , g \neq b\}$, and we call the estimation $\widehat{\tau g_o}^{S_b^c}$. Then, we obtain the residuals $\widehat{\tau u}_{it}^{S_b}= \tau x_{1it}-\widehat{\tau g_o^{S_b^c}(\tilde{x}_{it},z_{it})}$ for partition $S_b =\{i \in S_b\}$.
Finally, the estimator of $\theta_o$ for partition $S_b$ is the solution of the sample moment conditions such that:
where $\phi(O_i;\hat{\theta}_b,\widehat{\tau g_o}^{S_b^c}) = (\tilde{H}_{i}^{S_b})'V_i^{-1}M_{\tau_i}(y_i - \tilde{H_i}^{S_b}\hat{\theta}_b)$, $\tilde{ H}_{i}^{S_b}=[M_{\tau_i} x_{1i} \quad M_{\tau_i} \tilde{X}_{i} \quad \widehat{\tau u}_{i}^{S_b}]$, $\widehat{\tau u}_{i}^{S_b}$ is the vector stacking the residuals obtained in the first stage for partition $S_b$, $\tilde{X}_i$ a $T_i\times (K-1)$ vector collecting all time observations of $\tilde{X}_{it}$.
This is equivalent to estimate $\theta_o$ by performing generalized least squares regression of $\tau y_{it}$ on $\tau x_{1it}$, $\tau \tilde{x}_{it}$ and $\widehat{\tau u}_{it}^{S_b}$ or minimizing the following quadratic loss function: $$\hat{\theta}_b=\substack{argmin\\ \theta} Q_b(O_i;\theta ,\widehat{\tau g_o}^{S_b^c}),$$ with $Q_b(O_i;\theta ,\widehat{\tau g_o}^{S_b^c})=\frac{1}{n_{T,b}} \sum_{i \in S_b}(M_{\tau_i} y_{i}-M_{\tau_i} x_{1i}\beta_{1}-M_{\tau_i} \tilde{X}_{i}\beta_{2}-\widehat{\tau u_{i}}^{S_b}\rho)'V_i^{-1}(M_{\tau_i} y_{i}-M_{\tau_i} x_{1i}\beta_{1}-M_{\tau_i} \tilde{X}_{i}\beta_{2}- \widehat{\tau u_{i}}^{S_b}\rho)$. Importantly, we can just perform ordinary least squares estimation and use a sandwich type variance-covariance matrix estimator for robust inference. In addition, OLS is the most efficient if $M_{\tau_i}=D_i$, $T = 2$, and $\omega_{it}$ is homoskedastic and uncorrelated across both $i$ and $t$.
Second Step:
We average $\hat{\theta}_b$ to obtain the cross-fitting estimator of $\theta_o$ as:
$$\hat{\theta}_o=\frac{1}{B}\sum_{b=1}^B \hat{\theta}_b.$$
Thus, the cross-fitting estimator of our parameter of interest $\beta_{1o}$ is:
$$\hat{\beta}_{1o}=\frac{1}{B}\sum_{b=1}^B \hat{\beta}_b.$$
In order to deal with the dependency of the estimator $\hat{\beta}_{1o}$ on the particular sample split, we follow Chernozhukov2018EJ by proposing to repeate the procedure $ss$ times. The estimates are aggregated either by using the mean or the median. In addition, a correction term is added to the variance estimator of $\hat{\beta}_{1o}$ (See Subsection (ref)).
\begin{example_un}{First-difference transformation (continued):}
As described in the previous section, we can write the first-differenced structural equation for observation $i$ at period $t$ as follows:
\sloppy First, we partition the set $\{1, 2, ..., N\}$ in B subsets $S_1$, $S_2$, ..., $S_B$, and for each $b \in B$ we estimate $\mathbb{E}[\Delta x_{1it}|\tilde{x}_{it},\tilde{x}_{it-1},z_{it},z_{it-1}]$ using a super learner with partition $S_b^c$. Then, we get the residuals $\widehat{\Delta u}_{it}^{S_b}=\Delta x_{1it}-\widehat{\Delta g_o^{S_b^c}(\tilde{x}_{it},z_{it})}$ for partition $S_b$. We obtain the estimator of $\theta_o$ for partition $S_b$ as the solution of the sample moment conditions such that:
where $\phi(O_i;\hat{\theta}_b,\widehat{\Delta g}_o^{S_b^c}) = (\tilde{H}_i^{S_b})' V_i^{-1} (D_i y_i - \tilde{H}_i^{S_b} \hat{\theta}_b) $, $\tilde{H}_{i}^{S_b}=[D_i x_{1i} \quad D_i \tilde{x}_{i} \quad \widehat{\Delta u}_{i}^{S_b}]$, and $\widehat{\Delta u}_{i}^{S_b}$ are the residuals obtained in the first stage for partition $S_b$.
Finally, since $\hat{\beta}_{1o,b} \in \hat{\theta}_{o,b}$ we estimate our parameter of interest $\beta_{1o}$ using:
$$\hat{\beta}_{1o}=\frac{1}{B}\sum_{b=1}^B \hat{\beta}_{1o,b}.$$
\end{example_un}
In this subsection, we show that the SLCF estimator is a debiased Machine Learning (DML) estimator. The reason is that the moment condition (ref) has a Neyman orthogonal equivalent. Without loss of generality, we assume that there are not covariates such that the model (ref) simplifies to:
If we stack up the the observations in matrix form, we obtain:
Now, we premultiply (ref) by the projection matrix in the orthogonal space of $\tau u_i $ ($A_i = I_i - M_{\tau_i} u_i( u_i'M_{\tau_i}'M_{\tau_i} u_i)^{-1} u_i'M_{\tau_i}'$), we obtain:
Using the transformed model (ref), we know that the score for identification of $\beta_{1}$ is:
The proof is presented in the appendix (Subsection (ref)).
The proof is presented in the appendix (Subsection (ref)).
The proof is not presented as it follows a similar argument as the proof of Theorem (ref).
The estimator of the variance-covariance matrix of the SLCF estimator of $\theta_o$ is given by:
where:
In order to take into account the sample splitting uncertainty, we can repeat the cross-fitting process several times as explained in Section (ref) and algorithm (ref). We also need to add correction term to the variance-covariance estimator $\hat{\Sigma}_o$ given by:
with $\hat{\theta}_o = \frac{1}{SS} \sum_{ss = 1}^{SS}\hat{\theta}_{ss}$.
The Super Learner Control function estimator is not numerical equivalent to a 2SLS estimator that plugs in the estimated transformed nuisance parameter $\tau g_o(\tilde{x}_{it}, z_{it})$. Without loss of generality, we assume that exogenous covariates are not present (model (ref) to (ref)), and $V_{i} = I_{i}$. As described in Section (ref), one could be tempted to estimate the parameter of insterest using the naive 2SLS estimator:
The IV estimator is given by:
We test the performance of the proposed estimation methods using a Monte Carlo simulation experiment. As competing methods, we consider Within OLS, and Within 2SLS estimators including polynomial transformations of the covariates and instruments of degree 1 and 5.
We generate 100 samples from a triangular simultaneous model with a linear structural equation and a non-linear reduced-form equation with $N=1000$ and $T=2$. Both equations present additive disturbance terms, and individual-specific effects. The structural equation presents one endogenous regressor $x_{1it}$ and one exogenous covariate $x_{2it}$.
where:
The unknown function $g(x_{2it},z_{it})$ is similar to the ones proposed by GuoWP2022, and depends on the parameter $a$ that controls the nonlinearity in the functional form and the strength of the linear correlation between $x_{1it}$ and $z_{it}$. As in GuoWP2022, a larger $a$ causes higher nonlinearity and lower linear correlation between $x_{1it}$ and $z_{it}$. However, we cannot compare our results to theirs because they consider a model for cross-sectional data.
In this section, we present the results of the simulation experiment. In figure (ref), we present the average estimated parameter of interest $\beta_{1o}$ for different values of $a$. The methods are Within OLS (WOLS), Within 2SLS (W2SLS) using polynomial of degree 1 of the covariate and the instrumental variables, Within 2SLS using polynomial of degree 5 of the covariate and the instrumental variables (W2SLS_polynomial), the Super Learner Control Function estimator on the First-Differenced data (FDSLCF), and the Super Learner Control Function estimator on the Within-transformed data (WSLCF). We also plot the p-value of the estimated control functions using the first-differenced and the within transformed data. As $a$ is larger, the nonlinearity in the reduced form equation is more important. This translates in a lower correlation between the endogenous regressor and the instrumental variables. As a result, the Within-2SLS estimators present a decay in their performance. In contrast, the Super Learner Control Function estimators remain stable around the true value of the parameter of interest $\beta_{1o}=1$.
In panel (a) of figure (ref), we present the boxplot of the estimated parameter of interest $\beta_{1o}$ for a highly nonlinear reduced form equation (DGP 1 with $a=5$), and in panel (b) for a linear reduced form equation (DGP 1 with $a=1$). The estimated parameters with the Super Learner Control Function estimators are tightly concentrated around the true value of the parameter of interest. In contrast, the Within estimates present wide variation around the true value.
For the estimation of the nuisance parameters, we used a Super Learner. The base learners used for the Super Learner are a linear model, a neural network, and the mean. The architecture of the neural network is one hidden layer with two neurons, the activation function is the sigmoid one, and the output function is logistic, the regularization is 0, and the maximum number of iterations is 100.
As an illustration, we estimate the causal effect of air pollution on student performance in the United States. The data on air pollution are obtained from DeryuginaAER2019. This dataset presents daily data on air pollutants: PM 2.5, ozone ($O_3$), carbon monoxide (CO), sulfur dioxide ($SO_2$), nitrogen dioxide ($NO_2$). It also contains data on wind direction, wind speed, and weather conditions at county level. We aggregate the data to yearly level by averaging across days. The data on student performance are obtained from Stanford Educational Data Archive (SEDA, FahleSEDA2024) for county-yearly level. We merge both datasets by county-year indicator. The final dataset is a panel data set with 787 counties for the period 2009 to 2013. We transform the data using the first-difference and within transformations.
We specify the following linear regression for county $d$, subject $s$, cohort $c$ at year $t$:
where $y_{dst}$ represents the average grade of subject $s$, cohort $c$, county $d$, subject $s$ at year $t$. $PM2.5_{dt}$ represents the yearly average of $PM2.5$ concentration at county $d$ in year $t$, $X_{dt}$ is a vector of control variables at county-year level including average school characteristics, social characteristics, and weather variables, $\phi_c$ corresponds to cohort specific effects, $\tilde{\alpha}_g$ grade specific effects, and $\psi_s$ corresponds to subject specific effects. $\alpha_{d}$ represents county fixed effects.
As mentioned before, $PM2.5_{dt}$ is endogenous. Thus, following DeryuginaAER2019 we use wind direction as an instrument of air pollution along with a nonlinear first-stage equation as follows:
with $WD_{dt}$ is the yearly average wind direction in county $d$ in year $t$, $\eta_{d}$ are county specific effects.
Using a First-difference transformation, we estimate that an increase in the concentration of PM 2.5 causes a reduction in test scores of 0.976 standard deviations that is statistically significant at a 5% significance level. The base learners used for the Super Learner are the mean, a linear model, a neural network, and a random forest. The architecture of the neural network is one hidden layer with two neurons, the activation function is the sigmoid one, and the output function is linear, the regularization is 0, and the maximum number of iterations is 100. The random forest is composed of 100 trees with minimum leaf size of 5, and two variables to select at each split. The number of sample splits is equal to 10, and the number of folds is 5. We use $\tilde{V}_i = I_i$ which is the most efficient in the case of the Within transformation. Inference is still valid as we used a sandwich type estimator of the variance-covariance matrix of the coefficients of the structural model.
This paper proposes a triangular simultaneous equation model for panel data with additive separable individual-specific fixed effects composed of a linear structural equation with a nonlinear reduced form equation to deal with the issue of weak instrumental variables due to nonlinearities in panel data settings. The parameter of interest is the structural parameter of the endogenous variable. The identification of this parameter is obtained under the assumption of available exclusion restrictions and using a control function approach. Estimating the parameter of interest is done using an estimator that we call Super Learner Control Function estimator (SLCFE). The estimation procedure is composed of two main steps and cross-fitting. We estimate the control function with a super learner using sample splitting. In the following step, we use the estimated control function to control for endogeneity in the structural equation. Cross-fitting is done across the individual dimension. The estimator is consistent and asymptotically normal achieving a parametric rate of convergence. We perform a Monte Carlo simulation to test the performance of the estimators proposed. We conclude that the estimator perform well, provided that we can learn accurately the nuisance parameter in the first stage.