Mostly Harmless Machine Learning: Learning Optimal Instruments in Linear IV Models

Mostly Harmless Machine Learning: Learning Optimal Instruments in Linear IV Models

Introduction

Instrumental variable (IV) designs are a popular method in empirical economics. Over 30% of all NBER working papers and top journal publications considered by currie2020technology include some discussion of instrumental variables. The vast majority of IV designs used in practice are linear IV estimated via two-stage least squares (TSLS), a familiar technique in standard introductions to econometrics and causal inference angrist2008mostly. Standard TSLS, however, leaves on the table some variation provided by the instruments that may improve precision of estimates, as it only exploits variation that is linearly related to the endogenous regressors. In the event that the instrument has a low linear correlation with the endogenous variable, but nevertheless predicts the endogenous variable well through a nonlinear transformation, we should expect TSLS to perform poorly in terms of both estimation precision and inference robustness. In particular, in some cases, TSLS would provide spuriously precise but biased estimates andrews2018weak. Such nonlinear settings become increasingly plausible when exogenous variation includes high dimensional data or alternative data, such as text, images, or other complex attributes like weather. We show that off-the-shelf machine learning techniques provide a general-purpose toolbox for leveraging such complex variation, improving instrument strength and estimate quality.

Replacing the first stage linear regression with more flexible specifications does not come without cost in terms of stronger identifying assumptions. The validity of TSLS hinges only upon the restriction that the instrument is linearly uncorrelated with unobserved disturbances in the response variable. Relaxing the linearity requires that endogenous residuals are mean zero conditional on the exogenous instruments, which is stronger. However, it is rare that a researcher has a compelling reason to believe the weaker non-correlation assumption, but rejects the slightly stronger mean-independence assumption. Indeed, whenever researchers contemplate including higher order polynomials of the instruments, they are implicitly accepting stronger assumptions than TSLS allows. In fact, by not exploiting the nonlinearities, TSLS may accidentally make a strong instrument weak, and deliver spuriously precise inference: dieterle2016simple and references therein find that several applied microeconomics papers have conclusions that are sensitive to the specification (linear vs. quadratic) of the first-stage.

A more serious identification concern with leveraging machine learning in the first-stage comes from the parametric functional form in the second stage. When there are exogenous covariates that are included in the parametric structural specification, nonlinear transformations of these covariates could in principle be valid instruments, and provide variation that precisely estimates the parameter of interest. For example, in the standard IV setup of $Y = D^\t\tau + X^\t\beta + U$ where $X$ is an exogenous covariate, imposing $\E[U \mid X] = 0$ would formally result in $X^2, X^3$, etc. being valid “excluded” instruments. However, given that the researcher's stated source of identification comes from excluded instruments, such “identifying variation” provided by covariates is more of an artifact of parametric specification than any serious information from the data that relates to the researcher's scientific inquiry.

One principled response to the above concern is to make the second stage structural specification likewise nonparametric, thereby including an infinite dimensional parameter to estimate, making the empirical design a nonparametric instrumental variable (NPIV) design. Significant theoretical and computational progress have been made in this regard newey2003instrumental,ai2003efficient,ai2007estimation,horowitz2007nonparametric,severini2012efficiency,ai2012semiparametric,hartford2017deep,dikkala2020minimax,chen2012estimation,chen2015sieve,chernozhukov2018double,chernozhukov2016locally. However, regrettably, NPIV has received relatively little attention in applied work in economics, potentially due to theoretical complications, difficulty in interpretation and troubleshooting, and computational scalability. Moreover, in some cases parametric restrictions on structural functions come from theoretical considerations or techniques like log-linearization, where estimated parameters have intuitive theoretical interpretation and policy relevance. In these cases the researcher may have compelling reasons to stick with parametric specifications.

In the spirit of being user-friendly to practitioners, this paper considers estimation and inference in an instrumental variable model where the second stage structural relationship is linear, while allowing for as much nonlinearity in the instrumental variable as possible, without creating unintended and spurious identifying variation from included covariates. Our results provide intuition and justification for using machine learning methods in instrumental variable designs. We show that with sample-splitting, under weak consistency conditions, a simple estimator that uses the predicted values of endogenous and included regressors as technical instruments is consistent, asymptotically normal, and semiparametrically efficient. The constructed instrumental variable also readily provides weak instrument diagnostics and robust procedures. Moreover, standard diagnostics like out-of-sample prediction quality are directly related to quality of estimates. In the presence of exogenous covariates that are parametrically included in the second-stage structural function, adapting machine learning techniques requires caution to avoid spurious identification from functional forms of the included covariates. To that end, we formulate and analyze the problem as a sequential moment restriction, and develop estimators that utilize machine learning for extracting nonlinear variation from and only from instruments.

