Revisiting the Many Instruments Problem using Random Matrix Theory

Revisiting the Many Instruments Problem using Random Matrix Theory

Introduction

In this paper, we propose a new bias-adjusted two-stage least-squares (2SLS) Ridge estimator for a simultaneous equations model with a single endogenous regressor and many instrumental variables (IV). Our motivation comes from recent results in the theory of random matrices (RMT), which give a new understanding of the bias induced by many instruments.

It is well-known that IV estimators perform poorly under many and/or weak instruments. In these settings, IV estimators are biased towards the ordinary least-squares (OLS) estimate bound1995problems,staiger1997instrumental,stock2000gmm,chao2005consistent. The standard IV estimator is a 2SLS procedure, where we first use least-squares to predict the endogenous variable and use these predictions in a second least-squares regression to estimate the causal effect. Thus, the accuracy of 2SLS estimation depends crucially on the prediction quality in the first-stage regression. We would expect a method that reduces the predictive risk in the first stage to perform better in the second-stage regression. This reasoning motivated a series of papers to replace OLS in the first stage with statistical learning tools (such as belloni2012sparse, gold2020inference and angrist2022machine for Lasso regressions or okui2011instrumental and hansen2014instrumental for Ridge regressions). In this paper, we are interested in settings with dense first-stage parameters such that Ridge is the preferred choice in the first-stage regression.

Bias-adjusted 2SLS estimation can be motivated by the Silverstein equation. This link to random matrix theory allows us to improve IV estimation even in settings where no closed-form expression of the Silverstein equation exists. We first consider a standard many IV setting where the number of instruments is large relative to, but still smaller than, the sample size. In this setting, our bias approximation nests the well-known results for least-squares in the first stage hahn2002new,hahn2002notes,bun2011comparison, and a special case of our bias-adjustment coincides with the well-known bias-adjusted 2SLS estimator nagar1959bias,chao2005consistent,donald2001choosing. Bias-adjustments are costly in the sense that they reduce the overall signal strength of the instruments. We show that using Ridge in the first stage reduces the implicit price of bias-adjustments. This leads to a trade-off, which allows for less costly estimation of the causal effect parameter and comes along with improved asymptotic properties in the second-stage regression. Additionally, the results mentioned above generalize to settings where the number of instruments is larger than the number of observations, which has yet to be discussed in the IV literature. Finally, we derive the optimal tuning parameter of Ridge regressions in simultaneous equations models, which comprises the well-known results for single equation models as special case with uncorrelated error terms.

To derive our results, we make the following assumptions. We assume that each instrument has only a small, independent random effect on the endogenous variable. This can be considered an average-case analysis over dense parameters. chamberlain2004random also use a random effects assumption to derive a quasi-maximum likelihood estimator in an instrumental variables regression with many instruments. Moreover, we assume a high-dimensional asymptotic regime for the instrumental variables. In such a setting, we would expect Ridge to outperform OLS with respect to prediction risk dobriban_wager. While being standard in the RMT literature, these assumptions are restrictive. However, they allow us to derive insightful results that nest existing approaches and may serve as a basis for future conceptual developments.

