The RATS Software Forum

Code: Select all

*
* Replication File for Uhlig (2005), "What are the effects of monetary policy on
* output? Results from an agnostic identification procedure." Journal of Monetary
* Economics, 52, pp 381-419. Pure sign restriction approach
*
open data uhligdata.xls
calendar 1965 1 12
data(format=xls,org=columns) 1965:01 2003:12 gdpc1 gdpdef cprindex totresns bognonbr fedfunds
*
set gdpc1    = log(gdpc1)*100.0
set gdpdef   = log(gdpdef)*100.0
set cprindex = log(cprindex)*100.0
set totresns = log(totresns)*100.0
set bognonbr = log(bognonbr)*100.0
*
system(model=varmodel)
variables gdpc1 gdpdef cprindex fedfunds bognonbr totresns
lags 1 to 12
end(system)
estimate(noprint)
*
dec vect[strings] vl(6)
compute vl=||"real GDP","GDP price defl","Comm. Price Ind.","Fed Funds Rate","Nonborr. Reserv.","Total Reserves"||
*
* n1 is the number of draws from the posterior of the VAR
* n2 is the number of draws from the unit sphere for each draw for the VAR
* nvar is the number of variables
* nstep is the number of IRF steps to compute
* KMAX is the "K" value for the number of steps constrained
*
compute n1=200
compute n2=200
compute nkeep=1000
compute nvar=6
compute nstep=60
compute KMAX=5
*
* This is the standard setup for MC integration of an OLS VAR
*
compute sxx    =%decomp(%xx)
compute svt    =%decomp(inv(%nobs*%sigma))
compute betaols=%modelgetcoeffs(varmodel)
compute ncoef  =%rows(sxx)
compute wishdof=%nobs-ncoef
dec rect ranc(ncoef,nvar)
*
* Most draws are going to get rejected. We allow for up to nkeep good ones. The
* variable accept will count the number of accepted draws. GOODRESP will be a
* RECT(nsteps,nvar) at each accepted draw.
*
declare vect[rect] goodresp(nkeep) goodfevd(nkeep)
declare vector ik a(nvar) ones(nvar)
declare series[rect] irfsquared
compute ones=%const(1.0)
*
compute accept=0
infobox(action=define,progress,lower=1,upper=n1) "Monte Carlo Integration"
do draws=1,n1
   *
   * Make a draw from the posterior for the VAR and compute its impulse
   * responses.
   *
   compute sigmad  =%ranwisharti(svt,wishdof)
   compute swish   =%decomp(sigmad)
   compute ranc    =%ran(1.0)
   compute betau   =sxx*ranc*tr(swish)
   compute betadraw=betaols+betau
   compute %modelsetcoeffs(varmodel,betadraw)
   impulse(noprint,model=varmodel,decomp=swish,results=impulses,steps=nstep)
   gset irfsquared 1 1     = %xt(impulses,t).^2
   gset irfsquared 2 nstep = irfsquared{1}+%xt(impulses,t).^2
   *
   * Do the subdraws over the unit sphere. These give the weights on the
   * orthogonal components.
   *
   do subdraws=1,n2
      compute a=%ran(1.0),a=a/sqrt(%normsqr(a))
      *
      * Check that the responses have the correct signs for steps 1 to KMAX+1
      * (+1 because in the paper, the steps used 0-based subscripts, rather than
      * the 1 based subscripts used by RATS).
      *
      do k=1,KMAX+1
         compute ik=%xt(impulses,k)*a
         if ik(4)<0.or.ik(3)>0.or.ik(2)>0.or.ik(5)>0
            branch 105
      end do k
      *
      * This is an accepted draw. Copy the information out. If we have enough
      * good ones, drop out of the overall loop.
      *
      compute accept=accept+1
      dim goodresp(accept)(nstep,nvar) goodfevd(accept)(nstep,nvar)
      ewise goodresp(accept)(i,j)=ik=%xt(impulses,i)*a,ik(j)
      ewise goodfevd(accept)(i,j)=ik=(irfsquared(i)*(a.^2))./(irfsquared(i)*ones),ik(j)
      if accept>=nkeep
         break
   :105
   end do subdraws
   if accept>=nkeep
      break
   infobox(current=draws)
end do draws
infobox(action=remove)
*
* Post-processing. Graph the mean of the responses along with the 16% and 84%-iles
*
clear upper lower resp
*
spgraph(vfields=3,hfields=2,footer="Figure 6. Impulse Responses with Pure-Sign Approach")
do i=1,nvar
   smpl 1 accept
   do k=1,nstep
      set work = goodresp(t)(k,i)
      compute frac=%fractiles(work,||.16,.84||)
      compute lower(k)=frac(1)
      compute upper(k)=frac(2)
      compute resp(k)=%avg(work)
   end do k
*
   smpl 1 nstep
   graph(ticks,number=0,picture="##.##",header="Impulse Responses for "+vl(i)) 3
   # resp
   # upper / 2
   # lower / 2
end do i
*
spgraph(done)

spgraph(vfields=3,hfields=2,footer="Figure 8. Fraction of Variance Explained with Pure-Sign Approach")
do i=1,nvar
   compute minlower=maxupper=0.0
   smpl 1 accept
   do k=1,nstep
      set work = goodfevd(t)(k,i)
      compute frac=%fractiles(work,||.16,.50,.84||)
      compute lower(k)=frac(1)
      compute upper(k)=frac(3)
      compute resp(k)=frac(2)
   end do k
*
   smpl 1 nstep
   graph(ticks,number=0,min=0.0,picture="##.##",header="Fraction Explained for "+vl(i)) 3
   # resp
   # upper / 2
   # lower / 2
end do i
*
spgraph(done)