econ61001_asn1.R
2020-12-27
library(wooldridge)
## Warning: package ‘wooldridge’ was built under R version 4.0.3
library(MASS)
##
## Attaching package: ‘MASS’
## The following object is masked from ‘package:wooldridge’:
##
## cement
data(‘campus’)
##### Question 1
# calculate the sample mean of each explanatory variable, and round each value to 3 digits
round(apply(campus[, c(“enroll”, “police”, “priv”)], 2, mean), 3)
## enroll police priv
## 16076.351 20.495 0.124
# calculate the sample standard deviation of each explanatory variable, and round each value to 3 digits
round(apply(campus[, c(“enroll”, “police”, “priv”)], 2, sd), 3)
## enroll police priv
## 12298.986 15.631 0.331
##### Question 2
y <- campus$crime
x1 <- campus$enroll
x2 <- campus$police
# priv is a categorical variable
x3 <- as.factor(campus$priv)
# fitting the linear model using OLS
ols <- lm(y ~ x1+x2+x3)
summary(ols)
##
## Call:
## lm(formula = y ~ x1 + x2 + x3)
##
## Residuals:
## Min 1Q Median 3Q Max
## -617.88 -123.67 1.68 82.41 1144.99
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -174.07758 47.67878 -3.651 0.000431 ***
## x1 0.02486 0.00290 8.572 2.14e-13 ***
## x2 7.78662 2.26744 3.434 0.000890 ***
## x31 75.05962 78.57454 0.955 0.341920
## —
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1
##
## Residual standard error: 241.5 on 93 degrees of freedom
## Multiple R-squared: 0.7338, Adjusted R-squared: 0.7252
## F-statistic: 85.46 on 3 and 93 DF, p-value: < 2.2e-16
# extract the residuals
u <- resid(ols)
# x3^2 is omitted
f <- formula(u2 ~ (x1+x2+x3)^2+I(x1^2)+I(x2^2))
lm.fit <- lm(f, data=data.frame(u2=u^2, x1=x1, x2=x2, x3=x3))
summary(lm.fit)
##
## Call:
## lm(formula = f, data = data.frame(u2 = u^2, x1 = x1, x2 = x2,
## x3 = x3))
##
## Residuals:
## Min 1Q Median 3Q Max
## -230536 -55751 -10147 27014 1141347
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.867e+04 4.944e+04 -0.782 0.4362
## x1 1.084e+01 4.983e+00 2.176 0.0323 *
## x2 -4.039e+03 4.351e+03 -0.928 0.3559
## x31 4.177e+04 1.099e+05 0.380 0.7049
## I(x1^2) -6.415e-06 1.150e-04 -0.056 0.9556
## I(x2^2) 1.082e+02 7.753e+01 1.395 0.1664
## x1:x2 -1.468e-01 1.464e-01 -1.003 0.3185
## x1:x31 -6.174e+00 1.428e+01 -0.432 0.6665
## x2:x31 1.872e+03 7.555e+03 0.248 0.8049
## —
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ‘ 1
##
## Residual standard error: 147000 on 88 degrees of freedom
## Multiple R-squared: 0.1893, Adjusted R-squared: 0.1156
## F-statistic: 2.568 on 8 and 88 DF, p-value: 0.01446
n <- nrow(campus)
# the F-test statistic
F <- summary(lm.fit)$fstatistic
# the F-statistic
round(F[1], 3)
## value
## 2.568
k <- F[2]
# the p-value for F-test version
round(pf(F[1],F[2],F[3], lower.tail = FALSE), 3)
## value
## 0.014
# the R-squared of the regression
R2 <- summary(lm.fit)$r.squared
# the LM statistic
LM <- nrow(campus) * R2
# the LM statistic
round(LM, 3)
## [1] 18.358
# the p-value for LM-test version
round(pchisq(LM, df=k, lower.tail = FALSE), 3)
## [1] 0.019
##### Question 3
mc <- 1000
F.result <- vector('logical', length = mc)
LM.result <- vector('logical', length = mc)
for (i in 1:mc) {
set.seed(i)
u <- rnorm(n)
m <- lm(f,data=data.frame(u2=u^2, x1=x1, x2=x2, x3=x3))
R2 <- summary(m)$r.squared
F.result[i] <- (R2/k/((1-R2)/(n-k-1))) > F[1]
LM.result[i] <- n*R2 > LM
}
# size F
sum(F.result)/mc
## [1] 0.05
# size LM
sum(LM.result)/mc
## [1] 0.05
Sigma <- diag(x1, n, n)
for (i in 1:mc) {
set.seed(i)
u <- mvrnorm(1, rep(0, n), Sigma)
m <- lm(f,data=data.frame(u2=u^2, x1=x1, x2=x2, x3=x3))
R2 <- summary(m)$r.squared
F.result[i] <- (R2/k/((1-R2)/(n-k-1)))> F[1]
LM.result[i] <- n*R2 > LM
}
# Power F D1
sum(F.result)/mc
## [1] 0.718
# Power LM D1
sum(LM.result)/mc
## [1] 0.718
Sigma <- diag(x1+x2, n, n)
for (i in 1:mc) {
set.seed(i)
u <- mvrnorm(1, rep(0, n), Sigma)
m <- lm(f,data=data.frame(u2=u^2, x1=x1, x2=x2, x3=x3))
R2 <- summary(m)$r.squared
F.result[i] <- (R2/k/((1-R2)/(n-k-1)))> F[1]
LM.result[i] <- n*R2 > LM
}
# Power F D2
sum(F.result)/mc
## [1] 0.713
# Power LM D2
sum(LM.result)/mc
## [1] 0.713