\paragraph{Related Literature.} The core techniques that allow for the construction of our estimators follow from chamberlain1987asymptotic,chamberlain1992comment. The ideas in our proofs are also familiar in the double machine learning chernozhukov2018double,belloni2012sparse and semiparametrics literatures liu2020deep; our arguments, however, follow from elementary techniques that are accessible to graduate students and are self-contained. Our proposed estimator is similar to the split-sample IV or jackknife IV estimators in angrist1999jackknife, but we do not restrict ourselves to linear settings or linear smoothers. Using nonlinear or machine learning in the first stage of IV settings is considered by xu2021weakiv (for probit), hansen2014instrumental (for ridge), belloni2012sparse,chernozhukov2015post (for lasso), and bai2010instrumental (for boosting), among others; and our work can be viewed as providing a simple, unified analysis for practitioners, much in the spirit of chernozhukov2018double. To the best of our knowledge, we are the first to formally explore practical complications of making the first stage nonlinear in a context with exogenous covariates. Finally, we view our work as counterpoint to the recent work by angrist2019machine, which is more pessimistic about combining machine learning with instrumental variables---a point we explore in detail in (ref).

Main theoretical results

We consider the standard cross-sectional setup where the data $(R_i)_ {i=1}^N = (Y_i, D_i, X_i, W_i)_ {i=1}^N \iid P$ are sampled from some infinite population. $Y_i$ is some outcome variable, $D_i$ is a set of endogenous treatment variables, $X_i$ is a set of exogenous controls, and $W_i$ is a set of instrumental variables. The researcher is willing to argue that $W_i$ is exogenously or quasi-experimentally assigned. Moreover, the researcher believes that $W_i$ provides a source of variation that “identifies” an effect $\tau$ of $D_i$ on $Y_i$. We denote the endogenous variables and covariates as $T_i\equiv [1, D_i^\t, X_i^\t]^\t$ and the excluded instrument and covariates as the technical instruments $Z_i \equiv [1, W_i^\t, X_i^\t]^\t$.

A typical specification in empirical economics is the linear instrumental variables specification:

We believe that often the researcher is willing to assume more than that $U_i$ is uncorrelated with $(X_i, W_i)$. Common introductions of instrumental variables angrist2008mostly,angrist2001instrumental stress that instruments induce variation in $D_i$ and are otherwise unrelated to $U_i$, and that a common source of instruments is natural experiments. We argue that these narratives imply a stronger form of exogeneity than TSLS requires. After all, a symmetric mean-zero random variable $S$ is uncorrelated with $S^2$, but one would be hard pressed to say that $S^2$ is unrelated to $S$. Moreover, the condition $\E[W_i U_i] = 0$, strictly speaking, does not automatically make polynomial expansions of $W_i$ valid instruments, yet using higher order polynomials is common in empirical research, indicating that the conditional restriction $\E[U_i \mid W_i] = 0$ more accurately captures the assumptions imposed in many empirical projects. With this in mind, we will assume mean independence throughout the paper: $\E[U_i \mid W_i] = 0$. This stronger exogeneity assumption allows researcher to extract more identifying variation from instruments, but doing so calls for more flexible machinery for dealing with the first stage.\footnote{Moreover, under conditions such that our proposed estimator is efficient, we can test mean independence assuming $\E[W_iU_i] = 0$, since TSLS and our proposed estimator are two estimators that generate a Hausman test.}

No covariates

Let us first consider the case in which we do not have exogenous covariates $X_i$. Our mean-independence restrictions give rise to a conditional moment restriction, $ \E[Y_i - T_i^\t \theta \mid W_i] = 0 $, where $\theta = (\alpha, \tau^\t, \beta^\t)^\t$. The conditional moment restriction encodes an infinite set of unconditional moment restrictions: \[ \text{For all square integrable $\tilde \upsilon$}: \E[\tilde \upsilon(W_i)(Y_i - T_i^\t\theta)] = 0. \] chamberlain1987asymptotic finds that all relevant statistical information in a conditional moment restriction is contained in a single unconditional moment restriction involving an optimal instrument $\Upsilon^\star$, and the unconditional moment restriction with the optimal instrument delivers semiparametrically efficient estimation and inference. In our case, $\Upsilon^\star(W_i) = \frac{1}{\sigma^2(W_i)}[1, \mu(W_i)^\t]^\t$, where $ \mu(W_i) \equiv \E[D_i \mid W_i]$ and $\sigma^2(W_i) = \E[U_i^2 \mid W_i]$. We estimate $\Upsilon^\star$ with $\hat\Upsilon$ and form a plug-in estimator for $\theta$:

