# install.packages("leaps")
# install.packages("ivreg")
library("leaps")
library("ivreg")
data = read.csv("http://hedibert.org/wp-content/uploads/2024/11/card.csv",header=TRUE)
attach(data)
n = nrow(data)
OLS regression with
several exogenous variables
outcome = lwage
treatment = educ
W = cbind(exper,expersq,black,smsa,south,smsa66,
reg662,reg663,reg664,reg665,reg666,reg667,reg668,reg669)
X = cbind(treatment,W)
reg.ols = lm(outcome~X)
summary(reg.ols)
##
## Call:
## lm(formula = outcome ~ X)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.62326 -0.22141 0.02001 0.23932 1.33340
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.6208069 0.0742327 62.248 < 2e-16 ***
## Xtreatment 0.0746933 0.0034983 21.351 < 2e-16 ***
## Xexper 0.0848320 0.0066242 12.806 < 2e-16 ***
## Xexpersq -0.0022870 0.0003166 -7.223 6.41e-13 ***
## Xblack -0.1990123 0.0182483 -10.906 < 2e-16 ***
## Xsmsa 0.1363846 0.0201005 6.785 1.39e-11 ***
## Xsouth -0.1479550 0.0259799 -5.695 1.35e-08 ***
## Xsmsa66 0.0262417 0.0194477 1.349 0.17733
## Xreg662 0.0963671 0.0358979 2.684 0.00730 **
## Xreg663 0.1445400 0.0351244 4.115 3.97e-05 ***
## Xreg664 0.0550756 0.0416573 1.322 0.18623
## Xreg665 0.1280248 0.0418395 3.060 0.00223 **
## Xreg666 0.1405174 0.0452469 3.106 0.00192 **
## Xreg667 0.1179810 0.0448025 2.633 0.00850 **
## Xreg668 -0.0564361 0.0512579 -1.101 0.27098
## Xreg669 0.1185697 0.0388301 3.054 0.00228 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3723 on 2994 degrees of freedom
## Multiple R-squared: 0.2998, Adjusted R-squared: 0.2963
## F-statistic: 85.48 on 15 and 2994 DF, p-value: < 2.2e-16
IV regression with
several exogenous variables
iv = nearc4
Z = cbind(iv,W)
reg.iv = ivreg(outcome ~ X | Z)
summary(reg.iv)
##
## Call:
## ivreg(formula = outcome ~ X | Z)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.83164 -0.24075 0.02428 0.25208 1.42760
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.6661519 0.9248295 3.964 7.54e-05 ***
## Xtreatment 0.1315038 0.0549637 2.393 0.01679 *
## Xexper 0.1082711 0.0236586 4.576 4.92e-06 ***
## Xexpersq -0.0023349 0.0003335 -7.001 3.12e-12 ***
## Xblack -0.1467758 0.0538999 -2.723 0.00650 **
## Xsmsa 0.1118084 0.0316620 3.531 0.00042 ***
## Xsouth -0.1446715 0.0272846 -5.302 1.23e-07 ***
## Xsmsa66 0.0185311 0.0216086 0.858 0.39119
## Xreg662 0.1007678 0.0376857 2.674 0.00754 **
## Xreg663 0.1482588 0.0368141 4.027 5.78e-05 ***
## Xreg664 0.0498971 0.0437398 1.141 0.25406
## Xreg665 0.1462719 0.0470639 3.108 0.00190 **
## Xreg666 0.1629029 0.0519096 3.138 0.00172 **
## Xreg667 0.1345722 0.0494023 2.724 0.00649 **
## Xreg668 -0.0830770 0.0593313 -1.400 0.16155
## Xreg669 0.1078142 0.0418137 2.578 0.00997 **
##
## Diagnostic tests:
## df1 df2 statistic p-value
## Weak instruments 1 2994 13.256 0.000276 ***
## Wu-Hausman 1 2993 1.168 0.279973
## Sargan 0 NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3883 on 2994 degrees of freedom
## Multiple R-Squared: 0.2382, Adjusted R-squared: 0.2343
## Wald test: 51.01 on 15 and 2994 DF, p-value: < 2.2e-16
OLS regression with
selection of exogenous variables
X1 = X
data = data.frame(y=outcome,X=X1)
p = ncol(X1)
regs = regsubsets(y~X,data=data,method="exhaustive",nvmax=p)
bic = summary(regs)$bic
ii = 1:p
coef(regs,ii[bic==min(bic)])
## (Intercept) Xtreatment Xexper Xexpersq Xblack Xsmsa
## 4.713907129 0.074619220 0.084431565 -0.002270828 -0.190925636 0.157572489
## Xsouth Xreg663 Xreg668
## -0.116091147 0.054617876 -0.143638697
ols = lm(outcome ~ treatment+exper+expersq+black+smsa+south+reg663+reg668)
summary(ols)
##
## Call:
## lm(formula = outcome ~ treatment + exper + expersq + black +
## smsa + south + reg663 + reg668)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.58335 -0.22424 0.01746 0.24290 1.34410
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.7139071 0.0676384 69.693 < 2e-16 ***
## treatment 0.0746192 0.0034959 21.344 < 2e-16 ***
## exper 0.0844316 0.0066260 12.742 < 2e-16 ***
## expersq -0.0022708 0.0003168 -7.169 9.46e-13 ***
## black -0.1909256 0.0175808 -10.860 < 2e-16 ***
## smsa 0.1575725 0.0155371 10.142 < 2e-16 ***
## south -0.1160911 0.0159501 -7.278 4.29e-13 ***
## reg663 0.0546179 0.0183605 2.975 0.00296 **
## reg668 -0.1436387 0.0418657 -3.431 0.00061 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3728 on 3001 degrees of freedom
## Multiple R-squared: 0.2961, Adjusted R-squared: 0.2942
## F-statistic: 157.8 on 8 and 3001 DF, p-value: < 2.2e-16
IV regression with
selection of exogenous variables
data = data.frame(y=outcome,Z=Z)
p = ncol(data)
regs = regsubsets(y~Z,data=data,method="exhaustive",nvmax=p,force.in=Z[,1])
bic = summary(regs)$bic
ii = 1:p
coef(regs,ii[bic==min(bic)])
## (Intercept) Ziv Zexper Zexpersq Zblack Zsmsa
## 5.93718986 0.04594223 0.05358652 -0.00220199 -0.26286597 0.18248462
## Zsouth Zreg663
## -0.12741608 0.06092994
reg.iv = ivreg(outcome ~ treatment+exper+expersq+black+smsa+south+reg663 |
iv+exper+expersq+black+smsa+south+reg663)
summary(reg.iv)
##
## Call:
## ivreg(formula = outcome ~ treatment + exper + expersq + black +
## smsa + south + reg663 | iv + exper + expersq + black + smsa +
## south + reg663)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.83851 -0.23743 0.02038 0.25642 1.46049
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.6610748 0.8361146 4.379 1.23e-05 ***
## treatment 0.1365294 0.0495914 2.753 0.005939 **
## exper 0.1095925 0.0214619 5.106 3.49e-07 ***
## expersq -0.0023032 0.0003357 -6.860 8.33e-12 ***
## black -0.1254105 0.0532781 -2.354 0.018643 *
## smsa 0.1271510 0.0303902 4.184 2.95e-05 ***
## south -0.0862605 0.0241822 -3.567 0.000367 ***
## reg663 0.0661039 0.0193090 3.423 0.000627 ***
##
## Diagnostic tests:
## df1 df2 statistic p-value
## Weak instruments 1 3002 16.621 4.68e-05 ***
## Wu-Hausman 1 3001 1.762 0.185
## Sargan 0 NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3928 on 3002 degrees of freedom
## Multiple R-Squared: 0.2184, Adjusted R-squared: 0.2166
## Wald test: 104.1 on 7 and 3002 DF, p-value: < 2.2e-16
LS0tCnRpdGxlOiAiR2F1c3NpYW4gbGluZWFyIHJlZ3Jlc3Npb24gd2l0aCBJViIKc3VidGl0bGU6ICJvdXRjb21lPWxvZ3dhZ2UgLSB0cmVhdG1lbnQ9ZWR1YyAtIGl2PW5lYXJjNCIKYXV0aG9yOiAiSGVkaWJlcnQgRnJlaXRhcyBMb3BlcyIKZGF0ZTogIjE1LzExLzIwMjUiCm91dHB1dDoKICBodG1sX2RvY3VtZW50OgogICAgdG9jOiB0cnVlCiAgICB0b2NfZGVwdGg6IDMKICAgIHRvY19jb2xsYXBzZWQ6IHRydWUKICAgIGNvZGVfZG93bmxvYWQ6IHRydWUKICAgIG51bWJlcl9zZWN0aW9uczogdHJ1ZQogIHBkZl9kb2N1bWVudDoKICAgIHRvYzogdHJ1ZQogICAgdG9jX2RlcHRoOiAnMycKLS0tCgpgYGB7ciBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGVjaG8gPSBUUlVFKQpgYGAKCgpgYGB7ciBmaWcuYWxpZ249J2NlbnRlcicsIGZpZy53aWR0aD04LCBmaWcuaGVpZ2h0PTZ9CiMgaW5zdGFsbC5wYWNrYWdlcygibGVhcHMiKSAKIyBpbnN0YWxsLnBhY2thZ2VzKCJpdnJlZyIpCmxpYnJhcnkoImxlYXBzIikKbGlicmFyeSgiaXZyZWciKQoKZGF0YSA9IHJlYWQuY3N2KCJodHRwOi8vaGVkaWJlcnQub3JnL3dwLWNvbnRlbnQvdXBsb2Fkcy8yMDI0LzExL2NhcmQuY3N2IixoZWFkZXI9VFJVRSkgCgphdHRhY2goZGF0YSkgCgpuID0gbnJvdyhkYXRhKSAKYGBgCgoKIyBPTFMgcmVncmVzc2lvbiB3aXRoIHNldmVyYWwgZXhvZ2Vub3VzIHZhcmlhYmxlcwoKYGBge3J9Cm91dGNvbWUgICA9IGx3YWdlCnRyZWF0bWVudCA9IGVkdWMKVyAgICAgICAgID0gY2JpbmQoZXhwZXIsZXhwZXJzcSxibGFjayxzbXNhLHNvdXRoLHNtc2E2NiwKICAgICAgICAgICAgICAgICAgcmVnNjYyLHJlZzY2MyxyZWc2NjQscmVnNjY1LHJlZzY2NixyZWc2NjcscmVnNjY4LHJlZzY2OSkKWCAgICAgICAgID0gY2JpbmQodHJlYXRtZW50LFcpCnJlZy5vbHMgICA9IGxtKG91dGNvbWV+WCkgCnN1bW1hcnkocmVnLm9scykKYGBgCgojIElWIHJlZ3Jlc3Npb24gd2l0aCBzZXZlcmFsIGV4b2dlbm91cyB2YXJpYWJsZXMKCmBgYHtyfQppdiAgICAgPSBuZWFyYzQKWiAgICAgID0gY2JpbmQoaXYsVykKcmVnLml2ID0gaXZyZWcob3V0Y29tZSB+IFggfCBaKQpzdW1tYXJ5KHJlZy5pdikKYGBgCgoKIyBPTFMgcmVncmVzc2lvbiB3aXRoIHNlbGVjdGlvbiBvZiBleG9nZW5vdXMgdmFyaWFibGVzCgpgYGB7cn0KWDEgICA9IFgKZGF0YSA9IGRhdGEuZnJhbWUoeT1vdXRjb21lLFg9WDEpCnAgICAgPSBuY29sKFgxKQpyZWdzID0gcmVnc3Vic2V0cyh5flgsZGF0YT1kYXRhLG1ldGhvZD0iZXhoYXVzdGl2ZSIsbnZtYXg9cCkKYmljID0gc3VtbWFyeShyZWdzKSRiaWMKaWkgID0gMTpwCmNvZWYocmVncyxpaVtiaWM9PW1pbihiaWMpXSkKb2xzID0gbG0ob3V0Y29tZSB+IHRyZWF0bWVudCtleHBlcitleHBlcnNxK2JsYWNrK3Ntc2Erc291dGgrcmVnNjYzK3JlZzY2OCkKc3VtbWFyeShvbHMpCmBgYAoKIyBJViByZWdyZXNzaW9uIHdpdGggc2VsZWN0aW9uIG9mIGV4b2dlbm91cyB2YXJpYWJsZXMKCmBgYHtyfQpkYXRhID0gZGF0YS5mcmFtZSh5PW91dGNvbWUsWj1aKQpwICAgID0gbmNvbChkYXRhKQpyZWdzID0gcmVnc3Vic2V0cyh5flosZGF0YT1kYXRhLG1ldGhvZD0iZXhoYXVzdGl2ZSIsbnZtYXg9cCxmb3JjZS5pbj1aWywxXSkKYmljID0gc3VtbWFyeShyZWdzKSRiaWMKaWkgID0gMTpwCmNvZWYocmVncyxpaVtiaWM9PW1pbihiaWMpXSkKCnJlZy5pdiA9IGl2cmVnKG91dGNvbWUgfiB0cmVhdG1lbnQrZXhwZXIrZXhwZXJzcStibGFjaytzbXNhK3NvdXRoK3JlZzY2MyB8IAogICAgICAgICAgICAgICAgICAgICAgICAgaXYrZXhwZXIrZXhwZXJzcStibGFjaytzbXNhK3NvdXRoK3JlZzY2MykKc3VtbWFyeShyZWcuaXYpCmBgYAo=