Question 1.Obtain monthly time series for Crude Oil Prices: West Texas Intermediate (WTI), FRED/MCOILWTICO, and the monthly time series for US Regular Conventional Gas Price,FRED/GASREGCOVM. 2.Create a single time series plot with two log prices log for the sample 1995M1-2017M4. 3.Perform unit root tests for both logs. 4.Determine the number of lags to include in cointegration analysis using Schwarz information criterion. 5. estimate a bivariate VEC model for the variables log. 6.Are the adjustment parameters α1 and α2 in the estimated VEC model statistically significant? Are their signs consistent with error correction mechanism that moves the system back to the long run equilibrium, whenever there is a disruption ? 7.eestimate the VEC model with a restriction α2= 0. 8.What is the intuition for imposing the restriction in (7)? 9.Create and plot sequence of one month ahead forecasts (using either fixed, rolling or recursive scheme) for the period 2011M1-2017M4. Obtain the RMSE for the forecast of the gas price,p GAS t. 10.Plot the ACF and PACF . 11.Create a sequence of one month ahead forecast. 12.Compare the RMSE of the forecast. 13.Based on your results in (12), is there any advantage to using the VEC model to forecast the price of the gas, instead of using a simple univariate AR/MA/ARMA model?
library(Quandl)
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(zoo)
library(xts)
library(dygraphs)
library(knitr)
library(forecast)
library(urca)
library(vars)
## Loading required package: MASS
## Loading required package: strucchange
## Loading required package: sandwich
## Loading required package: lmtest
library(tseries)
Quandl.api_key('rxVQhZ8_nxxeo2yy4Uz4')
COP <- Quandl("FRED/MCOILWTICO", type="ts")
RCGP <- Quandl("FRED/GASREGCOVM", type="ts")
ly1=100*log(COP)
ly2=100*log(RCGP)
ly1 <- window(ly1, start=c(1995,1), end=c(2017,4))
## Warning in window.default(x, ...): 'end' value not changed
ly2 <- window(ly2, start=c(1995,1), end=c(2017,4))
plot(cbind(ly1,ly2), col=c(2,4), plot.type = "single", xlab="Years 1995-2017")
legend("right", c(expression(Log_Crude_Oil_Prices),expression(Log_Regular_Conventional_Gas_Price)), col=c(2,4), lty=1, bty = "n")
dl.y1 <- diff(ly1,lag=1)
dl.y2 <- diff(ly2,lag=1)
par(mfrow=c(2,1))
plot(dl.y1, xlab=" Years 1995-2017", ylab="", main="Diff Log Crude Oil Prices")
plot(dl.y2, xlab=" Years 1995-2017", ylab="", main="Diff Log Regular Conventional Gas Price")
adf.test(ly1)
##
## Augmented Dickey-Fuller Test
##
## data: ly1
## Dickey-Fuller = -1.6586, Lag order = 6, p-value = 0.7197
## alternative hypothesis: stationary
kpss.test(ly1, null="Trend")
## Warning in kpss.test(ly1, null = "Trend"): p-value smaller than printed p-
## value
##
## KPSS Test for Trend Stationarity
##
## data: ly1
## KPSS Trend = 0.80651, Truncation lag parameter = 3, p-value = 0.01
y1.urkpss <- ur.kpss(ly1, type="tau", lags="short")
summary(y1.urkpss)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.5563
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
y1.ers1 <- ur.ers(ly1, type="P-test", model="trend")
summary(y1.ers1)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 10.1644
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
All tests (ADF, KPSS, ERS) show that y1 is not stationary.
adf.test(dl.y1)
## Warning in adf.test(dl.y1): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: dl.y1
## Dickey-Fuller = -6.63, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationary
kpss.test(dl.y1, null="Trend")
## Warning in kpss.test(dl.y1, null = "Trend"): p-value greater than printed
## p-value
##
## KPSS Test for Trend Stationarity
##
## data: dl.y1
## KPSS Trend = 0.052193, Truncation lag parameter = 3, p-value = 0.1
dl.y1.urkpss <- ur.kpss(dl.y1, type="tau", lags="short")
summary(dl.y1.urkpss)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.0514
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
dl.y1.ers1 <- ur.ers(dl.y1, type="P-test", model="trend")
summary(dl.y1.ers1)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 0.5881
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
All tests Y1 is stationary.
adf.test(ly2)
##
## Augmented Dickey-Fuller Test
##
## data: ly2
## Dickey-Fuller = -1.6181, Lag order = 6, p-value = 0.7368
## alternative hypothesis: stationary
kpss.test(ly2, null="Trend")
## Warning in kpss.test(ly2, null = "Trend"): p-value smaller than printed p-
## value
##
## KPSS Test for Trend Stationarity
##
## data: ly2
## KPSS Trend = 0.7202, Truncation lag parameter = 3, p-value = 0.01
y2.urkpss <- ur.kpss(ly2, type="tau", lags="short")
summary(y2.urkpss)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.5035
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
y2.ers1 <- ur.ers(ly2, type="P-test", model="trend")
summary(y2.ers1)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 8.1303
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
All three performed tests (ADF, KPSS, ERS) sugeest that y2 is not stationary.
adf.test(dl.y2)
## Warning in adf.test(dl.y2): p-value smaller than printed p-value
##
## Augmented Dickey-Fuller Test
##
## data: dl.y2
## Dickey-Fuller = -7.9022, Lag order = 6, p-value = 0.01
## alternative hypothesis: stationary
kpss.test(dl.y2, null="Trend")
## Warning in kpss.test(dl.y2, null = "Trend"): p-value greater than printed
## p-value
##
## KPSS Test for Trend Stationarity
##
## data: dl.y2
## KPSS Trend = 0.037877, Truncation lag parameter = 3, p-value = 0.1
dl.y2.urkpss <- ur.kpss(dl.y2, type="tau", lags="short")
summary(dl.y2.urkpss)
##
## #######################
## # KPSS Unit Root Test #
## #######################
##
## Test is of type: tau with 5 lags.
##
## Value of test-statistic is: 0.0403
##
## Critical value for a significance level of:
## 10pct 5pct 2.5pct 1pct
## critical values 0.119 0.146 0.176 0.216
dl.y2.ers1 <- ur.ers(dl.y2, type="P-test", model="trend")
summary(dl.y2.ers1)
##
## ###############################################
## # Elliot, Rothenberg and Stock Unit Root Test #
## ###############################################
##
## Test of type P-test
## detrending of series with intercept and trend
##
## Value of test-statistic is: 0.0533
##
## Critical values of P-test are:
## 1pct 5pct 10pct
## critical values 3.96 5.62 6.89
first difference logarithmic of y2 is stationary.
y <- cbind(ly1, ly2)
y <- window(y, start=c(1995,1), end=c(2010,12))
VARselect(y, type="const", lag.max=12)
## $selection
## AIC(n) HQ(n) SC(n) FPE(n)
## 6 6 2 6
##
## $criteria
## 1 2 3 4 5
## AIC(n) 7.353484 7.050363 7.000602 7.015203 6.968809
## HQ(n) 7.396638 7.122286 7.101294 7.144664 7.127039
## SC(n) 7.459916 7.227750 7.248944 7.334499 7.359059
## FPE(n) 1561.637976 1153.310834 1097.380193 1113.619041 1063.280429
## 6 7 8 9 10
## AIC(n) 6.908009 6.931080 6.921235 6.926159 6.945912
## HQ(n) 7.095008 7.146848 7.165772 7.199464 7.247987
## SC(n) 7.369214 7.463239 7.524349 7.600227 7.690935
## FPE(n) 1000.758328 1024.392235 1014.716198 1020.183056 1041.114353
## 11 12
## AIC(n) 6.956168 6.983254
## HQ(n) 7.287012 7.342867
## SC(n) 7.772146 7.870186
## FPE(n) 1052.556306 1082.326306
Schwarz information criterion suggests using 2 lags.
y.CA <- ca.jo(y, ecdet="const", type="trace", K=2, spec="transitory", season=12)
summary(y.CA)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , without linear trend and constant in cointegration
##
## Eigenvalues (lambda):
## [1] 1.980704e-01 9.957529e-03 -3.903128e-18
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 1.90 7.52 9.24 12.97
## r = 0 | 43.84 17.85 19.96 24.60
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## ly1.l1 ly2.l1 constant
## ly1.l1 1.000000 1.0000000 1.0000000
## ly2.l1 -1.587125 -0.1377993 0.2504144
## constant -275.337070 -429.0868103 -322.0681760
##
## Weights W:
## (This is the loading matrix)
##
## ly1.l1 ly2.l1 constant
## ly1.d 0.04824429 -0.008302568 1.915537e-16
## ly2.d 0.24802732 -0.003053918 9.713652e-16
y.CA <- ca.jo(y, ecdet="const", type="eigen", K=2, spec="transitory", season=12)
summary(y.CA)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: maximal eigenvalue statistic (lambda max) , without linear trend and constant in cointegration
##
## Eigenvalues (lambda):
## [1] 1.980704e-01 9.957529e-03 -3.903128e-18
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 1.90 7.52 9.24 12.97
## r = 0 | 41.94 13.75 15.67 20.20
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## ly1.l1 ly2.l1 constant
## ly1.l1 1.000000 1.0000000 1.0000000
## ly2.l1 -1.587125 -0.1377993 0.2504144
## constant -275.337070 -429.0868103 -322.0681760
##
## Weights W:
## (This is the loading matrix)
##
## ly1.l1 ly2.l1 constant
## ly1.d 0.04824429 -0.008302568 1.915537e-16
## ly2.d 0.24802732 -0.003053918 9.713652e-16
y.VEC <- cajorls(y.CA, r=1)
y.VEC
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## ly1.d ly2.d
## ect1 0.04824 0.24803
## sd1 6.08302 4.57391
## sd2 2.31774 2.22917
## sd3 8.03228 6.21726
## sd4 4.09514 5.92019
## sd5 5.27266 6.26258
## sd6 4.71141 4.85558
## sd7 4.94215 3.27466
## sd8 5.02513 4.23428
## sd9 3.40821 4.01838
## sd10 1.76899 -0.01183
## sd11 -0.20628 -0.33635
## ly1.dl1 0.27002 0.16276
## ly2.dl1 -0.11007 0.18754
##
##
## $beta
## ect1
## ly1.l1 1.000000
## ly2.l1 -1.587125
## constant -275.337070
summary(y.VEC$rlm)
## Response ly1.d :
##
## Call:
## lm(formula = ly1.d ~ ect1 + sd1 + sd2 + sd3 + sd4 + sd5 + sd6 +
## sd7 + sd8 + sd9 + sd10 + sd11 + ly1.dl1 + ly2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -26.112 -4.685 1.140 6.214 19.117
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.04824 0.08132 0.593 0.55375
## sd1 6.08302 3.00078 2.027 0.04416 *
## sd2 2.31774 3.05064 0.760 0.44842
## sd3 8.03228 2.99738 2.680 0.00807 **
## sd4 4.09514 3.07874 1.330 0.18520
## sd5 5.27266 3.14608 1.676 0.09552 .
## sd6 4.71141 3.17557 1.484 0.13969
## sd7 4.94215 3.14817 1.570 0.11825
## sd8 5.02513 3.07894 1.632 0.10445
## sd9 3.40821 3.04442 1.119 0.26446
## sd10 1.76899 3.03716 0.582 0.56101
## sd11 -0.20628 2.96477 -0.070 0.94461
## ly1.dl1 0.27002 0.11374 2.374 0.01867 *
## ly2.dl1 -0.11007 0.14766 -0.745 0.45702
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.325 on 176 degrees of freedom
## Multiple R-squared: 0.134, Adjusted R-squared: 0.06508
## F-statistic: 1.945 on 14 and 176 DF, p-value: 0.02469
##
##
## Response ly2.d :
##
## Call:
## lm(formula = ly2.d ~ ect1 + sd1 + sd2 + sd3 + sd4 + sd5 + sd6 +
## sd7 + sd8 + sd9 + sd10 + sd11 + ly1.dl1 + ly2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -19.3100 -2.3394 0.3539 3.5036 11.8312
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.24803 0.04799 5.169 6.35e-07 ***
## sd1 4.57391 1.77075 2.583 0.010606 *
## sd2 2.22917 1.80017 1.238 0.217249
## sd3 6.21726 1.76874 3.515 0.000559 ***
## sd4 5.92019 1.81675 3.259 0.001343 **
## sd5 6.26258 1.85649 3.373 0.000913 ***
## sd6 4.85558 1.87389 2.591 0.010367 *
## sd7 3.27466 1.85772 1.763 0.079681 .
## sd8 4.23428 1.81687 2.331 0.020912 *
## sd9 4.01838 1.79650 2.237 0.026556 *
## sd10 -0.01183 1.79221 -0.007 0.994740
## sd11 -0.33635 1.74950 -0.192 0.847763
## ly1.dl1 0.16276 0.06712 2.425 0.016318 *
## ly2.dl1 0.18754 0.08713 2.152 0.032732 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.913 on 176 degrees of freedom
## Multiple R-squared: 0.4769, Adjusted R-squared: 0.4353
## F-statistic: 11.46 on 14 and 176 DF, p-value: < 2.2e-16
adjustment parameters are estimated as: α1 is not significant and α2 is significant at the level of 1%.and the signs are not consistent with error correction mechanism.
rest.alpha <- matrix(c(1,0), c(2,1))
y.CA.ralpha <- alrtest(y.CA, A=rest.alpha, r=1)
summary(y.CA.ralpha)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Estimation and testing under linear restrictions on beta
##
## The VECM has been estimated subject to:
## beta=H*phi and/or alpha=A*psi
##
## [,1]
## [1,] 1
## [2,] 0
##
## Eigenvalues of restricted VAR (lambda):
## [1] 0.082 0.000 0.000
##
## The value of the likelihood ratio test statistic:
## 25.68 distributed as chi square with 1 df.
## The p-value of the test statistic is: 0
##
## Eigenvectors, normalised to first column
## of the restricted VAR:
##
## [,1]
## RK.ly1.l1 1.0000
## RK.ly2.l1 -1.5589
## RK.constant -278.3270
##
## Weights W of the restricted VAR:
##
## [,1]
## [1,] -0.2446
## [2,] 0.0000
VEC with restriction on alpha_2
y.VEC.rest1 <- cajorls(y.CA.ralpha, r=1)
y.VEC.rest1
## $rlm
##
## Call:
## lm(formula = substitute(form1), data = data.mat)
##
## Coefficients:
## ly1.d ly2.d
## ect1 0.01989 0.22771
## sd1 6.16105 4.65599
## sd2 2.35313 2.25799
## sd3 8.05228 6.27209
## sd4 4.01057 5.86157
## sd5 5.08712 6.11692
## sd6 4.39434 4.58973
## sd7 4.58845 2.97139
## sd8 4.73235 3.97373
## sd9 3.18491 3.80830
## sd10 1.57007 -0.20095
## sd11 -0.32171 -0.44939
## ly1.dl1 0.29063 0.18202
## ly2.dl1 -0.12807 0.17142
##
##
## $beta
## ect1
## ly1.l1 1.00000
## ly2.l1 -1.55894
## constant -278.32704
summary(y.VEC.rest1$rlm)
## Response ly1.d :
##
## Call:
## lm(formula = ly1.d ~ ect1 + sd1 + sd2 + sd3 + sd4 + sd5 + sd6 +
## sd7 + sd8 + sd9 + sd10 + sd11 + ly1.dl1 + ly2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -25.946 -4.719 1.290 6.271 18.953
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.01989 0.08040 0.247 0.80487
## sd1 6.16105 3.00238 2.052 0.04165 *
## sd2 2.35313 3.05307 0.771 0.44189
## sd3 8.05228 2.99965 2.684 0.00796 **
## sd4 4.01057 3.08111 1.302 0.19473
## sd5 5.08712 3.14682 1.617 0.10776
## sd6 4.39434 3.17109 1.386 0.16758
## sd7 4.58845 3.14109 1.461 0.14586
## sd8 4.73235 3.07381 1.540 0.12546
## sd9 3.18491 3.04160 1.047 0.29648
## sd10 1.57007 3.03530 0.517 0.60562
## sd11 -0.32171 2.96558 -0.108 0.91374
## ly1.dl1 0.29063 0.11263 2.580 0.01069 *
## ly2.dl1 -0.12807 0.14716 -0.870 0.38534
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.332 on 176 degrees of freedom
## Multiple R-squared: 0.1325, Adjusted R-squared: 0.06354
## F-statistic: 1.921 on 14 and 176 DF, p-value: 0.02696
##
##
## Response ly2.d :
##
## Call:
## lm(formula = ly2.d ~ ect1 + sd1 + sd2 + sd3 + sd4 + sd5 + sd6 +
## sd7 + sd8 + sd9 + sd10 + sd11 + ly1.dl1 + ly2.dl1 - 1, data = data.mat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -19.4165 -1.8241 0.8595 3.8376 11.9833
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## ect1 0.22771 0.04790 4.754 4.12e-06 ***
## sd1 4.65599 1.78847 2.603 0.010020 *
## sd2 2.25799 1.81866 1.242 0.216049
## sd3 6.27209 1.78684 3.510 0.000569 ***
## sd4 5.86157 1.83537 3.194 0.001665 **
## sd5 6.11692 1.87451 3.263 0.001323 **
## sd6 4.58973 1.88896 2.430 0.016114 *
## sd7 2.97139 1.87110 1.588 0.114070
## sd8 3.97373 1.83102 2.170 0.031329 *
## sd9 3.80830 1.81183 2.102 0.036984 *
## sd10 -0.20095 1.80807 -0.111 0.911631
## sd11 -0.44939 1.76655 -0.254 0.799493
## ly1.dl1 0.18202 0.06709 2.713 0.007331 **
## ly2.dl1 0.17142 0.08766 1.955 0.052119 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.963 on 176 degrees of freedom
## Multiple R-squared: 0.4661, Adjusted R-squared: 0.4236
## F-statistic: 10.97 on 14 and 176 DF, p-value: < 2.2e-16
It is not possible to achieve the long run equilibrium, if we do not restrict one of the alphas to be equal to 0.
Rolling Scheme one-month ahead forecast
ly1=100*log(COP)
ly2=100*log(RCGP)
y.fcst <- cbind(ly1, ly2)
y.fcst <- na.trim(y.fcst)
y.fcst <- window(y.fcst, start=1995+0, end=2017+3/12)
## Warning in window.default(x, ...): 'end' value not changed
fstM <- 1995+0
lstM <- 2010+11/12
y.1 <- window(y.fcst, end=lstM)
y.2 <- window(y.fcst, start=lstM+1/12)
y.VAR.f1 <-data.frame()
y.VAR.f2 <-data.frame()
for(i in 1:length(y.2[,2]))
{
y <- window(y.fcst, end = lstM + (i-1)/12)
y.CA <- ca.jo(y, ecdet="const", type="eigen", K=3, spec="transitory")
y.VAR <- vec2var(y.CA, r=1)
y.VAR.updt <- predict(y.VAR, n.ahead=1)
y.VAR.f1 <-rbind(y.VAR.f1, as.ts(y.VAR.updt$fcst$ly1))
y.VAR.f2 <-rbind(y.VAR.f2, as.ts(y.VAR.updt$fcst$ly2))
}
y.VAR.f1 <- ts(y.VAR.f1, start=2011, frequency=12)
y.VAR.f2 <- ts(y.VAR.f2, start=2011, frequency=12)
y.pgas <-as.ts(y.fcst[,2])
y.pgas1 <-as.ts(y.1[,2])
y.pgas2 <-as.ts(y.2[,2])
par(mfrow=c(1,1))
plot(y.pgas, type="l", pch=16, lty="solid", main="One month ahead forecast for gas price", xlim=c(2010,2017), ylim=c(45,140))
lines(y.pgas1, type="o", pch=16, lty="solid")
lines(y.VAR.f2[,1], type="o", pch=16, lty="dashed", col="blue")
VEC.f.rol <- as.ts(y.pgas2)
accuracy(VEC.f.rol, y.VAR.f2[,1])
## ME RMSE MAE MPE MAPE ACF1
## Test set -1.41565 4.930924 3.547289 -1.667412 3.922441 0.05211672
## Theil's U
## Test set 0.5858428
RMSE is 4.930924.
par(mfrow=c(1,2))
acf(dl.y2, type="correlation", lag=48, xlab="Lag", ylab="ACF", main="Gas Price")
acf(dl.y2, type="partial", lag=48, xlab="Lag", ylab="PACF", main="Gas Price")
results of ACF and PACF, we can suggest different ARMA models.
arma22 <- Arima(y.pgas1, order=c(2, 0, 2))
tsdiag(arma22, gof.lag=36)
#(11)
rol.f <-zoo()
for(i in 1: length(y.2[,2]))
{
y <- window(dl.y2, end = lstM+(i-1)/12)
rol.updt <- arima(y, order=c(2,0,2))
rol.f <- c(rol.f, forecast(rol.updt, 1)$mean)
}
rol.f <-as.ts(rol.f)
par(mfrow=c(1,1))
plot(rol.f, type="o", pch=16, xlim=c(2010, 2017), main="One month ahead forecast for ARMA (2,2)")
lines(rol.f, type="p", pch=16, lty="solid", col="blue")
lines(dl.y2, type="l", pch=16, lty="solid")
lines(ly1, type="o", pch=16, lty="dashed")
accuracy(rol.f, dl.y2)
## ME RMSE MAE MPE MAPE ACF1 Theil's U
## Test set -0.5683734 5.221744 3.814195 Inf Inf -0.06078359 NaN
VEC model preforms a better forcast.