Chapter 7 Multivariate TS Analysis

# preamble setting the directories and loading packages
rm(list = ls())
setwd("C:/Users/Tejendra/Desktop/FoldersOnDesktop/UdemyCourse/Section7")
require(tidyverse)
require(tidymodels)
require(data.table)
require(tidyposterior)
require(tsibble)  #tsibble for time series based on tidy principles
require(fable)  #for forecasting based on tidy principles
require(ggfortify)  #for plotting timeseries
require(forecast)  #for forecast function
require(tseries)
require(chron)
require(lubridate)
require(directlabels)
require(zoo)
require(lmtest)
require(TTR)  #for smoothing the time series
require(MTS)
require(vars)
require(fUnitRoots)
require(lattice)
require(grid)

### Multivariate Time Series Datasets
# Generating a random dataframe
set.seed(40)
x = rnorm(100, 1)
y = rnorm(100, 30)
z = rnorm(100, 500)

xyz = data.frame(x, y, z)

class(xyz)

## [1] "data.frame"

# Converting a data.frame into mts
mymts = ts(xyz,
           frequency = 12,
           start = c(1940, 4))
mymts

##                    x        y        z
## Apr 1940  1.47773904 29.83632 499.5348
## May 1940  1.49618282 29.57833 500.8470
## Jun 1940  0.14041570 30.18632 500.6993
## Jul 1940  0.17094004 33.41376 500.5783
## Aug 1940  0.67842692 29.88003 498.8215
## Sep 1940 -0.30377040 29.39793 498.4307
## Oct 1940 -0.42148660 29.97486 499.4152
## Nov 1940  2.74491495 29.10156 498.9308
## Dec 1940  0.71172064 29.65150 500.1091
## Jan 1941 -0.30886572 29.66069 500.5443
## Feb 1941  0.93054781 30.34087 499.6077
## Mar 1941 -0.22492668 32.01089 498.4770
## Apr 1941  1.80899626 28.50112 501.2086
## May 1941  0.50784966 31.43644 501.1260
## Jun 1941  1.45269393 29.73437 500.6066
## Jul 1941  1.99963310 29.35371 498.6896
## Aug 1941  1.46702960 29.57435 499.6374
## Sep 1941  1.37605199 29.68409 500.0454
## Oct 1941  2.70349095 29.96825 499.6607
## Nov 1941 -0.03546152 29.04966 501.9248
## Dec 1941  2.32812210 29.84856 499.7123
## Jan 1942  0.40571289 29.49766 498.6045
## Feb 1942  2.61128458 28.72180 500.4350
## Mar 1942 -0.11267383 29.41831 500.5698
## Apr 1942 -0.46018298 29.60005 499.7704
## May 1942  1.73215598 30.67877 500.7278
## Jun 1942 -0.61033943 29.86472 500.5900
## Jul 1942  1.33206736 30.32859 499.6907
## Aug 1942  1.76085616 31.30714 499.6516
## Sep 1942 -0.85366955 31.60610 501.2060
## Oct 1942  1.79115157 29.95238 499.9182
## Nov 1942 -0.28174039 31.88079 499.5227
## Dec 1942  0.22012266 30.16986 500.9506
## Jan 1943  2.43834228 30.46903 500.2272
## Feb 1943  0.76517046 29.30462 501.4385
## Mar 1943  1.61219714 30.54121 499.8786
## Apr 1943 -0.45847063 28.51256 499.2938
## May 1943  1.24276907 32.63799 499.2861
## Jun 1943  0.68260789 30.06999 500.0912
## Jul 1943  1.85937333 28.61832 499.6628
## Aug 1943  2.34415507 30.33621 499.5449
## Sep 1943  1.31555720 29.74085 499.1711
## Oct 1943  3.04544486 30.51940 500.2317
## Nov 1943  0.81752947 30.47371 499.0529
## Dec 1943 -0.48393823 30.77232 501.0480
## Jan 1944  1.39716563 31.01997 498.3103
## Feb 1944  2.56965755 30.66090 501.7181
## Mar 1944  2.27023296 29.67146 500.5579
## Apr 1944  0.12855665 29.49982 499.1805
## May 1944  1.36088326 29.09545 500.3088
## Jun 1944  1.91588753 29.84240 500.4543
## Jul 1944  1.57451813 28.72539 499.2390
## Aug 1944  0.51459597 29.62820 499.4080
## Sep 1944  1.01574267 29.96911 501.3104
## Oct 1944  2.01421468 30.80206 501.3224
## Nov 1944  1.05814926 29.62449 499.8477
## Dec 1944  3.26829549 29.19358 500.2824
## Jan 1945  1.70865319 29.56288 499.7287
## Feb 1945  2.57247878 28.81038 498.0113
## Mar 1945  1.88720230 31.18393 500.0002
## Apr 1945 -0.41166914 29.01367 498.7472
## May 1945  1.80377513 29.01178 500.8901
## Jun 1945  0.44761960 32.77873 500.5178
## Jul 1945  1.42573830 30.32139 500.1941
## Aug 1945  0.22952536 29.30319 501.1462
## Sep 1945  0.37141318 30.91817 498.6691
## Oct 1945  0.98161035 28.43434 499.0585
## Nov 1945  0.34002587 30.74192 500.0088
## Dec 1945  0.48957134 28.47928 499.2432
## Jan 1946  1.64637043 30.25345 499.1858
## Feb 1946  0.53402838 29.99529 501.5672
## Mar 1946  2.37546031 30.09232 500.1563
## Apr 1946  2.19759957 30.22118 501.4116
## May 1946  0.38367022 28.16100 500.9055
## Jun 1946 -0.03894626 29.26891 499.0285
## Jul 1946  1.39954949 29.76422 500.0708
## Aug 1946  1.22058087 31.97935 499.8144
## Sep 1946  1.39729109 29.56818 500.6330
## Oct 1946  1.02702970 29.09356 501.1575
## Nov 1946  2.23407288 29.47902 501.3850
## Dec 1946  1.03060886 31.77104 501.4232
## Jan 1947 -0.36391409 31.91089 499.9674
## Feb 1947  1.90495376 28.87947 499.6365
## Mar 1947 -1.06942831 29.30213 497.3460
## Apr 1947  1.21278838 30.20043 501.0113
## May 1947  1.12402274 28.39637 500.9910
## Jun 1947  1.05353083 31.21420 499.5156
## Jul 1947  1.33005175 27.35457 501.7041
## Aug 1947  2.66493296 30.53251 501.5721
## Sep 1947  0.69573063 29.79494 499.8651
## Oct 1947 -0.41182405 28.44126 499.6831
## Nov 1947  2.42143748 28.23419 501.0438
## Dec 1947  0.34725532 29.96483 501.3890
## Jan 1948 -0.67167091 30.39152 499.5824
## Feb 1948  0.12753867 28.55380 500.9327
## Mar 1948  1.83922250 31.28148 499.0765
## Apr 1948 -0.16800389 30.63600 499.7599
## May 1948  2.08739965 29.04089 501.1305
## Jun 1948  1.85917397 30.20759 501.3492
## Jul 1948  0.89464946 30.01724 499.7520

