# Lewbel (1997)'s Higher Moments IV Approach for Multiplicative Model

I wonder how I can introduce Lewbel (1997)'s higher moments IV approach in a multiplicative / log-log model.

Assume the following linear model: We know that e.g. $$Y_t=5+1 X_{1t}+1 X_{2t}+1X_{3t}+\epsilon_t$$. But we do not observe $$X_{1t}$$ and $$X_{2t}$$ directly but $$Z_{1t}$$ and $$Z_{2t}$$ which are $$X_{1t}$$ and $$X_{2t}$$ with a measurement error $$Z_{1t}=X_{1t}+\phi_t$$ resp. $$Z_{2t}=X_{2t}+\eta_t$$.

If we now run the regression of $$Y$$ on $$Z_{1}$$, $$Z_{2}$$ and $$X_{3}$$ we will get biased estimates for the coefficients of $$Z_{1}$$ and $$Z_{2}$$ since these are now correlated with our error.

This can be easily shown in R:

N <- 1000000 # I set N to 1 Million just to check the consistency of the estimates
X1 <- rlnorm(N,0,1) + 10
X2 <- rlnorm(N,0,1) + 10
X3 <- rlnorm(N,0,1) + 10
Z1 <- X1 + rnorm(N)
Z2 <- X2 + rnorm(N)

Y <- 5 + 1 * X1 + 1 * X2 + 1 * X3 + rnorm(N)

summary(lm(Y ~ Z1 + Z2 + X3))

Call:
lm(formula = Y ~ Z1 + Z2 + X3)

Residuals:
Min     1Q Median     3Q    Max
-7.355 -1.094 -0.022  1.066 26.798

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 9.1508216  0.0144121   634.9   <2e-16 ***
Z1          0.8238412  0.0006830  1206.3   <2e-16 ***
Z2          0.8201544  0.0006892  1190.0   <2e-16 ***
X3          0.9994776  0.0007529  1327.5   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.626 on 999996 degrees of freedom
Multiple R-squared:  0.8221,    Adjusted R-squared:  0.8221
F-statistic: 1.541e+06 on 3 and 999996 DF,  p-value: < 2.2e-16


Lewbel shows that without any assumptions regarding the measurement error distribution we get consistent estimates if we use one or more of the following instruments:

1. $$(G_t - \bar G)$$ where $$G_t$$ is any function of the exogenous regressors so e.g. $$G_t=X_{3t}^3$$ resulting in $$(X_{3t}^3 - \overline {X_{3t}^3})$$
2. $$(G_t - \bar G)(Z_t - \bar Z)$$,
3. $$(G_t - \bar G)(Y_t - \bar Y)$$,
4. $$(Y_t - \bar Y)(Z_t - \bar Z)$$

Since we have two endogenous regressors we can now use 2SLS with e.g. $$(Y_t - \bar Y)(Z_{1t} - \bar Z_1)$$ and $$(Y_t - \bar Y)(Z_{2t} - \bar Z_2)$$ and the exogenous regressor $$X_3$$ as our predictors for $$Z_{1}$$ and $$Z_{2}$$.

instrument1 <- (Y - mean(Y))*(Z1 - mean(Z1))
instrument2 <- (Y - mean(Y))*(Z2 - mean(Z2))

library(ivreg)
summary(ivreg(Y ~ Z1+Z2+X3| instrument1+instrument2+X3))

Call:
ivreg(formula = Y ~ Z1 + Z2 + X3 | instrument1 + instrument2 +
X3)

Residuals:
Min         1Q     Median         3Q        Max
-7.9941360 -1.1703872  0.0004406  1.1662657  8.7512704

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 5.0200059  0.0240095   209.1   <2e-16 ***
Z1          0.9985308  0.0013820   722.5   <2e-16 ***
Z2          0.9995701  0.0012804   780.7   <2e-16 ***
X3          1.0000704  0.0008015  1247.8   <2e-16 ***

Diagnostic tests:
df1   df2 statistic p-value
Weak instruments (Z1) 2e+00 1e+06    191302  <2e-16 ***
Weak instruments (Z2) 2e+00 1e+06    244406  <2e-16 ***
Wu-Hausman            2e+00 1e+06     30562  <2e-16 ***
Sargan                0e+00    NA        NA      NA
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.731 on 999996 degrees of freedom
Multiple R-Squared: 0.7985, Adjusted R-squared: 0.7985
Wald test: 8.901e+05 on 3 and 999996 DF,  p-value: < 2.2e-16


