Vector Autoregression or VAR can be estimated using several packages. We will use the “vars” and “tsDyn” packages to illustrate the application of VAR in R. In addition, we will apply some goodness of fit statistics and plots to assess the fit of the VAR model.
In this post, we will apply the VAR model in level as well as VAR in differences. For the goodness of fit, we will discuss the test and plot for stability check of the model and the Portmanteau test to check for serial correlation. After the estimation and goodness of fit, we will estimate the Impulse Response Functions (IRFs) and FEVDs in R.
Related Posts
Estimation, Goodness of Fit, Impulse responses after VAR in R
knitr::opts_chunk$set(echo = TRUE)
library(tidyverse)
library(urca)
library(stargazer)
library(vars)
library(haven)
library(tsDyn)
library(tsutils)
data(denmark)
data_var <- denmark[,c("LRM","LRY","IBO","IDE")]
view(data_var)
var_matrix <- cbind(denmark$LRM, denmark$LRY, denmark$IBO, denmark$IDE)
colnames(var_matrix, do.NULL = FALSE)## [1] "col1" "col2" "col3" "col4"colnames(var_matrix) <- c("LRM","LRY","IBO","IDE")Basic Plots
It is always a good option to plot the time series before moving on to any kind of analysis. Even basic plots can give a lot of information about the behavior of variables. Here, we will use simple line graphs using ts.plot to plot each time series.
The “par” function can be used to set several parameters for the graphs. We have 4 time series in this dataset and we can use the “mfrow” option to choose the layout of the 4 plots. The “mar” option can be used to set margins between the plots. We also have options to choose colors for the different regions of the plots such as the axis, lines, labels and so on.
par(col = "skyblue", col.axis = "mediumpurple", mfrow = c(2,2), mar = c(2,2,2,2))
ts.plot(data_var$LRM)
ts.plot(data_var$LRY)
ts.plot(data_var$IBO)
ts.plot(data_var$IDE)
Lag-length Selection
The first and one of the most important step in the application of VAR or Vector Autoregressive Model is selecting the lag-length. The “vars” package provides a functions for this purpose: “VARselect”. This functions uses various Information Criteria to choose appropriate lag-length for the VAR model. It reports the values of AIC (Akaike Information Criteria), SIC (Schwarz or Bayesian), HQIC (Hannan-Quinn) and FPE (Final Prediction Error) and the lags that minimize these criteria.
The “VARselect” function provides various useful options: 1) lag.max: the maximum number of lags to consider for the model.
type: the type of determininstic terms to include in the model. We can use “const” to include only the constant in the model, “trend” to include only trend, “both” to include constant as well as trend, and “none” to estimate without constant and trend.
season: allows the inclusion of seasonal dummies in the VAR model and its lag selection process.
exogen: this option can be used to include exogenous variables in the VAR model. You can construct a separate matrix of exogenous variables and use this option to include those in the model estimation and lag selection.
VARselect(var_matrix, lag.max = 5, type = "const", season = 4)## $selection
## AIC(n) HQ(n) SC(n) FPE(n)
## 2 1 1 2
##
## $criteria
## 1 2 3 4 5
## AIC(n) -3.511335e+01 -3.518424e+01 -3.500713e+01 -3.488959e+01 -3.483300e+01
## HQ(n) -3.464736e+01 -3.448526e+01 -3.407515e+01 -3.372462e+01 -3.343503e+01
## SC(n) -3.388966e+01 -3.334870e+01 -3.255974e+01 -3.183036e+01 -3.116192e+01
## FPE(n) 5.692226e-16 5.448356e-16 6.871711e-16 8.507607e-16 1.050901e-15VAR Estimation
The VAR function available in the “vars” package can be used to estimate the VAR model. Both the VARselect and VAR function come from the same package. The options available for both the functions are similar. The only difference is that the VAR function provides an option “p” to specify the lags of VAR model. The other options of type, season and exogen are similar to that in VARselect.
To customize the output table in a specific format, we can also use stargazer library and function on the results of the VAR model.
var1 <- VAR(var_matrix, p = 1, type = "const", season = 4)
summary(var1)##
## VAR Estimation Results:
## =========================
## Endogenous variables: LRM, LRY, IBO, IDE
## Deterministic variables: const
## Sample size: 54
## Log Likelihood: 664.799
## Roots of the characteristic polynomial:
## 0.9755 0.8495 0.8495 0.6655
## Call:
## VAR(y = var_matrix, p = 1, type = "const", season = 4L)
##
##
## Estimation results for equation LRM:
## ====================================
## LRM = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.820257 0.080604 10.176 2.33e-13 ***
## LRY.l1 0.081758 0.115447 0.708 0.482399
## IBO.l1 -1.186809 0.280062 -4.238 0.000107 ***
## IDE.l1 0.667722 0.390332 1.711 0.093883 .
## const 1.758377 0.484526 3.629 0.000711 ***
## sd1 -0.050303 0.009360 -5.374 2.48e-06 ***
## sd2 -0.023512 0.009202 -2.555 0.013991 *
## sd3 -0.035384 0.009054 -3.908 0.000304 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.02342 on 46 degrees of freedom
## Multiple R-Squared: 0.9796, Adjusted R-squared: 0.9765
## F-statistic: 316 on 7 and 46 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation LRY:
## ====================================
## LRY = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.092422 0.084742 1.091 0.281
## LRY.l1 0.756673 0.121373 6.234 1.29e-07 ***
## IBO.l1 0.029088 0.294439 0.099 0.922
## IDE.l1 -0.410897 0.410370 -1.001 0.322
## const 0.397467 0.509399 0.780 0.439
## sd1 -0.006643 0.009841 -0.675 0.503
## sd2 -0.005054 0.009675 -0.522 0.604
## sd3 -0.003697 0.009519 -0.388 0.700
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.02462 on 46 degrees of freedom
## Multiple R-Squared: 0.9016, Adjusted R-squared: 0.8867
## F-statistic: 60.23 on 7 and 46 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation IBO:
## ====================================
## IBO = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.0004092 0.0348901 0.012 0.9907
## LRY.l1 0.0161188 0.0499719 0.323 0.7485
## IBO.l1 1.1044825 0.1212271 9.111 7.22e-12 ***
## IDE.l1 -0.3259861 0.1689584 -1.929 0.0599 .
## const -0.0883127 0.2097310 -0.421 0.6757
## sd1 0.0015179 0.0040516 0.375 0.7097
## sd2 0.0062114 0.0039834 1.559 0.1258
## sd3 0.0047736 0.0039193 1.218 0.2294
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.01014 on 46 degrees of freedom
## Multiple R-Squared: 0.9077, Adjusted R-squared: 0.8937
## F-statistic: 64.64 on 7 and 46 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation IDE:
## ====================================
## IDE = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.012734 0.020553 0.620 0.53862
## LRY.l1 0.009339 0.029438 0.317 0.75249
## IBO.l1 0.215441 0.071414 3.017 0.00415 **
## IDE.l1 0.632235 0.099532 6.352 8.57e-08 ***
## const -0.205992 0.123551 -1.667 0.10226
## sd1 -0.002543 0.002387 -1.066 0.29220
## sd2 -0.002671 0.002347 -1.138 0.26091
## sd3 -0.001994 0.002309 -0.864 0.39222
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.005971 on 46 degrees of freedom
## Multiple R-Squared: 0.8478, Adjusted R-squared: 0.8247
## F-statistic: 36.62 on 7 and 46 DF, p-value: < 2.2e-16
##
##
##
## Covariance matrix of residuals:
## LRM LRY IBO IDE
## LRM 5.484e-04 3.350e-04 -1.188e-04 -5.001e-05
## LRY 3.350e-04 6.061e-04 -2.583e-05 -3.479e-05
## IBO -1.188e-04 -2.583e-05 1.027e-04 2.334e-05
## IDE -5.001e-05 -3.479e-05 2.334e-05 3.566e-05
##
## Correlation matrix of residuals:
## LRM LRY IBO IDE
## LRM 1.0000 0.5810 -0.5003 -0.3576
## LRY 0.5810 1.0000 -0.1035 -0.2366
## IBO -0.5003 -0.1035 1.0000 0.3856
## IDE -0.3576 -0.2366 0.3856 1.0000stargazer(var1$varresult, type = "text", style = "aer", title = "VAR(1)", column.labels = c("LRM","LRY","IBO","IDE"))##
## VAR(1)
## =======================================================================
## y
## LRM LRY IBO IDE
## (1) (2) (3) (4)
## -----------------------------------------------------------------------
## LRM.l1 0.820*** 0.092 0.000 0.013
## (0.081) (0.085) (0.035) (0.021)
##
## LRY.l1 0.082 0.757*** 0.016 0.009
## (0.115) (0.121) (0.050) (0.029)
##
## IBO.l1 -1.187*** 0.029 1.104*** 0.215***
## (0.280) (0.294) (0.121) (0.071)
##
## IDE.l1 0.668* -0.411 -0.326* 0.632***
## (0.390) (0.410) (0.169) (0.100)
##
## const 1.758*** 0.397 -0.088 -0.206
## (0.485) (0.509) (0.210) (0.124)
##
## sd1 -0.050*** -0.007 0.002 -0.003
## (0.009) (0.010) (0.004) (0.002)
##
## sd2 -0.024** -0.005 0.006 -0.003
## (0.009) (0.010) (0.004) (0.002)
##
## sd3 -0.035*** -0.004 0.005 -0.002
## (0.009) (0.010) (0.004) (0.002)
##
## Observations 54 54 54 54
## R2 0.980 0.902 0.908 0.848
## Adjusted R2 0.977 0.887 0.894 0.825
## Residual Std. Error (df = 46) 0.023 0.025 0.010 0.006
## F Statistic (df = 7; 46) 316.039*** 60.230*** 64.645*** 36.615***
## -----------------------------------------------------------------------
## Notes: ***Significant at the 1 percent level.
## **Significant at the 5 percent level.
## *Significant at the 10 percent level.VAR Goodness of fit
We can check the stability of the model with the help of “stability” function. The results of this function can also be plotted to check the stability condition. The “sctest” function can be used to estimate the p-value of stability test.
Testing for autocorrelation or serial correlation is essential in a VAR model. If the model suffers from serial correlation, we must increase the number of lags in the model. The function “serial.test” can be used to apply the Portmanteau Test. The option lags.pt is used to specify the number of lags for the Portmanteau Test.
stability_var1 <- stability(var1)
par(col = "skyblue", col.axis = "mediumpurple", mfrow = c(2,2))
plot(stability_var1$stability$LRM)
plot(stability_var1$stability$LRY)
plot(stability_var1$stability$IBO)
plot(stability_var1$stability$IDE)
sctest(var1)##
## M-fluctuation test
##
## data: var1
## f(efp) = 1.4194, p-value = 0.6862serial.test(var1, lags.pt = 8)##
## Portmanteau Test (asymptotic)
##
## data: Residuals of VAR object var1
## Chi-squared = 145.74, df = 112, p-value = 0.01767VAR(2) Estimation
The VAR(1) model or the VAR with 1 lag suffered from the problem of serial correlation. Therefore, we have to increase the number of lags in order to remove serial correlation.
We can estimate VAR(2) by specifying p = 2 in the VAR function. Everything else can stay the same because we only need to increase the lags.
The stability test and the Portmanteau test indicate that the model is stable and there is no serial correlation.
var2 <- VAR(var_matrix, p = 2, type = "const", season = 4)
summary(var2)##
## VAR Estimation Results:
## =========================
## Endogenous variables: LRM, LRY, IBO, IDE
## Deterministic variables: const
## Sample size: 53
## Log Likelihood: 678.644
## Roots of the characteristic polynomial:
## 0.9725 0.7713 0.7713 0.6734 0.6734 0.6051 0.2716 0.2716
## Call:
## VAR(y = var_matrix, p = 2, type = "const", season = 4L)
##
##
## Estimation results for equation LRM:
## ====================================
## LRM = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + LRM.l2 + LRY.l2 + IBO.l2 + IDE.l2 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 1.014228 0.201655 5.030 1.02e-05 ***
## LRY.l1 0.013753 0.166549 0.083 0.93459
## IBO.l1 -1.180148 0.393173 -3.002 0.00456 **
## IDE.l1 0.176409 0.598347 0.295 0.76961
## LRM.l2 -0.194958 0.176615 -1.104 0.27609
## LRY.l2 0.096016 0.157827 0.608 0.54630
## IBO.l2 0.138489 0.436917 0.317 0.75288
## IDE.l2 0.461713 0.576711 0.801 0.42798
## const 1.582925 0.547678 2.890 0.00613 **
## sd1 -0.055917 0.010563 -5.294 4.34e-06 ***
## sd2 -0.016458 0.009426 -1.746 0.08831 .
## sd3 -0.039480 0.008961 -4.406 7.40e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.02165 on 41 degrees of freedom
## Multiple R-Squared: 0.9842, Adjusted R-squared: 0.9799
## F-statistic: 231.9 on 11 and 41 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation LRY:
## ====================================
## LRY = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + LRM.l2 + LRY.l2 + IBO.l2 + IDE.l2 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.689838 0.206656 3.338 0.00180 **
## LRY.l1 0.646384 0.170679 3.787 0.00049 ***
## IBO.l1 0.280519 0.402923 0.696 0.49023
## IDE.l1 -0.587402 0.613185 -0.958 0.34370
## LRM.l2 -0.504019 0.180994 -2.785 0.00807 **
## LRY.l2 0.044561 0.161741 0.276 0.78431
## IBO.l2 0.377122 0.447752 0.842 0.40453
## IDE.l2 -0.060277 0.591012 -0.102 0.91926
## const -0.389553 0.561260 -0.694 0.49155
## sd1 -0.025121 0.010825 -2.321 0.02536 *
## sd2 0.007339 0.009660 0.760 0.45179
## sd3 -0.011369 0.009183 -1.238 0.22276
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.02218 on 41 degrees of freedom
## Multiple R-Squared: 0.9272, Adjusted R-squared: 0.9076
## F-statistic: 47.44 on 11 and 41 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation IBO:
## ====================================
## IBO = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + LRM.l2 + LRY.l2 + IBO.l2 + IDE.l2 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.0654218 0.0802292 0.815 0.4195
## LRY.l1 0.1179269 0.0662622 1.780 0.0825 .
## IBO.l1 1.3825683 0.1564254 8.839 4.83e-11 ***
## IDE.l1 0.0858928 0.2380548 0.361 0.7201
## LRM.l2 -0.0509340 0.0702669 -0.725 0.4727
## LRY.l2 -0.1356369 0.0627923 -2.160 0.0367 *
## IBO.l2 -0.3009861 0.1738292 -1.732 0.0909 .
## IDE.l2 -0.2532475 0.2294467 -1.104 0.2761
## const -0.0641545 0.2178960 -0.294 0.7699
## sd1 -0.0000689 0.0042027 -0.016 0.9870
## sd2 0.0073995 0.0037503 1.973 0.0553 .
## sd3 0.0048269 0.0035652 1.354 0.1832
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.008613 on 41 degrees of freedom
## Multiple R-Squared: 0.9401, Adjusted R-squared: 0.924
## F-statistic: 58.46 on 11 and 41 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation IDE:
## ====================================
## IDE = LRM.l1 + LRY.l1 + IBO.l1 + IDE.l1 + LRM.l2 + LRY.l2 + IBO.l2 + IDE.l2 + const + sd1 + sd2 + sd3
##
## Estimate Std. Error t value Pr(>|t|)
## LRM.l1 0.065001 0.050881 1.278 0.208607
## LRY.l1 -0.001606 0.042023 -0.038 0.969697
## IBO.l1 0.370309 0.099203 3.733 0.000576 ***
## IDE.l1 0.950624 0.150972 6.297 1.64e-07 ***
## LRM.l2 -0.068678 0.044563 -1.541 0.130962
## LRY.l2 0.021744 0.039822 0.546 0.588003
## IBO.l2 -0.227189 0.110241 -2.061 0.045699 *
## IDE.l2 -0.264860 0.145513 -1.820 0.076037 .
## const -0.071218 0.138187 -0.515 0.609062
## sd1 -0.004189 0.002665 -1.572 0.123709
## sd2 -0.001087 0.002378 -0.457 0.650085
## sd3 -0.002730 0.002261 -1.208 0.234120
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.005462 on 41 degrees of freedom
## Multiple R-Squared: 0.8862, Adjusted R-squared: 0.8557
## F-statistic: 29.04 on 11 and 41 DF, p-value: 7.092e-16
##
##
##
## Covariance matrix of residuals:
## LRM LRY IBO IDE
## LRM 4.686e-04 2.541e-04 -8.290e-05 -3.636e-05
## LRY 2.541e-04 4.921e-04 -1.586e-05 -2.883e-05
## IBO -8.290e-05 -1.586e-05 7.418e-05 1.178e-05
## IDE -3.636e-05 -2.883e-05 1.178e-05 2.983e-05
##
## Correlation matrix of residuals:
## LRM LRY IBO IDE
## LRM 1.0000 0.5292 -0.4447 -0.3075
## LRY 0.5292 1.0000 -0.0830 -0.2379
## IBO -0.4447 -0.0830 1.0000 0.2503
## IDE -0.3075 -0.2379 0.2503 1.0000stargazer(var2$varresult, type = "text", style = "aer", title = "VAR(2)", column.labels = c("LRM","LRY","IBO","IDE"))##
## VAR(2)
## =======================================================================
## y
## LRM LRY IBO IDE
## (1) (2) (3) (4)
## -----------------------------------------------------------------------
## LRM.l1 1.014*** 0.690*** 0.065 0.065
## (0.202) (0.207) (0.080) (0.051)
##
## LRY.l1 0.014 0.646*** 0.118* -0.002
## (0.167) (0.171) (0.066) (0.042)
##
## IBO.l1 -1.180*** 0.281 1.383*** 0.370***
## (0.393) (0.403) (0.156) (0.099)
##
## IDE.l1 0.176 -0.587 0.086 0.951***
## (0.598) (0.613) (0.238) (0.151)
##
## LRM.l2 -0.195 -0.504*** -0.051 -0.069
## (0.177) (0.181) (0.070) (0.045)
##
## LRY.l2 0.096 0.045 -0.136** 0.022
## (0.158) (0.162) (0.063) (0.040)
##
## IBO.l2 0.138 0.377 -0.301* -0.227**
## (0.437) (0.448) (0.174) (0.110)
##
## IDE.l2 0.462 -0.060 -0.253 -0.265*
## (0.577) (0.591) (0.229) (0.146)
##
## const 1.583*** -0.390 -0.064 -0.071
## (0.548) (0.561) (0.218) (0.138)
##
## sd1 -0.056*** -0.025** -0.000 -0.004
## (0.011) (0.011) (0.004) (0.003)
##
## sd2 -0.016* 0.007 0.007* -0.001
## (0.009) (0.010) (0.004) (0.002)
##
## sd3 -0.039*** -0.011 0.005 -0.003
## (0.009) (0.009) (0.004) (0.002)
##
## Observations 53 53 53 53
## R2 0.984 0.927 0.940 0.886
## Adjusted R2 0.980 0.908 0.924 0.856
## Residual Std. Error (df = 41) 0.022 0.022 0.009 0.005
## F Statistic (df = 11; 41) 231.913*** 47.439*** 58.465*** 29.037***
## -----------------------------------------------------------------------
## Notes: ***Significant at the 1 percent level.
## **Significant at the 5 percent level.
## *Significant at the 10 percent level.stability_var2 <- stability(var2)
par(col = "skyblue", col.axis = "mediumpurple", mfrow = c(2,2))
plot(stability_var2$stability$LRM)
plot(stability_var2$stability$LRY)
plot(stability_var2$stability$IBO)
plot(stability_var2$stability$IDE)
sctest(var2)##
## M-fluctuation test
##
## data: var2
## f(efp) = 1.4956, p-value = 0.6697serial.test(var2, lags.pt = 8)##
## Portmanteau Test (asymptotic)
##
## data: Residuals of VAR object var2
## Chi-squared = 110.87, df = 96, p-value = 0.1423Impulse Response Functions and FEVDs
One of the ,most important application of VAR models is the ability to estimate Impulse Response Functions and the Forecast Error Variance Decomposition. These provide valuable insights into the behavior of variables to shocks and can be useful in policy analysis.
irf_var2 <- irf(var2, n.ahead = 16, ortho = FALSE,impulse = c("LRM","LRY","IBO","IDE"), response = c("LRM","LRY"))
plot(irf_var2)
#irf_var2$irf$LRM
fevd_var2 <- fevd(var2, n.ahead = 8)
par(col = "mediumpurple", col.axis = "mediumpurple", mfrow = c(4,1), mar = c(2,2,2,2))
plot(fevd_var2)
fevd_var2$LRM## LRM LRY IBO IDE
## [1,] 1.0000000 0.000000000 0.0000000 0.000000000
## [2,] 0.9351761 0.002192469 0.0619773 0.000654151
## [3,] 0.7881634 0.011124917 0.1948010 0.005910660
## [4,] 0.6293382 0.022255308 0.3325729 0.015833585
## [5,] 0.5022989 0.029609745 0.4402478 0.027843600
## [6,] 0.4138380 0.032309199 0.5144322 0.039420631
## [7,] 0.3557796 0.031726551 0.5633801 0.049113727
## [8,] 0.3187657 0.029447462 0.5953646 0.056422190Forecasts
VAR models are often used for forecasting purposes and have been observed to provide good results. The predict function can be used to estimate the forecasts after VAR estimation. The option n.ahead is used to specify the number of time periods to forecast into the future.
for_var2 <- predict(var2, n.ahead = 8)
par(col = "skyblue", col.axis = "mediumpurple", mfrow = c(4,1), mar = c(2,2,2,2))
plot(for_var2)
fanchart(for_var2,
colors = c("skyblue4","skyblue3","skyblue2","skyblue1","skyblue",
"purple4","purple3","purple2","purple1"),
col.y = "skyblue")
VAR in differences
The model applied above is mis-specified because the variables are non-stationary at levels. In such a case, we can either use VAR in differences to ensure stationarity or apply VECM if the variables are cointegrated.
We can use the vars package to estimate VAR in differences by differencing the dataset using the “diff” function. The estimation can be carried out with VAR function using the differenced dataset. However, this method does not provide level forecasts when we use predict function after estimation. Instead, the output is forecasts in differenced variables and level forecasts must be calculated manually.
To overcome this, we can use the “lineVar” function from the tsDyn package. This function provides an additional option I = “diff to specify VAR in differences. This function can also be used to estimate a VECM. The disadvantage of this function is that there is no option to specify seasonal dummies like season option in vars package. Therefore, we have to include the seasonal dummies as exogenous variables with the help of exogen option.
time <- ts(denmark$ENTRY, deltat = 1/4, start = c(1974,1), end = c(1987, 3))
season_dummy <- seasdummy(nrow(data_var), y = time)
colnames(season_dummy) <- c("Q1", "Q2", "Q3")
view(season_dummy)
# var_exo_dummies <- VAR(var_matrix, p = 2, type = "const", exogen = season_dummy)
diff_var <- diff(as.matrix(data_var), lag = 1)
#view(diff_var)
VARselect(diff_var, lag.max = 5, type = "const", season = 4)## $selection
## AIC(n) HQ(n) SC(n) FPE(n)
## 1 1 1 1
##
## $criteria
## 1 2 3 4 5
## AIC(n) -3.471471e+01 -3.449955e+01 -3.412300e+01 -3.406793e+01 -3.384134e+01
## HQ(n) -3.424597e+01 -3.379645e+01 -3.318552e+01 -3.289609e+01 -3.243512e+01
## SC(n) -3.347923e+01 -3.264634e+01 -3.165205e+01 -3.097925e+01 -3.013491e+01
## FPE(n) 8.486279e-16 1.083141e-15 1.673727e-15 1.959679e-15 2.902021e-15var1_diff <- lineVar(var_matrix, lag = 1, include = "const", model = "VAR", I = "diff", exogen = season_dummy)
summary(var1_diff)## #############
## ###Model VAR
## #############
## Full sample size: 55 End sample size: 53
## Number of variables: 4 Number of estimated slope parameters 32
## AIC -1849.251 BIC -1786.202 SSR 0.05289925
##
## Intercept LRM -1 LRY -1
## Equation LRM 0.0344(0.0065)*** 0.3967(0.1664)* -0.1366(0.1561)
## Equation LRY 0.0072(0.0065) 0.4989(0.1648)** -0.1417(0.1546)
## Equation IBO -0.0047(0.0023). 0.0574(0.0596) 0.1410(0.0559)*
## Equation IDE 0.0016(0.0016) 0.0399(0.0419) 0.0172(0.0393)
## IBO -1 IDE -1 Q1
## Equation LRM -0.7524(0.3836). -0.3169(0.5768) -0.0660(0.0108)***
## Equation LRY 0.0864(0.3798) -0.4091(0.5712) -0.0213(0.0107).
## Equation IBO 0.3918(0.1374)** 0.1817(0.2067) 7.8e-06(0.0039)
## Equation IDE 0.3626(0.0965)*** 0.1599(0.1451) -0.0036(0.0027)
## Q2 Q3
## Equation LRM -0.0105(0.0098) -0.0422(0.0096)***
## Equation LRY 0.0041(0.0097) -0.0120(0.0095)
## Equation IBO 0.0073(0.0035)* 0.0046(0.0034)
## Equation IDE -0.0020(0.0025) -0.0027(0.0024)for_season <- seasdummy(nrow(data_var)+12, y = time)
for_hor_seas <- for_season[56:67,1:3]
var1_diff_fcst <- predict(var1_diff, n.ahead = 12, exoPred = for_hor_seas)
forecasts <- rbind(var_matrix, var1_diff_fcst)
par(mfrow = c(4,1), mar = c(2,2,2,2))
names <- c("LRM","LRY","IBO","IDE")
for(i in 1:4)
{
matplot(forecasts[,i], type=c("l"), col = c("skyblue"),
main = names[i] )
abline(v= nrow(var_matrix), col="mediumpurple")
}