library(readxl)
library(Metrics)
## Warning: package 'Metrics' was built under R version 4.4.1
library(tseries)
## Warning: package 'tseries' was built under R version 4.4.1
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(vars)
## Warning: package 'vars' was built under R version 4.4.1
## Loading required package: MASS
## Loading required package: strucchange
## Warning: package 'strucchange' was built under R version 4.4.1
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
## Loading required package: sandwich
## Warning: package 'sandwich' was built under R version 4.4.1
## Loading required package: urca
## Loading required package: lmtest
currencies <- read_excel("C:/Users/Admin/Downloads/currencies.xls")
currencies
## # A tibble: 7,142 × 4
## Date EUR GBP JPY
## <dttm> <dbl> <dbl> <dbl>
## 1 1998-12-14 00:00:00 0.847 0.593 117.
## 2 1998-12-15 00:00:00 0.849 0.595 116.
## 3 1998-12-16 00:00:00 0.851 0.597 116.
## 4 1998-12-17 00:00:00 0.848 0.594 115.
## 5 1998-12-18 00:00:00 0.850 0.594 116.
## 6 1998-12-19 00:00:00 0.849 0.594 116.
## 7 1998-12-20 00:00:00 0.851 0.594 115.
## 8 1998-12-21 00:00:00 0.855 0.595 117.
## 9 1998-12-22 00:00:00 0.854 0.596 116.
## 10 1998-12-23 00:00:00 0.858 0.597 116.
## # ℹ 7,132 more rows
adf.test(currencies$EUR)
##
## Augmented Dickey-Fuller Test
##
## data: currencies$EUR
## Dickey-Fuller = -1.6103, Lag order = 19, p-value = 0.7434
## alternative hypothesis: stationary
adf.test(currencies$GBP)
##
## Augmented Dickey-Fuller Test
##
## data: currencies$GBP
## Dickey-Fuller = -1.6029, Lag order = 19, p-value = 0.7466
## alternative hypothesis: stationary
adf.test(currencies$JPY)
##
## Augmented Dickey-Fuller Test
##
## data: currencies$JPY
## Dickey-Fuller = -1.7792, Lag order = 19, p-value = 0.6719
## alternative hypothesis: stationary
currencies
## # A tibble: 7,142 × 4
## Date EUR GBP JPY
## <dttm> <dbl> <dbl> <dbl>
## 1 1998-12-14 00:00:00 0.847 0.593 117.
## 2 1998-12-15 00:00:00 0.849 0.595 116.
## 3 1998-12-16 00:00:00 0.851 0.597 116.
## 4 1998-12-17 00:00:00 0.848 0.594 115.
## 5 1998-12-18 00:00:00 0.850 0.594 116.
## 6 1998-12-19 00:00:00 0.849 0.594 116.
## 7 1998-12-20 00:00:00 0.851 0.594 115.
## 8 1998-12-21 00:00:00 0.855 0.595 117.
## 9 1998-12-22 00:00:00 0.854 0.596 116.
## 10 1998-12-23 00:00:00 0.858 0.597 116.
## # ℹ 7,132 more rows
currencies$reur = c (NA ,100* diff ( log ( currencies$EUR ) ) )
currencies$rgbp = c (NA ,100* diff ( log ( currencies$GBP ) ) )
currencies$rjpy = c (NA ,100* diff ( log ( currencies$JPY ) ) )
currencies = currencies [ -1 ,]
adf.test(currencies$reur)
## Warning in adf.test(currencies$reur): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: currencies$reur
## Dickey-Fuller = -18.535, Lag order = 19, p-value = 0.01
## alternative hypothesis: stationary
adf.test(currencies$rgbp)
## Warning in adf.test(currencies$rgbp): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: currencies$rgbp
## Dickey-Fuller = -19.246, Lag order = 19, p-value = 0.01
## alternative hypothesis: stationary
adf.test(currencies$rjpy)
## Warning in adf.test(currencies$rjpy): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: currencies$rjpy
## Dickey-Fuller = -18.927, Lag order = 19, p-value = 0.01
## alternative hypothesis: stationary
VAR(currencies[c("reur","rgbp","rjpy") ] , p = 2)
##
## VAR Estimation Results:
## =======================
##
## Estimated coefficients for equation reur:
## =========================================
## Call:
## reur = reur.l1 + rgbp.l1 + rjpy.l1 + reur.l2 + rgbp.l2 + rjpy.l2 + const
##
## reur.l1 rgbp.l1 rjpy.l1 reur.l2 rgbp.l2
## 0.1474965506 -0.0183557359 -0.0070978299 -0.0118082131 0.0066228786
## rjpy.l2 const
## -0.0054272893 0.0001368541
##
##
## Estimated coefficients for equation rgbp:
## =========================================
## Call:
## rgbp = reur.l1 + rgbp.l1 + rjpy.l1 + reur.l2 + rgbp.l2 + rjpy.l2 + const
##
## reur.l1 rgbp.l1 rjpy.l1 reur.l2 rgbp.l2 rjpy.l2
## -0.025270542 0.221361804 -0.039016040 0.046926637 -0.067794083 0.003286917
## const
## 0.002825689
##
##
## Estimated coefficients for equation rjpy:
## =========================================
## Call:
## rjpy = reur.l1 + rgbp.l1 + rjpy.l1 + reur.l2 + rgbp.l2 + rjpy.l2 + const
##
## reur.l1 rgbp.l1 rjpy.l1 reur.l2 rgbp.l2
## 0.0410612623 -0.0708456993 0.1324572148 -0.0188915075 0.0249075437
## rjpy.l2 const
## 0.0149565492 -0.0004125908
VARselect(currencies[c("reur","rgbp","rjpy")], lag.max = 10)
## $selection
## AIC(n) HQ(n) SC(n) FPE(n)
## 4 2 1 4
##
## $criteria
## 1 2 3 4 5
## AIC(n) -5.433453428 -5.43726462 -5.437666647 -5.438909241 -5.438464346
## HQ(n) -5.429472196 -5.43029747 -5.427713567 -5.425970237 -5.422539418
## SC(n) -5.421888924 -5.41702674 -5.408755386 -5.401324601 -5.392206329
## FPE(n) 0.004367985 0.00435137 0.004349621 0.004344219 0.004346152
## 6 7 8 9 10
## AIC(n) -5.437389700 -5.436794113 -5.435808057 -5.435335655 -5.433668221
## HQ(n) -5.418478849 -5.414897337 -5.410925357 -5.407467031 -5.402813673
## SC(n) -5.382458304 -5.373189339 -5.363529904 -5.354384124 -5.344043311
## FPE(n) 0.004350826 0.004353418 0.004357713 0.004359772 0.004367048
model = VAR(currencies[c("reur","rgbp","rjpy") ] , p = 1) #Given name to call it latter
model
##
## VAR Estimation Results:
## =======================
##
## Estimated coefficients for equation reur:
## =========================================
## Call:
## reur = reur.l1 + rgbp.l1 + rjpy.l1 + const
##
## reur.l1 rgbp.l1 rjpy.l1 const
## 0.1456393628 -0.0168336370 -0.0081023309 0.0001747216
##
##
## Estimated coefficients for equation rgbp:
## =========================================
## Call:
## rgbp = reur.l1 + rgbp.l1 + rjpy.l1 + const
##
## reur.l1 rgbp.l1 rjpy.l1 const
## -0.015871266 0.206799579 -0.036538233 0.002664532
##
##
## Estimated coefficients for equation rjpy:
## =========================================
## Call:
## rjpy = reur.l1 + rgbp.l1 + rjpy.l1 + const
##
## reur.l1 rgbp.l1 rjpy.l1 const
## 0.0378371246 -0.0663624959 0.1336151098 -0.0003202048
causality(model, cause=c("reur", "rgbp"))$Granger
##
## Granger causality H0: reur rgbp do not Granger-cause rjpy
##
## data: VAR object model
## F-Test = 7.7511, df1 = 2, df2 = 21408, p-value = 0.0004315
causality(model, cause=c("reur"))$Granger
##
## Granger causality H0: reur do not Granger-cause rgbp rjpy
##
## data: VAR object model
## F-Test = 4.1141, df1 = 2, df2 = 21408, p-value = 0.01635
grangertest(reur ~ rgbp, order = 1, data=currencies)
## Granger causality test
##
## Model 1: reur ~ Lags(reur, 1:1) + Lags(rgbp, 1:1)
## Model 2: reur ~ Lags(reur, 1:1)
## Res.Df Df F Pr(>F)
## 1 7137
## 2 7138 -1 0.9912 0.3195
When test restiction: using F-test –> model 1 có restrict, model 2 không có restrict Q: y_1 Granger cause y_2 equivalent test H0: y_1 does not Granger cause y_2 –> test H0: alpha_12 = beta_12 = 0 using p-value = 0.3 > 0.05 –> not reject H0 –> y_1 does not Granger cause y_2 (phần 13 trong slide) –> 3! granger test
grangertest(rjpy ~ rgbp, order = 1, data=currencies)
## Granger causality test
##
## Model 1: rjpy ~ Lags(rjpy, 1:1) + Lags(rgbp, 1:1)
## Model 2: rjpy ~ Lags(rjpy, 1:1)
## Res.Df Df F Pr(>F)
## 1 7137
## 2 7138 -1 9.69 0.00186 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(reur ~ rjpy, order = 1, data=currencies)
## Granger causality test
##
## Model 1: reur ~ Lags(reur, 1:1) + Lags(rjpy, 1:1)
## Model 2: reur ~ Lags(reur, 1:1)
## Res.Df Df F Pr(>F)
## 1 7137
## 2 7138 -1 0.4202 0.5169
grangertest(rjpy ~ reur, order = 1, data=currencies)
## Granger causality test
##
## Model 1: rjpy ~ Lags(rjpy, 1:1) + Lags(reur, 1:1)
## Model 2: rjpy ~ Lags(rjpy, 1:1)
## Res.Df Df F Pr(>F)
## 1 7137
## 2 7138 -1 5e-04 0.9815
grangertest(reur ~ rgbp, order = 1, data=currencies)
## Granger causality test
##
## Model 1: reur ~ Lags(reur, 1:1) + Lags(rgbp, 1:1)
## Model 2: reur ~ Lags(reur, 1:1)
## Res.Df Df F Pr(>F)
## 1 7137
## 2 7138 -1 0.9912 0.3195
ir = irf(model, n.ahead = 20)
plot(ir)
vd = fevd(model, n.ahead = 20)
par(mar=c(1,1,1,1))
plot(vd)
forecasts1 <- predict(model, n.ahead = 12)
forecasts1[1]
## $fcst
## $fcst$reur
## fcst lower upper CI
## [1,] -1.659409e-02 -0.9178834 0.8846953 0.9012893
## [2,] -1.823675e-03 -0.9112383 0.9075909 0.9094146
## [3,] -6.306415e-05 -0.9096242 0.9094981 0.9095611
## [4,] 1.282602e-04 -0.9094356 0.9096921 0.9095639
## [5,] 1.443687e-04 -0.9094196 0.9097083 0.9095639
## [6,] 1.444208e-04 -0.9094195 0.9097084 0.9095639
## [7,] 1.439514e-04 -0.9094200 0.9097079 0.9095639
## [8,] 1.437791e-04 -0.9094202 0.9097077 0.9095639
## [9,] 1.437307e-04 -0.9094202 0.9097077 0.9095639
## [10,] 1.437183e-04 -0.9094202 0.9097077 0.9095639
## [11,] 1.437152e-04 -0.9094202 0.9097076 0.9095639
## [12,] 1.437145e-04 -0.9094202 0.9097076 0.9095639
##
## $fcst$rgbp
## fcst lower upper CI
## [1,] -0.022958915 -0.8322457 0.7863279 0.8092868
## [2,] -0.001676255 -0.8261020 0.8227494 0.8244257
## [3,] 0.002345004 -0.8227657 0.8274557 0.8251107
## [4,] 0.003160391 -0.8219847 0.8283055 0.8251451
## [5,] 0.003334861 -0.8218120 0.8284818 0.8251469
## [6,] 0.003373586 -0.8217734 0.8285206 0.8251470
## [7,] 0.003382381 -0.8217646 0.8285294 0.8251470
## [8,] 0.003384407 -0.8217626 0.8285314 0.8251470
## [9,] 0.003384877 -0.8217621 0.8285319 0.8251470
## [10,] 0.003384987 -0.8217620 0.8285320 0.8251470
## [11,] 0.003385012 -0.8217620 0.8285320 0.8251470
## [12,] 0.003385018 -0.8217620 0.8285320 0.8251470
##
## $fcst$rjpy
## fcst lower upper CI
## [1,] -3.933914e-03 -0.9186254 0.9107576 0.9146915
## [2,] 4.990315e-05 -0.9237681 0.9238679 0.9238180
## [3,] -2.712992e-04 -0.9243809 0.9238383 0.9241096
## [4,] -5.144610e-04 -0.9246376 0.9236087 0.9241232
## [5,] -5.938230e-04 -0.9247177 0.9235301 0.9241239
## [6,] -6.153958e-04 -0.9247393 0.9235085 0.9241239
## [7,] -6.208461e-04 -0.9247448 0.9235031 0.9241239
## [8,] -6.221758e-04 -0.9247461 0.9235018 0.9241239
## [9,] -6.224944e-04 -0.9247464 0.9235014 0.9241239
## [10,] -6.225700e-04 -0.9247465 0.9235014 0.9241239
## [11,] -6.225879e-04 -0.9247465 0.9235014 0.9241239
## [12,] -6.225921e-04 -0.9247465 0.9235014 0.9241239
forecasts1$fcst$reur[,1]
## [1] -1.659409e-02 -1.823675e-03 -6.306415e-05 1.282602e-04 1.443687e-04
## [6] 1.444208e-04 1.439514e-04 1.437791e-04 1.437307e-04 1.437183e-04
## [11] 1.437152e-04 1.437145e-04
mae_reur <- mae(currencies$reur, forecasts1$fcst$reur[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
rmse_reur <- rmse(currencies$reur, forecasts1$fcst$reur[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
mae_rgbp <- mae(currencies$rgbp, forecasts1$fcst$rgbp[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
rmse_rgbp <- rmse(currencies$rgbp, forecasts1$fcst$rgbp[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
mae_rjpy <- mae(currencies$rjpy, forecasts1$fcst$rjpy[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
rmse_rjpy <- rmse(currencies$rjpy, forecasts1$fcst$rjpy[, 1])
## Warning in actual - predicted: longer object length is not a multiple of
## shorter object length
print(paste("MAE (reur):", mae_reur))
## [1] "MAE (reur): 0.311525918482539"
print(paste("RMSE (reur):", rmse_reur))
## [1] "RMSE (reur): 0.463932234934771"
print(paste("MAE (rgbp):", mae_rgbp))
## [1] "MAE (rgbp): 0.276267294860245"
print(paste("RMSE (rgbp):", rmse_rgbp))
## [1] "RMSE (rgbp): 0.420874426417948"
print(paste("MAE (rjpy):", mae_rjpy))
## [1] "MAE (rjpy): 0.314258624083522"
print(paste("RMSE (rjpy):", rmse_rjpy))
## [1] "RMSE (rjpy): 0.471385335908988"