Content
- Simulate a VAR(1) model
- Financial data
- Modeling
- Residual checking
- Prediction
- Co-integration
1. Simulate a VAR(1) model
require(MTS)
## Loading required package: MTS
require(ggplot2)
## Loading required package: ggplot2
(p1=matrix(c(0.5,-0.25,-1,0.5),2,2))
## [,1] [,2]
## [1,] 0.50 -1.0
## [2,] -0.25 0.5
(sig=matrix(c(1,0,0,1),2,2))
## [,1] [,2]
## [1,] 1 0
## [2,] 0 1
m1 <- VARMAsim(400,arlags=c(1),phi=p1,sigma=sig)
qplot(x = 1:length(m1$series[,1]), y = m1$series[,1],
geom = "line", xlab = "time index",
ylab = "column1")
qplot(x = 1:length(m1$series[,2]), y = m1$series[,2],
geom = "line", xlab = "time index",
ylab = "column2")
acf(m1$series[,1])
acf(m1$series[,2])
2. RiskMetrics-IGARCH
- Data loading and visualization
da = read.table("q-gdp-ukcaus.txt",header=T)
dim(da)
## [1] 126 5
head(da)
## year mon uk ca us
## 1 1980 1 172436 624794 5908500
## 2 1980 4 169359 623433 5787400
## 3 1980 7 169038 623215 5776600
## 4 1980 10 167180 630215 5883500
## 5 1981 1 166052 645957 6005700
## 6 1981 4 166393 651954 5957800
gdp = log(da[,3:5])
qplot(1:length(gdp[,1]), gdp[,1],type="l",
xlab = "time index", ylab = "UK(log)")
qplot(1:length(gdp[,2]), gdp[,2],type="l",
xlab = "time index", ylab = "CA(log)")
qplot(1:length(gdp[,3]), gdp[,3],type="l",
xlab = "time index", ylab = "US(log)")
- Computes the multivariate Ljung-Box statistics for cross-correlation matrices
mq(gdp,10)
## Ljung-Box Statistics:
## m Q(m) df p-value
## [1,] 1 331 9 0
## [2,] 2 613 18 0
## [3,] 3 853 27 0
## [4,] 4 1064 36 0
## [5,] 5 1255 45 0
## [6,] 6 1434 54 0
## [7,] 7 1605 63 0
## [8,] 8 1768 72 0
## [9,] 9 1920 81 0
## [10,] 10 2063 90 0
- Compute the growth rate of gdp. ‘2’ indicates columns
zt=apply(gdp,2,diff)
mq(zt,10) ##The results show that the 3-dimensional series is strongly serially correlated
## Ljung-Box Statistics:
## m Q(m) df p-value
## [1,] 1.0 79.1 9.0 0
## [2,] 2.0 118.1 18.0 0
## [3,] 3.0 143.3 27.0 0
## [4,] 4.0 163.9 36.0 0
## [5,] 5.0 174.7 45.0 0
## [6,] 6.0 185.0 54.0 0
## [7,] 7.0 199.8 63.0 0
## [8,] 8.0 214.1 72.0 0
## [9,] 9.0 219.7 81.0 0
## [10,] 10.0 226.7 90.0 0
3. Modeling
m1 = VARorderI(zt)
## selected order: aic = 8
## selected order: bic = 1
## selected order: hq = 4
## M statistic and its p-value
## Mstat pv
## [1,] 115.806 0.0000000
## [2,] 19.481 0.0213966
## [3,] 23.054 0.0060759
## [4,] 31.546 0.0002385
## [5,] 14.840 0.0954378
## [6,] 11.959 0.2156223
## [7,] 13.921 0.1251771
## [8,] 27.994 0.0009562
## [9,] 12.384 0.1925391
## [10,] 6.985 0.6386805
## [11,] 16.819 0.0516252
## [12,] 11.173 0.2640476
## [13,] 14.962 0.0919874
## Summary table:
## p AIC BIC HQ M(p) p-value
## [1,] 0 -30.158 -30.158 -30.158 0.000 0.00000000
## [2,] 1 -31.090 -30.885 -31.007 115.806 0.00000000
## [3,] 2 -31.129 -30.717 -30.961 19.481 0.02139658
## [4,] 3 -31.206 -30.585 -30.954 23.054 0.00607588
## [5,] 4 -31.374 -30.542 -31.036 31.546 0.00023853
## [6,] 5 -31.374 -30.329 -30.950 14.840 0.09543785
## [7,] 6 -31.346 -30.085 -30.834 11.959 0.21562235
## [8,] 7 -31.341 -29.862 -30.740 13.921 0.12517713
## [9,] 8 -31.502 -29.802 -30.812 27.994 0.00095624
## [10,] 9 -31.484 -29.562 -30.704 12.384 0.19253910
## [11,] 10 -31.402 -29.254 -30.530 6.985 0.63868053
## [12,] 11 -31.448 -29.072 -30.483 16.819 0.05162521
## [13,] 12 -31.423 -28.816 -30.365 11.173 0.26404756
## [14,] 13 -31.455 -28.615 -30.302 14.962 0.09198740
m2 = VAR(zt,1) ## based on BIC, estimate a VAR(1) model
## Constant term:
## Estimates: 0.001713324 0.001182869 0.002785892
## Std.Error: 0.0006790162 0.0007193106 0.0007877173
## AR coefficient matrix
## AR( 1 )-matrix
## [,1] [,2] [,3]
## [1,] 0.434 0.189 0.0373
## [2,] 0.185 0.245 0.3917
## [3,] 0.322 0.182 0.1674
## standard error
## [,1] [,2] [,3]
## [1,] 0.0811 0.0827 0.0872
## [2,] 0.0859 0.0877 0.0923
## [3,] 0.0940 0.0960 0.1011
##
## Residuals cov-mtx:
## [,1] [,2] [,3]
## [1,] 2.893347e-05 1.965508e-06 6.619853e-06
## [2,] 1.965508e-06 3.246932e-05 1.686272e-05
## [3,] 6.619853e-06 1.686272e-05 3.893867e-05
##
## det(SSE) = 2.721916e-14
## AIC = -31.09086
## BIC = -30.88722
## HQ = -31.00813
value = c(1,1,1,1,1,0,0,1,0)
value1 = matrix(value,3,3,byrow=T)
(value2 = rbind(c(1,0,1),value1)) ##For constant term
## [,1] [,2] [,3]
## [1,] 1 0 1
## [2,] 1 1 1
## [3,] 1 1 0
## [4,] 0 1 0
m22 = VAR(zt,1,fixed=value2) ##fixed: removing insignificant estimates.
## Constant term:
## Estimates: 0.001769493 0 0.004079231
## Std.Error: 0.0006639398 0 0.0007440541
## AR coefficient matrix
## AR( 1 )-matrix
## [,1] [,2] [,3]
## [1,] 0.446 0.209 0.000
## [2,] 0.219 0.276 0.421
## [3,] 0.498 0.000 0.000
## standard error
## [,1] [,2] [,3]
## [1,] 0.0757 0.0685 0.0000
## [2,] 0.0839 0.0862 0.0912
## [3,] 0.0844 0.0000 0.0000
##
## Residuals cov-mtx:
## [,1] [,2] [,3]
## [1,] 2.897757e-05 1.930763e-06 6.817902e-06
## [2,] 1.930763e-06 3.320102e-05 1.606269e-05
## [3,] 6.817902e-06 1.606269e-05 4.356240e-05
##
## det(SSE) = 3.315141e-14
## AIC = -30.94169
## BIC = -30.80593
## HQ = -30.88654
4. Residual checking
mq(m2$residuals,10) ##The results show that the residuals are not serially correlated
## Ljung-Box Statistics:
## m Q(m) df p-value
## [1,] 1.00 9.66 9.00 0.38
## [2,] 2.00 17.53 18.00 0.49
## [3,] 3.00 26.88 27.00 0.47
## [4,] 4.00 45.07 36.00 0.14
## [5,] 5.00 52.91 45.00 0.20
## [6,] 6.00 58.52 54.00 0.31
## [7,] 7.00 66.50 63.00 0.36
## [8,] 8.00 81.90 72.00 0.20
## [9,] 9.00 92.83 81.00 0.17
## [10,] 10.00 103.90 90.00 0.15
mq(m2$residuals,24,adj=9) ##The results show that the residuals are not serially correlated (1% significance level)
## Ljung-Box Statistics:
## m Q(m) df p-value
## [1,] 1.00 9.66 0.00 1.00
## [2,] 2.00 17.53 9.00 0.04
## [3,] 3.00 26.88 18.00 0.08
## [4,] 4.00 45.07 27.00 0.02
## [5,] 5.00 52.91 36.00 0.03
## [6,] 6.00 58.52 45.00 0.09
## [7,] 7.00 66.50 54.00 0.12
## [8,] 8.00 81.90 63.00 0.06
## [9,] 9.00 92.83 72.00 0.05
## [10,] 10.00 103.90 81.00 0.04
## [11,] 11.00 107.82 90.00 0.10
## [12,] 12.00 119.23 99.00 0.08
## [13,] 13.00 132.59 108.00 0.05
## [14,] 14.00 142.52 117.00 0.05
## [15,] 15.00 153.51 126.00 0.05
## [16,] 16.00 158.83 135.00 0.08
## [17,] 17.00 165.14 144.00 0.11
## [18,] 18.00 171.03 153.00 0.15
## [19,] 19.00 184.52 162.00 0.11
## [20,] 20.00 193.41 171.00 0.12
## [21,] 21.00 198.35 180.00 0.17
## [22,] 22.00 206.70 189.00 0.18
## [23,] 23.00 211.62 198.00 0.24
## [24,] 24.00 223.23 207.00 0.21
MTSdiag(m2, gof = 24, adj = 9, level = F)
## [1] "Covariance matrix:"
## uk ca us
## uk 2.92e-05 1.98e-06 6.67e-06
## ca 1.98e-06 3.27e-05 1.70e-05
## us 6.67e-06 1.70e-05 3.93e-05
## CCM at lag: 0
## [,1] [,2] [,3]
## [1,] 1.0000 0.0641 0.197
## [2,] 0.0641 1.0000 0.474
## [3,] 0.1972 0.4742 1.000
## Simplified matrix:
## CCM at lag: 1
## . . .
## . . .
## . . .
## CCM at lag: 2
## . . .
## . . .
## . . .
## CCM at lag: 3
## . . .
## . . .
## . . .
## CCM at lag: 4
## . . -
## . . .
## . - .
## CCM at lag: 5
## . . .
## . . .
## . . .
## CCM at lag: 6
## . . .
## . . .
## . - .
## CCM at lag: 7
## . . .
## . . .
## . . .
## CCM at lag: 8
## . . .
## . . .
## . . .
## CCM at lag: 9
## . . .
## . . .
## . . .
## CCM at lag: 10
## . . .
## . + .
## . . .
## CCM at lag: 11
## . . .
## . . .
## . . .
## CCM at lag: 12
## . . .
## . . .
## . . .
## CCM at lag: 13
## . - .
## . . .
## . . .
## CCM at lag: 14
## . . .
## . . .
## . . .
## CCM at lag: 15
## . . .
## . . .
## . . .
## CCM at lag: 16
## . . .
## . . .
## . . .
## CCM at lag: 17
## . . .
## . . .
## . . .
## CCM at lag: 18
## . . .
## . . .
## . . .
## CCM at lag: 19
## . . .
## . . +
## . . .
## CCM at lag: 20
## . . .
## . . .
## . . .
## CCM at lag: 21
## . . .
## . . .
## . . .
## CCM at lag: 22
## . . .
## . . .
## . . .
## CCM at lag: 23
## . . .
## . . .
## . . .
## CCM at lag: 24
## . . .
## . . .
## . . .
## Hit Enter for p-value plot of individual ccm:
## Hit Enter to compute MQ-statistics:
##
## Ljung-Box Statistics:
## m Q(m) df p-value
## [1,] 1.00 9.66 0.00 1.00
## [2,] 2.00 17.53 9.00 0.04
## [3,] 3.00 26.88 18.00 0.08
## [4,] 4.00 45.07 27.00 0.02
## [5,] 5.00 52.91 36.00 0.03
## [6,] 6.00 58.52 45.00 0.09
## [7,] 7.00 66.50 54.00 0.12
## [8,] 8.00 81.90 63.00 0.06
## [9,] 9.00 92.83 72.00 0.05
## [10,] 10.00 103.90 81.00 0.04
## [11,] 11.00 107.82 90.00 0.10
## [12,] 12.00 119.23 99.00 0.08
## [13,] 13.00 132.59 108.00 0.05
## [14,] 14.00 142.52 117.00 0.05
## [15,] 15.00 153.51 126.00 0.05
## [16,] 16.00 158.83 135.00 0.08
## [17,] 17.00 165.14 144.00 0.11
## [18,] 18.00 171.03 153.00 0.15
## [19,] 19.00 184.52 162.00 0.11
## [20,] 20.00 193.41 171.00 0.12
## [21,] 21.00 198.35 180.00 0.17
## [22,] 22.00 206.70 189.00 0.18
## [23,] 23.00 211.62 198.00 0.24
## [24,] 24.00 223.23 207.00 0.21
## Hit Enter to obtain residual plots:
5. Prediction
VARpred(m2, h=6) ## R command to compute prediction
## orig 125
## Forecasts at origin: 125
## uk ca us
## [1,] 0.002301 0.002467 0.003630
## [2,] 0.003314 0.003634 0.004582
## [3,] 0.004010 0.004480 0.005280
## [4,] 0.004498 0.005089 0.005774
## [5,] 0.004843 0.005522 0.006125
## [6,] 0.005088 0.005829 0.006373
## Standard Errors of predictions:
## [,1] [,2] [,3]
## [1,] 0.005379 0.005698 0.006240
## [2,] 0.006031 0.006764 0.006787
## [3,] 0.006309 0.007159 0.007043
## [4,] 0.006444 0.007346 0.007168
## [5,] 0.006511 0.007439 0.007230
## [6,] 0.006545 0.007485 0.007261
## Root mean square errors of predictions:
## [,1] [,2] [,3]
## [1,] 0.005464 0.005789 0.006339
## [2,] 0.040704 0.054184 0.039966
## [3,] 0.028021 0.035325 0.028630
## [4,] 0.020427 0.025433 0.020945
## [5,] 0.015215 0.018788 0.015718
## [6,] 0.011744 0.014341 0.012275
6. Co-integration
require(quantmod)
## Loading required package: quantmod
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
## Loading required package: TTR
##
## Attaching package: 'TTR'
## The following object is masked from 'package:MTS':
##
## VMA
## Version 0.4-0 included new data defaults. See ?getSymbols.
getSymbols("BHP")
## 'getSymbols' currently uses auto.assign=TRUE by default, but will
## use auto.assign=FALSE in 0.5-0. You will still be able to use
## 'loadSymbols' to automatically load data. getOption("getSymbols.env")
## and getOption("getSymbols.auto.assign") will still be checked for
## alternate defaults.
##
## This message is shown once per session and may be disabled by setting
## options("getSymbols.warning4.0"=FALSE). See ?getSymbols for details.
## [1] "BHP"
getSymbols("VALE")
## [1] "VALE"
bhp = log(BHP$BHP.Adjusted)
vale = log(VALE$VALE.Adjusted)
(m1 = lm(bhp~vale))
##
## Call:
## lm(formula = bhp ~ vale)
##
## Coefficients:
## (Intercept) vale
## 2.6567 0.4261
wt = bhp-0.4261*vale
autoplot(wt,colour = "red")
- Co-integration tests: Johansen’s co-integration test
require(urca)
## Loading required package: urca
head(BHP)
## BHP.Open BHP.High BHP.Low BHP.Close BHP.Volume BHP.Adjusted
## 2007-01-03 39.67 39.69 38.62 38.87 6285400 24.54328
## 2007-01-04 37.90 38.19 37.48 37.74 5121000 23.82978
## 2007-01-05 37.64 37.69 36.88 37.16 4412500 23.46355
## 2007-01-08 37.55 37.66 36.93 37.33 4674400 23.57089
## 2007-01-09 37.77 37.77 36.37 37.65 8331600 23.77294
## 2007-01-10 37.30 38.05 37.02 37.84 4646300 23.89291
head(VALE)
## VALE.Open VALE.High VALE.Low VALE.Close VALE.Volume
## 2007-01-03 15.070 15.610 14.305 14.410 19475000
## 2007-01-04 14.320 14.320 14.015 14.235 20641000
## 2007-01-05 14.235 14.240 13.555 13.755 20187600
## 2007-01-08 13.845 14.285 13.795 14.235 21989000
## 2007-01-09 14.205 14.205 13.530 13.820 21337600
## 2007-01-10 13.625 14.355 13.565 14.190 18589000
## VALE.Adjusted
## 2007-01-03 10.263877
## 2007-01-04 10.139231
## 2007-01-05 9.797338
## 2007-01-08 10.139231
## 2007-01-09 9.843638
## 2007-01-10 10.107178
tail(BHP)
## BHP.Open BHP.High BHP.Low BHP.Close BHP.Volume BHP.Adjusted
## 2019-04-11 56.39 56.63 56.12 56.35 851200 56.35
## 2019-04-12 57.03 57.17 56.73 56.83 1440000 56.83
## 2019-04-15 56.19 56.40 55.96 56.02 1222700 56.02
## 2019-04-16 56.07 56.18 55.07 55.37 2713700 55.37
## 2019-04-17 54.95 55.19 54.48 55.12 1964600 55.12
## 2019-04-18 54.83 55.00 54.67 54.95 1137400 54.95
dim(BHP)
## [1] 3095 6
par(mfcol=c(2,1))
m1 = lm(bhp~vale)
summary(m1)
##
## Call:
## lm(formula = bhp ~ vale)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.67076 -0.10064 -0.00332 0.13101 0.31926
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.656659 0.014569 182.35 <2e-16 ***
## vale 0.426078 0.005637 75.59 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1727 on 3093 degrees of freedom
## Multiple R-squared: 0.6488, Adjusted R-squared: 0.6487
## F-statistic: 5713 on 1 and 3093 DF, p-value: < 2.2e-16
x = cbind(bhp,vale) ##Take a sequence of vector and combine by columns
m1 = ar(x)
m1$order
## [1] 6
m2 = ca.jo(x, K=6)#maximum eigenvalue test
summary(m2)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: maximal eigenvalue statistic (lambda max) , with linear trend
##
## Eigenvalues (lambda):
## [1] 0.0041244184 0.0008515475
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 2.63 6.50 8.18 11.65
## r = 0 | 12.77 12.91 14.90 19.19
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## BHP.Adjusted.l6 VALE.Adjusted.l6
## BHP.Adjusted.l6 1.0000000 1.000000
## VALE.Adjusted.l6 -0.3622562 -1.965875
##
## Weights W:
## (This is the loading matrix)
##
## BHP.Adjusted.l6 VALE.Adjusted.l6
## BHP.Adjusted.d -0.009210101 -7.222733e-05
## VALE.Adjusted.d -0.010014343 6.121659e-04
m3 = ca.jo(x, K=6, type=c("trace")) #trace co-integration test
summary(m3)
##
## ######################
## # Johansen-Procedure #
## ######################
##
## Test type: trace statistic , with linear trend
##
## Eigenvalues (lambda):
## [1] 0.0041244184 0.0008515475
##
## Values of teststatistic and critical values of test:
##
## test 10pct 5pct 1pct
## r <= 1 | 2.63 6.50 8.18 11.65
## r = 0 | 15.40 15.66 17.95 23.52
##
## Eigenvectors, normalised to first column:
## (These are the cointegration relations)
##
## BHP.Adjusted.l6 VALE.Adjusted.l6
## BHP.Adjusted.l6 1.0000000 1.000000
## VALE.Adjusted.l6 -0.3622562 -1.965875
##
## Weights W:
## (This is the loading matrix)
##
## BHP.Adjusted.l6 VALE.Adjusted.l6
## BHP.Adjusted.d -0.009210101 -7.222733e-05
## VALE.Adjusted.d -0.010014343 6.121659e-04