This is numerically equivalent to estimating (ref) with two-stage weighted least-squares with $\hat \mu(W_i)$ as an instrument and weighting with $1/\sigma^2(W_i)$. In particular, if $U_i$ is homoskedastic, the optimal instrument is simply $[1, \mu(W_i)^\t]^\t$, and two-stage least-squares with an estimate of $\mu(W_i)$ returns an estimate $\hat\theta_N$.\footnote{This approach should not be confused with what many applied researchers think of when they think of two-stage least squares, namely directly regressing $Y_i$ on the estimated instrument $\hat\Upsilon(W_i)$ by OLS - i.e. $\hat\theta_N = \pr{\frac{1}{N} \sum_{i=1}^N \hat\Upsilon(W_i) \hat\Upsilon(W_i)^\t}^{-1} \pr{\frac{1}{N} \sum_{i=1}^N \hat\Upsilon(W_i) Y_i}.$ This is what angrist2008mostly term the “forbidden regression”, and it will not generally return consistent estimates of $\theta$.} Under heteroskedasticity, this instrument is no longer optimal (in the sense of semiparametric efficiency) but remains valid. Therefore, we shall refer to the instrument with the weighting $1/\sigma^2 (W_i)$ as the optimal instrument under efficient weighting and the instrument without $1/\sigma^2(W_i)$ as the optimal instrument under identity weighting.

Under identity-weighting, estimating $\Upsilon^\star$ amounts to learning $\mu (W_i) \equiv \E[D_i \mid W_i]$, which is well-suited to machine learning techniques; this is only slightly complicated by the estimation of $\sigma^2(W_i) \equiv \E[U_i^2 \mid W_i]$ under efficient weighting. One might worry that the preliminary estimation of $\Upsilon^\star$ complicates asymptotic analysis of $\hat \theta_N$. Under a simple sampling-splitting scheme, however, we state a high-level condition for consistency, normality, and efficiency of $\hat \theta_N$. Though it simplifies the proof and potentially weakens regularity conditions, sample-splitting does reduce the amount of data used to estimate the optimal instrument $\Upsilon^\star$, but such problems can be effectively mitigated by $k$-fold sample-splitting: 20-fold sample-splitting, for instance, limits the loss of data to $5\%$ at the cost of 20 computations that can be effectively parallelized. Such concerns notwithstanding, we focus our exposition to two-fold sample-splitting.

Specifically, assume $N=2n$ for simplicity and let $S_1, S_2 \subset [N]$ be the two subsamples with size $n$. Under identity-weighting, for $j \in \{1,2\}$, form $\hat\Upsilon^{(j)}$ by estimating $\mu(W_i)$ with data from the other sample, $S_{-j}$. An estimator for $\mu$ may be a neural network or a random forest trained via empirical risk minimization, or a penalized linear regression such as elastic net.\footnote{With $k$-fold sample-splitting, $S_{-j}$ is the union of all sample-split folds other than the $j$-th one.} The estimated instrument $\hat\Upsilon_i$ is then formed by evaluating $\hat\Upsilon^{(j)} (W_i)$ for all $ i \in S_j$. We may then use (ref) to form an (identity-weighted) estimator of $\theta$ by plugging in $\hat\Upsilon$. Under efficient weighting, on each $S_{-j}$, we would use the identity-weighted estimator of $\theta$ as an initial estimator to obtain an estimate of $U_i$, and similarly predict $U_i^2$ with $W_i$ to form an estimate of $\sigma^2(W_i)$. The resulting estimated optimal instrument under efficient weighting may then be plugged into (ref) and form an efficient-weighted estimator. We term such estimators the machine learning split-sample (MLSS) estimators. The pseudocode for major procedures considered in this paper is collected in (ref).

(ref) shows that the MLSS estimator is consistent and asymptotically normal when the first-stage estimator $\hat\Upsilon^{(j)}$ converges to a strong instrument. Moreover, it is semiparametrically efficient when $\hat\Upsilon^{(j)}$ is consistent for the optimal instrument $\Upsilon^\star(W_i) \equiv [1, \mu (W_i)^\t]/\sigma^2(W_i)$ in $L^2(W)$ norm. The $L^2$ consistency condition\footnote{In many-instrument settings under bekker1994alternative-type asymptotic sequences, there may be no consistent estimator of the optimal instrument in the absence of sparsity assumptions raskutti2011minimax.} is not strong---in particular, it is weaker than the $L^2$ consistency at $o(N^ {-1/4})$-rate conditions commonly required in the double machine learning and semiparametrics literatures chernozhukov2018double,\footnote{We are not claiming that the MLSS procedure has any advantage over the double machine learning literature, but simply that the statistical problem here is sufficiently well-behaved such that we enjoy weaker conditions than is typically required.} where such conditions are considered mild.\footnote{The nuisance parameter $\E[D_i \mid W_i]$ in this setting enjoys higher-order orthogonality property described in mackey2018orthogonal. In particular, it is infinite-order orthogonal, thereby requiring no rate condition to work. Intuitively, estimation error in $\Upsilon(\cdot)$ has no effect on the moment condition $\E[\Upsilon (Y_i - \alpha - T_i^\t \theta)] = 0$ holding, and this feature of the problem makes the estimation robust to estimation of $\Upsilon$.}

