#----------------------------------------
# CAU 6
#(i) Hãy hồi quy log(wage) (log của tiền lương) theo số năm học (educ). Hồi quy IV (biến
#công cụ): Vì ta nghi ngờ rằng educ có thể là biến nội sinh (endogenous), nên hãy sử
#dụng sibs (số anh/chị/em ruột) làm biến công cụ (instrument) cho educ trong mô hình
#IV regression. Dùng lệnh “stargazer” để so sánh kết quả của hai mô hình hồi quy
#(OLS và IV).
#-----------------------------------------
library(wooldridge)
library(AER)
## Loading required package: car
## Loading required package: carData
## Loading required package: lmtest
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
## Loading required package: sandwich
## Loading required package: survival
library(stargazer)
##
## Please cite as:
## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.
## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer
data(wage2, package="wooldridge")
str(wage2)
## 'data.frame': 935 obs. of 17 variables:
## $ wage : int 769 808 825 650 562 1400 600 1081 1154 1000 ...
## $ hours : int 40 50 40 40 40 40 40 40 45 40 ...
## $ IQ : int 93 119 108 96 74 116 91 114 111 95 ...
## $ KWW : int 35 41 46 32 27 43 24 50 37 44 ...
## $ educ : int 12 18 14 12 11 16 10 18 15 12 ...
## $ exper : int 11 11 11 13 14 14 13 8 13 16 ...
## $ tenure : int 2 16 9 7 5 2 0 14 1 16 ...
## $ age : int 31 37 33 32 34 35 30 38 36 36 ...
## $ married: int 1 1 1 1 1 1 0 1 1 1 ...
## $ black : int 0 0 0 0 0 1 0 0 0 0 ...
## $ south : int 0 0 0 0 0 0 0 0 0 0 ...
## $ urban : int 1 1 1 1 1 1 1 1 0 1 ...
## $ sibs : int 1 1 1 4 10 1 1 2 2 1 ...
## $ brthord: int 2 NA 2 3 6 2 2 3 3 1 ...
## $ meduc : int 8 14 14 12 6 8 8 8 14 12 ...
## $ feduc : int 8 14 14 12 11 NA 8 NA 5 11 ...
## $ lwage : num 6.65 6.69 6.72 6.48 6.33 ...
## - attr(*, "time.stamp")= chr "25 Jun 2011 23:03"
ols1 <- lm(lwage ~ educ, data=wage2)
iv1 <- ivreg(lwage ~ educ | sibs, data=wage2)
stargazer(ols1, iv1,
type = "text",
title = "So sánh OLS và IV (sibs làm công cụ cho educ)",
keep = c("educ"),
covariate.labels = c("Số năm học ở trường"),
dep.var.labels = "log(wage)")
##
## So sánh OLS và IV (sibs làm công cụ cho educ)
## ====================================================================
## Dependent variable:
## -------------------------------------
## log(wage)
## OLS instrumental
## variable
## (1) (2)
## --------------------------------------------------------------------
## Số năm học ở trường 0.060*** 0.122***
## (0.006) (0.026)
##
## --------------------------------------------------------------------
## Observations 935 935
## R2 0.097 -0.009
## Adjusted R2 0.096 -0.010
## Residual Std. Error (df = 933) 0.400 0.423
## F Statistic 100.700*** (df = 1; 933)
## ====================================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
print("TRA LOI (i): OLS (6%) có thể bị đánh giá thấp tác động thật của giáo dục. IV, dùng sibs là biến công cụ, cho hệ số lớn hơn (12,2%).")
## [1] "TRA LOI (i): OLS (6%) có thể bị đánh giá thấp tác động thật của giáo dục. IV, dùng sibs là biến công cụ, cho hệ số lớn hơn (12,2%)."
#----------------------------------
#(ii)
#Biến brthord thể hiện thứ tự sinh trong gia đình. brthord = 1 nếu là con cả, = 2 nếu là
#con thứ hai, v.v.). Giải thích mối tương quan giữa educ và brthord. Hãy giải thích tại
#sao educ và brthord có thể tương quan âm (negatively correlated).
cor(wage2$educ, wage2$brthord, use = "complete.obs")
## [1] -0.2049925
print("TRA LOI (ii): nguồn lực bị trải đều cho nhiều con, người sinh sau thiệt thòi hơn.Nói cách khác, cha mẹ thường đặt kỳ vọng cao vào con đầu, các con sau có thể chịu ảnh hưởng từ hoàn cảnh kinh tế hoặc kỳ vọng giảm dần. Con sau sinh khi cha mẹ đã lớn tuổi hơn, có thể có thu nhập thấp hơn hoặc ít thời gian đầu tư cho học hành...giả thích tương quan âm.")
## [1] "TRA LOI (ii): nguồn lực bị trải đều cho nhiều con, người sinh sau thiệt thòi hơn.Nói cách khác, cha mẹ thường đặt kỳ vọng cao vào con đầu, các con sau có thể chịu ảnh hưởng từ hoàn cảnh kinh tế hoặc kỳ vọng giảm dần. Con sau sinh khi cha mẹ đã lớn tuổi hơn, có thể có thu nhập thấp hơn hoặc ít thời gian đầu tư cho học hành...giả thích tương quan âm."
#-------------------------------------
#(iii)
#Dùng brthord làm biến công cụ: Hãy sử dụng brthord làm biến công cụ cho educ
#trong phương trình hồi quy
#log(wage) = 𝛽0 + 𝛽1 educ + u (1)
#IV với brthord làm công cụ:
iv2 <- ivreg(lwage ~ educ | brthord, data = wage2)
summary(iv2)
##
## Call:
## ivreg(formula = lwage ~ educ | brthord, data = wage2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8532 -0.2557 0.0435 0.2970 1.3033
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.03040 0.43295 11.619 < 2e-16 ***
## educ 0.13064 0.03204 4.078 4.97e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4215 on 850 degrees of freedom
## Multiple R-Squared: -0.02862, Adjusted R-squared: -0.02983
## Wald test: 16.63 on 1 and 850 DF, p-value: 4.975e-05
#---------------------------
#(iv) OLS có sibs
#log(wage) = 𝛽0 + 𝛽1 educ + 𝛽2 sibs + u (2)
ols2 <- lm(lwage ~ educ + sibs, data = wage2)
summary(ols2)
##
## Call:
## lm(formula = lwage ~ educ + sibs, data = wage2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.9206 -0.2446 0.0423 0.2699 1.2650
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.068773 0.089116 68.099 < 2e-16 ***
## educ 0.056037 0.006123 9.152 < 2e-16 ***
## sibs -0.015133 0.005832 -2.595 0.00961 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3991 on 932 degrees of freedom
## Multiple R-squared: 0.1039, Adjusted R-squared: 0.102
## F-statistic: 54.03 on 2 and 932 DF, p-value: < 2.2e-16
print("TRA LOI: Khi số năm học tăng thêm 1 năm, thu nhập tăng trung bình 5.6%. Khi số anh/chị/em tăng thêm 1 người, thu nhập giảm khoảng 1.5%")
## [1] "TRA LOI: Khi số năm học tăng thêm 1 năm, thu nhập tăng trung bình 5.6%. Khi số anh/chị/em tăng thêm 1 người, thu nhập giảm khoảng 1.5%"
#--------------------------------
#(v)
print("Ước lượng phương trình (2) bằng hồi quy IV. Biến educ bị nội sinh, nên dùng biến công cụ brthord.")
## [1] "Ước lượng phương trình (2) bằng hồi quy IV. Biến educ bị nội sinh, nên dùng biến công cụ brthord."
iv_educ_brthord <- ivreg(lwage ~ educ + sibs | brthord + sibs, data = wage2)
summary(iv_educ_brthord)
##
## Call:
## ivreg(formula = lwage ~ educ + sibs | brthord + sibs, data = wage2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.84808 -0.26227 0.03841 0.29901 1.30836
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.938527 1.055690 4.678 3.37e-06 ***
## educ 0.136994 0.074681 1.834 0.0669 .
## sibs 0.002111 0.017372 0.122 0.9033
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.427 on 849 degrees of freedom
## Multiple R-Squared: -0.05428, Adjusted R-squared: -0.05676
## Wald test: 10.9 on 2 and 849 DF, p-value: 2.124e-05
#SE(educ) = 0.0747; SE(sibs) = 0.017 nhỏ, nhưng hệ số gần bằng 0
#Khi thêm sibs vào phương trình, công cụ brthord và sibs có tương quan, làm công cụ yếu hơn. Do đó, phần biến thiên của educ được giải thích bởi brthord sau khi trừ ảnh hưởng sibs rất nhỏ -> Bước 1 yếu -> SE tăng mạnh.
# Kiem tra
summary(iv_educ_brthord, diagnostics = TRUE)
##
## Call:
## ivreg(formula = lwage ~ educ + sibs | brthord + sibs, data = wage2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.84808 -0.26227 0.03841 0.29901 1.30836
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.938527 1.055690 4.678 3.37e-06 ***
## educ 0.136994 0.074681 1.834 0.0669 .
## sibs 0.002111 0.017372 0.122 0.9033
##
## Diagnostic tests:
## df1 df2 statistic p-value
## Weak instruments 1 849 7.155 0.00762 **
## Wu-Hausman 1 848 1.343 0.24683
## Sargan 0 NA NA NA
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.427 on 849 degrees of freedom
## Multiple R-Squared: -0.05428, Adjusted R-squared: -0.05676
## Wald test: 10.9 on 2 and 849 DF, p-value: 2.124e-05
#F-stat bước 1 thấp (<10)
print("Khi thêm sibs vào mô hình, công cụ brthord trở nên yếu hơn, dẫn đến sai số chuẩn của hệ số ước lượng educ tăng lên đáng kể.")
## [1] "Khi thêm sibs vào mô hình, công cụ brthord trở nên yếu hơn, dẫn đến sai số chuẩn của hệ số ước lượng educ tăng lên đáng kể."
#---------------------------------
#(vi)
w2 <- subset(wage2, !is.na(educ) & !is.na(brthord) & !is.na(sibs))
first_stage <- lm(educ ~ brthord + sibs, data = w2)
w2$educ_hat <- fitted(first_stage)
cor(w2$educ_hat, w2$sibs)
## [1] -0.9294818
print("TRA LOI (vi): Kết quả cho thấy hệ số tương quan giữa educ_hat và sibs là −0.93, tức là tương quan âm rất mạnh. Điều này có nghĩa rằng phần dự đoán của giáo dục (từ công cụ brthord) gần như phụ thuộc hoàn toàn vào số lượng ace. Khi hai biến educ_hat và sibs được đưa đồng thời vào phương trình bước hai, hiện tượng đa cộng tuyến xảy ra, làm cho sai số chuẩn của giá trị ước lượng (educ) trong mô hình IV tăng mạnh, đúng như kết quả đã quan sát ở phần (v).")
## [1] "TRA LOI (vi): Kết quả cho thấy hệ số tương quan giữa educ_hat và sibs là −0.93, tức là tương quan âm rất mạnh. Điều này có nghĩa rằng phần dự đoán của giáo dục (từ công cụ brthord) gần như phụ thuộc hoàn toàn vào số lượng ace. Khi hai biến educ_hat và sibs được đưa đồng thời vào phương trình bước hai, hiện tượng đa cộng tuyến xảy ra, làm cho sai số chuẩn của giá trị ước lượng (educ) trong mô hình IV tăng mạnh, đúng như kết quả đã quan sát ở phần (v)."