#####################################
#Lab 3 Instrumental variables Stat 209 2/15/19
Rogosa R-Session
#####################################
In this session output I did not duplicate all the explanation
in the lab text, but did all the basic analyses and added some notes.
So these two versions should be looked at side by side
First run through I use the traditional tsls function from the sem package.
Then at the end I repeat with the (Stata-phonic) ivreg command from the
newer AER package.
Plus a bit of Dylan Small's ivpack (robust s.e., improved CI).
##### Examples from Wooldridge, Introductory Econometrics
####### Chapters 15 and 16
## stata results available from
http://fmwww.bc.edu/gstat/examples/wooldridge/wooldridge15.html
http://fmwww.bc.edu/gstat/examples/wooldridge/wooldridge16.html
> # Lab 3
> mroz87 = read.table( "http://statweb.stanford.edu/~rag/stat209/Mroz87.dat", header = T)
> names(mroz87)
[1] "lfp" "hours" "kids5" "kids618" "age" "educ"
[7] "wage" "repwage" "hushrs" "husage" "huseduc" "huswage"
[13] "faminc" "mtr" "motheduc" "fatheduc" "unem" "city"
[19] "exper" "nwifeinc" "wifecoll" "huscoll"
> # Variable definition in lab script also linked, main page
# lfp Dummy variable for labor-force participation.
# hours Wife's hours of work in 1975.
# kids5 Number of children 5 years old or younger.
# kids618 Number of children 6 to 18 years old.
# age Wife's age.
# educ Wife's educational attainment, in years.
# wage Wife's average hourly earnings, in 1975 dollars.
# repwage Wife's wage reported at the time of the 1976 interview.
# hushrs Husband's hours worked in 1975.
# husage Husband's age.
# huseduc Husband's educational attainment, in years.
# huswage Husband's wage, in 1975 dollars.
# faminc Family income, in 1975 dollars.
# mtr Marginal tax rate facing the wife.
# motheduc Wife's mother's educational attainment, in years.
# fatheduc Wife's father's educational attainment, in years.
# unem Unemployment rate in county of residence, in percentage points.
# city Dummy variable = 1 if live in large city, else 0.
# exper Actual years of wife's previous labor market experience.
# nwifeinc Non-wife income.
# wifecoll Dummy variable for wife's college attendance.
# huscoll Dummy variable for husband's college attendance.
> library(sem)
> # I already had sem installed; if not need to install.packages("sem") or "AER" to do this with ivreg; see Appendix
> help(tsls)
> # tsls is traditional tool for IV; ivreg from AER the newer alternative (see appendix at end of lab)
> # outcome variable is log-wage for the 428 working women
> attach(mroz87)
> table(lfp)
lfp
0 1
325 428
> table(wage > 0)
FALSE TRUE
325 428
> detach(mroz87)
> poswage = subset(mroz87, wage > 0) # my new data subset
> poswage$logWage = log( poswage$wage ) # adding logWage to the data set for session
> attach(poswage)
> names(poswage)
[1] "lfp" "hours" "kids5" "kids618" "age" "educ"
[7] "wage" "repwage" "hushrs" "husage" "huseduc" "huswage"
[13] "faminc" "mtr" "motheduc" "fatheduc" "unem" "city"
[19] "exper" "nwifeinc" "wifecoll" "huscoll" "logWage"
> length(logWage)
[1] 428
> table(lfp) #all the women in this subset are in the workforce
lfp
1
428
> # so the poswage data subset is the 428 working women, and appended
> # to that subset is logWage
# see what you have by take log(wage)
> logWage[1:10]
[1] 1.21015366 0.32851207 1.51413773 0.09212329 1.52427210 1.55648008
[7] 2.12025954 2.05963416 0.75433635 1.54489939
> wage[1:10]
[1] 3.3540 1.3889 4.5455 1.0965 4.5918 4.7421 8.3333 7.8431 2.1262 4.6875
> log(3.354)
[1] 1.210154
> min(logWage)
[1] -2.054164
> min(wage) #not paid so well
[1] 0.1282
> quantile(wage)
0% 25% 50% 75% 100%
0.12820 2.26260 3.48190 4.97075 25.00000
> table(wage < 1)
FALSE TRUE
411 17
> # so we have a few very poorly paid women in this dataset it seems < 1$/hr
> # onto fitting regression (predictor educ)
> lm.posWage = lm(logWage ~ educ)
> summary(lm.posWage)
Call:
lm(formula = logWage ~ educ)
Residuals:
Min 1Q Median 3Q Max
-3.10256 -0.31473 0.06434 0.40081 2.10029
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.1852 0.1852 -1.000 0.318
educ 0.1086 0.0144 7.545 2.76e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.68 on 426 degrees of freedom
Multiple R-Squared: 0.1179, Adjusted R-squared: 0.1158
F-statistic: 56.93 on 1 and 426 DF, p-value: 2.761e-13
> # we get a highly significant slope (but not big Rsq),
> # year increase in educ, fit increases 11 cents hourly wage
> cor(logWage,educ)^2 # R-squared for the OLS equation
[1] 0.1178826
> # now the IV fit using fatheduc as instrument (omitted vars concern)
> cor.test(educ, fatheduc)
Pearson's product-moment correlation
data: educ and fatheduc
t = 9.4255, df = 426, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.3337579 0.4908623
sample estimates:
cor
0.4154030
> # significant, fairly strong (moderate at least) correlation educ and instrument (IV s.e. inflated by ~2.4)
> instr.posWage = tsls(logWage ~ educ, ~ fatheduc)
> summary(instr.posWage)
2SLS Estimates
Model Formula: logWage ~ educ
Instruments: ~fatheduc
Residuals:
Min. 1st Qu. Median Mean 3rd Qu. Max.
-3.09e+00 -3.39e-01 5.25e-02 -1.99e-14 4.04e-01 2.07e+00
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.44110 0.44610 0.9888 0.32332
educ 0.05917 0.03514 1.6839 0.09294
Residual standard error: 0.6894 on 426 degrees of freedom
> # educ slope reduced by half, se increased, no longer signif
> # you can show by math (with the regression recursion) and I mentioned in class when I showed the
> # Woolridge textbook results for these data, that if the omitted variable(s)are
> # positively related to predictor, the OLS regression will overestimate the slope
> # easy to see from the week 1, regression recursion relations sheet, example is linked week 6
> cov(logWage, fatheduc)/cov(educ, fatheduc) # IV estimate of slope "by hand"
[1] 0.05917348
> .0144/.4154 # inflation of OLS standard error for the IV regression
[1] 0.03466538
> #matches IV s.e.
> # do two stage least squares as 2 stages
> lm.1 = lm(educ ~ fatheduc)
> lm.2 = lm(logWage ~ fitted(lm.1))
> summary(lm.2)
Call:
lm(formula = logWage ~ fitted(lm.1))
Residuals:
Min 1Q Median 3Q Max
-3.21264 -0.37631 0.05632 0.41728 2.06040
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.44110 0.46711 0.944 0.346
fitted(lm.1) 0.05917 0.03680 1.608 0.109
Residual standard error: 0.7219 on 426 degrees of freedom
Multiple R-Squared: 0.006034, Adjusted R-squared: 0.003701
F-statistic: 2.586 on 1 and 426 DF, p-value: 0.1086
> # matches the output from tsls pretty closely
#### more on R-square. as we saw in the Woolridge stata output and you can see
# from ivreg in appendix reported R-square from the IV fit is .0934, smaller
# than the OLS R-square. You can compute that by 1 - SSresiduals/SSoutcome
# for the IV fit. I do that by hand below (to answer a prior class question).
# Intuitively, the instrument is not as good a predictor of outcome as
# the (flawed) OLS predictor. The ivreg function does give this R-square (see Appendix).
# > residuals(instr.posWage)
# [1] 0.058968499 -0.822673098 0.362952568 -1.059061876 0.254739979 0.405294911 0.732380450 0.908449000
# [9] -0.396848816 0.393714227 0.250736453 0.432260420 -0.417231990 -0.332816019 0.269993024 -0.793983319
# .........
# [417] 0.464990512 -0.360203911 -0.234894433 0.408927460 -0.206501358 0.090083425 0.117931765 -0.194811762
# [425] 0.517671938 0.559070004 0.075263162 0.255303904
# > mr = mean(residuals(instr.posWage))
# > mr
# [1] 3.956327e-14
# > ssres = sum((residuals(instr.posWage) - mr)^2)
# > ssres
# [1] 202.4601
# > ssposw = sum((logWage - mean(logWage))^2)
# > 1 - ssres/ssposw #matches ivreg output
# [1] 0.09343841
>
> # move on to Task 2, adding exper to the prediction of logWage (use exper and exper^2)
> # exper is clearly a possible important variable in predicting wage
> lm.posWage2 = lm(logWage ~ educ + exper + I(exper^2))
# I made more typing for myself by not creating a separate variable for exper^2
# here's the OLS estimation, adding experience to the prediction of wage
> summary(lm.posWage2)
Call:
lm(formula = logWage ~ educ + exper + I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-3.08404 -0.30627 0.04952 0.37498 2.37115
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.5220406 0.1986321 -2.628 0.00890 **
educ 0.1074896 0.0141465 7.598 1.94e-13 ***
exper 0.0415665 0.0131752 3.155 0.00172 **
I(exper^2) -0.0008112 0.0003932 -2.063 0.03974 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.6664 on 424 degrees of freedom
Multiple R-Squared: 0.1568, Adjusted R-squared: 0.1509
F-statistic: 26.29 on 3 and 424 DF, p-value: 1.302e-15
> # R-sq a little higher using exper; educ coeff almost identical to single predictor eq
> # but the regression may also have omitted variable bias
> # textbooks use use motheduc and fatheduc as exogenous instruments for this regression
# since exper is also regarded as exogenous, exper and exper^2 are included as instruments
by Woolridge in his stata output, the main IV emphasis is on fatheduc, motheduc as
instruments for educ, exper is just put in seeking better prediction (Rsq still measily .13)
> instr.posWage2 = tsls(logWage ~ educ + exper + I(exper^2), ~fatheduc + motheduc + exper + I(exper^2))
> summary(instr.posWage2)
2SLS Estimates
Model Formula: logWage ~ educ + exper + I(exper^2)
Instruments: ~fatheduc + motheduc + exper + I(exper^2)
Residuals:
Min. 1st Qu. Median Mean 3rd Qu. Max.
-3.10e+00 -3.20e-01 5.51e-02 1.87e-15 3.69e-01 2.35e+00
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.048100 0.4003281 0.1202 0.904419
educ 0.061397 0.0314367 1.9530 0.051474
exper 0.044170 0.0134325 3.2883 0.001092
I(exper^2) -0.000899 0.0004017 -2.2380 0.025740
Residual standard error: 0.6747 on 424 degrees of freedom
> # once again using IV for omitted vars cuts down the return to educ estimate, also more than doubles se of coeff
> # you can check these results against textbook (Wooldridge) versions using stata
> # go to http://fmwww.bc.edu/gstat/examples/wooldridge/wooldridge15.html Ex 15.5 (also 15.1)
>
> ########################################################################################
> ########################################################################################
> # onto Task 3 simultaneous eqs: more complex modelling, we have a supply and demand eq
> # supply, outcome hours with logWage educ age kidslt6 nwifeinc as predictors
> # demand, logwage outcome, with hours educ and exper (Task 2) as predictors
> # these are regarded as linked eqs as outcome from supply is predictor in demand
## eqs written out in week 6 class handout on simultaneity and reciprocal effects
> # so here we have linked eqs predictor in one is the outcome in the other, hours wages
> # IV TSLS allows estimation of both, you can mess around with reduced form etc,
> # A guide to reduced form, estimation etc at http://elsa.berkeley.edu/eml/ra_reader/18-simul.pdf
> # Preliminaries, test for distinct elements in each equations set of instruments
# scan of Wooldridge p.562 at http://www-stat.stanford.edu/~rag/stat209/woolp562.pdf
> # do the nested F-tests for identifiability [also see Sargan test in ivreg]
> lm.hours = lm( hours ~ educ + age + kids5 + nwifeinc + exper + I(exper^2))
> summary(lm.hours)
Call:
lm(formula = hours ~ educ + age + kids5 + nwifeinc + exper +
I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-1521.23 -556.81 63.86 448.71 3310.63
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 1691.0312 312.4039 5.413 1.04e-07 ***
educ -17.7573 16.3887 -1.084 0.279205
age -15.7096 5.6871 -2.762 0.005991 **
kids5 -294.8382 96.8761 -3.043 0.002485 **
nwifeinc -0.1066 3.6261 -0.029 0.976570
exper 51.7395 14.5003 3.568 0.000401 ***
I(exper^2) -0.5758 0.4390 -1.312 0.190374
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 729.8 on 421 degrees of freedom
Multiple R-squared: 0.1287, Adjusted R-squared: 0.1163
F-statistic: 10.36 on 6 and 421 DF, p-value: 1.014e-10
> lm.hoursSub = lm( hours ~ educ + age + kids5 + nwifeinc )
# I wrote out the submodel, Karen's use of update is more elegant
> anova(lm.hoursSub, lm.hours)
Analysis of Variance Table
Model 1: hours ~ educ + age + kids5 + nwifeinc
Model 2: hours ~ educ + age + kids5 + nwifeinc + exper + I(exper^2)
Res.Df RSS Df Sum of Sq F Pr(>F)
1 423 248021309
2 421 224200428 2 23820881 22.365 5.875e-10 ***
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
> # Not surprisingly, since exper was signif in the full model, we can reject the composite hypoth that both coeff = 0
> # many ways to test the composite hypoth (from elementary regression class, increment-to-Rsq is equivalent)
> # same process for second eq identifiability
# Now go to the direct IV estimation, using all exogenous vars as instruments
# Results match Wooldridge text, see also
# http://fmwww.bc.edu/gstat/examples/wooldridge/wooldridge16.html, Ex 16.5
> sem.hours = tsls( hours ~ logWage + educ + age + kids5 + nwifeinc, instruments =~ educ + age + kids5 + nwifeinc + exper + I(exper^2) )
> summary(sem.hours)
2SLS Estimates
Model Formula: hours ~ logWage + educ + age + kids5 + nwifeinc
Instruments: ~educ + age + kids5 + nwifeinc + exper + I(exper^2)
Residuals:
Min. 1st Qu. Median Mean 3rd Qu. Max.
-4.57e+03 -6.54e+02 -3.69e+01 -4.58e-12 5.70e+02 8.37e+03
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2225.662 574.564 3.8737 0.0001242
logWage 1639.556 470.576 3.4841 0.0005454
educ -183.751 59.100 -3.1092 0.0020032
age -7.806 9.378 -0.8324 0.4056640
kids5 -198.154 182.929 -1.0832 0.2793249
nwifeinc -10.170 6.615 -1.5374 0.1249417
Residual standard error: 1354.2045 on 422 degrees of freedom
> sem.logWage = tsls( logWage ~ hours + educ + exper + I(exper^2), instruments =~ educ + exper + I(exper^2) + age + kids5 + nwifeinc )
> summary(sem.logWage)
2SLS Estimates
Model Formula: logWage ~ hours + educ + exper + I(exper^2)
Instruments: ~educ + exper + I(exper^2) + age + kids5 + nwifeinc
Residuals:
Min. 1st Qu. Median Mean 3rd Qu. Max.
-3.50e+00 -2.93e-01 3.21e-02 1.48e-14 3.65e-01 2.46e+00
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.6557254 0.3377883 -1.9412 5.289e-02
hours 0.0001259 0.0002546 0.4945 6.212e-01
educ 0.1103300 0.0155244 7.1069 5.077e-12
exper 0.0345824 0.0194916 1.7742 7.675e-02
I(exper^2) -0.0007058 0.0004541 -1.5543 1.209e-01
Residual standard error: 0.6794 on 423 degrees of freedom
> # now isn't it interesting that with the more complex model, coeff for educ in Wage prediction is back up to
> # the value and signif given by the original, naive single predictor eq estimated by OLS !!!
> # Being an economist isn't as easy as it looks?
>
==================================================================================================
==================================================================================================
Appendix: Using ivreg command in package AER
Replicate tsls fits in Lab3; AER package linked in week 6
> install.packages("AER")
--- Please select a CRAN mirror for use in this session ---
also installing the dependencies car, lmtest, sandwich, strucchange, zoo
package 'car' successfully unpacked and MD5 sums checked
package 'lmtest' successfully unpacked and MD5 sums checked
package 'sandwich' successfully unpacked and MD5 sums checked
package 'strucchange' successfully unpacked and MD5 sums checked
package 'zoo' successfully unpacked and MD5 sums checked
package 'AER' successfully unpacked and MD5 sums checked
The downloaded packages are in
C:\Documents and Settings\Administrator\Local Settings\Temp\RtmpfPNXMz\downloaded_packages
updating HTML package descriptions
> library(AER)
Loading required package: car
Loading required package: lmtest
Loading required package: zoo
Attaching package: 'zoo'
The following object(s) are masked from package:base :
as.Date.numeric
Loading required package: sandwich
Loading required package: strucchange
Loading required package: survival
Loading required package: splines
> help(ivreg)
ivreg package:AER R Documentation
Instrumental-Variable Regression
Description:
Fit instrumental-variable regression by two-stage least squares.
This is equivalent to direct instrumental-variables estimation
when the number of instruments is equal to the number of
predictors.
Usage:
ivreg(formula, instruments, data, subset, na.action, weights, offset,
contrasts = NULL, model = TRUE, y = TRUE, x = FALSE, ...)
Arguments:
formula, instruments: formula specification(s) of the regression
relationship and the instruments. Either 'instruments' is
missing and 'formula' has three parts as in 'y ~ x1 + x2 | z1
+ z2 + z3' (recommended) or 'formula' is 'y ~ x1 + x2' and
'instruments' is a one-sided formula '~ z1 + z2 + z3' (only
for backward compatibility).
data: an optional data frame containing the variables in the model.
By default the variables are taken from the environment from
which 'ivreg' is called.
subset: an optional vector specifying a subset of observations to be
used in fitting the model.
na.action: a function that indicates what should happen when the data
contain 'NA's. The default is set by the 'na.action' option.
weights: an optional vector of weights to be used in the fitting
process.
offset: an optional offset that can be used to specify an a priori
known component to be included during fitting.
contrasts: an optional list. See the 'contrasts.arg' of
'model.matrix.default'.
model, x, y: logicals. If 'TRUE' the corresponding components of the
fit (the model frame, the model matrices , the response) are
returned.
...: further arguments passed to 'ivreg.fit'.
#Redo Lab3 IV fits
> mroz87 = read.table( "http://www-stat.stanford.edu/~rag/stat209/Mroz87.dat", header = T)
> poswage = subset(mroz87, wage > 0) # my new data subset
> poswage$logWage = log( poswage$wage ) # adding logWage to the data set for session
> attach(poswage)
> ivreg1 = ivreg(logWage ~ educ| fatheduc)
> summary(ivreg1)
Call:
ivreg(formula = logWage ~ educ | fatheduc)
Residuals:
Min 1Q Median 3Q Max
-3.0870 -0.3393 0.0525 0.4042 2.0677
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.44110 0.44610 0.989 0.323
educ 0.05917 0.03514 1.684 0.093 .
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6894 on 426 degrees of freedom
Multiple R-Squared: 0.09344, Adjusted R-squared: 0.09131
Wald test: 2.835 on 1 and 426 DF, p-value: 0.09294
> # matches tsls Task1
> confint(ivreg1) # parameter confidence intervals as usual
2.5 % 97.5 %
(Intercept) -0.433239996 1.3154468
educ -0.009703131 0.1280501
####### use the diagnostics option of ivreg (recent)
> ivreg1 = ivreg(logWage ~ educ| fatheduc)
> summary(ivreg1, diagnostics = TRUE)
Call:
ivreg(formula = logWage ~ educ | fatheduc)
Residuals:
Min 1Q Median 3Q Max
-3.0870 -0.3393 0.0525 0.4042 2.0677
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.44110 0.44610 0.989 0.3233
educ 0.05917 0.03514 1.684 0.0929 .
Diagnostic tests:
df1 df2 statistic p-value
Weak instruments 1 426 88.84 <2e-16 ***
Wu-Hausman 1 425 2.47 0.117
Sargan 0 NA NA NA
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6894 on 426 degrees of freedom
Multiple R-Squared: 0.09344, Adjusted R-squared: 0.09131
Wald test: 2.835 on 1 and 426 DF, p-value: 0.09294
> confint(ivreg1) # parameter confidence intervals as usual
2.5 % 97.5 %
(Intercept) -0.433239996 1.3154468
educ -0.009703131 0.1280501
##Investigate the 'Weak instruments' entry
> cor(educ, fatheduc)
[1] 0.415403
> cor.test(educ, fatheduc)
Pearson's product-moment correlation
data: educ and fatheduc
t = 9.4255, df = 426, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.3337579 0.4908623
sample estimates:
cor
0.415403
> 9.4255^2
[1] 88.84005
## just a test (like in the main lab of correlation between predictor and instrument
# the Hausman test is often used to show there's a difference between OLS and IV, but what does that tell you?
# ivreg diagnostics finds non-signicance even though IV est cuts the OLS estimate in half and IV estimate
# for comparison OLS fit
> lm1 = lm(logWage ~ educ)
> summary(lm1)
Call:
lm(formula = logWage ~ educ)
Residuals:
Min 1Q Median 3Q Max
-3.10256 -0.31473 0.06434 0.40081 2.10029
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.1852 0.1852 -1.000 0.318
educ 0.1086 0.0144 7.545 2.76e-13 ***
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.68 on 426 degrees of freedom
Multiple R-squared: 0.1179, Adjusted R-squared: 0.1158
F-statistic: 56.93 on 1 and 426 DF, p-value: 2.761e-13
> confint(lm1)
2.5 % 97.5 %
(Intercept) -0.54926726 0.1788736
educ 0.08034506 0.1369522
## to see the details of the Hausman test (and do it by hand) see pdf page 46 onward from
http://personal.rhul.ac.uk/uhte/006/ec2203/Lecture%2015_IVestimation.pdf
also the Basel class notes linked in week 6
## The Sargan test is a statistical test used for testing over-identifying restrictions
## relevant to simultaneous eqs, which we investigated in identifiability in lab script
> ivreg2 = ivreg(logWage ~ educ + exper + I(exper^2) | fatheduc + motheduc + exper + I(exper^2))
> summary(ivreg2)
Call:
ivreg(formula = logWage ~ educ + exper + I(exper^2) | fatheduc +
motheduc + exper + I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-3.0986 -0.3196 0.0551 0.3689 2.3493
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0481003 0.4003281 0.120 0.90442
educ 0.0613966 0.0314367 1.953 0.05147 .
exper 0.0441704 0.0134325 3.288 0.00109 **
I(exper^2) -0.0008990 0.0004017 -2.238 0.02574 *
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6747 on 424 degrees of freedom
Multiple R-Squared: 0.1357, Adjusted R-squared: 0.1296
Wald test: 8.141 on 3 and 424 DF, p-value: 2.787e-05
> # matches tsls Task 2 exactly
### summary with diagnostics
> summary(ivreg2, diagnostics = TRUE)
Call:
ivreg(formula = logWage ~ educ + exper + I(exper^2) | fatheduc +
motheduc + exper + I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-3.0986 -0.3196 0.0551 0.3689 2.3493
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0481003 0.4003281 0.120 0.90442
educ 0.0613966 0.0314367 1.953 0.05147 .
exper 0.0441704 0.0134325 3.288 0.00109 **
I(exper^2) -0.0008990 0.0004017 -2.238 0.02574 *
Diagnostic tests:
df1 df2 statistic p-value
Weak instruments 2 423 55.400 <2e-16 ***
Wu-Hausman 1 423 2.793 0.0954 .
Sargan 1 NA 0.378 0.5386
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6747 on 424 degrees of freedom
Multiple R-Squared: 0.1357, Adjusted R-squared: 0.1296
Wald test: 8.141 on 3 and 424 DF, p-value: 2.787e-05
> confint(ivreg2)
2.5 % 97.5 %
(Intercept) -0.7365283088 0.8327289181
educ -0.0002181633 0.1230114191
exper 0.0178432261 0.0704975626
I(exper^2) -0.0016862590 -0.0001116803
> ivreg3 = ivreg(hours ~ logWage + educ + age + kids5 + nwifeinc | educ + age + kids5 + nwifeinc + exper + I(exper^2) )
> summary(ivreg3)
Call:
ivreg(formula = hours ~ logWage + educ + age + kids5 + nwifeinc |
educ + age + kids5 + nwifeinc + exper + I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-4570.13 -654.08 -36.94 569.86 8372.91
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2225.662 574.564 3.874 0.000124 ***
logWage 1639.556 470.576 3.484 0.000545 ***
educ -183.751 59.100 -3.109 0.002003 **
age -7.806 9.378 -0.832 0.405664
kids5 -198.154 182.929 -1.083 0.279325
nwifeinc -10.170 6.615 -1.537 0.124942
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 1354 on 422 degrees of freedom
Multiple R-Squared: -2.008, Adjusted R-squared: -2.043
Wald test: 3.441 on 5 and 422 DF, p-value: 0.004648
> ivreg4 = ivreg(logWage ~ hours + educ + exper + I(exper^2) | educ + age + kids5 + nwifeinc + exper + I(exper^2) )
> summary(ivreg4)
Call:
ivreg(formula = logWage ~ hours + educ + exper + I(exper^2) |
educ + age + kids5 + nwifeinc + exper + I(exper^2))
Residuals:
Min 1Q Median 3Q Max
-3.49800 -0.29307 0.03208 0.36486 2.45912
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.6557254 0.3377883 -1.941 0.0529 .
hours 0.0001259 0.0002546 0.494 0.6212
educ 0.1103300 0.0155244 7.107 5.08e-12 ***
exper 0.0345824 0.0194916 1.774 0.0767 .
I(exper^2) -0.0007058 0.0004541 -1.554 0.1209
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6794 on 423 degrees of freedom
Multiple R-Squared: 0.1257, Adjusted R-squared: 0.1174
Wald test: 19.03 on 4 and 423 DF, p-value: 2.108e-14
>
> confint(ivreg3)
2.5 % 97.5 %
(Intercept) 1099.53686 3351.786795
logWage 717.24422 2561.866954
educ -299.58477 -67.917788
age -26.18666 10.574476
kids5 -556.68882 160.380224
nwifeinc -23.13425 2.795067
> confint(ivreg4)
2.5 % 97.5 %
(Intercept) -1.3177783230 0.0063274473
hours -0.0003731274 0.0006249278
educ 0.0799028211 0.1407571860
exper -0.0036203894 0.0727851026
I(exper^2) -0.0015957505 0.0001842115
####################### Another bolt-on, try some of the tools in Dylan Small's ivpack ########################
> install.packages("AER")
> install.packages("ivpack")
> library(AER)
> library(ivpack)
> mroz87 = read.table( "http://statweb.stanford.edu/~rag/stat209/Mroz87.dat", header = T)
> poswage = subset(mroz87, wage > 0) # my new data subset
> poswage$logWage = log( poswage$wage ) # adding logWage to the data set for session
> attach(poswage)
> ivreg1 = ivreg(logWage ~ educ| fatheduc, x = TRUE)
> summary(ivreg1)
Call:
ivreg(formula = logWage ~ educ | fatheduc, x = TRUE)
Residuals:
Min 1Q Median 3Q Max
-3.0870 -0.3393 0.0525 0.4042 2.0677
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.44110 0.44610 0.989 0.3233
educ 0.05917 0.03514 1.684 0.0929 .
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6894 on 426 degrees of freedom
Multiple R-Squared: 0.09344, Adjusted R-squared: 0.09131
Wald test: 2.835 on 1 and 426 DF, p-value: 0.09294
> confint(ivreg1)
2.5 % 97.5 %
(Intercept) -0.433239996 1.3154468
educ -0.009703131 0.1280501
> anderson.rubin.ci(ivreg1) # from ivpack (primarily for Wald estimator)
$confidence.interval
[1] "[ -0.0142227779599083 , 0.127132921928464 ]"
# Compute robust to heteroskedasticity standard errors for an instrumental variables analysis. These
# are the Huber-White standard errors for an instrumental variable analysis as described in White (1982).
> robust.se(ivreg1)
[1] "Robust Standard Errors"
t test of coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.441103 0.464287 0.9501 0.3426
educ 0.059173 0.036943 1.6017 0.1100
> ivreg2 = ivreg(logWage ~ educ + exper + I(exper^2) | fatheduc + motheduc + exper + I(exper^2), x = TRUE)
> summary(ivreg2)
Call:
ivreg(formula = logWage ~ educ + exper + I(exper^2) | fatheduc +
motheduc + exper + I(exper^2), x = TRUE)
Residuals:
Min 1Q Median 3Q Max
-3.0986 -0.3196 0.0551 0.3689 2.3493
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.0481003 0.4003281 0.120 0.90442
educ 0.0613966 0.0314367 1.953 0.05147 .
exper 0.0441704 0.0134325 3.288 0.00109 **
I(exper^2) -0.0008990 0.0004017 -2.238 0.02574 *
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
Residual standard error: 0.6747 on 424 degrees of freedom
Multiple R-Squared: 0.1357, Adjusted R-squared: 0.1296
Wald test: 8.141 on 3 and 424 DF, p-value: 2.787e-05
> confint(ivreg2)
2.5 % 97.5 %
(Intercept) -0.7365283088 0.8327289181
educ -0.0002181633 0.1230114191
exper 0.0178432261 0.0704975626
I(exper^2) -0.0016862590 -0.0001116803
> robust.se(ivreg2) # a little bigger than ivreg
[1] "Robust Standard Errors"
t test of coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.04810030 0.42778460 0.1124 0.910527
educ 0.06139663 0.03318243 1.8503 0.064969 .
exper 0.04417039 0.01547356 2.8546 0.004521 **
I(exper^2) -0.00089897 0.00042807 -2.1001 0.036314 *
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
>
END
Lab3 Rogosa session 2019
==================================================