Formally, regularity conditions are stated in (ref). The first condition simply states that the nuisance estimation attains some limit as sample size tends to infinity, which is a similar requirement as sampling stability in lei2018distribution. The second condition states that the limit is a strong instrument. The third condition assumes bounded moments so as to ensure a central limit theorem. The last condition, which is only required for semiparametric efficiency, states that the nuisance estimation is consistent for the optimal instrument in $L^2$ norm. For consistency of standard error estimates, we assume more bounded moments in (ref).

Exogenous covariates

The presence of covariates $X_i$ complicates the analysis considerably. Under the researcher's model, both $W_i$ and $X_i$ are considered exogenous, and thus we may assume $\E[U_i \mid Z_i] = 0$ and use it as a conditional moment restriction, under which the efficient instrument is $\var(U_i \mid Z_i)^ {-1} \E[T_i \mid Z_i]$ and our analysis from the previous section continues to apply mutatis mutandis. However, if the researcher maintains a linear specification $Y_i = T_i^\t \theta + U_i$, estimating $\theta$ based on the conditional moment restriction $\E[U_i \mid Z_i] = 0$ may inadvertently “identify” $\theta$ through nonlinear behavior in $X_i$ rather than the variation in $W_i$. Such a specification may allow the researcher to precisely estimate $\theta$ even when the instrument $W_i$ is completely irrelevant, when, say, higher-order polynomial terms in the scalar $X_i$, $X_i^2, X_i^3$, are strongly correlated with $D_i$, perhaps due to misspecification of the linear moment condition. There may well be compelling reasons why these nonlinear terms in $X_i$ allow for identification of $\tau$ under an economic or causal model; however, they are likely not the researcher's stated source of identification, and allowing their influence to leak into the estimation procedure undermines credibility of the statistical exercise.

One idea to resolve such a conundrum is to make the structural function nonparametric as well, and convert the model to a nonparametric instrumental variable regression newey2003instrumental,ai2003efficient,ai2007estimation,ai2012semiparametric,chen2012estimation (See (ref) for discussion).\footnote{More recently, chernozhukov2018double derive Neyman-orthogonal moment conditions assuming a partially linear second stage in Section 4.2 of their paper.} Another idea, which we undertake in this paper, is to weaken the moment condition and rule nonlinearities in $X_i$ as inadmissible for inference.

To that end, we analyze the statistical restrictions implied by the model and consider relaxations. The conditional moment restriction $\E[U_i \mid Z_i] = 0$ is equivalent to the following orthgonality constraint

Condition (ref) is too strong, since it allows nonlinear transforms of $X_i$ to be valid instruments. A natural idea is to restrict the class of allowable instruments $\Upsilon (W_i, X_i)$ to those that are partially linear in $X_i$, $\Upsilon(W_i, X_i) = h(W_i) + X_i^\t\ell$, thereby deliberately discarding information from nonlinear transformations of $X_i$. Doing so yields the following family of orthogonality constraints:

We may view (ref) as imposing an orthogonality condition on the structural errors $U_i$ that is intermediate between that of TSLS and that of (ref). In particular, if we define $\E_{\mathrm{(PL)}}[\cdot \mid X_i, W_i]$ as a projection operator that projects onto partially linear functions of $(X_i, W_i)$: \[ \E_{\mathrm{(PL)}}[U_i \mid X_i, W_i]\equiv \argmin_{\substack{h(X_i, W_i) \\ h (X_i, W_i) = X_i^\t \ell + g(W_i)}} \E\pr{U_i - h(X_i, W_i)}^2, \] then requiring (ref) is equivalent to requiring orthogonality under this partially linear projection operator:

In contrast, the $\cov(U_i, Z_i) = 0$ orthogonality requirement of TSLS can be written as $\E_{\mathrm{(L)}} [\cdot \mid Z_i] = 0$, where $\E_{\mathrm{(L)}}[\cdot \mid Z_i]$ is analogously defined as a projection operator onto linear functions of $Z_i$. We see that (ref) is a natural interpolation between the respective orthogonality structures on the errors $U_i$ induced by the TSLS and the conditional moment restrictions.

The moment restriction corresponding to (ref) is the following sequential moment restriction

We see that (ref) is a natural interpolation between the usual unconditional moment condition, $\E[Z_i U_i] = 0$, and the conditional moment restriction that may be spurious $\E[U_i \mid Z_i] =0$, by only allowing nonlinear information in $W_i$ to be used for estimation and inference.

Having set up the estimation problem as (equivalently) characterized by (ref), (ref), or (ref), efficient estimation is discussed by chamberlain1992comment. In particular, the optimal instrument under identity weighting takes the convenient form