# Standard exploratory tools
plot(mymts)

theme_set(theme_bw())
autoplot(mymts) +
  ggtitle("Time Series Plot of the `mymts' Time-Series") +
  theme(plot.title = element_text(hjust = 0.5)) #for centering the text

head(mymts)

##                   x        y        z
## Apr 1940  1.4777390 29.83632 499.5348
## May 1940  1.4961828 29.57833 500.8470
## Jun 1940  0.1404157 30.18632 500.6993
## Jul 1940  0.1709400 33.41376 500.5783
## Aug 1940  0.6784269 29.88003 498.8215
## Sep 1940 -0.3037704 29.39793 498.4307

class(mymts)

## [1] "mts"    "ts"     "matrix"

# Our further exercise dataset
class(EuStockMarkets)

## [1] "mts"    "ts"     "matrix"

head(EuStockMarkets)

## Time Series:
## Start = c(1991, 130) 
## End = c(1991, 135) 
## Frequency = 260 
##              DAX    SMI    CAC   FTSE
## 1991.496 1628.75 1678.1 1772.8 2443.6
## 1991.500 1613.63 1688.5 1750.5 2460.2
## 1991.504 1606.51 1678.6 1718.0 2448.2
## 1991.508 1621.04 1684.1 1708.1 2470.4
## 1991.512 1618.16 1686.6 1723.1 2484.7
## 1991.515 1610.61 1671.6 1714.3 2466.8

plot(EuStockMarkets)

autoplot(EuStockMarkets) +
  ggtitle("Time Series Plot of the `EuStockMarkets' Time-Series") +
  theme(plot.title = element_text(hjust = 0.5))

# Main packages - problem: both have different functions VAR
## Testing for stationarity
### tseries - standard test adt.test
apply(EuStockMarkets, 2, adf.test)

## Warning in FUN(newX[, i], ...): p-value greater than printed p-value

## $DAX
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -0.82073, Lag order = 12, p-value = 0.9598
## alternative hypothesis: stationary
## 
## 
## $SMI
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -0.522, Lag order = 12, p-value = 0.9808
## alternative hypothesis: stationary
## 
## 
## $CAC
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -0.24897, Lag order = 12, p-value = 0.99
## alternative hypothesis: stationary
## 
## 
## $FTSE
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -1.9736, Lag order = 12, p-value = 0.5895
## alternative hypothesis: stationary

