Serial Correlation and Heteroskedasticity in Time Series Regressions (Chapter 12)

Comments will follow…


library(foreign)
library(lmtest)
library(sandwich)
library(tseries)

download.file('http://fmwww.bc.edu/ec-p/data/wooldridge/phillips.dta','phillips.dta',mode="wb")
download.file('http://fmwww.bc.edu/ec-p/data/wooldridge/prminwge.dta','prminwge.dta',mode="wb")
download.file('http://fmwww.bc.edu/ec-p/data/wooldridge/barium.dta','barium.dta',mode="wb")
download.file('http://fmwww.bc.edu/ec-p/data/wooldridge/intdef.dta','intdef.dta',mode="wb")
download.file('http://fmwww.bc.edu/ec-p/data/wooldridge/nyse.dta','nyse.dta',mode="wb")

# Example 12.1
phillips<-read.dta('phillips.dta')

# Static Phillips curve
lm.12.1.1<-lm(inf ~ unem, data=phillips)
res.static<-lm.12.1.1$res
res.static_1<-c(NA,lm.12.1.1$res[1:48])
lm.12.1.1.test<-lm(res.static ~ res.static_1)
summary(lm.12.1.1.test)

# Augmented Phillips curve
lm.12.1.2<-lm((inf-inf_1) ~ unem, data=phillips)
res.aug<-lm.12.1.2$res
res.aug_1<-c(NA,lm.12.1.2$res[1:47])
lm.12.1.2.test<-lm(res.aug ~ res.aug_1)
summary(lm.12.1.2.test)

# Heteroskedasticity robust t statistic
coeftest(lm.12.1.1.test,vcov=vcovHC(lm.12.1.1.test,type="HC0"))
coeftest(lm.12.1.2.test,vcov=vcovHC(lm.12.1.2.test,type="HC0"))

# Durbin-Watson Test
# from package lmtest
dwtest(inf ~ unem, data=phillips)
dwtest((inf-inf_1) ~ unem, data=phillips)

# Example 12.2
prminwge<-read.dta('prminwge.dta')

lm.12.2<-lm(lprepop ~ lmincov + lusgnp + lprgnp + t,data=prminwge)
lm.12.2.res<-lm.12.2$res
lm.12.2.res_1<-c(NA,lm.12.2$res[1:(length(lm.12.2.res)-1)])
lm.12.2.test<-lm(lm.12.2.res ~ lmincov + lusgnp + lprgnp + t + lm.12.2.res_1,data=prminwge)
summary(lm.12.2.test)

# For comparison the test with strict exogenous regressors
summary(lm(lm.12.2.res ~ lm.12.2.res_1))

# Example 12.3
barium<-read.dta('barium.dta')

lm.12.3<-lm(lchnimp ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6, data=barium)

lm.12.3.res<-lm.12.3$res
lm.12.3.res_1<-c(NA,lm.12.3$res[1:(length(lm.12.3.res)-1)])
lm.12.3.res_2<-c(NA,NA,lm.12.3$res[1:(length(lm.12.3.res)-2)])
lm.12.3.res_3<-c(NA,NA,NA,lm.12.3$res[1:(length(lm.12.3.res)-3)])

data<-cbind(barium,lm.12.3.res)
data<-cbind(data,lm.12.3.res_1)
data<-cbind(data,lm.12.3.res_2)
data<-cbind(data,lm.12.3.res_3)

lm.12.3.test<-lm(lm.12.3.res ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6 +
 lm.12.3.res_1 + lm.12.3.res_2 +lm.12.3.res_3, data=data)

lm.12.3.test.res<-lm(lm.12.3.res ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6,
 data=data[4:length(data[,1]),])

anova(lm.12.3.test,lm.12.3.test.res)

# Seasonal autocorrelation
data<-cbind(data,lm.12.3.res_12=c(rep(NA,12),lm.12.3$res[1:(length(lm.12.3.res)-12)]))
summary(lm(lm.12.3.res ~ lm.12.3.res_12,data=data))

# Include regressors
summary(lm(lm.12.3.res ~ lchempi + lgas + lrtwex + befile6 + affile6 + afdec6 +
 lm.12.3.res_12,data=data))

# Example 12.4
barium<-read.dta('barium.dta')

# Download my prais package
library(prais)
pw.1<-prais.winsten(lchnimp ~ lchempi + lgas + lrtwex + befile6 + 
 affile6 + afdec6, data=barium)
pw.1

# Example 12.5
phillips<-read.dta('phillips.dta')

# Column 1
lm.12.5<-lm(inf ~ unem, data=phillips)
summary(lm.12.5)

# Column 2
pw.2<-prais.winsten(inf ~ unem, data=phillips)
pw.2

# Example 12.6
intdef<-read.dta('intdef.dta')
summary(lm.12.6.1<-lm(i3 ~ inf + def, data=intdef))
res<-lm.12.6.1$residuals
res_1<-c(NA,res[-length(res)])
summary(lm(res ~ -1 + res_1))

lm.12.6.2<-lm(diff(i3) ~ diff(inf) + diff(def),data=intdef)
summary(lm.12.6.2)

# Correlation between i3t and i3t-1
with(intdef,cor(i3,c(NA,i3[-length(i3)]),use="pairwise.complete.obs"))

# Regression of et on et-1
res<-lm.12.6.2$residuals
res_1<-res_1<-c(NA,res[-length(res)])
summary(lm(res ~ -1 + res_1))

# Example 12.7
prminwge<-read.dta('prminwge.dta')

lm.12.7.1<-lm(lprepop ~ lmincov + lusgnp + lprgnp + t,data=prminwge)
summary(lm.12.7.1)

# Still looking for the SC/heteroskedasticiy robust se.

pw.12.7.1<-prais.winsten(lprepop ~ lmincov + lusgnp + lprgnp + t,data=prminwge)
pw.12.7.1

# Example 12.8
nyse<-read.dta('nyse.dta')

lm.12.8<-lm(return ~ return_1,data=nyse)
summary(lm.12.8)

res<-lm.12.8$residuals^2
lm.12.8.bp<-lm(res ~ return_1,data=na.omit(nyse))
summary(lm.12.8.bp)

# Example 12.9
nyse<-read.dta('nyse.dta')

lm.12.9<-lm(return ~ return_1,data=nyse)
summary(lm.12.9)

u<-lm.12.8$residuals
u_1<-c(NA,u[-length(u)])
usq<-lm.12.8$residuals^2
usq_1<-c(NA,usq[-length(usq)])

lm.12.9.arch<-lm(usq ~ usq_1)
summary(lm.12.9.arch)

lm.12.9.3<-lm(u ~ u_1)
summary(lm.12.9.3)
Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s