I used your varma-bekk code to turn it into a varmax bekk(2,2) code. could you please check that i have done it correctly. I also added the commands to compute the standerdized residuals and the rho but i am not sure that this is correct because you define the varma residuals(residv) differently. and how do i get the conditional variances for each equation?
if i want to add another ar lag do i just creat a new artwo(n,n) matrix and add to to the mean nonlin command line and add artwo*yv{2} term to the frml residv equation?
last question, if i want the bekk to be lower triangular i can add the restrictions manually using the nonlin command and set each component (1,2) in the arch and garch term to zero?
VARMAX(1,1)-BEKK(2,2)
Code: Select all
smpl 2002:12:30 2010:10:5
set x1 = 100.0*usa
set x2 = 100.0*arg
set xoil = 100.0*oil
set xcus3t = 100.0*cus3t
*
compute n=2
compute gstart=2003:1:20,gend=2010:10:5
*
dec series[vect] yv eg
dec frml[vect] residv
dec vect[series] u(n)
*
* The paths of the covariance matrices and uu' are saved in the
* SERIES[SYMM] named H and UU. UX and HX are used for the current values
* of the residual vector and H matrices
*
declare series[symm] h uu
*
* ux is used when extracting a u vector
*
declare symm hx(n,n)
declare vect ux(n)
*
* These will be the parameters for the mean equations. These are adjusted to add
* variance or covariance terms as needed.
*
dec vect b(n)
dec rect ar(n,n) ma(n,n) exo(n,n)
nonlin(parmset=meanparms) b ar ma exo
*
* Mean model = VARMA(1,1)
*
frml residv = yv-b-ar*yv{1}-ma*%xt(u,t-1)-exo*eg
*
gset yv = ||x1,x2||
gset eg = ||xoil,xcus3t||
*
* Run preliminary VAR(4) to get estimates of residuals
*
linreg x1 / u(1)
# constant x1{1 to 4} x2{1 to 4} xoil xcus3t
linreg x2 / u(2)
# constant x1{1 to 4} x2{1 to 4} xoil xcus3t
*
linreg x1
# constant x1{1} x2{1} u(1){1} u(2){1} xoil xcus3t
compute b(1)=%beta(1),ar(1,1)=%beta(2),ar(1,2)=%beta(3),$
ma(1,1)=%beta(4),ma(1,2)=%beta(5),$
exo(1,1)=%beta(6),exo(1,2)=%beta(7)
linreg x2
# constant x1{1} x2{1} u(1){1} u(2){1} xoil xcus3t
compute b(2)=%beta(1),ar(2,1)=%beta(2),ar(2,2)=%beta(3),$
ma(2,1)=%beta(4),ma(2,2)=%beta(5),$
exo(2,1)=%beta(6),exo(2,2)=%beta(7)
vcv(matrix=rr,noprint)
# u
*
* These are used to initialize pre-sample variances.
*
gset h * gend = rr
gset uu * gend = rr
set u(1) = 0.0
set u(2) = 0.0
*
declare frml[symm] hf
*
frml logl = $
hx = hf(t) , $
h(t) = hx, $
ux = residv , $
uu(t) = %outerxx(ux), $
%pt(u,t,ux),$
%logdensity(hx,ux)
*
* This does a simple GARCH(1,1) model for the variance.
*
dec symm vcs(n,n) vas(n,n) vbs(n,n) vastwo(n,n) vbstwo(n,n)
compute vcs=rr,vbs=%const(0.50),vas=%const(0.05),vastwo=%const(0.05),vbstwo=%const(0.5)
nonlin(parmset=garchparms) vcs vas vbs vastwo vabtwo
frml hf = vcs+vbs.*h{1}+vas.*uu{1}+vbstwo*h{2}+vastwo*uu{2}
maximize(trace,parmset=meanparms+garchparms,pmethod=simplex,piters=10,method=bfgs,iters=400) logl gstart gend
*
* This does a BEKK model for the variance
*
dec rect vab(n,n) vbb(n,n) vabtwo(n,n) vbbtwo(n,n)
compute vab=.05*%identity(n),vbb=.50*%identity(n),vabtwo=0.05*%identity(n),vbbtwo=0.05*%identity(n)
nonlin(parmset=garchparms) vcs vab vbb vabtwo vbbtwo
compute vcs=%decomp(rr)
frml hf = vcs*tr(vcs)+vbb*h{1}*tr(vbb)+vab*uu{1}*tr(vab)+vbbtwo*h{2}*tr(vbbtwo)+vabtwo*uu{2}*tr(vabtwo)
maximize(trace,parmset=meanparms+garchparms,pmethod=simplex,piters=10,method=bfgs,iters=400) logl gstart gend
set x1sd = u(1)/sqrt(h(t)(1,1))
set x2sd = u(2)/sqrt(h(t)(2,2))
set x1sdsq = x1sd**2
set x2sdsq = x2sd**2
set rho12 = h(t)(1,2)/sqrt(h(t)(1,1)*h(t)(2,2))
set rho21 = h(t)(2,1)/sqrt(h(t)(1,1)*h(t)(2,2))
graph
# rho12
@mvqstat(lags=8,dfc=4)
# x1sd x2sd
@mvqstat(lags=8,dfc=6)
# x1sdsq x2sdsqHashem