# Alternative: lib fUnitRoots, function
apply(EuStockMarkets, 2, adfTest, 
      lags=0, #maximum number of lags used for error term correction
      type="c", #type of unit root regression
      title = "ADF Test for EuStockMarkets Data") #title of the project

## Warning in FUN(newX[, i], ...): p-value greater than printed p-value

## Warning in FUN(newX[, i], ...): p-value greater than printed p-value

## Warning in FUN(newX[, i], ...): p-value greater than printed p-value

## $DAX
## 
## Title:
##  ADF Test for EuStockMarkets Data
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 0
##   STATISTIC:
##     Dickey-Fuller: 1.9429
##   P VALUE:
##     0.99 
## 
## Description:
##  Mon Aug 19 23:37:21 2019 by user: Tejendra
## 
## 
## $SMI
## 
## Title:
##  ADF Test for EuStockMarkets Data
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 0
##   STATISTIC:
##     Dickey-Fuller: 2.2138
##   P VALUE:
##     0.99 
## 
## Description:
##  Mon Aug 19 23:37:21 2019 by user: Tejendra
## 
## 
## $CAC
## 
## Title:
##  ADF Test for EuStockMarkets Data
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 0
##   STATISTIC:
##     Dickey-Fuller: 1.2494
##   P VALUE:
##     0.99 
## 
## Description:
##  Mon Aug 19 23:37:21 2019 by user: Tejendra
## 
## 
## $FTSE
## 
## Title:
##  ADF Test for EuStockMarkets Data
## 
## Test Results:
##   PARAMETER:
##     Lag Order: 0
##   STATISTIC:
##     Dickey-Fuller: 0.2207
##   P VALUE:
##     0.9735 
## 
## Description:
##  Mon Aug 19 23:37:21 2019 by user: Tejendra

# Differencing the whole mts
stnry = diffM(EuStockMarkets) #difference operation on a vector of time series. Default order of differencing is 1.

# Retest
apply(stnry, 2, adf.test)

## Warning in FUN(newX[, i], ...): p-value smaller than printed p-value

## Warning in FUN(newX[, i], ...): p-value smaller than printed p-value

## Warning in FUN(newX[, i], ...): p-value smaller than printed p-value

## Warning in FUN(newX[, i], ...): p-value smaller than printed p-value

## $DAX
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -9.9997, Lag order = 12, p-value = 0.01
## alternative hypothesis: stationary
## 
## 
## $SMI
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -10.769, Lag order = 12, p-value = 0.01
## alternative hypothesis: stationary
## 
## 
## $CAC
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -11.447, Lag order = 12, p-value = 0.01
## alternative hypothesis: stationary
## 
## 
## $FTSE
## 
##  Augmented Dickey-Fuller Test
## 
## data:  newX[, i]
## Dickey-Fuller = -10.838, Lag order = 12, p-value = 0.01
## alternative hypothesis: stationary

## VAR modeling
plot.ts(stnry)

autoplot(ts(stnry,
            start = c(1990,130),
            frequency = 260)) +
  ggtitle("Time Series Plot of the stationary `EuStockMarkets' Time-Series")

# Lag order identification
#We will use two different functions, from two different packages to identify the lag order for the VAR model. Both functions are quite similar to each other but differ in the output they produce. vars::VAR is a more powerful and convinient function to identify the correct lag order. 
VARselect(stnry, 
          type = "none", #type of deterministic regressors to include. We use none becasue the time series was made stationary using differencing above. 
          lag.max = 10) #highest lag order

## $selection
## AIC(n)  HQ(n)  SC(n) FPE(n) 
##      9      1      1      9 
## 
## $criteria
##                   1            2            3            4            5
## AIC(n) 2.527062e+01 2.527564e+01 2.526566e+01 2.525844e+01 2.525725e+01
## HQ(n)  2.528823e+01 2.531088e+01 2.531850e+01 2.532891e+01 2.534534e+01
## SC(n)  2.531840e+01 2.537122e+01 2.540902e+01 2.544959e+01 2.549619e+01
## FPE(n) 9.438206e+10 9.485771e+10 9.391500e+10 9.324010e+10 9.312964e+10
##                   6            7            8            9           10
## AIC(n) 2.525408e+01 2.525692e+01 2.525696e+01 2.525073e+01 2.525455e+01
## HQ(n)  2.535978e+01 2.538023e+01 2.539789e+01 2.540927e+01 2.543071e+01
## SC(n)  2.554080e+01 2.559143e+01 2.563926e+01 2.568081e+01 2.573242e+01
## FPE(n) 9.283467e+10 9.309877e+10 9.310329e+10 9.252533e+10 9.288047e+10

# Creating a VAR model with vars
var.a <- vars::VAR(stnry,
                   lag.max = 10, #highest lag order for lag length selection according to the choosen ic
                   ic = "AIC", #information criterion
                   type = "none") #type of deterministic regressors to include
