The RATS Software Forum

Code: Select all

system(model=var)
variables dgs3mo spret curve
lags 1 to 5
det constant
end(system)
estimate(resids=rvar)
*
compute gend   = %regend()
compute gstart = %regstart()
*
dec series[symm] h uu
*
gset h  * gend = %sigma
gset uu * gend = %sigma
*
* The covariance matrix is modelled in vec form as
*
* vec(H)=BL x Psi + BL x Gamma x B**-2 diag(H(t-1)) +
*     BL x Lambda x B**-2 diag (UU'(t-1))
*
* where Psi, Gamma, Lambda, and B are to be estimated, with B being the
* inverse of a matrix A with unit diagonal. BL is a matrix of squared
* interactions of members of B. The free coefficients are listed as
*
*    psi    n vector
*    gamma  n x n matrix
*    lambda n x n matrix
*    ap     n x n-1 matrix with the n-1 elements in each row
*            giving the non-diagonal elements of A.
*
compute n=3
compute nstack=n*(n+1)/2
dec vect psi(n)
dec rect gamma(n,n) lambda(n,n)
dec rect ap(n,n-1)
dec vect gammav(n) lambdav(n)
*
dec vect garchc
dec rect garcha garchb
*
*
* The function initrf is called as an initialization function for
* MAXIMIZE to compute the time-invariant matrices in the formula for
* vectorized covariance. Those matrices are called GARCHA, GARCHB and
* GARCHC.
*
function initrf
local integer i j k fill
dec rect a b bl
dim a(n,n)
dim bl(nstack,n)
*
ewise a(i,j)=%if(i==j,1.0,ap(i,j-(j>i)))
compute b=inv(a)
do k=1,n
   compute fill=0
   do i=1,n
      do j=1,i
         compute fill=fill+1
         compute bl(fill,k)=b(i,k)*b(j,k)
      end do j
   end do i
end do k
compute garchc=bl*psi
compute garcha=bl*gamma*a*a
compute garchb=bl*lambda*a*a
end
*
*
dec frml[symm] hf
*
frml hf = %vectosymm(garchc+garcha*%xdiag(h{1})+garchb*%xdiag(uu{1}),n)
frml logl = $
    hx = hf(t) , $
    ux = %xt(rvar,t), $
    h(t)=hx, uu(t)=%outerxx(ux), $
    %logdensity(hx,ux)
*
* Initial guess values
*
compute psi=%xdiag(%sigma)*.10
compute ap=%const(0.0)
compute gammav=%const(.80)
compute lambdav=%const(.10)
compute gamma=%diag(gammav)
compute lambda=%diag(lambdav)
*
* Estimate the restricted model with only diagonal elements in gamma and
* lambda.
*
nonlin(parmset=garchparms) psi gammav lambdav ap
*
maximize(trace,start=%(gamma=%diag(gammav),lambda=%diag(lambdav),initrf()),$
  parmset=garchparms,pmethod=simplex,piters=5,method=bfgs,iters=200) logl gstart gend
*
* Free up lambda
*
nonlin(parmset=garchparms) psi gammav lambda ap
maximize(trace,start=%(gamma=%diag(gammav),initrf()),$
  parmset=garchparms,pmethod=simplex,piters=5,method=bfgs,iters=200) logl gstart gend
*
* Convert from conditional variances of the non-orthogonalized
* residuals, to the structural ones.
*
set hi gstart gend = aha=(a*h*tr(a)),aha(1,1)
set hs gstart gend = aha=(a*h*tr(a)),aha(2,2)
set hy gstart gend = aha=(a*h*tr(a)),aha(3,3)
spgraph(footer="Conditional Variances of Structural Innovations",vfields=3)
graph(header="Innovation to Short Rate")
# hi
graph(header="Innovation to Stock Prices")
# hs
graph(header="Innovation to Yield Curve Slope")
# hy
spgraph(done)


* conditional correlations
*
set rho12 = h(t)(1,2)/sqrt(h(t)(1,1)*h(t)(2,2))
set rho13 = h(t)(1,3)/sqrt(h(t)(1,1)*h(t)(3,3))
set rho23 = h(t)(2,3)/sqrt(h(t)(2,2)*h(t)(3,3))
spgraph(footer="Conditional Correlations between financial variables",vfields=3)
graph(header="Short Rate and Stock Prices")
# rho12
graph(header="Short Rate and Yield Curve Slope")
# rho13
graph(header="Yield Curve Slope and Stock Prices")
# rho23
spgraph(done)
*
*