平稳VAR模型一

knitr::opts_chunk$set(echo = TRUE, message = F, warning = F)

Tsay书的第一章通过上两篇得到总结，主要目的是清楚多元模型的数据结构和模型目标。第二章用以了解 Stationary VAR 时间序列的基本知识，主要包括:

• 模型成分和结构
• 平稳性特征/条件
• 模型测量
• 模型筛选
• 预测

总体的，一个VAR(p)表达与AR(p)类似:

zt=ϕ0+∑i=1pϕizt−i+at

参见上篇详细表达。 值得注意的是，VAR(p) 的 p 取值于最大lag， 如一个观测对象数量为2的VAR(2)模型， 可以由一个AR(1)和一个AR(2)模型组成， 如下, 例子中z2,t−1 如果不能为z_{1t} 的预测提供有效信息, 依照Granger
Causality ϕ1,12=0, 原因是它不能提高预测的准确性。

z1t=ϕ10+ϕ1,11z1,t−1+ϕ1,12z2,t−1+a1tz2t=ϕ20+ϕ1,21z1,t−1+ϕ1,22z2,t−1+a2t

另外，现实中P的值都很小，以下在尽量不失总体性的情况下以如下三元VAR(2)模型为例子。 对 VAR(2) 模型结构和成分分析主要考虑moment分析，平稳性， MA变形以及单元序列变形三方面。

⎛⎝⎜z1tz2tz3t⎞⎠⎟=⎛⎝⎜530⎞⎠⎟+⎡⎣⎢0.47−0.350.470.210.340.2300.470.23⎤⎦⎥⎛⎝⎜z1,t−1z2,t−1z3,t−1⎞⎠⎟+⎡⎣⎢0−0.190.300.180000⎤⎦⎥⎛⎝⎜z1,t−1z2,t−1z3,t−1⎞⎠⎟+⎛⎝⎜a1ta2ta3t⎞⎠⎟,withΣa=⎡⎣⎢0.2850.0260.0690.0260.2870.1370.0690.1370.357⎤⎦⎥.

本篇仅考虑前两个Monment，即均值和协方差，公式和计算过程如下。

• 均值

μ=(Ik−ϕ1−ϕ2)−1ϕ0

• 方差和前p-1个lag 的协方差矩阵:

vec(Γ∗0)=[I(kp)2−Φ⊗Φ]−1vec(Γb)with:Φ=⎡⎣⎢⎢⎢⎢⎢⎢⎢⎢ϕ1I0⋮0ϕ20I⋮⋯⋯⋯⋱⋯ϕp−100⋮Iϕp00⋮0⎤⎦⎥⎥⎥⎥⎥⎥⎥⎥,andΓb=[Γa0k0k0k]generally,Γℓ=∑i=1pϕiΓℓ−i∀ℓ>0

使用上述VAR(2)数值进行如下计算

library(MTS)p0 = matrix(c(5, 3, 0), 3, 1)p1 = matrix(c(.47, -.35, .47, .21, .34, .23, 0, .47, .23), 3, 3)p2 = matrix(c(0,-.19, .3, 0, .18, 0, 0,0,0),3,3)sig = matrix(c(.285, .026, .069,.026,.287,.137,.069,.137,.357),3,3)k =  nrow(sig)Ik= diag(k)mu = solve(Ik - p1 - p2) %*% p0rownames(mu) = c("z1","z2","z3")colnames(mu) = c(mean)muzt = VARMAsim(30000, arlags = c(1,2), phi = cbind(p1,p2), cnst = c(p0),              sigma = sig, set.seed(100))meansm = matrix(apply(zt$series, 2, mean), nrow =1) "mean throung simulation"meansmP1 = diag(ncol(p1))P2 = P1 * 0PR1 = cbind(p1, p2)PR2 = cbind(P1, P2)P = rbind(PR1, PR2)PG1 = matrix(0, nrow(sig),ncol(sig))G2 = G1G3 = G1GR1 = cbind(sig, G1)GR2 = cbind(G2, G3)G = rbind(GR1, GR2)Gk = 3p = 2nVI = (k*p)^2VI = diag(nVI)VKP = kronecker(P, P)G0 = solve(VI - VKP) %*% matrix(as.numeric(G), ncol = 1)nG = length(G)^0.5GM = matrix(G0, nrow = nG)g0 = GM[1:nrow(sig), 1:ncol(sig)]g1 = GM[1:nrow(sig), (ncol(g0)+1):ncol(GM)]g2 = p1 %*% g1 + p2 %*% g0GamaM = cbind(g0, g1, g2)rownames(GamaM) = c("z1","z2","z3")colnames(GamaM) = paste("L", rep(c(1,2,3),  3), rep(c(0,1,2),  each =3), sep = "")GamaMsm = diag(diag(g0))se = sm^0.5s = solve(se)rho0 = s %*% g0  %*% srho1 = s %*% g1  %*% srho2 = s %*% g2  %*% sRohM = cbind(rho0, rho1, rho2)dimnames(RohM) = dimnames(GamaM)RohMzt = zt$seriesnr = nrow(zt)zt1 = zt[1: (nr - 2),]zt2 = zt[2:(nr -1),]zt3 = zt[3:nr,]Rohsim = cbind(cor(zt), cor(zt2,zt1), cor(zt3, zt1))"check throung simulation, outcomes differ somehow from the parameters's one if T is small"Rohsim 

Quelle

Rmarkdwon版本
Tsay, Ruey S(2014): *Multivariate Time Series Analysis
With R and Financial Applications*, Chapter 2.