# install.packages('lmtest')
library(readr) # read csv
library(lmtest) # stat tests
library(sandwich) # cov matrixes
library(readxl) # csv and excel files
library(psych) # descriptive statistics
library(knitr) # for tables
library(xtable) # for latex tables
library(ggplot2) # for graphs
library(car) # for linear hypothesis testing
earn_data <- read_csv('C:/Users/Popov/Documents/Studies/NES_studies/R/Econometrics_1/HA1/cps99_ps1.csv',
show_col_types = FALSE)
ahe_yrseduc <- lm(ahe ~ yrseduc, data = earn_data)
rob_se_ahe_yrseduc<- sqrt(diag(vcovHC(ahe_yrseduc, type = "HC1"))) #for plots
print(coeftest(ahe_yrseduc, vcov = vcovHC, type = 'HC1'))
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.609050 0.647177 -4.0314 5.653e-05 ***
## yrseduc 1.321859 0.050043 26.4145 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ahe_2 <- lm(ahe ~ yrseduc + female, data = earn_data)
rob_se_ahe_2<- sqrt(diag(vcovHC(ahe_2, type = "HC1"))) #for plots
print(coeftest(ahe_2, vcov = vcovHC, type = 'HC1'))
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.383958 0.630403 -2.1954 0.0282 *
## yrseduc 1.341470 0.048854 27.4587 <2e-16 ***
## female -3.395995 0.228019 -14.8935 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
H_{0}: \(female = 0\)
H_{1}: \(female \neq 0\)
\(\alpha\) = 0.05
f_H0 <- c("female = 0")
linearHypothesis(ahe_2, f_H0, vcov = vcovHC(ahe_2, type = "HC1"))
## Linear hypothesis test
##
## Hypothesis:
## female = 0
##
## Model 1: restricted model
## Model 2: ahe ~ yrseduc + female
##
## Note: Coefficient covariance matrix supplied.
##
## Res.Df Df F Pr(>F)
## 1 3779
## 2 3778 1 221.82 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Based of the test, null hypothesis is rejected.
cor_t <- cor.test(earn_data$yrseduc, earn_data$female)
# Display the correlation and its significance
print(cor_t)
##
## Pearson's product-moment correlation
##
## data: earn_data$yrseduc and earn_data$female
## t = 1.7746, df = 3779, p-value = 0.07605
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.003023847 0.060676025
## sample estimates:
## cor
## 0.02885538
P.S. The t-statistic for the correlation coefficient in the cor.test function is computed using the formula: \(t = \frac{1 - r^2}{r \sqrt{N - 2}}\),
where \(r\) is the sample correlation coefficient, \(N\) is the number of pairs of observations.
ahe_3 <- lm(ahe ~ yrseduc + ba, data = earn_data)
rob_se_ahe_3<- sqrt(diag(vcovHC(ahe_3, type = "HC1"))) #for plots
print(coeftest(ahe_3, vcov = vcovHC, type = 'HC1'))
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.245578 0.894054 -1.3932 0.16365
## yrseduc 1.204144 0.074538 16.1547 < 2e-16 ***
## ba 0.786030 0.465545 1.6884 0.09142 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
c_int_ba <- confint(ahe_3, coef = "ba", vcov. = vcovHC(ahe_yrseduc, type = "HC1"))
cat("\n95% Confidence Interval for β3 (slope):\n", c_int_ba["ba",])
##
## 95% Confidence Interval for β3 (slope):
## -0.1332442 1.705304
```