Unfortunately it looks like this approach does not work with a multiplicative model. If we change our data generating process to $$Y_t=e^5X_{1t}X_{2t}X_{3t}e^{\epsilon_t}$$ and still assume that $$Z_{1t}=X_{1t}+\phi_t$$ and $$Z_{2t}=X_{2t}+\eta_t$$ we get biased estimates in the OLS estimation again (of course)

N <- 1000000
X1 <- rlnorm(N,0,1) + 10
X2 <- rlnorm(N,0,1) + 10
X3 <- rlnorm(N,0,1) + 10
Z1 <- X1 + rnorm(N)
Z2 <- X2 + rnorm(N)

Y <- exp(5)* X1 * X2 * X3 * exp(rnorm(N))

Call:
lm(formula = log(Y) ~ log(Z1) + log(Z2) + log(X3))

Residuals:
Min      1Q  Median      3Q     Max
-4.7435 -0.6778  0.0011  0.6795  4.5768

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.342661   0.026737   237.2   <2e-16 ***
log(Z1)     0.725010   0.005936   122.1   <2e-16 ***
log(Z2)     0.724755   0.005942   122.0   <2e-16 ***
log(X3)     1.003215   0.007018   142.9   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.005 on 999996 degrees of freedom
Multiple R-squared:  0.04782,   Adjusted R-squared:  0.04781
F-statistic: 1.674e+04 on 3 and 999996 DF,  p-value: < 2.2e-16


If I now create the instruments according to our new model as $$(ln(Y_t) - \overline {ln(Y)})(ln(Z_{1t}) - \overline {ln(Z_1)})$$ and $$(ln(Y_t) - \overline {ln(Y)})(ln(Z_{2t}) - \overline {ln(Z_2)})$$ and run the 2SLS regression I get indeed a lot lot better but still somehow biased estimates.

instrument1 <- (log(Y) - mean(log(Y)))*(log(Z1) - mean(log(Z1)))
instrument2 <- (log(Y) - mean(log(Y)))*(log(Z2) - mean(log(Z2)))
summary(ivreg(log(Y) ~ log(Z1)+log(Z2)+log(X3)| instrument1+instrument2+log(X3)))

Call:
ivreg(formula = log(Y) ~ log(Z1) + log(Z2) + log(X3) | instrument1 +
instrument2 + log(X3))

Residuals:
Min         1Q     Median         3Q        Max
-4.7072431 -0.6797044  0.0001246  0.6815530  4.7018362

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept)  4.73916    0.08389   56.49   <2e-16 ***
log(Z1)      1.04275    0.02347   44.43   <2e-16 ***
log(Z2)      1.06425    0.02398   44.38   <2e-16 ***
log(X3)      1.00338    0.00704  142.53   <2e-16 ***

Diagnostic tests:
df1   df2 statistic p-value
Weak instruments (log(Z1)) 2e+00 1e+06   34397.2  <2e-16 ***
Weak instruments (log(Z2)) 2e+00 1e+06   32928.5  <2e-16 ***
Wu-Hausman                 2e+00 1e+06     204.8  <2e-16 ***
Sargan                     0e+00    NA        NA      NA
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.008 on 999996 degrees of freedom
Multiple R-Squared: 0.04198,    Adjusted R-squared: 0.04198
Wald test:  8073 on 3 and 999996 DF,  p-value: < 2.2e-16


Compared to the OLS estimates with the true regressors $$X_1$$ and $$X_2$$ it looks like there is still a bias of around .05 for the 2SLS estimates for $$Z_1$$ and $$Z_2$$.

Call:
lm(formula = log(Y) ~ log(X1) + log(X2) + log(X3))

Residuals:
Min      1Q  Median      3Q     Max
-4.7467 -0.6739  0.0009  0.6754  4.5881

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.964183   0.029565   167.9   <2e-16 ***
log(X1)     1.008105   0.006988   144.3   <2e-16 ***
log(X2)     1.004447   0.006997   143.6   <2e-16 ***
log(X3)     1.002263   0.006979   143.6   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9997 on 999996 degrees of freedom
Multiple R-squared:  0.05851,   Adjusted R-squared:  0.0585
F-statistic: 2.071e+04 on 3 and 999996 DF,  p-value: < 2.2e-16


Is there anything I am doing wrong here? Thanks!

I actually found the problem. The higher moments IV approach requires the endogenous regressors to be non-normal. In my example above, the logs of $$Z_1$$ and $$Z_2$$ are almost normal. If I change rlnorm(N,0,1) + 10 in the data generating process to e.g. rlnorm(N,0,1.5) + 10 also the logs of $$Z_1$$ and $$Z_2$$ remain non-normal and the IVs for the log-log model work well.