The RATS Software Forum

Code: Select all

compute eps0=rv(ecpolicy).^2 
compute sigma0=%xdiag(hf(blackwed))

Code: Select all

*
* Example of a CC-spillover model with a single dummy break applied to the mean,
* spillover variances and CC matrix as in Sun, Tong and Yan(2009), "Market
* liberalization within a country," Journal of Empirical Finance, vol. 16, no 1,
* 18-41.
*
* This is applied to the data set used in the GARCHMV.RPF example.
*
all 3720
open data xrates.xls
data(format=xls,org=columns) / DATE USXFRA USXSUI

*
set xfra = 100.0*log(USXFRA/USXFRA{1})
set xsui = 100.0*log(USXSUI/USXSUI{1})
*
* The base model is a CC with VARIANCES=SPILLOVER.
*
garch(mv=cc,variances=spillover,rvectors=rv,derives=dd) / xfra xsui
*
* Somewhat interestingly (given that it's just one of our standard datasets), the
* fluctuation statistics show that a distinct break in the CC coefficient (at
* around entry 2725) *is* a distinct possibility.
*
@flux
# dd
acc dd(11) / ddcc
graph(footer="Accumulated gradients of CC parameter")
# ddcc
*
set dummy = t>=2725
*
* Estimation using MAXIMIZE
* The initial few lines of this set the estimation range, which needs to
* be done explicitly, and the number of variables.
*
compute gstart=2,gend=3720
compute n=2
*
dec vect[series] u(n)
dec series[symm] h uu
*
dec vect[equation] eqn(n)
dec vect[vect] bmean(n)
nonlin(parmset=meanparms) bmean
*
* Mean equations
*
equation eqn(1) xfra
# constant dummy
equation eqn(2) xsui
# constant dummy
*
* Estimate by OLS to get guess values
*
do i=1,n
   linreg(equation=eqn(i))
   compute bmean(i)=%beta
   set u(i) = %resids
end do i
*
* Get the covariance matrix of the OLS residuals to get presample values for the H
* matrices and to use for guess values for the GARCH parameters.
*
vcv
# u
compute rr=%sigma
*
* 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)=%eqnrvalue(eqn(i),t,bmean(i))), $
  ux = %xt(u,t),$
  h(t)=hx, uu(t)=%outerxx(ux),$
  %logdensity(hx,ux)
*****************************************************
*
* VCV is the VECTOR of variance constants
* VBV is the VECTOR of "GARCH" parameters
* VAR is the RECTANGULAR of standard spillover parameters
* VARD is the RECTANGULAR of dummied spillover parameters
*
dec vect vcv(n) vbv(n)
dec rect var(n,n) vard(n,n)
*
* QC is the subdiagonal of the standard CC
* QCD is the subdiagonal of the dummied CC
*
dec symm qc(n-1,n-1) qcd(n-1,n-1)
*
* This computes the VECTOR of variances allowing for the dummies in the spillover
* terms.
*
function varspilloverdummy time
type vector varspilloverdummy
type integer time
*
local integer i j
local vect    hv(n)
*
do i=1,n
   compute hv(i)=vcv(i)+vbv(i)*h(time-1)(i,i)
   do j=1,n
      compute hv(i)+=var(i,j)*uu(time-1)(j,j)
      compute hv(i)+=vard(i,j)*uu(time-1)(j,j)*dummy(time)
   end do j
end do i
*
compute varspilloverdummy=hv
end

* This converts the subdiagonal correlation to a full correlation matrix.
*
function %%corrfromsubd corrsubd
type symmetric %%corrfromsubd corrsubd
*
dim %%corrfromsubd(%rows(corrsubd)+1,%rows(corrsubd)+1)
ewise %%corrfromsubd(i,j)=%if(i==j,1.0,corrsubd(i-1,j))
end
*
* Compute the covariance matrix from the time-varying correlation matrix
* subdiagonal, and the spillover variance.
*
frml hf = %corrtocv(%%corrfromsubd(qc+qcd*dummy),varspilloverdummy(t))

* Guess values for GARCH parameters
*
ewise vbv(i)=.45
ewise var(i,j)=.05*(i==j)
ewise vcv(i)=rr(i,i)*(1-vbv(i)-var(i,i))
*
compute vard=%zeros(n,n)
ewise qc(i,j)=.5
compute qcd=%zeros(n-1,n-1)
*
nonlin(parmset=garchparms) vcv vbv var vard qc qcd
maximize(parmset=meanparms+garchparms,$
   pmethod=simplex,piters=10,method=bfgs,iters=400) logl gstart gend

* VIRF with historical incidents
*
sstats(max) / %if(date==800201,t,0)>>xblackwed %if(date==930802,t,0)>>xecpolicy
compute blackwed=fix(xblackwed),ecpolicy=fix(xecpolicy)

compute nstep=400
spgraph(vfields=2,hfields=2,footer="Figure 1",$
  xlabels=||"1992 Sept 16","1993 Aug 2"||,$
  ylabels=||"FRA/USD Variance","SUI/USD Variance"||)
*
* Black Wednesday shocks. These are computed using a baseline of the
* estimated volatility state, so they are excess over the predicted
* covariance.
*
compute eps0=rv(blackwed).^2
compute sigma0=%xdiag(hf(blackwed))
compute shock=eps0-sigma0
*
* This generates responses for the 2 elements of the covariance matrix,
* which will be (in order) the (1,1) and (2,2).
*
dec vect[series] sept92virf(2)
do i=1,2
   set sept92virf(i) 1 nstep = 0.0
end do i
*
*
do step=1,nstep
   if step==1
      compute hvec=%diag(vbv)*shock
   else
      compute hvec=(%diag(vbv)+%diag(%xdiag(var)))*hvec
   compute %pt(sept92virf,step,hvec)
end do step

do i=1,2
   graph(min=0.0,column=1,row=i,picture="*.###",vticks=5)
   # sept92virf(i) 1 nstep
end do i


* EC Policy Change shock. These use a baseline of zero.
*
compute eps0=rv(ecpolicy).^2
compute sigma0=%xdiag(hf(ecpolicy))
compute shock=eps0-sigma0
*
dec vect[series] aug93virf(2)
do i=1,2
   set aug93virf(i) 1 nstep = 0.0
end do i
*
do step=1,nstep
   if step==1
      compute hvec=%diag(vbv)*shock
   else
      compute hvec=(%diag(vbv)+%diag(%xdiag(vard)))*hvec
   compute %pt(aug93virf,step,hvec)
end do step
*
do i=1,2
   graph(min=0.0,column=2,row=i,picture="*.###",vticks=5)
   # aug93virf(i) 1 nstep
end do i
spgraph(done)