summary(var.a)

## 
## VAR Estimation Results:
## ========================= 
## Endogenous variables: DAX, SMI, CAC, FTSE 
## Deterministic variables: none 
## Sample size: 1850 
## Log Likelihood: -33712.408 
## Roots of the characteristic polynomial:
## 0.817 0.817 0.8116 0.8116 0.7915 0.7915 0.7864 0.7864 0.7784 0.7784 0.7579 0.7579 0.7541 0.7541 0.7537 0.7537 0.7473 0.7421 0.7421 0.7295 0.7295 0.7285 0.7153 0.7153 0.6723 0.6723 0.6696 0.6696 0.6616 0.6616 0.6551 0.577 0.577 0.4544 0.289 0.1213
## Call:
## vars::VAR(y = stnry, type = "none", lag.max = 10, ic = "AIC")
## 
## 
## Estimation results for equation DAX: 
## ==================================== 
## DAX = DAX.l1 + SMI.l1 + CAC.l1 + FTSE.l1 + DAX.l2 + SMI.l2 + CAC.l2 + FTSE.l2 + DAX.l3 + SMI.l3 + CAC.l3 + FTSE.l3 + DAX.l4 + SMI.l4 + CAC.l4 + FTSE.l4 + DAX.l5 + SMI.l5 + CAC.l5 + FTSE.l5 + DAX.l6 + SMI.l6 + CAC.l6 + FTSE.l6 + DAX.l7 + SMI.l7 + CAC.l7 + FTSE.l7 + DAX.l8 + SMI.l8 + CAC.l8 + FTSE.l8 + DAX.l9 + SMI.l9 + CAC.l9 + FTSE.l9 
## 
##           Estimate Std. Error t value Pr(>|t|)    
## DAX.l1   0.0096570  0.0424492   0.227 0.820065    
## SMI.l1  -0.1008170  0.0297641  -3.387 0.000721 ***
## CAC.l1   0.0752689  0.0465795   1.616 0.106285    
## FTSE.l1  0.0730055  0.0366170   1.994 0.046328 *  
## DAX.l2   0.0190453  0.0423265   0.450 0.652792    
## SMI.l2  -0.0172409  0.0298939  -0.577 0.564188    
## CAC.l2   0.0687124  0.0465965   1.475 0.140487    
## FTSE.l2 -0.0804753  0.0369389  -2.179 0.029489 *  
## DAX.l3  -0.0676359  0.0423179  -1.598 0.110155    
## SMI.l3   0.0135412  0.0299928   0.451 0.651696    
## CAC.l3   0.0484694  0.0466586   1.039 0.299032    
## FTSE.l3  0.0409793  0.0369675   1.109 0.267783    
## DAX.l4  -0.0501669  0.0422723  -1.187 0.235480    
## SMI.l4   0.0162536  0.0300860   0.540 0.589099    
## CAC.l4   0.1001510  0.0469324   2.134 0.032981 *  
## FTSE.l4 -0.0451988  0.0369319  -1.224 0.221170    
## DAX.l5   0.0109497  0.0424940   0.258 0.796687    
## SMI.l5  -0.0978623  0.0303192  -3.228 0.001270 ** 
## CAC.l5   0.0731622  0.0469140   1.559 0.119054    
## FTSE.l5 -0.0254787  0.0369942  -0.689 0.491086    
## DAX.l6  -0.0121897  0.0424062  -0.287 0.773800    
## SMI.l6   0.0246183  0.0303677   0.811 0.417660    
## CAC.l6   0.0871855  0.0468724   1.860 0.063039 .  
## FTSE.l6  0.0007736  0.0369967   0.021 0.983320    
## DAX.l7   0.0786601  0.0425103   1.850 0.064421 .  
## SMI.l7  -0.0050826  0.0302543  -0.168 0.866604    
## CAC.l7  -0.0691098  0.0466880  -1.480 0.138981    
## FTSE.l7 -0.0418380  0.0370855  -1.128 0.259406    
## DAX.l8  -0.0336346  0.0425857  -0.790 0.429743    
## SMI.l8   0.0963209  0.0304325   3.165 0.001576 ** 
## CAC.l8  -0.1180253  0.0466645  -2.529 0.011515 *  
## FTSE.l8  0.0517022  0.0371047   1.393 0.163666    
## DAX.l9  -0.0262047  0.0423049  -0.619 0.535714    
## SMI.l9   0.0052002  0.0303603   0.171 0.864022    
## CAC.l9   0.1359369  0.0467408   2.908 0.003678 ** 
## FTSE.l9 -0.0068091  0.0369929  -0.184 0.853983    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 32.09 on 1814 degrees of freedom
## Multiple R-Squared: 0.05129, Adjusted R-squared: 0.03246 
## F-statistic: 2.724 on 36 and 1814 DF,  p-value: 2.04e-07 
## 
## 
## Estimation results for equation SMI: 
## ==================================== 
## SMI = DAX.l1 + SMI.l1 + CAC.l1 + FTSE.l1 + DAX.l2 + SMI.l2 + CAC.l2 + FTSE.l2 + DAX.l3 + SMI.l3 + CAC.l3 + FTSE.l3 + DAX.l4 + SMI.l4 + CAC.l4 + FTSE.l4 + DAX.l5 + SMI.l5 + CAC.l5 + FTSE.l5 + DAX.l6 + SMI.l6 + CAC.l6 + FTSE.l6 + DAX.l7 + SMI.l7 + CAC.l7 + FTSE.l7 + DAX.l8 + SMI.l8 + CAC.l8 + FTSE.l8 + DAX.l9 + SMI.l9 + CAC.l9 + FTSE.l9 
## 
##          Estimate Std. Error t value Pr(>|t|)    
## DAX.l1   0.035132   0.052025   0.675 0.499578    
## SMI.l1  -0.038299   0.036478  -1.050 0.293895    
## CAC.l1   0.046647   0.057087   0.817 0.413969    
## FTSE.l1  0.127516   0.044877   2.841 0.004541 ** 
## DAX.l2   0.006278   0.051875   0.121 0.903691    
## SMI.l2   0.018350   0.036637   0.501 0.616532    
## CAC.l2   0.104672   0.057108   1.833 0.066983 .  
## FTSE.l2 -0.096675   0.045272  -2.135 0.032859 *  
## DAX.l3  -0.148622   0.051864  -2.866 0.004210 ** 
## SMI.l3   0.004229   0.036759   0.115 0.908422    
## CAC.l3   0.094768   0.057184   1.657 0.097644 .  
## FTSE.l3  0.131679   0.045307   2.906 0.003701 ** 
## DAX.l4  -0.175243   0.051808  -3.383 0.000733 ***
## SMI.l4   0.029175   0.036873   0.791 0.428904    
## CAC.l4   0.124249   0.057520   2.160 0.030895 *  
## FTSE.l4  0.011514   0.045263   0.254 0.799239    
## DAX.l5   0.007207   0.052080   0.138 0.889954    
## SMI.l5  -0.089506   0.037159  -2.409 0.016106 *  
## CAC.l5   0.070892   0.057497   1.233 0.217751    
## FTSE.l5 -0.037913   0.045339  -0.836 0.403156    
## DAX.l6  -0.072106   0.051972  -1.387 0.165490    
## SMI.l6   0.011650   0.037218   0.313 0.754308    
## CAC.l6   0.102452   0.057446   1.783 0.074681 .  
## FTSE.l6 -0.001026   0.045343  -0.023 0.981944    
## DAX.l7   0.147987   0.052100   2.840 0.004555 ** 
## SMI.l7  -0.012999   0.037079  -0.351 0.725941    
## CAC.l7  -0.123208   0.057220  -2.153 0.031432 *  
## FTSE.l7 -0.049168   0.045451  -1.082 0.279502    
## DAX.l8   0.008599   0.052192   0.165 0.869153    
## SMI.l8   0.089777   0.037298   2.407 0.016182 *  
## CAC.l8  -0.099393   0.057191  -1.738 0.082397 .  
## FTSE.l8 -0.019262   0.045475  -0.424 0.671921    
## DAX.l9   0.072664   0.051848   1.401 0.161245    
## SMI.l9  -0.091853   0.037209  -2.469 0.013657 *  
## CAC.l9   0.081425   0.057285   1.421 0.155371    
## FTSE.l9  0.068442   0.045338   1.510 0.131322    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 39.33 on 1814 degrees of freedom
## Multiple R-Squared: 0.05981, Adjusted R-squared: 0.04115 
## F-statistic: 3.206 on 36 and 1814 DF,  p-value: 7.227e-10 
## 
## 
## Estimation results for equation CAC: 
## ==================================== 
## CAC = DAX.l1 + SMI.l1 + CAC.l1 + FTSE.l1 + DAX.l2 + SMI.l2 + CAC.l2 + FTSE.l2 + DAX.l3 + SMI.l3 + CAC.l3 + FTSE.l3 + DAX.l4 + SMI.l4 + CAC.l4 + FTSE.l4 + DAX.l5 + SMI.l5 + CAC.l5 + FTSE.l5 + DAX.l6 + SMI.l6 + CAC.l6 + FTSE.l6 + DAX.l7 + SMI.l7 + CAC.l7 + FTSE.l7 + DAX.l8 + SMI.l8 + CAC.l8 + FTSE.l8 + DAX.l9 + SMI.l9 + CAC.l9 + FTSE.l9 
## 
##          Estimate Std. Error t value Pr(>|t|)   
## DAX.l1  -0.001041   0.034340  -0.030  0.97582   
## SMI.l1  -0.071184   0.024078  -2.956  0.00315 **
## CAC.l1   0.043492   0.037681   1.154  0.24857   
## FTSE.l1  0.082268   0.029622   2.777  0.00554 **
## DAX.l2   0.014488   0.034241   0.423  0.67225   
## SMI.l2  -0.027912   0.024183  -1.154  0.24857   
## CAC.l2   0.083842   0.037695   2.224  0.02626 * 
## FTSE.l2 -0.063578   0.029882  -2.128  0.03350 * 
## DAX.l3  -0.032437   0.034234  -0.948  0.34350   
## SMI.l3   0.031449   0.024263   1.296  0.19508   
## CAC.l3  -0.059480   0.037745  -1.576  0.11524   
## FTSE.l3  0.023769   0.029905   0.795  0.42683   
## DAX.l4  -0.112680   0.034197  -3.295  0.00100 **
## SMI.l4   0.045902   0.024338   1.886  0.05946 . 
## CAC.l4   0.071056   0.037967   1.872  0.06143 . 
## FTSE.l4 -0.020521   0.029877  -0.687  0.49225   
## DAX.l5  -0.040047   0.034376  -1.165  0.24418   
## SMI.l5  -0.040002   0.024527  -1.631  0.10308   
## CAC.l5   0.044130   0.037952   1.163  0.24507   
## FTSE.l5 -0.011466   0.029927  -0.383  0.70166   
## DAX.l6  -0.010487   0.034305  -0.306  0.75987   
## SMI.l6   0.017464   0.024566   0.711  0.47724   
## CAC.l6   0.046108   0.037918   1.216  0.22415   
## FTSE.l6 -0.002253   0.029929  -0.075  0.94000   
## DAX.l7   0.093443   0.034389   2.717  0.00665 **
## SMI.l7  -0.011696   0.024475  -0.478  0.63280   
## CAC.l7  -0.058576   0.037769  -1.551  0.12110   
## FTSE.l7 -0.059667   0.030001  -1.989  0.04687 * 
## DAX.l8   0.012292   0.034450   0.357  0.72128   
## SMI.l8   0.026246   0.024619   1.066  0.28653   
## CAC.l8  -0.102523   0.037750  -2.716  0.00667 **
## FTSE.l8  0.048842   0.030016   1.627  0.10387   
## DAX.l9   0.019936   0.034223   0.583  0.56028   
## SMI.l9  -0.025465   0.024560  -1.037  0.29995   
## CAC.l9   0.048093   0.037812   1.272  0.20357   
## FTSE.l9 -0.003901   0.029926  -0.130  0.89630   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 25.96 on 1814 degrees of freedom
## Multiple R-Squared: 0.04515, Adjusted R-squared: 0.0262 
## F-statistic: 2.382 on 36 and 1814 DF,  p-value: 8.732e-06 
## 
## 
## Estimation results for equation FTSE: 
## ===================================== 
## FTSE = DAX.l1 + SMI.l1 + CAC.l1 + FTSE.l1 + DAX.l2 + SMI.l2 + CAC.l2 + FTSE.l2 + DAX.l3 + SMI.l3 + CAC.l3 + FTSE.l3 + DAX.l4 + SMI.l4 + CAC.l4 + FTSE.l4 + DAX.l5 + SMI.l5 + CAC.l5 + FTSE.l5 + DAX.l6 + SMI.l6 + CAC.l6 + FTSE.l6 + DAX.l7 + SMI.l7 + CAC.l7 + FTSE.l7 + DAX.l8 + SMI.l8 + CAC.l8 + FTSE.l8 + DAX.l9 + SMI.l9 + CAC.l9 + FTSE.l9 
## 
##          Estimate Std. Error t value Pr(>|t|)    
## DAX.l1   0.025600   0.039796   0.643  0.52012    
## SMI.l1  -0.085874   0.027904  -3.078  0.00212 ** 
## CAC.l1  -0.003879   0.043668  -0.089  0.92923    
## FTSE.l1  0.165616   0.034328   4.824 1.52e-06 ***
## DAX.l2   0.023464   0.039681   0.591  0.55438    
## SMI.l2  -0.023240   0.028025  -0.829  0.40708    
## CAC.l2   0.028324   0.043684   0.648  0.51683    
## FTSE.l2 -0.031301   0.034630  -0.904  0.36619    
## DAX.l3  -0.052914   0.039673  -1.334  0.18245    
## SMI.l3   0.012312   0.028118   0.438  0.66154    
## CAC.l3   0.057729   0.043742   1.320  0.18709    
## FTSE.l3  0.007780   0.034657   0.224  0.82240    
## DAX.l4  -0.054187   0.039630  -1.367  0.17170    
## SMI.l4   0.043084   0.028206   1.527  0.12681    
## CAC.l4   0.078160   0.043999   1.776  0.07583 .  
## FTSE.l4 -0.083589   0.034624  -2.414  0.01587 *  
## DAX.l5   0.001615   0.039838   0.041  0.96767    
## SMI.l5  -0.042176   0.028424  -1.484  0.13803    
## CAC.l5   0.102931   0.043982   2.340  0.01938 *  
## FTSE.l5 -0.069017   0.034682  -1.990  0.04674 *  
## DAX.l6  -0.027039   0.039756  -0.680  0.49651    
## SMI.l6   0.058310   0.028470   2.048  0.04069 *  
## CAC.l6   0.094202   0.043943   2.144  0.03219 *  
## FTSE.l6 -0.088315   0.034684  -2.546  0.01097 *  
## DAX.l7   0.054056   0.039853   1.356  0.17514    
## SMI.l7   0.052084   0.028363   1.836  0.06648 .  
## CAC.l7  -0.065521   0.043770  -1.497  0.13458    
## FTSE.l7 -0.055592   0.034768  -1.599  0.11000    
## DAX.l8  -0.004926   0.039924  -0.123  0.90181    
## SMI.l8   0.057267   0.028530   2.007  0.04488 *  
## CAC.l8  -0.080907   0.043748  -1.849  0.06456 .  
## FTSE.l8  0.017291   0.034786   0.497  0.61919    
## DAX.l9   0.003630   0.039661   0.092  0.92709    
## SMI.l9  -0.017471   0.028463  -0.614  0.53942    
## CAC.l9   0.069078   0.043819   1.576  0.11510    
## FTSE.l9  0.012475   0.034681   0.360  0.71911    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 30.08 on 1814 degrees of freedom
## Multiple R-Squared: 0.05941, Adjusted R-squared: 0.04075 
## F-statistic: 3.183 on 36 and 1814 DF,  p-value: 9.501e-10 
## 
## 
## 
## Covariance matrix of residuals:
##         DAX    SMI   CAC  FTSE
## DAX  1025.7  940.4 616.8 647.6
## SMI   940.4 1537.5 650.9 722.0
## CAC   616.8  650.9 672.0 522.2
## FTSE  647.6  722.0 522.2 903.2
## 
## Correlation matrix of residuals:
##         DAX    SMI    CAC   FTSE
## DAX  1.0000 0.7488 0.7430 0.6729
## SMI  0.7488 1.0000 0.6403 0.6127
## CAC  0.7430 0.6403 1.0000 0.6703
## FTSE 0.6729 0.6127 0.6703 1.0000

