I would like to implement (in R) an instrumental variable (IV) estimator, that takes the most general form (here not 2SLS or GMM!): $$ \beta_{IV} = (Z'X)^{-1}Z'Y $$

I could code this in the naive way, inverting $Z'X$, and computing the rest, but I remember that this is very inefficient in term of computing, and that in the similar case of OLS $(X'X)^{-1}X'Y$, one should use rather a QR decomposition rather than inverting $X'X$.

However, in the case of my IV estimator, I guess the standard QR won't apply, so I am not sure whether there is any similar technique/decomposition that could be applied in this case?

Any idea of an efficient technique to compute $(Z'X)^{-1}Z'Y$ without inverting $(Z'X)^{-1}$?


  • $\begingroup$ Please expand your use of the abbreviation "IV" at first use. It appears to mean instrumental variable in this case - you should say so explicitly rather than leaving people to guess. $\endgroup$
    – Glen_b
    Nov 27 '13 at 21:23
  • $\begingroup$ It looks to me like QR decomposition would be a good way to proceed. Alternatively, if that doesn't work for some reason, singular value decomposition might do. $\endgroup$
    – Glen_b
    Nov 27 '13 at 21:47
  • $\begingroup$ Indeed, Part 5 of this document seems to outline just such a calculation. $\endgroup$
    – Glen_b
    Nov 27 '13 at 21:59
  • $\begingroup$ Great, exactly what I was looking for, thanks a lot Glen_b!!! And indeed, seems to confirm that this can be done with a QR! $\endgroup$
    – Matifou
    Nov 28 '13 at 18:15

OLS and IV estimators are both solutions to systems of equations, which can be computed using QR decompositions.

Linear regression model

The linear regression model is $$ \begin{align} \boldsymbol{Y} &= \mathbf{X}\boldsymbol{\beta} + \boldsymbol{\varepsilon}\\ \mathbb{E}(\boldsymbol{\varepsilon} \mid \mathbf{X}) &= \boldsymbol{0} \end{align} $$ The OLS estimator ($\hat{\boldsymbol{\beta}}_{OLS}$) is the solution to the system of equations: $$ \left(\mathbf{X}'\mathbf{X}\right)\hat{\boldsymbol{\beta}}_{OLS} = \mathbf{X}'\boldsymbol{Y} $$

This computation can be easily performed using a call to qr.solve. I have added formula handling because that makes the interface neat.

# download and pre-process the data

# download the Mroz dataset
download.file(url = 'http://people.stern.nyu.edu/wgreene/Text/tables/TableF4-1.txt', 
              destfile = 'mroz.txt')

# read in the file
dfMroz = read.table('mroz.txt', header = TRUE, skip = 36)
names(dfMroz) = tolower(names(dfMroz))

# linear regression of wife's wages
# using lm()
lmMroz = lm(log(ww) ~ we + ax. + I(ax.^2), 
             data = dfMroz[dfMroz$ww > 0, ])

# computation by hand
## define the formula on which the OLS estimate is based
formulaMrozLM = log(ww) ~ we + ax. + I(ax.^2)

## get the model matrices
mfMrozLM = model.frame(formulaMrozLM, data = dfMroz, subset = ww > 0)
vYMrozLM = model.response(mfMrozLM)
mXMrozLM = model.matrix(formulaMrozLM, data = mfMrozLM)

# solve using solve() and qr.solve()
# solve(crossprod(mXMrozLM), crossprod(mXMrozLM, vYMrozLM))
vBetaLMHat = qr.solve(mXMrozLM, vYMrozLM)

Instrumental variables model

The instrumental variables model is $$ \begin{align} \boldsymbol{Y} &= \mathbf{X}\boldsymbol{\beta} + \boldsymbol{\varepsilon}\\ \mathbf{X} &= \mathbf{Z}\mathbf{\Pi} + \mathbf{V}\\ \mathbb{E}(\boldsymbol{\varepsilon} \mid \mathbf{Z}) &= \boldsymbol{0} \end{align} $$ where I have included the instrument exogeneity condition, but excluded the relevance condition.

The IV estimator for the case that there are as many instruments as there are endogenous variables, is defined as the solution to the system of equations

$$ \left(\mathbf{Z}'\mathbf{X}\right)\hat{\boldsymbol{\beta}}_{IV} = \mathbf{Z}'\boldsymbol{Y} $$

In the general case, where there are more instrumental variables than endogenous variables in the system, we have the 2SLS estimator which is an IV estimator, for $\mathbf{Z} = \hat{\mathbf{X}}$, where $\hat{\mathbf{X}}$ is the linear projection of the matrix $\mathbf{X}$ onto the columns of the matrix $\mathbf{Z}$, that is, $\hat{\mathbf{X}} = \mathbf{Z}\hat{\mathbf{\Pi}}$. This projection is computed by way of a solution of the auxiliary system of equations: $$ \left(\mathbf{Z}'\mathbf{Z}\right)\hat{\mathbf{\Pi}} = \mathbf{Z}'\mathbf{X} $$

# intrumental variables regression of wife's wages
# computing using ivreg() in AER
ivMroz = ivreg(log(ww) ~ ax. + I(ax.^2) + kl6 + k618 + we |
                 ax. + I(ax.^2) + kl6 + k618 + wmed + wfed + he, 
               data = dfMroz[dfMroz$ww > 0, ])

# NOTE this is provided by the 'Formula' package of Achim Zeileis
#   and this is used for brevity only, and these matrices can 
#   easily created by hand
formulaMrozIV = as.Formula(log(ww) ~ ax. + I(ax.^2) + kl6 + k618 + we |
                        ax. + I(ax.^2) + kl6 + k618 + wmed + wfed + he)
mfMrozIV = model.frame(formulaMrozIV, data = dfMroz, subset = ww > 0)
vYMrozIV = model.response(mfMrozIV)
mXMrozIV = model.matrix(formulaMrozIV, data = mfMrozIV, rhs = 1)
mZMrozIV = model.matrix(formulaMrozIV, data = mfMrozIV, rhs = 2)

# compute the mXHat
# NOTE in the case of the pure IV estimator, where the number 
#   of endogenous variables is the same as the number of exogenous
#   variables, this computation is unnecessary
mXHat = mZMrozIV %*% qr.solve(mZMrozIV, mXMrozIV)

# compute the 2SLS (IV) estimator
vBetaIVHat = qr.solve(mXHat, vYMrozIV)

The IV estimator in your equation seems hard to understand. The length of the instrument vector $Z$ may not equal that of the regressor vector $X$. If they are not equal, the matrix $Z'X$ is not a square matrix; its inversion would be ill defined.

One possible alternative is to introduce a positive-definite weighting matrix $S$ such that $\beta_{IV} = ((X'Z)S(Z'X))^{-1} (X'Z)S(Z'y)$. This returns to the GMM framework. The weighting matrix $S$ could be an identity matrix or a more efficient alternative. In this case, the inversion of the matrix could be done through the usual QR decomposition. I would suggest the textbook "Econometrics" by Fumio Hayashi, which has a more detailed description of the IV estimators. Good luck!


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.