which is simply $(1, X_i)$, along with the best partially linear prediction of the endogenous treatment $D_i$ from $W_i, X_i$. Observe that the only difference between (ref) and chamberlain1987asymptotic's optimal instrument under homoskedasticity is modifying $\E$ into $\E_{\mathrm{(PL)}}$. Implementing (ref) is straightforward, as by robinson1988root,\footnote{Moreover, it is easy to impose partial linear structure on certain estimators, including series regression and feedforward neural networks, and in those cases we may minimize squared error directly without Robinson's transformation.} partially linear regression is reducible to fitting two nonparametric regressions $\E[D_i \mid W_i]$ and $\E[X_i \mid W_i]$, and forming the following prediction function \[ \E_{\mathrm{(PL)}}\bk{D_i \mid X_i, W_i} = \E[D_i \mid W_i] + \E_{\mathrm{(L)}}[D_i \mid \pr{X_i - \E[X_i \mid W_i]}]. \]

Efficient estimation in the heteroskedastic case is more complex. The optimal instrument is the vector \[ \Upsilon^\star(Z_i) = \frac{\E[T_i \mid W_i]}{\sigma^2(W_i)} + \E\bk{T_i \tilde X_i^\t} \E \bk{U_i^2 \tilde X_i \tilde X_i^\t}^{-1}\tilde X_i , \quad \tilde X_i \equiv X_i - \frac{\E[X_i U_i^2 \mid W_i]} {\sigma^2(W_i)} \] and the associated set of unconditional moment restrictions are

The intuition for (ref) is the following: the two moment conditions $\E[U_i X_i] = \E[U_i \mid W_i] = 0$ provide non-orthogonal information for $\theta$ that prevents us from applying the optimal instrument on each moment condition. However, we may orthogonalize one against the other.\footnote{Orthogonal here does not refer to Neyman orthogonality chernozhukov2018double, but simply means that the two moments are uncorrelated.} In particular, the moment condition $\E[\tilde X_i U_i] = 0$ is orthogonal to $\E[U_i \mid W_i]$ in the sense that $\E[\tilde X_i U_i \cdot U_i \mid W_i] = 0$. Indeed, the term $\frac{\E[X_i U_i^2 \mid W_i]} {\sigma^2(W_i)}U_i \equiv \frac{\ip{X_iU_i, U_i}}{\ip{U_i, U_i}} U_i $ is constructed to be the projection of $X_i U_i$ onto $U_i$ under the inner product $\ip{A,B} = \E[AB \mid W_i]$.

As before, complications in nonparametric estimation can be avoided by sample splitting, where nuisance parameters are estimated on $K-1$ folds of the data and the moment condition is evaluated on the remaining fold. As a summary across our settings, we collect the nuisance parameters that require a first-step estimation in (ref). The estimator \[ \hat\theta^{\mathrm{MLSS}}_{N} = \pr{\frac{1}{N} \sum_{i=1}^N \hat\Upsilon(Z_i) T_i^\t }^{-1} \frac{1}{N} \sum_{i=1}^N \hat\Upsilon(Z_i) Y_i \] remains the same as (ref) and is subjected to the same analysis in (ref)---under conditions 1--3 in (ref), $\hat\theta^{\mathrm{MLSS}}_{N}$ is consistent and asymptotically normal, and additionally, it is semiparametrically efficient if the $L^2$-limit of $\hat\Upsilon$ coincides with the optimal instrument $\Upsilon^\star$ in their respective settings.

Lastly, we make two remarks about the case with exogenous covariates, assuming identity weighting for tractability. First, it is possible for $\E[D_i \mid W_i] = 0$ and for (ref) to generate precise estimates of the coefficient $\tau$. The reason is that it is possible for the partially linear specification $D_i = h(W_i) + X_i^\t \ell + V_i$ to generate nonzero $h (W_i)$ but zero conditional expectation, in much the same way that some regression coefficients may be zero without adjusting for $X_i$, but nonzero when adjusted for $X_i$. Whether or not this makes $W_i$ a plausibly exogenous and strong instrument is likely to be context specific. A robustness check may be generated by replacing $\E_{\mathrm{(PL)}}[D_i \mid X_i, W_i]$ with $\E[D_i \mid W_i]$, which delivers consistent and asymptotically normal estimates (assuming strong instrument) at the cost of efficiency. Second, a Frisch-Waugh-Lovell- or double machine learning-like procedure of first partialling out $X_i$ from $Y_i, D_i$ and then treating

as a conditional moment restriction also delivers consistent and asymptotically normal estimates. However, using the “optimal” instrument for (ref)---the predicted residual $\E[D_i - \E_{\mathrm{(L)}}[D_i \mid X_i] \mid W_i]$---does not achieve semiparametric efficiency, since it uses the information in the sequential moment restriction (ref) separately, without considering them jointly and orthogonalizing one against the other, resulting in efficiency loss.

Discussion

“Forbidden regression.” Nonlinearities in the first stage are often discouraged due to a “forbidden regression,” where the researcher regresses $Y$ on \smash{$\hat D$} estimated via nonlinear methods, motivated by a heuristic explanation for TSLS. As angrist2001instrumental point out, this regression is inconsistent, and consistent estimation follows from using $\hat D$ as an instrument for $D$, as we do, rather than replacing $D$ with $\hat D$---in the case where the first-stage is linear, the two estimates are numerically equivalent, but not in general.

