?>

 

library(dynlm);library(stargazer)
data(hseinv, package=’wooldridge’) ###從hseinv資料庫裡取資料
str(hseinv)
hseinv

windows()
plot((invpc)~price,data=hseinv,
pch=19,col=’blue’
)
abline(lm(invpc~price,data=hseinv),col=2)

windows()
plot(log(invpc)~price,data=hseinv)
abline(lm(log(invpc)~price,data=hseinv),col=2)

# Define Yearly time series beginning in 1947
tsdata <- ts(hseinv[,c(1,4,9)], start=1947)
tsdata2 <- (hseinv[,c(1,4,9)])

windows()
plot(tsdata) ###時間序列資料

 

windows()
plot(tsdata2)

 

# Linear regression of model with lags:
res1 <- dynlm(log(invpc) ~ log(price) , data=tsdata) ###不考慮時間
res2 <- dynlm(log(invpc) ~ log(price) + trend(tsdata), data=tsdata) ###考慮時間因素 ###trend是把資料庫從1開始排起

# Pretty regression table
stargazer(res1,res2, type=”text”)

 

# Linear regression of model with lags:
res1_1 <- lm(log(invpc) ~ log(price) , data=hseinv) ###不考慮時間
res2_1 <- lm(log(invpc) ~ log(price) + year , data=hseinv) ###考慮時間因素

# Pretty regression table
stargazer(res1,res2,res1_1,res2_1, type=”text”)

 

hseinv$deyear1= hseinv$year-mean(hseinv[,1]) ###deyear是去掉時間的意思###去掉時間平均
hseinv

res2_2 <- lm(log(invpc) ~ log(price) + deyear1, data=hseinv)

hseinv$deyear2= hseinv$year-1947 #mean(hseinv[,1]) ###去掉1947年
hseinv

res2_3 <- lm(log(invpc) ~ log(price) + deyear2, data=hseinv)

# Pretty regression table
stargazer(res1,res2,res1_1,res2_1,res2_2,res2_3, type=”text”)

# how to get the -0.913 ???
hseinv$deyear4= hseinv$year-1946
res2_4 <- lm(log(invpc) ~ log(price) + deyear4, data=hseinv)

# Pretty regression table
stargazer(res1,res2,res2_4, type=”text”)

 

windows()
plot(hseinv$year,res2$residuals)

windows()
pacf <- pacf(res2$residuals)
pacf

   

 

windows()
acf <- acf(res2$residuals)
acf

write.csv(file=’res.csv’, res2$residuals) ###write把東西寫出去

res

res2 <- dynlm(log(invpc) ~ log(price) + trend(tsdata), data=tsdata)
res3 <- dynlm(log(invpc) ~ log(price) + trend(tsdata)+lag(log(price),1), data=tsdata)
stargazer(res2,res3, type=”text”)

 

windows()
ar1 <- pacf(res3$residuals)
ar1

 

res4 <- dynlm(log(invpc) ~ log(price) + trend(tsdata)
+lag(log(invpc),1)+lag(log(invpc),3), data=tsdata)
stargazer(res2,res3,res4, type=”text”)

 

windows()
pacf_13 <- pacf(res4$residuals)

 

windows()
acf_13 <- acf(res4$residuals)

 

windows()
ar13 <- pacf(res4$residuals)
ar13

 

plot(res4,1)


library(lmtest)
bptest(res1)


bptest(res2)


bptest(res1_1)


bptest(res2_1)


bptest(res2_2)


bptest(res2_3)

bptest(res2_4)

bptest(res3)

bptest(res4)

221201 計量經濟學_1
?>