*
* Example of two-step estimates of various DCC models. This is the
* technique described in Cappiello, Engle & Sheppard(2006),  "Asymmetric
* Dynamics in the Correlations of Global Equity and Bond Returns,"
* Journal of Financial Econometrics, vol. 4, no 4, pages 537-572,
* applied to a different data set.
*
open data garch.asc
data(format=free,org=columns) 1 1867 bp cd dm jy sf
*
set dlogbp = 100*log(bp/bp{1})
set dlogcd = 100*log(cd/cd{1})
set dlogdm = 100*log(dm/dm{1})
set dlogjy = 100*log(jy/jy{1})
set dlogsf = 100*log(sf/sf{1})
*
compute n=5
dec vect[series] y(n)
set y(1) = dlogbp
set y(2) = dlogcd
set y(3) = dlogdm
set y(4) = dlogjy
set y(5) = dlogsf
*
* Estimate univariate GARCH models, saving the residuals and variances.
* Create VECT[SERIES] with the standardized residuals and the negative
* part of the standardized residuals.
*
dec vect[series] u(n) h(n) eps(n) eta(n)
do i=1,n
   garch(p=1,q=1,hseries=h(i),resids=u(i),asymm) / y(i)
   set eps(i) = u(i)/sqrt(h(i))
   set eta(i) = %min(eps(i),0.0)
end do i
*
* The recursions in the DCC models need a lagged value, so this makes
* sure there is space for them when setting GSTART (start parameter used
* in the MAXIMIZE instructions). This isn't a problem in this data set,
* since there's a data point lost in converting to returns, but if the
* original data were already returns, GARCH would be able to start at 1.
* This makes sure GSTART is at least 2.
*
compute gstart=%imax(%regstart(),2),gend=%regend()
*
* Keep track of the number of parameters across all the univariate
* models.
*
compute uniparms=n*%nreg
*
* These are empirical estimates
*
vcv(nocenter,matrix=qbar)
# eps
vcv(nocenter,matrix=nbar)
# eta
*************************************************************************
*************************************************************************
*
* The same overall setup will be used for all models, with the
* restrictions imposed through the parameter sets.
*
* These are the general parameters
*
dec vect aq(n) bq(n) gq(n)
compute aq=%const(sqrt(.1))
compute bq=%const(sqrt(.85))
compute gq=%const(0.0)
*
dec frml[symm] qf
dec series[symm] q ee nn hh
dec symm qcmatrix
*
gset q  = qbar
gset ee = qbar
gset nn = nbar
gset hh = qbar
gset ee gstart gend = %outerxx(%xt(eps,t))
gset nn gstart gend = %outerxx(%xt(eta,t))
*
* This computes the "intercept" in the evolution matrix for Q.
*
function StartQC
compute qcmatrix=qbar-%diag(aq)*qbar*%diag(aq)-%diag(bq)*qbar*%diag(bq)-%diag(gq)*nbar*%diag(gq)
eigen qcmatrix qceigen
end StartQC
*
* This computes the Q matrix at t, given the lags of Q and the
* standardized residuals.
*
frml qf   = qcmatrix+%diag(aq)*ee{1}*%diag(aq)+%diag(bq)*q{1}*%diag(bq)+%diag(gq)*nn{1}*%diag(gq)
*
* This computes the log likelihood at t given the values of the
* (unstandardized) residuals.
*
frml logl = (q=qf),hh=%corrtocv(%cvtocorr(q),%xt(h,t)),%logdensity(hh,%xt(u,t))
*************************************************************************
*
* DCC - this won't be quite the same as the DCC you would get with the
* GARCH instruction  because it uses the two-step likelihood.
*
* Note that AFIX and BFIX enter the calculation of the QC as squares, so
* the stability condition is AFIX^2+BFIX^2<1
*
compute bfix=sqrt(.85),afix=sqrt(.10)
nonlin afix bfix
compute gq=%const(0.0)
maximize(start=%(bq=%const(bfix),aq=%const(afix),StartQC()),method=bfgs,title="Two-step DCC") logl gstart gend
disp "DCC BIC" -2.0*%logl+(%nreg+uniparms)*log(%nobs)
*
* ADCC
*
compute bfix=sqrt(.85),afix=sqrt(.10),gfix=0.0
nonlin afix bfix gfix
maximize(start=%(bq=%const(bfix),aq=%const(afix),gq=%const(gfix),StartQC()),method=bfgs,title="Two-step ADCC") logl gstart gend
disp "ADCC BIC" -2.0*%logl+(%nreg+uniparms)*log(%nobs)
*
* GDDCC (Generalized Diagonal DCC)
*
nonlin aq bq
compute aq=%const(sqrt(.1))
compute bq=%const(sqrt(.85))
compute gq=%const(0.0)
maximize(start=%(StartQC()),pmethod=simplex,piters=5,method=bfgs,title="GDDCC") logl gstart gend
disp "GDDCC BIC" -2.0*%logl+(%nreg+uniparms)*log(%nobs)
*
* AGDDCC (Asymmetric Generalized Diagonal DCC)
*
compute gq=%const(0.0)
nonlin aq bq gq
maximize(start=%(StartQC()),pmethod=simplex,piters=5,method=bfgs,title="AGDDCC") logl gstart gend
disp "AGDDCC BIC" -2.0*%logl+(%nreg+uniparms)*log(%nobs)