Made by 1A182901-2 Yen Cheng Hsuan
Library Packages
library(readr)
library(dplyr)
library(astsa)
library(forecast)
library(ggplot2)
library(SciViews)
Import/Clean data
data <- read.csv("~/Desktop/早稲田/2018 Spring Course/TSA/TSA2018_PS1_Data.csv")
ffr <- as.numeric(data$FFR)
date <- as.character(as.Date(data$Date))
res <- as.data.frame(cbind(date,ffr))
res$date <- as.Date(res$date)
res$ffr <- as.numeric(as.character(res$ffr))
3.a
first_graph <- ggplot(res,aes(res$date,res$ffr))+
xlab("Date")+
ylab("FFR(%)")+
geom_line()+
scale_x_date(date_breaks="5 year",date_labels = "%Y")
first_graph
#ggsave("3.a.png", width = 9, height = 6)
3.b
data$IRS <- data$r5 - data$Tbill
irs <- data$IRS
res2 <- as.data.frame(cbind(date,irs))
res2$date <- as.Date(res2$date)
res2$irs <- as.numeric(as.character(res2$irs))
second_graph <- ggplot(res2,aes(res2$date,res2$irs))+
xlab("Date")+
ylab("Interest Rate Spread")+
geom_line()+
scale_x_date(date_breaks="5 year",date_labels = "%Y")
second_graph
#ggsave("3.b.png", width = 9, height = 6)
3.c
t <- data$X
res2 <- cbind(t,res2)
res3 <- as.data.frame(acf2(res2$irs))
res3
## ACF PACF
## 1 0.86 0.86
## 2 0.68 -0.21
## 3 0.55 0.11
## 4 0.41 -0.18
## 5 0.28 -0.01
## 6 0.15 -0.14
## 7 0.07 0.14
## 8 0.04 0.01
## 9 -0.03 -0.18
## 10 -0.13 -0.12
## 11 -0.20 -0.04
## 12 -0.22 0.08
## 13 -0.24 -0.08
## 14 -0.26 0.03
## 15 -0.25 0.00
## 16 -0.22 -0.02
## 17 -0.19 -0.01
## 18 -0.13 0.17
## 19 -0.05 0.05
## 20 0.00 -0.06
## 21 0.04 -0.03
## 22 0.06 -0.07
## 23 0.07 0.07
## 24 0.10 0.01
## 25 0.10 0.01
#write.csv(res3,file = "3.c.csv")
3.d
res4 <- arima(res2$irs,c(7,0,0))
res5 <- arima(res2$irs,c(1,0,1))
res4
##
## Call:
## arima(x = res2$irs, order = c(7, 0, 0))
##
## Coefficients:
## ar1 ar2 ar3 ar4 ar5 ar6 ar7 intercept
## 1.1057 -0.4441 0.3891 -0.2899 0.2156 -0.2928 0.1328 1.1970
## s.e. 0.0680 0.1003 0.1038 0.1050 0.1031 0.0994 0.0675 0.1689
##
## sigma^2 estimated as 0.2087: log likelihood = -135.66, aic = 289.32
res5
##
## Call:
## arima(x = res2$irs, order = c(1, 0, 1))
##
## Coefficients:
## ar1 ma1 intercept
## 0.7572 0.3771 1.1934
## s.e. 0.0541 0.0952 0.1813
##
## sigma^2 estimated as 0.2235: log likelihood = -142.75, aic = 293.5
3.e
MSPE
res6.1 <- filter(res2,t<163)
res6.2 <- filter(res2,t>162)
res7 <- arima(res6.1$irs,c(7,0,0))
res8 <- arima(res6.1$irs,c(1,0,1))
res9 <- predict(res7, n.ahead = 50,newxreg = NULL,se.fit = TRUE)
res10 <- predict(res8, n.ahead = 50,newxreg = NULL,se.fit = TRUE)
res_pr_ar7 <- as.data.frame(res9$pred[1:50])
res_pr_ar7 <- as.data.frame(cbind(res_pr_ar7$`res9$pred[1:50]`,res6.2$irs))
colnames(res_pr_ar7) <- c("pre","real")
res_pr_ar1ma1 <- as.data.frame(res10$pred[1:50])
res_pr_ar1ma1 <- as.data.frame(cbind(res_pr_ar1ma1$`res10$pred[1:50]`,res6.2$irs))
colnames(res_pr_ar1ma1) <- c("pre","real")
res_pr_ar7$er <- ((res_pr_ar7$pre) - (res_pr_ar7$real))^2
res_pr_ar1ma1$er <- ((res_pr_ar1ma1$pre) - (res_pr_ar1ma1$real))^2
MSPE_ar7 <- sum(res_pr_ar7$er)/50
MSPE_ar1ma1 <- sum(res_pr_ar1ma1$er)/50
cat("MSPE of ar7 = ",MSPE_ar7)
## MSPE of ar7 = 0.7461955
cat("\nMSPE of ar1ma1 = ",MSPE_ar1ma1)
##
## MSPE of ar1ma1 = 0.7415018
GrangerNewbold
gnt1 <- list()
i=1
for (i in 1:50){
res_mo <- filter(res2,t<162+i)
res_ar7 <- arima(res_mo$irs,c(7,0,0))
res_ar1ma1 <- arima(res_mo$irs,c(1,0,1))
res_p_ar7 <- predict(res_ar7, n.ahead = 1,newxreg = NULL,se.fit = TRUE)
res_p_ar1ma1 <- predict(res_ar1ma1, n.ahead = 1,newxreg = NULL,se.fit = TRUE)
gnt1$e_ar7[i] <- res2$irs[162+i] - (res_p_ar7$pred)
gnt1$e_ar1ma1[i] <- res2$irs[162+i] - (res_p_ar1ma1$pred)
}
gntdf <- as.data.frame(gnt1)
gntdf$xi <- gntdf$e_ar7 + gntdf$e_ar1ma1
gntdf$zi <- gntdf$e_ar7 - gntdf$e_ar1ma1
gn <- cor.test(gntdf$xi,gntdf$zi)
gn
##
## Pearson's product-moment correlation
##
## data: gntdf$xi and gntdf$zi
## t = 1.2276, df = 48, p-value = 0.2256
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.1091767 0.4318484
## sample estimates:
## cor
## 0.1744731
3.f
cpinf <- list()
i=1
for (i in 1:211){
cpinf[i] <- log(data$CPI[i+1]/data$CPI[i])*100
}
cpinf <- as.numeric(cpinf)
acf2(cpinf)
## ACF PACF
## [1,] 0.76 0.76
## [2,] 0.65 0.17
## [3,] 0.65 0.29
## [4,] 0.57 -0.06
## [5,] 0.54 0.11
## [6,] 0.52 0.01
## [7,] 0.44 -0.07
## [8,] 0.37 -0.07
## [9,] 0.34 -0.03
## [10,] 0.33 0.08
## [11,] 0.30 -0.01
## [12,] 0.29 0.08
## [13,] 0.28 0.03
## [14,] 0.30 0.12
## [15,] 0.27 -0.06
## [16,] 0.27 0.04
## [17,] 0.28 0.01
## [18,] 0.27 0.00
## [19,] 0.27 0.01
## [20,] 0.30 0.09
## [21,] 0.30 0.03
## [22,] 0.30 0.06
## [23,] 0.26 -0.13
## [24,] 0.24 -0.01
## [25,] 0.19 -0.16
res11 <- arima(cpinf, c(3,0,0))
cat("\nResult of ar(3)\n\n")
##
## Result of ar(3)
res11
##
## Call:
## arima(x = cpinf, order = c(3, 0, 0))
##
## Coefficients:
## ar1 ar2 ar3 intercept
## 0.5767 -0.0152 0.2953 0.9344
## s.e. 0.0655 0.0769 0.0657 0.2150
##
## sigma^2 estimated as 0.218: log likelihood = -139.27, aic = 288.54
cat("BIC",BIC(res11),"\n")
## BIC 305.3004
cat("\nResult of auto.arima\n\n")
##
## Result of auto.arima
auto.arima(cpinf,max.d = 0)
## Series: cpinf
## ARIMA(2,0,2) with non-zero mean
##
## Coefficients:
## ar1 ar2 ma1 ma2 mean
## 0.4793 0.4234 0.1325 -0.3878 0.9228
## s.e. 0.1559 0.1377 0.1457 0.0694 0.2339
##
## sigma^2 estimated as 0.2234: log likelihood=-139.33
## AIC=290.65 AICc=291.07 BIC=310.77