# Residual diagnostics
#serial.test function takes the VAR model as the input.  
serial.test(var.a)

## 
##  Portmanteau Test (asymptotic)
## 
## data:  Residuals of VAR object var.a
## Chi-squared = 183.36, df = 112, p-value = 2.444e-05

#selecting the variables
# Granger test for causality
#for causality function to give reliable results we need all the variables of the multivariate time series to be stationary. 
causality(var.a, #VAR model
          cause = c("DAX")) #cause variable. If not specified then first column of x is used. Multiple variables can be used. 

## $Granger
## 
##  Granger causality H0: DAX do not Granger-cause SMI CAC FTSE
## 
## data:  VAR object var.a
## F-Test = 1.7314, df1 = 27, df2 = 7256, p-value = 0.01074
## 
## 
## $Instant
## 
##  H0: No instantaneous causality between: DAX and SMI CAC FTSE
## 
## data:  VAR object var.a
## Chi-squared = 759.19, df = 3, p-value < 2.2e-16

## Forecasting VAR models
fcast = predict(var.a, n.ahead = 25) # we forecast over a short horizon because beyond short horizon prediction becomes unreliable or uniform
par(mar = c(2.5,2.5,2.5,2.5))
plot(fcast)

# Forecasting the DAX index
DAX = fcast$fcst[1]; DAX # type list