\paragraph{Interpretation under heterogeneous treatment effects.} Assume $D_i$ is binary and suppose $Y_i = D_i Y_i(1) + (1-D_i)Y_i(0)$. Suppose the treatment follows a Roy model, $D_i = \one \pr{\mu(W_i) \ge V_i}$, where $V_i \sim \Unif[0,1]$. In this setting, the conditional moment restriction (ref) is misspecified, since it assumes constant treatment effects, and different choices of the instrument would generate estimators that converge to different population quantities. Nevertheless, the results of heckman2005structural (Section 4) show that different choices of the instrument generate estimators that estimate different weightings of marginal treatment effects (MTEs); moreover, optimal instruments, whether under identity weighting or optimal weighting, correspond to convex averages of MTEs, whereas no such guarantees are available for linear IV estimators with $W_i$ as the instrument, without assuming that $\E[D_i \mid W_i]$ is linear. The weights on the MTEs are explicitly stated in (ref).

\paragraph{Weak IV detection and robust inference.}

A major practical motivation for our work, following bai2010instrumental, is to use machine learning to rescue otherwise weak instruments due to a lack of linear correlation; nonetheless, the instrument may be irredeemably weak, and providing weak-instrument robust inference is important in practice. Relatedly, xu2021weakiv and antoine2019identification also consider weak IV inference with nonlinear first-stages; the benefits of split-sampling in the presence of many or weak instruments are recently exploited by mikusheva2020inference and date to dufour2003identification, angrist1999jackknife, staiger1994instrumental, and references therein; kaji2019theory proposes a general theory of weak identification in semiparametric settings.

On weak IV detection, our procedure produces estimated optimal instruments, which result in just-identified moment conditions. As a result, in models with a single endogenous treatment variable, the stock2005testing $F$-statistic rule-of-thumb has its exact interpretation\footnote{Namely, the worst-case bias of TSLS exceeds 10% of the worst-case bias of OLS andrews2018weak.} regardless of homo- or heteroskedasticity andrews2018weak, and the first stage $F$-statistic may be used as a tool for detecting weak instruments.

Pre-testing for weak instruments distorts downstream inferences. Alternatively, weak IV robust inferences, which are inferences of $\tau$ that are valid regardless of instrument strength, are often preferred. The procedure we propose, under identity-weighting, is readily compatible with simple robust procedures. In particular, on each subsample $S_j$, we may perform the Anderson--Rubin test anderson1949estimation and combining the results across subsamples via Bonferroni correction. For a confidence interval at the 95% nominal level with two-fold sample-splitting, this amounts to intersecting two 97.5%-nominal AR confidence intervals on each subsample, and these confidence intervals may be computed using off-the-shelf software implementations of AR confidence intervals.

More formally, consider the null hypothesis $H_0: \tau = \tau_0$. Consider a Frisch--Waugh--Lovell procedure that partials out the covariates $X_i$. Let the residuals be $U_i (\tau_0) \equiv Y_i -D_i^\t \tau_0$, and let $\tilde U_i(\tau_0) \equiv U_i(\tau_0) - \tilde\delta^\t \bar X_i$ be the residual $U_i (\tau_0)$ after partialling out $\bar X_i \equiv [1, X_i^\t]^\t$. Suppose the estimated instrument takes the form $\hat\Upsilon(Z_i) = [1, \hat\upsilon(Z_i)^\t, X_i^\t]^\t$, where $\dim\hat\upsilon(Z_i) = \dim D_i$; this requirement is satisfied under identity weighting. Similarly, partial out $\bar X_i$ from the estimated instrument to obtain $\tilde \upsilon^{(j)} \equiv \hat\upsilon^{(j)} - \tilde \lambda^\t \bar X_i$. Finally, consider the covariance $V_{n,j} (\tau_0)$ between the residual and the instrument, after partialling out the covariates $\bar X_i$, and let $\Omega_{n,j}$ be an estimate of $V_{n,j}$'s variance matrix: i.e., \[ V_{n,j}(\tau_0) \equiv \frac{1}{\sqrt{n}} \sum_{i\inS_j} \hat\upsilon^{(j)} (Z_i)\tilde U_i(\tau_0) \text{ and } \Omega_{n,j}(\tau_0) \equiv \frac{1} {n}\sum_{i\in S_j} \tilde U_i (\tau_0)^{2} \tilde \upsilon^{(j)} (Z_i)\tilde \upsilon^{(j)} (Z_i)^\t. \] The Anderson--Rubin statistic on the $j$-th subsample is then defined as the normalized magnitude of $V_{n,j}$: $\operatorname{AR}_j(\tau_0) \equiv V_{n,j}^\t \Omega_{n,j}^{-1} V_{n,j}.$ Under $H_0$, by virtue of the exclusion restriction, we should expect $V_{n,j}$ be mean-zero Gaussian, and thus $\operatorname{AR}_j$ should be $\chi^2$. Indeed, (ref) shows that on each subsample, under mild bounded moment conditions that ensure convergence ((ref)), $\operatorname{AR}_j(\tau_0)$ attains a limiting $\chi^2$ distribution. Under weak IV asymptotics, it is not necessarily the case that the AR statistics are asymptotically uncorrelated across subsamples, and so we resort to the Bonferroni procedure in outputting a single confidence interval.