Closest to our study in the RMT literature are dobriban_wager and hastie2022surprises who analyze the performance of Ridge in a single equation model. In the econometrics literature, our study is close to hahn2002new, hahn2002notes and bun2011comparison, who analyze the bias of 2SLS under many and/or weak instruments. Further, our study contributes to several econometric approaches that adjust for the 2SLS bias nagar1959bias,chamberlain2004random,chao2005consistent,hansen2014instrumental,spiess2017bias. In economics, our study is closely related to angrist2022machine, who compare bias-adjusted IV estimators with estimators based on Lasso in the first stage regression. They conclude that statistical learning techniques are less encouraging in IV settings. In contrast to this conclusion, we show that bias-adjusted 2SLS estimators can indeed be motivated by results from the statistical learning literature. These allow us to improve the bias-adjustment, which leads to a better finite-sample performance. Our computational results can be replicated using software available from \href{https://github.com/farbmacher/RMT-many-instruments}{https://github.com/farbmacher/RMT-many-instruments}

Model and 2SLS Estimation

We have a random i.i.d. sample $\left\{ y_i,x_i,z_i' \right\} _{i=1}^{n}$ drawn from the following simultaneous equations model

where $y_i$ and $x_i$ are the entries of the $n\times 1$ vectors $y$ and $x$, and $z_i$ are arranged as the rows of the $n \times k$ matrix $Z$. The aim is to estimate the causal effect $\beta$. Throughout, we assume that the following model conditions hold.

Assumption (ref) imposes homoskedasticity and joint normality of $\epsilon$ and $\nu$. If $\sigma_{\varepsilon\nu}\neq 0 $, $x$ is endogenous, and the OLS estimator of $\beta$ is inconsistent. The instrumental variables $Z$ can be used to overcome the endogeneity. The standard IV estimator is a 2SLS procedure. First, we use least-squares to predict $\widehat{x}=P_0 x$ with $P_0=Z(Z'Z)^{-1}Z' \, $. In a second least-squares regression, we use $\widehat{x} \, $ to estimate the causal effect $\beta$. A first-stage prediction based on Ridge is defined as $\widehat{x}_\lambda=P_\lambda x$ with $P_{\lambda}=Z(Z'Z/n+\lambda I_k)^{-1}Z'/n \, $ being the Ridge smoother and $\lambda$ being a tuning parameter. When we set $\lambda\downarrow 0$ in settings with fewer instruments than observations, we get 2SLS. If we set $\lambda\downarrow 0$ in settings with more instruments than observations, we obtain a ridgeless 2SLS estimator.\footnote{For a detailed discussion of ridgeless least-squares and its interpretation as min-norm least-squares estimator see, for example, hastie2022surprises.}

Following the RMT literature, we consider a high-dimensional asymptotic regime for the instrumental variables.

See dobriban_wager for a detailed discussion of covariance matrices $\Sigma$ that are aligned with this setting. For example, the AR-1 model with $\Sigma_{ij}=\rho^{\vert i-j \vert}$, which is often used in econometrics, satisfies Assumption (ref).

Moreover, following the RMT literature, we assume that each instrument has only a small and independent random effect on the endogenous variable.

This assumption has been used before in the literature on many instruments, e.g. in a theoretical analysis without RMT by chamberlain2004random, and can be regarded as an average-case analysis over dense parameters. The concentration parameter, $\mu^2\equiv n\alpha^2/\sigma_\nu^2 \, $ is a standard measure of instrument strength. In a weak instrument setting, the concentration parameter is small relative to the sample size.

Throughout the paper, we illustrate our theoretical results using Monte-Carlo simulations based on a standard many IV setting

with $\varepsilon=\rho \nu +\sqrt{1-\rho^2} w$ where $\nu\sim N(0,1)$ and $w\sim N(0,1)$ are independent, and where we generally set $\rho=0.6$. Moreover, $\pi\sim N(0,\sigma^2_\pi)$ with $\sigma^2_\pi=\alpha^2/k$ and $\alpha^2=\gamma F \, $. This allows us to control the first-stage $F$ statistic, which we generally set equal to 5. We show simulation results for $Z\sim N(0,\Sigma)$ and a) isotropic $\Sigma=I_k$ or b) AR-1 model $\Sigma_{ij}=0.5^{\vert i-j\vert}$. Additionally, we set $n=200$, $k=150$ or $k=250$, and $\beta=0$.

Results from Random Matrix Theory

We start with a brief review of key results from the theory of random matrices, which help us to describe and adjust for the bias induced by Ridge estimation in the first-stage regression. From marchenko1967distribution and silverstein1995strong we know that the spectral distribution of the sample covariance matrix $\widehat\Sigma=Z'Z/n$ converges weakly, with probability 1, to a limiting distribution. This limiting distribution is called the empirical spectral distribution (ESD). Let $Y$ be a random variable distributed according to the ESD of $\widehat\Sigma$, then the well-known Stieltjes transform $m(-\lambda)$ is

The Stieltjes transform of the spectral measure of $\widehat\Sigma$ satisfies

Additionally, we need the companion Stieltjes transform $v(-\lambda)$, which is the Stieltjes transform of the ESD of the $n \times n$ matrix $\underline{\widehat\Sigma}=ZZ'/n\,$. Let $\underline{Y}$ be a random variable distributed according to the ESD of $\underline{\widehat\Sigma}\, $, then

and the Stieltjes transform of the spectral measure of $\underline{\widehat\Sigma}$ satisfies

To calculate the population quantities $m(-\lambda)$ and $v(-\lambda)$ as a benchmark in the simulations, we use the procedure described in dobriban_wager_supplement. For estimation, we follow dobriban2020wonder and use $\widehat{m}(-\lambda)$ and $\widehat{v}(-\lambda) \, $.

In the following, we repeatedly use well-known results in the RMT literature (the proofs are replicated in Appendix (ref) for convenience). First, we mention a useful representation of the Ridge smoother.

Second, we mention a useful link between the Stieltjes transform and its companion.

Third, we mention limits involving the Stieltjes transforms that help us to work out the links to the existing literature on many weak instruments.

We, moreover, make use of the Silverstein equation silverstein1995strong, which is

where $\frac{dH(t)}{dt}$ is the pdf of the PSD. The Silverstein equation links the limit PSD to the limit ESD. If instruments are equi-correlated with correlation coefficient $\rho_z$, $\Sigma$ has $k-1$ eigenvalues that equal $1-\rho_z$, thus, the pdf of the PSD is a point mass at $1-\rho_z$. Therefore, the Silverstein equation has a closed-form expression

In the special case of isotropic instruments ($\rho_z=0$), this simplifies to

hastie2022supplement.

2SLS Estimation and RMT

Bias Approximation

We derive a bias approximation for the causal effect $\beta$ when the first-stage regression is based on Ridge. Consider the first-stage Ridge estimates $\widehat\pi_\lambda$ from a regression of the endogenous variable $x$ on the instruments $Z$ given by

or, equivalently written as $\widehat\pi_\lambda=\big( Z'Z/n+\lambda I_k \big)^{-1}Z'x/n\, $. In the second stage, we use the predicted values $\widehat{x}_{\lambda}=Z\widehat\pi_\lambda$ in a just-identified 2SLS regression to obtain

In the following, we call $\widetilde\beta_\lambda$ the 2SLS-Ridge estimator to emphasize that it is a two-stage procedure with Ridge in the first stage and just-identified 2SLS in the second stage regression. For the limit $\lambda\downarrow 0$, we get the (ridgeless) 2SLS estimator, which corresponds to standard 2SLS if $\gamma<1$. As already mentioned in hansen2014instrumental, 2SLS-Ridge is biased. Theorem (ref) gives precise information about the bias induced by using Ridge in the first-stage regression. While the numerator only contains a noise term, the denominator is the sum of a signal term and a noise term. For convenience, we use the Stieltjes transform $m(-\lambda)$ for all terms related to the identification strength $\alpha^2$ and the companion Stieltjes transform $v(-\lambda)$ for all terms related to $\sigma^2_\nu$ or $\sigma_{\varepsilon\nu}$.

Before discussing ways to get rid of the bias term later in Section (ref), we briefly address the links of our bias approximation to the existing literature and the advantage of using Ridge compared to least-squares in the first-stage regression. The first part of the next corollary replicates the well-established 2SLS bias approximation using Theorem (ref) and Lemma (ref) (compare, for example, hahn2002new,hahn2002notes and chao2005consistent). The second part of the corollary extends these results to ridgeless 2SLS estimation.

Comparing these results with the bias of OLS hahn2002notes, which is

we find that the (ridgeless) 2SLS bias is never larger than the OLS bias (even if $\gamma>1$). The equivalence of the ridgeless 2SLS and OLS bias for $\gamma>1$ intuitively follows from the fact that the ridgeless 2SLS estimator interpolates the data in the first stage if $\gamma>1$, i.e. $x=Z\widehat\pi_{(\lambda\downarrow 0)}=\widehat{x}_{(\lambda\downarrow 0)}$. More interestingly, we can also show that increasing $\lambda$ reduces the bias relative to the (ridgeless) 2SLS bias (even if $\gamma>1$). To see this, we can rearrange the bias of 2SLS-Ridge (when we rule out the case $\lambda v(-\lambda)=1$, which occurs if and only if $\lambda\to\infty$)

$F \equiv \mu^2/k$ is the population first-stage $F$-statistic and $a(-\lambda)\equiv \gamma\frac{1-\lambda p(-\lambda)}{1-\lambda v(-\lambda)}$ can be interpreted as a Ridge amplifier of the first-stage $F$-statistic. We show that this amplifier equals unity for the standard 2SLS (i.e., $\lambda=0$). Additionally, we show that its first derivative close to the ridgeless limit is positive, which suggests that there is a value of $\lambda$ for which the bias of 2SLS-Ridge is smaller compared to (ridgeless) 2SLS.

Figure (ref) depicts the bias approximation from Theorem (ref) and the simulated bias of $\widetilde\beta_\lambda$ for different values of $\lambda$. In the left part of Figure (ref), we see the standard bias approximation for 2SLS if $\lambda\downarrow 0$ and $\gamma<1$, while in the right part we can observe that the ridgeless 2SLS bias equals the OLS bias if $\gamma>1$, which illustrates the results in Corollary (ref). Independent of $\gamma$, we see that increasing $\lambda$ reduces the bias of 2SLS-Ridge. However, while increasing $\lambda$ reduces the bias of $\widetilde\beta_\lambda$, we cannot expect to get rid of it fully. Therefore, we propose a bias-adjustment for the 2SLS-Ridge estimator in the following.

Bias-Adjustment

The standard bias-adjusted 2SLS (aka Nagar's estimator, nagar1959bias; donald2001choosing; chao2005consistent) has several representations

Analogous to the bias-adjusted OLS projection matrix, $P_0 - \gamma I_n \, $, we propose a bias-adjusted Ridge smoother

For the first term, we know $\operatorname{tr}(P_\lambda)/n\to_{a.s.} 1-\lambda v(-\lambda)$. Therefore, setting

we know from the Silverstein equation that $\operatorname{tr}\big(S_\lambda(s_1)\big)/n \to_{a.s.} 0$. Unfortunately, the integral depends on unknown population values and cannot be derived in general. There are, however, several approaches to proceed. First, we may consider the case with $\lambda\downarrow 0$ and $\gamma<1$. Rearranging $s_1$ and using the dominated convergence theorem together with Lemma (ref), we can derive that

because $\lim_{\lambda\downarrow 0} \lambda v(-\lambda) = 1-\gamma$ and hence $\lim_{\lambda\downarrow 0} \lambda / \lambda v(-\lambda) = 0$. Again invoking the Silverstein equation gives us

which, for $\gamma<1$, can be used to construct the standard bias-adjusted 2SLS. Second, the availability of a closed-form expression for the Silverstein equation can be used to debias the 2SLS-Ridge estimator for values of $\lambda$ different from zero and/or values of $\gamma$ larger than unity. We illustrate this for isotropic instruments ($\Sigma=I_k$). When we consider

with the limit Stieltjes transform for such instruments depending only on $\gamma$ and $\lambda$

then we know from the Silverstein equation that $\operatorname{tr}\big( S_\lambda(s_2)\big)/n\to_{a.s.} 0$.

Third, we can use an estimate of $v(-\lambda)$ to construct an $s$ for general matrices $\Sigma$ that fulfill Assumption (ref). Consider

with

and

such that by construction, $\operatorname{tr}\big(S_\lambda(s_3)\big)/n=0\, $. We focus on $S_\lambda(s_3)$, in the following simply referred to as $S_\lambda \, $, which gives us a bias-adjusted 2SLS-Ridge estimator that works for general values of $\lambda \, $,

The bias-adjustment utilizes that $P_\lambda$ acts differently on the signal and noise terms. By construction, the noise term is successfully eliminated while the signal is preserved. Theorem (ref) shows that $\widehat\beta_{\lambda}$ estimates the causal effect consistently.

Using Lemma (ref), we show in the following corollary that Theorem (ref) also replicates the standard bias-adjusted 2SLS estimator as a special case (if $\gamma<1$ and $\lambda\downarrow 0$). Additionally, it extends this estimator to ridgeless 2SLS (if $\gamma>1$ and $\lambda\downarrow 0$).

Figure (ref) shows that $\widehat\beta_\lambda$ is (by construction) unbiased for all values of $\lambda$ (including $\lambda\downarrow 0$). The bias-adjustment works also if $\gamma>1\, $ (given that $\lambda\neq 0$).

It is well known that the standard bias-adjusted 2SLS estimator (i.e., $\lambda=0$ and $\gamma<1$) can perform poorly in finite samples. A reason for this is that $\gamma$ in Corollary (ref) can be interpreted as costs for the bias-adjustment, which reduce the signal to $\alpha^2(1-\gamma)$. By inspecting the limit of the denominator of bias-adjusted 2SLS-Ridge derived in the proof of Theorem (ref) (see Appendix (ref)), we see two negative terms:

The first term ($-\lambda p(-\lambda)$), which appears also in the asymptotic bias approximation of $\widetilde\beta_\lambda$ (see Theorem (ref)), stems from using Ridge and the second term ($-\gamma p(-\lambda)$) stems from the bias-adjustment. In both cases, $p(-\lambda)$ can be interpreted as a price per unit, which is unity if $\lambda=0$ and declines with increasing $\lambda$. This decline leads to a trade-off, allowing for less costly bias-adjustments (with respect to the signal term).

The following corollary discusses this trade-off close to the ridgeless limit and shows under what conditions we can expect a stronger signal term if $\lambda$ increases. Compared with standard bias-adjusted 2SLS, the trade-off turns out to be positive if $\gamma \textrm{E}\left[ Y^{-1} \right]\geq 1$. This condition simplifies to $\gamma>1/2 \, $ for isotropic instruments (i.e., when $\textrm{E}\left[ T^{-1} \right]=1$) because from the Marchenko-Pastur equation, we know $\textrm{E}(Y^{-1})=\textrm{E}(T^{-1})/(1-\gamma)$.\footnote{In the next section, we additionally show that the asymptotic variance of $\widehat\beta_{\lambda}$ with $\lambda\neq 0$ turns out to be smaller than the asymptotic variance of the standard bias-adjusted 2SLS ($\lambda=0$) under the weaker condition $\textrm{E}[Y^{-1}]\geq1$, which, for the isotropic case, is always fulfilled ($\gamma>0$).}

Figure (ref) illustrates the trade-off described in Corollary (ref) for different values of $\gamma$ and $\lambda$. If $\gamma>1/2$, the signal can be improved compared to standard bias-adjusted 2SLS. If $\gamma>1$, we know for standard bias-adjusted 2SLS that $\lim_{\lambda\downarrow 0} f(-\lambda) = 0$ (using Lemma (ref)) implying that standard bias-adjusted 2SLS is infeasible. However, a positive value of $\lambda$ can reestablish the signal, which enables us to identify the causal effect even if $\gamma>1$. We observe similar patterns for the AR-1 covariance matrix (see Figure (ref) in Appendix (ref)).

Variance and Choice of the Tuning Parameter

Variance Estimation and Asymptotic Variance

For inference, we consider bekker1994alternative standard errors. Let $\widehat\varepsilon=y-x\widehat\beta_\lambda$, $\widehat\sigma_\varepsilon^2=\widehat\varepsilon'\widehat\varepsilon/n$, and $\widetilde{x} = x + \widehat\varepsilon \, (\widehat\varepsilon'x)/(\widehat\varepsilon'\widehat\varepsilon)$. Inspired by hansen2008estimation, we can estimate the variance of the standard bias-adjusted 2SLS estimator by

Following this structure, we propose to estimate the variance of the bias-adjusted 2SLS-Ridge estimator by

The corresponding asymptotic variance of $\widehat\beta_\lambda$ is (see Appendix (ref) for a derivation)

which for $\lambda\downarrow 0$ coincides with the well-known asymptotic variance of the standard bias-adjusted 2SLS estimator under normality of the errors bekker1994alternative, anatolyev2013instrumental,kolesar2015identification,anatolyev2019many

In Appendix (ref), we show that for any positive value of $\lambda$, the asymptotic variance of $\widehat\beta_\lambda$ is always smaller than the asymptotic variance of $\widehat\beta_0$. In Figure (ref), we present results from Monte-Carlo simulations over a grid of different $\lambda$ values. The results illustrate that the observed variation of $\widehat\beta_\lambda$, and the estimated Bekker-type variance are close to the asymptotic variance of the bias-adjusted 2SLS-Ridge estimator. To showcase convergence of the observed variance to the asymptotic variance, we use $n=1000$ (for results using $n=200$, see Figure (ref) in Appendix (ref)). For settings where $\gamma < 1$, we further display the asymptotic variances of the standard bias-adjusted 2SLS and the limited information maximum likelihood estimator (Liml) as benchmarks. As suggested by our asymptotic results, we see considerable efficiency improvements compared to the standard bias-adjusted 2SLS estimator, and, in some settings also to the Liml estimator. In Section (ref), we further illustrate for which settings bias-adjusted 2SLS-Ridge outperforms Liml regarding estimation accuracy.

To further illustrate the validity of the variance estimates, Figure (ref) shows for different settings that under the null-hypothesis, observed rejection frequencies for testing at the 5%-level are always close to their nominal value.

Optimal Choice of the Tuning Parameter

Since $\widehat\beta_\lambda$ is consistent for all values of $\lambda$, we aim to find the optimal value of $\lambda$ that minimizes the variance of $\widehat\beta_\lambda$. To find this value theoretically, we take the derivative of the asymptotic variance (Equation (ref)) and set it to zero. This gives us the following expression (for a complete derivation, see Appendix (ref)):

This expression is fulfilled for $\lambda^* = \frac{1+\rho_{\varepsilon\nu}^2}{F}$. Interestingly, for $\rho_{\epsilon \nu}=0$, this result coincides with the optimal $\lambda^*$ in a single equation model which is equal to $\frac{1}{F}$ dobriban_wager.

To obtain $\lambda^*$ empirically, we estimate $\rho_{\varepsilon\nu}$ using the formula for the OLS bias hahn2002notes:

with bias-adjusted 2SLS ($\hat{\beta}_{\downarrow 0}$) as consistent but inefficient preliminary estimate.

Finally, we can compute the optimal $\lambda$ using:

Figure (ref) in Appendix (ref) illustrates for different settings that estimating $\hat\lambda^*$ via Equations (ref) and (ref) gives a reliable estimate of $\lambda^* \, $.

Next, we compare the asymptotic variance formula of bias-adjusted 2SLS-Ridge using the optimal tuning parameter $\lambda^*$ (see Appendix (ref) for a derivation) to that of Liml:

In Figure (ref), we plot the relative difference of the asymptotic variances of the two estimators over a grid of $\gamma$ and $\rho_z$ values for two distinct values of $\rho_{\epsilon \nu}$. We can see that the improvement in efficiency over Liml can be substantial, especially when the number of instruments approaches the number of observations.

Given Figure (ref), we construct three settings for two different values of $\rho_{\epsilon \nu}$ each. The settings correspond to different relative performances of the two estimators. In setting A, our estimator is expected to perform worse than Liml with respect to estimation accuracy. Both estimators are expected to show similar performance in setting B, whereas in setting C, bias-adjusted 2SLS-Ridge should outperform Liml. In Table (ref), we report estimation results for 2SLS, standard bias-adjusted 2SLS, Liml, and bias-adjusted 2SLS-Ridge with $\hat{\lambda}^*$ for each of these settings. While all bias-adjusted estimators show little bias, Liml and bias-adjusted 2SLS-Ridge outperform the standard bias-adjusted 2SLS estimator in terms of accuracy in every setting. The relative difference in observed variation between Liml and bias-adjusted 2SLS-Ridge reflects the asymptotic results, although bias-adjusted 2SLS-Ridge performs slightly worse than expected for the case when $\rho_{\epsilon \nu} = 0.6$.

Conclusion

Instrumental variables estimation performs poorly when the number of instruments increases. In these settings, standard IV estimators are biased towards the OLS estimate. Traditional bias-adjustments can be motivated by the Silverstein equation. This link to random matrix theory allows us to improve IV estimation even in settings where no closed-form expression of the Silverstein equation exists. Bias-adjustments are costly in the sense that they reduce the overall signal strength of the instruments. We show that Ridge estimation of the first-stage parameters reduces the implicit price of bias-adjustments. This leads to a trade-off, allowing for less costly estimation of the causal effect which comes along with improved asymptotic properties. Our theoretical results nest existing methods for bias approximation and adjustment with OLS in the first-stage regression. Moreover, the results generalize to settings where the number of instruments is larger than the number of observations, which has not been discussed in the IV literature so far. As a byproduct, we derive the optimal tuning parameter of Ridge regressions in simultaneous equations models, which comprises the well-known results for single equation models as special case with uncorrelated error terms.