## $DAX
##               fcst     lower    upper       CI
##  [1,]   5.60084304 -57.29384 68.49552 62.89468
##  [2,]  -5.91965377 -69.10567 57.26636 63.18601
##  [3,] -13.00325572 -76.29843 50.29191 63.29517
##  [4,]  13.72698427 -49.63068 77.08465 63.35767
##  [5,] -36.35436670 -99.79645 27.08772 63.44208
##  [6,]  -0.64933493 -64.37759 63.07892 63.72826
##  [7,]   5.96623225 -57.93847 69.87093 63.90470
##  [8,]   8.82501176 -55.17502 72.82505 64.00004
##  [9,]   2.15548901 -62.14270 66.45368 64.29819
## [10,]   2.12089142 -62.42840 66.67018 64.54929
## [11,]  -0.80619025 -65.36242 63.75003 64.55622
## [12,]  -4.24776524 -68.81260 60.31707 64.56484
## [13,]   0.31282719 -64.25708 64.88274 64.56991
## [14,]  -2.10767272 -66.68072 62.46537 64.57304
## [15,]   1.30500015 -63.28275 65.89275 64.58775
## [16,]   0.23742105 -64.35338 64.82822 64.59080
## [17,]   1.78179652 -62.81456 66.37815 64.59635
## [18,]  -0.36003725 -64.96020 64.24012 64.60016
## [19,]   0.72242303 -63.87835 65.32319 64.60077
## [20,]  -0.69067282 -65.29176 63.91041 64.60109
## [21,]  -0.43392018 -65.03547 64.16762 64.60155
## [22,]  -0.39358579 -64.99520 64.20802 64.60161
## [23,]  -0.08492216 -64.68679 64.51694 64.60186
## [24,]   0.19373467 -64.40828 64.79575 64.60202
## [25,]   0.22812144 -64.37395 64.83020 64.60207