\paragraph{Connection between first-stage fitting and estimate quality.} By using machine learning in the first stage, one may be able to improve the quality of the first-stage fit, as measured by out-of-sample $R^2$. We now offer an argument as to why improving that fit may improve the mean squared error of the estimator.

Consider a setting with no covariates\footnote{In the presence of covariates under identity weighting, we may, without loss of generality, partial out the covariates after estimating the optimal instrument.} and i.i.d. $U_i$ and an estimated instrument $\hat\upsilon(Z_i)$,\footnote{i.e. $\hat\Upsilon(Z_i) = [1, \hat\upsilon(Z_i)]^\t.$} meant to approximate $\E[D_i \mid Z_i]$. In linear IV, $\hat\upsilon(Z_i)$ is the linear projection of $D_i$ onto $Z_i$. Define the extra-sample error\footnote{See, for instance, friedman2001elements Section 7.4.} of an estimator $\hat\tau$ based on $\hat\upsilon (Z_i)$ to be the random quantity

where $(Y_i,D_i,Z_i)_ {i=1}^n$ is a new and independent sample unrelated to the estimate $\hat\upsilon$, and we hold $\hat\upsilon$ fixed. The subscript $n$ denotes sample quantities such as sample variances and covariances. The quantity $\mathrm{Err} (\hat\upsilon)$ is an optimistic measure of the quality of using $\hat\upsilon$ as the instrument, as construction of $\hat\upsilon$ without sample-splitting typically introduces a bias term since $\cov(\hat\upsilon, U_i)$ cannot be assumed to be small. The following calculation shows that $\mathrm{Err}(\hat\upsilon)$ scales with the inverse out-of-sample $R^2$ of $\hat\upsilon$ as a predictor of $D_i$:

where $R^2_n$ is the out-of-sample $R^2$ of predicting $D_i$ with $\hat\upsilon(Z_i)$, and $R^2(\cdot,\cdot)$ is its limit in probability.\footnote{The $F$-statistic is a monotone transformation of the $R^2$, which also serves as an indication of estimation quality.}

The out-of-sample $R^2$, which can be readily computed from a split-sample procedure, therefore offers a useful indicator for quality of estimation. In particular, if one is comfortable with the strengthened identification assumptions, there is little reason not to use the model that achieves the best out-of-sample prediction performance on the split-sample. In some settings, this best-performing model will be linear regression, but in many settings it may not be, and attempting more complex tools may deliver considerable benefits.

Moreover, much of the discussion on using machine learning for instrumental variables analysis focuses on selecting (relevant or valid) instruments belloni2012sparse,angrist2019machine assuming some level of sparsity, motivated by statistical difficulties encountered when the number of instruments is high. In light of the heuristic above, a more precise framing is perhaps combining instruments to form a better prediction of the endogenous regressor, as noted by hansen2014instrumental.

\paragraph{(When) is machine learning useful?} We conclude this section by discussing our work relative to angrist2019machine, who note that using lasso and random forest methods in the first stage does not seem to provide large performance benefits in practice, on a simulation design based on the data of angrist1991does. We note that, per our discussion above in the connection between first-stage fitting and estimate quality, a good heuristic summary for the estimation precision is the $R^2$ between the fitted instrument and the true optimal instrument---$\E[D_i\mid W_i]$ in the homoskedastic case. It is quite possible that in some settings, the conditional expectation $\E [D_i \mid W_i]$ is estimated well with linear regression, and lasso or random forest do not provide large benefits in terms of out-of-sample prediction quality. Since angrist1991does's instruments are quarter-of-birth interactions and are hence binary, it is in fact likely that predicting $D$ with linear regression performs well relative to nonlinear or complex methods\footnote{Indeed, in some of our experiments calibrated to the angrist1991does data, using a simple gradient boosting method (lightgbm) does not outperform linear regression in terms of out-of-sample $R^2$ (0.039% vs. 0.05%).} in the setting. Whether or not machine learning methods work well relative to linear methods is something that the researcher may verify in practice, via evaluating performance on a hold-out set, which is standard machine learning practice but is not yet widely adopted in empirical economics. Indeed, we observe that in both real ((ref)) and Monte Carlo ((ref)) settings where the out-of-sample prediction quality of more complex machine learning methods out-perform linear regression, MLSS estimators perform better than TSLS.

Empirical Application

