*
* Example of bootstrapping for MV GARCH models
*
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})
***********************************************************************
*
* MV GARCH with DCC
*
garch(p=1,q=1,mv=dcc,rvectors=rv,hmatrices=h,mvhseries=hhs) / xjpn xfra xsui
compute gstart=%regstart(),gend=%regend()
set cor %regstart() %regend() = hhs(1,2)/sqrt(h(t)(1,1)*h(t)(2,2))
gra 1; # cor

* Transform residuals into standardized form. For a DCC model, there are
* two types of standardizations necessary - one set transforms the
* residuals down to fully uncorrelated with unit variances (RVSTD) and
* one which merely standardizes each variance but leaves correlations
* intact (RVSTDC). Also save the vectors of variances and squared
* residuals.
*
dec series[vect] rvstd rvstdc
dec series[vect] hv uuv
gset rvstd  gstart gend = %solve(%decomp(h(t)),rv(t))
gset rvstdc gstart gend = rv(t)./%xdiag(h(t))
gset hv     gstart gend = %xdiag(h(t))
gset uuv    gstart gend = rv(t).*rv(t)
*
* Set up the individual variance recursions
*
dec vect %%g_c(%nvar) %%g_a(%nvar) %%g_b(%nvar)
ewise %%g_c(i)=%beta(%nregmean+0*%nvar+i)
ewise %%g_a(i)=%beta(%nregmean+1*%nvar+i)
ewise %%g_b(i)=%beta(%nregmean+2*%nvar+i)
*
dec frml[vect] hvf
frml hvf = %%g_c+%%g_b.*hv{1}+%%g_a.*uuv{1}
*
* Pull out the DCC recursion coefficients
*
compute %%dcc_a=%beta(%nreg-1),%%dcc_b=%beta(%nreg)
*
* The "constant" term in the DCC recursion is based upon the
* unconditional correlation matrix of the data.
*
vcv(center)
# xjpn xfra xsui
compute qbar=%cvtocorr(%sigma)
*
* Set up the DCC recursion
*
dec series[symm] q
dec frml[symm] qf
frml qf = (1-%%dcc_a-%%dcc_b)*qbar+%%dcc_b*q{1}+%%dcc_a*%outerxx(rvstdc{1})
*
* The Q sequence can't be backed out of the information available from
* GARCH. This re-generates it.
*
gset q gstart   gstart = qbar
gset q gstart+1 gend   = qf(t)
*
* Bootstrap model to produce new data set.
*
compute span=6237
compute ndraws=20
compute bstart=gend+1,bend=gend+span
*
* Extend various objects out-of-sample. (Values don't really matter at
* this point).
*
gset uuv    bstart bend = uuv{1}
gset hv     bstart bend = hv{1}
gset h      bstart bend = h{1}
gset q      bstart bend = q{1}
gset rvstd  bstart bend = rvstd{1}
gset rvstdc bstart bend = rvstdc{1}
*
*
set cor 1 ndraws = 0.
*
compute [vector] ones=%fill(%nvar,1,1.0)
*
dec vect[series] correlations(6237)
do i = 1, 6237
 set correlations(i) 1 ndraws = 0.
end do i

*
infobox(action=define,progress,lower=1,upper=ndraws) "Progress"
do draw=1,ndraws
   boot entries bstart bend gstart gend
   *
   * Do the DCC recursion and the recursion for the univariate variances,
   * convert the pair into the covariance matrix, post-multiply the
   * shuffled standardized residuals by a square root of the variance,
   * save the squares for the next GARCH recusion, save the standardized
   * residuals for the next DCC recursion, and return the re-flated
   * residuals.
   *
   * To do the bootstrap, update Q(t) and HV(t). Convert Q to a
   * correlation matrix and create H(t) from the correlations and the
   * variances in HV. Use that to reflate the shuffled fully
   * standardized residuals.
   *
   gset rv bstart bend = q=qf,hv=hvf,h=%corrtocv(%cvtocorr(q),hv),$
     uboot=%decomp(h)*rvstd(entries(t)),uuv=uboot.*uboot,rvstdc=uboot./hv,uboot
   *
   * Re-estimate the GARCH model
   *
	garch(p=1,q=1,mv=dcc,rvectors=rv,hmatrices=h,mvhseries=hhs) / xjpn xfra xsui
	set cor %regstart() %regend() = hhs(1,2)/sqrt(h(t)(1,1)*h(t)(2,2))
	do j = 1 , 6237
		com correlations(j)(draw) = cor(j)
	end do j

   infobox(current=draw)
end do draw
infoxbox(action=remove)
*
*
stats(fractiles) cor
tab
pri / correlations(300)
gra 1; # cor