# Extracting the forecast column
x = DAX$DAX[,1]; x

##  [1]   5.60084304  -5.91965377 -13.00325572  13.72698427 -36.35436670
##  [6]  -0.64933493   5.96623225   8.82501176   2.15548901   2.12089142
## [11]  -0.80619025  -4.24776524   0.31282719  -2.10767272   1.30500015
## [16]   0.23742105   1.78179652  -0.36003725   0.72242303  -0.69067282
## [21]  -0.43392018  -0.39358579  -0.08492216   0.19373467   0.22812144

tail(EuStockMarkets)

## Time Series:
## Start = c(1998, 164) 
## End = c(1998, 169) 
## Frequency = 260 
##              DAX    SMI    CAC   FTSE
## 1998.627 5598.32 7952.9 4041.9 5680.4
## 1998.631 5460.43 7721.3 3939.5 5587.6
## 1998.635 5285.78 7447.9 3846.0 5432.8
## 1998.638 5386.94 7607.5 3945.7 5462.2
## 1998.642 5355.03 7552.6 3951.7 5399.5
## 1998.646 5473.72 7676.3 3995.0 5455.0

# Inverting the differencing
#To get the data to the original scale we invert the time series
#since the values are just difference from the previous value, to get the values on the original scale we add the last value from the DAX time series to the predicted values.
#the plot of the predicted values will also show that over longer horizon the predicted values are not reliable
x = cumsum(x) + 5473.72
par(mar = c(2.5,2.5,1,2.5)) #bottom, left, top, and right
plot.ts(x)