We consider an empirical application in the criminal justice setting of deepivinlaw, where we consider the causal effect of appellate court decisions at the U.S. circuit court level on lengths of criminal sentences at the U.S. district courts under the jurisdiction of the circuit court. deepivinlaw exploit the fact that appellate judges are randomly assigned, and use the characteristics of appellate judges---including age, party affiliation, education, and career backgrounds---as instrumental variables. In criminal justice cases, plaintiffs rarely appeal, as it would involve trying the defendant twice for the same offense---generally not permitted in the United States; therefore, an appellate court reversal would typically be in favor of defendants, and we may posit a causal channel in which such reversals affect sentencing; for instance, the district court may be more lenient as a result of a reversal, as would be naturally predicted if the reversal sets up a precedent in favor of the defendant.

To connect the empirical setting with our notation, the outcome variable $Y$ is the change in sentencing length before and after an appellate decision, measured in months, where positive values of $Y$ indicates that sentences become longer after the appellate court decision. The endogenous treatment variable $D$ is whether or not an appellate decision reverses a district court ruling. The instruments $W$ are the characteristics of the randomly assigned circuit judge presiding over the appellate case in question, and covariates $X$ contain textual features from the circuit case, represented by Doc2Vec embeddings le2014distributed.

We compute two estimators of the optimal instrument under identity weighting based on flexible methods that are often characterized as machine learning (random forest and LightGBM, for light gradient boosting machine), as well as a variety of linear or polynomial regression estimators, with or without sample-splitting. We present our results in (ref) (see the notes to the figure for precise definitions of each estimator). For all of the split-sample estimators we have three sets of point estimates and confidence intervals, corresponding to the different splits. We also have specifications that exclude (panel (a)) and include covariates (panel (b)); except where discussed below the results are consistent across both.

The machine learning estimators perform similarly across splits, reporting mildly negative point estimates of between $1$ and $2$ months reduction in sentencing, with confidence intervals that are reasonably tight but include a zero effect. Moreover, the Wald and Anderson--Rubin confidence intervals are similar, suggesting that the instrument constructed by the ML methods is sufficiently strong to result in inferences that are not distorted.

A natural benchmark to compare these results to is TSLS, with the instruments entering either linearly or quadratically (but not interacted). Linear TSLS estimates a slightly positive effect without covariates and a slightly negative one with them. Though it has a tight Wald confidence interval, the AR interval is quite large ($[-3.3, 23.8]$). This indicates a weak instruments problem, and indeed the first-stage $F$-statistic is only 1.5. Quadratic TSLS returns a point estimate that is close to the MLSS estimates, with a very tight Wald confidence interval; however, the corresponding Anderson--Rubin interval is empty. Anderson--Rubin tests test model misspecification jointly with point nulls of the structural coefficient, and may report an empty interval if the model is misspecified.\footnote{In our case, the endogenous treatment is binary, and so the only source of model misspecification is heterogeneous treatment effects. In that case, TSLS continues to estimate population objects that are (possibly nonconvex) averages of marginal treatment effects, and arguably researchers would nonetheless like non-empty confidence sets. One benefit of split-sample approaches is that the power of the Anderson--Rubin test is wholly directed to testing the structural parameter rather than testing overidentification, since the estimated instrument always results in a just-identified system. As a result, Anderson--Rubin intervals under split-sample approaches will never be empty. Another benefit of our approach, which we discuss in (ref) and (ref), is that MLSS consistently estimates a convex average of marginal treatment effects assuming the first-stage is consistent for the conditional mean of the endogenous treatment on the exogenous instruments.} Moreover, we should still be concerned with weak instruments: The first-stage $F$-statistic is only $2.2$ (excluding covariates) and $2.3$ (including covariates).

The sample-splitting estimators based on traditional polynomial expansions rather than machine learning all perform poorly, with out-of-sample $R^2$ close to zero and consequently huge confidence intervals (the point estimates also vary wildly across splits). Overall, the MLSS estimators successfully extract more variation from the instruments than the alternatives, and consequently deliver more statistical precision.

Conclusion

In this paper, we provide a simple and user-friendly analysis of incorporating flexible prediction into instrumental variable analysis in a manner that is familiar to applied researchers. In particular, we document via elementary techniques that a split-sample IV estimator with machine learning methods as the first stage inherits classical asymptotic and optimality properties of usual instrumental regression, requiring only weak conditions governing the consistency of the first stage prediction. In the presence of covariates, we also formalize moment conditions for instrumental regression that continues to leverage nonlinearities in the excluded instrument without creating spurious identification from the nonlinearities in the included covariates. Leveraging such nonlinearity in the first stage allows the user to extract more identifying variation from the instrumental variables and can have the potential of rescuing seemingly weak instruments into strong ones, as we demonstrate with simulated data and real data from a criminal justice context. Conventional components of an instrumental variable analysis, such as identification-robust confidence sets, extend seamlessly in the presence of a machine learning first stage. We believe that machine learning in IV settings is a mostly harmless addition to the empiricist's toolbox.