*
* GARCHMV.PRG
* Manual Example 12.2
*
open data g10xrate.xls
data(format=xls,org=columns) 1 6237 usxjpn usxfra usxsui
*
set xjpn = 100.0*log(usxjpn/usxjpn{1})
set xfra = 100.0*log(usxfra/usxfra{1})
set xsui = 100.0*log(usxsui/usxsui{1})
*
garch(p=1,q=1,iters=200,hmatrices=hh) / xjpn xfra xsui
garch(p=1,q=1,mv=bek,method=bfgs,iters=200,pmethod=simplex,piters=10) / xjpn xfra xsui
garch(p=1,q=1,mv=diag,hmatrices=hd,rvectors=rd) / xjpn xfra xsui
garch(p=1,q=1,mv=cc) / xjpn xfra xsui
garch(p=1,q=1,mv=dcc,method=bfgs) / xjpn xfra xsui
*
* Compute the covariance matrix of the standardized residuals from
* the diagonal GARCH
*
set z1 = rd(t)(1)/sqrt(hd(t)(1,1))
set z2 = rd(t)(2)/sqrt(hd(t)(2,2))
set z3 = rd(t)(3)/sqrt(hd(t)(3,3))
vcv(matrix=cc)
# z1 z2 z3
*
* Compute the correlations from the multivariate GARCH
*
set rho12 = hh(t)(1,2)/sqrt(hh(t)(1,1)*hh(t)(2,2))
set rho13 = hh(t)(1,3)/sqrt(hh(t)(1,1)*hh(t)(3,3))
set rho23 = hh(t)(2,3)/sqrt(hh(t)(2,2)*hh(t)(3,3))
graph(header="Correlation of JPN with FRA",vgrid=||cc(1,2)||)
# rho12
graph(header="Correlation of JPN with SUI",vgrid=||cc(1,3)||)
# rho13
graph(header="Correlation of FRA with SUI",vgrid=||cc(2,3)||)
# rho23
*
* AR(1) models for each
*
equation(constant) jpneq xjpn 1
equation(constant) fraeq xfra 1
equation(constant) suieq xsui 1
group ar1 jpneq fraeq suieq
garch(p=1,q=1,model=ar1,mv=dcc,pmethod=simplex,piter=20,method=bfgs,iters=200,trace) / xjpn xfra xsui
*
* Compute correlations into the forecast period, and graph them along with some of the final
* values from the actual sample. The GRID option puts a vertical line at the separation
* between actual data and forecasts.
*
garch(p=1,q=1,iters=200,hmatrices=hh,rvectors=us) / xjpn xfra xsui
@MVGarchFore(steps=100) hh us
set rho12 6238 6337 = hh(t)(1,2)/sqrt(hh(t)(1,1)*hh(t)(2,2))
set rho13 6238 6337 = hh(t)(1,3)/sqrt(hh(t)(1,1)*hh(t)(3,3))
set rho23 6238 6337 = hh(t)(2,3)/sqrt(hh(t)(2,2)*hh(t)(3,3))
graph(header="Correlation of JPN with FRA",vgrid=||cc(1,2)||,grid=(t==6237))
# rho12 6100 6337
graph(header="Correlation of JPN with SUI",vgrid=||cc(1,3)||,grid=(t==6237))
# rho13 6100 6337
graph(header="Correlation of FRA with SUI",vgrid=||cc(2,3)||,grid=(t==6237))
# rho23 6100 6337
*
* Estimation using MAXIMIZE
* The initial few lines of this set the estimation range, which needs to be done explicitly,
* and the number of variables. Then, vectors for the dependent variables, residuals and
* residuals formulas are set up. The SET instructions copy the dependent variables over into
* the slots in the vector of series.
*
compute gstart=2,gend=6237
compute n=3
dec vect[series] y(n) u(n)
dec vect[frml] resid(n)
set y(1) = xjpn
set y(2) = xfra
set y(3) = xsui
*
* This is specific to a mean-only model. It sets up the formulas (the &i are needed in the
* formula definitions when the FRML is defined in a loop), and estimates them using NLSYSTEM.
* This both initializes the mean parameters, and computes the unconditional covariance matrix.
* If you want more general mean equations, the simplest way to do that would be to define each
* FRML separately.
*
dec vect b(n)
nonlin(parmset=meanparms) b
do i=1,n
   frml resid(i) = (y(&i)-b(&i))
end do i
nlsystem(parmset=meanparms,resids=u) gstart gend resid
compute rr=%sigma
*
* The paths of the covariance matrices and uu' are saved in the SERIES[SYMM] names H and UU.
* UX and HX are used to pull in residuals 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 are used to initialize pre-sample variances.
*
gset h  * gend = rr
gset uu * gend = rr
*
* This is a standard (normal) log likelihood formula for any multivariate GARCH model.
* The difference among these will be in the definitions of HF and RESID. The function
* %XT pulls information out of a matrix of SERIES.
*
declare frml[symm] hf
*
frml logl = $
    hx = hf(t) , $
    %do(i,1,n,u(i)=resid(i)) , $
    ux = %xt(u,t), $
    h(t)=hx, uu(t)=%outerxx(ux), $
    %logdensity(hx,ux)
*****************************************************
*
* Standard GARCH(1,1)
*
dec symm vcs(n,n) vas(n,n) vbs(n,n)
compute vcs=rr,vbs=%const(0.05),vas=%const(0.05)
nonlin(parmset=garchparms) vcs vas vbs
frml hf = vcs+vbs.*h{1}+vas.*uu{1}
maximize(parmset=meanparms+garchparms,pmethod=simplex,piters=10,method=bfgs,iters=400) logl gstart gend
*****************************************************
*
* CCC
* The correlations are parameterized using an (n-1)x(n-1) matrix for
* the subdiagonal. The (i,j) element of this will actually be the
* correlation between i+1 and j.
*
dec symm qc(n-1,n-1)
dec vect vcv(n) vbv(n) vav(n)
*
function hfcccgarch time
type symm hfcccgarch
type integer time
do i=1,n
   compute hx(i,i)=vcv(i)+vav(i)*h(time-1)(i,i)+vbv(i)*uu(time-1)(i,i)
   do j=1,i-1
     compute hx(i,j)=qc(i-1,j)*sqrt(hx(j,j)*hx(i,i))
   end do j
end do i
compute hfcccgarch=hx
end
*
frml hf = hfcccgarch(t)
nonlin(parmset=garchparms) vcv vbv vav qc
compute vcv=%xdiag(rr),vbv=%const(0.05),vav=%const(0.05),qc=%const(0.0)
maximize(parmset=meanparms+garchparms,pmethod=simplex,piters=10,method=bfgs) logl gstart gend