# Adding data and forecast to one time series
DAXinv =ts(c(EuStockMarkets[,1], x),
           start = c(1991,130), 
           frequency = 260)
plot(DAXinv)

plot.ts(DAXinv[1786:1885])

DAXinv_datframe <- as.data.frame(DAXinv[1786:1885]) 
colnames(DAXinv_datframe) <- c("x")
head(DAXinv_datframe)

##         x
## 1 5337.75
## 2 5226.20
## 3 5264.62
## 4 5164.89
## 5 5270.61
## 6 5348.75

ggplot() + 
  geom_line(data = as.data.frame(DAXinv_datframe[1:75,]), aes(y = get("DAXinv_datframe[1:75, ]"), x = seq(1, 75)), color = "green") +
  geom_line(data = as.data.frame(DAXinv_datframe[76:100,]), aes(y = get("DAXinv_datframe[76:100, ]"), x = seq(76, 100)), color = "red") +
  ggtitle("Plot of forecast of the VAR model on `EuStockMarkets''s DAX time series") +
  theme(plot.title = element_text(hjust = 0.5)) +
  xlab("Time") + ylab("Value")

## Creating an advanced plot with visual separation
# Converting to object zoo
x = zoo(DAXinv[1786:1885])

# Advanced xyplot from lattice
xyplot(x, grid=TRUE, panel = function(x, y, ...){
  panel.xyplot(x, y, col="red", ...)
  grid.clip(x = unit(76, "native"), just=c("right"))
  panel.xyplot(x, y, col="green", ...) })

#we repeat the plots from above using the ggplot2 package 
# Inverting the differencing
x_dat_frame <- as.data.frame(x)
ggplot(x_dat_frame, aes(y = x, x = seq(1, length(x_dat_frame$x)))) + 
  geom_line() + 
  ggtitle("Plot of forecast of the VAR model on `EuStockMarkets''s DAX time series") +
  theme(plot.title = element_text(hjust = 0.5)) +
  xlab("Time") + ylab("Value")

# Adding data and forecast to one time series
autoplot(DAXinv) +
  ggtitle("Time Series Plot of the `DAXinv' Time-Series") +
  theme(plot.title = element_text(hjust = 0.5))