The following program illustrates the two forms of imposing the identifying
restrictions in structural decompositions. For illustration purposes and to
check that the restrictions are correctly imposed, we impose restrictions
that replicate the Cholesky factorization. (If all you want is the standard
impulse response based on the Cholesky factorization, there is no need to
use the structrual decomposition.) The program checks whether the factorization
matrix from the three methods match by computing the L-infinity norm
(maximum of absolute column sums) of the matrix difference. Theoretically,
the difference should be zero but because of floating errors the numerical
difference won't be exactly zero.
' replicate cholesky factorization using SVAR
' for version 4
' 11/5/99 h
' last checked 10/27/2000 h
' include subroutine
include sub_rmaxdiff.prg
'change path to program path
%path = @runpath
cd "{%path}"
' create workfile
workfile cholsvar q 1948:1 1979:3
' fetch data from database
fetch(d=data_svar) rgnp rinv m1
' take log of each series
series lgnp = log(rgnp)
series linv = log(rinv)
series lm1 = log(m1)
' estimate unrestricted VAR
var var1.ls 1 4 lgnp linv lm1 @ c
'-------------------------------------------------------------------
' method 1: short-run restrictions in text form
'-------------------------------------------------------------------
var1.cleartext(svar)
var1.append(svar) @e1 = c(1)*@u1
var1.append(svar) @e2 = -c(2)*@e1 + c(3)*@u2
var1.append(svar) @e3 = -c(4)*@e1 - c(5)*@e2 + c(6)*@u3
freeze(tab1) var1.svar(rtype=text,conv=1e-5)
'show tab1
' store estimated A and B matrices from method 1
matrix mata1 = var1.@svaramat
matrix matb1 = var1.@svarbmat
' compute factorization matrix
matrix fact1 = @inverse(mata1)*matb1
'-------------------------------------------------------------------
' method 2: short-run restrictions in pattern matrix
'-------------------------------------------------------------------
' create pattern matrices
matrix(3,3) pata
' fill matrix in row major order
pata.fill(by=r) 1,0,0, na,1,0, na,na,1
matrix(3,3) patb
patb.fill(by=r) na,0,0, 0,na,0, 0,0,na
freeze(tab2) var1.svar(rtype=patsr,conv=1e-5,namea=pata,nameb=patb)
'show tab2
' store estimated A and B matrices from method 2
matrix mata2 = var1.@svaramat
matrix matb2 = var1.@svarbmat
' compute factorization matrix
matrix fact2 = @inverse(mata2)*matb2
'-------------------------------------------------------------------
' method 3: built-in cholesky
'-------------------------------------------------------------------
' do impulse response with default Cholesky ordering
var1.impulse(10)
' store cholesky factor
matrix fact3 = var1.@impfact
'-------------------------------------------------------------------
' check whether results match
'-------------------------------------------------------------------
' compute difference in factorization matrices
scalar diff13
call sub_rmaxdiff(fact1, fact3, diff13)
scalar diff23
call sub_rmaxdiff(fact2, fact3, diff23)
' display L-infinity norm of matrix difference in table
table(2,2) check
setcolwidth(check,1,20)
check(1,1) = "method1 - method3:"
check(1,2) = diff13
check(2,1) = "method2 - method3:"
check(2,2) = diff23
show check
The following program replicates the long-run restriction model by
Blanchard-Quah (1989). (The data are not the same and the results won't
exactly match.) The DY (output growth) and U (unemployment) series are already
demeaned following the procedure by Blanchard-Quah (1989). The program
estimates the factorization matrix by three different methods which should
all yield the same result.
' Blanchard-Quah long-run restriction
' verify estimates using three methods
' for version 4
' 11/5/99 h
' last checked 10/27/2000 h
' include subroutine
include sub_rmaxdiff.prg
' create workfile
workfile blanquah q 1948:1 1987:4
' fetch data from database (already demeaned)
fetch(d=data_svar) dy u
' estimate VAR (no constant)
var var1.ls 1 8 dy u @
'-------------------------------------------------------------------
' method 1: text form long-run restrictions
'-------------------------------------------------------------------
var1.cleartext(svar)
var1.append(svar) @lr1(@u1)=0
freeze(tab1) var1.svar(rtype=text,conv=1e-4)
show tab1
' store estimated A and B matrices from method 1
matrix mata1 = var1.@svaramat ' A should be identity
matrix matb1 = var1.@svarbmat
'-------------------------------------------------------------------
' method 2: long-run restrictions in short-run form
' *not recommended; only for checking*
'-------------------------------------------------------------------
' get unit (cumulated) long-run response
var1.impulse(imp=u)
matrix clr = var1.@lrrsp
var1.cleartext(svar)
var1.append(svar) @e1 = c(1)*@u1 + c(2)*@u2
var1.append(svar) @e2 = -clr(1,1)*c(1)/clr(1,2)*@u1 + c(4)*@u2
freeze(tab2) var1.svar(rtype=text,conv=1e-5)
show tab2
' store estimated A and B matrices from method 2
matrix mata2 = var1.@svaramat ' A should be identity
matrix matb2 = var1.@svarbmat
'-------------------------------------------------------------------
' method 3: moment matching (used in Blanchard-Quah (1989) paper)
' *only works for just-identified models*
'-------------------------------------------------------------------
' get residual covariance
sym rcov = var1.@residcov
' and do cholesky factorization
matrix chol = @cholesky(rcov)
' factorize unit long-run response
matrix p = clr * chol
' factorized B matrix
matrix(2,2) q
' impose sign normalization by taking positive root
q(1,1) = @sqrt( 1/(1+p(1,1)*p(1,1)/p(1,2)/p(1,2)) )
q(2,1) = -p(1,1)*q(1,1)/p(1,2)
q(1,2) = -q(2,1)
q(2,2) = q(1,1)
' get implied B matrix
matrix matb3 = chol * q
'-------------------------------------------------------------------
' check whether results match
'-------------------------------------------------------------------
' difference in A matrices
matrix eye2 = @identity(2) ' truth
scalar adiff1
call sub_rmaxdiff(mata1, eye2, adiff1)
scalar adiff2
call sub_rmaxdiff(mata2, eye2, adiff2)
' difference in B matrices
scalar bdiff13
call sub_rmaxdiff(matb1, matb3, bdiff13)
scalar bdiff23
call sub_rmaxdiff(matb2, matb3, bdiff23)
' display L-infinity norm of matrix difference in table
table(2,3) check
setcolwidth(check,1,10)
setcolwidth(check,2,15)
setcolwidth(check,3,15)
check(1,2) = "A - eye2"
check(1,3) = "B - method3"
check(2,1) = "method1:"
check(3,1) = "method2:"
check(2,2) = adiff1
check(3,2) = adiff2
check(2,3) = bdiff13
check(3,3) = bdiff23
show check
The following program replicates the short-run restriction exercises in
Sims (1986). (The data are not the same and the results do not match exactly.)
Note that while Sims assumes a diagonal covariance matrix for the structural
innovations, EViews assumes an identity covariance matrix. For this reason, we
need to estimate the standard deviation of the structural shocks as elements of
the B matrix.
' replicate Sims (1986) with similar data
' for version 4
' 1/12/2000 h
' last checked 10/27/2000 h
'change path to program path
%path = @runpath
cd "{%path}"
' create workfile
workfile sims1 q 1948:1 1979:3
' fetch data from database
fetch(d=data_svar) rgnp rinv pidgnp m1 tb3sec unemp
' transform series
series lgnp = log(rgnp)
series linv = log(rinv)
series lprc = log(pidgnp)
series lm1 = log(m1)
'----------------------------------------------------------------------
' replicate Chart 1 (p.11) *no Bayesian priors imposed*
'----------------------------------------------------------------------
var var1.ls 1 4 lgnp linv lprc lm1 unemp tb3sec @ c
freeze(chart1) var1.impulse(32,m,imp=chol,se=a)
show chart1
' plot only some of money innovation responses
freeze(chart1a) var1.impulse(32,m,imp=chol,se=a) lgnp lprc tb3sec @ lm1
'----------------------------------------------------------------------
' replicate Chart 2 (p.13)
'----------------------------------------------------------------------
' reorder VAR to confirm to Sims' (1986) ordering
var var1.ls 1 4 tb3sec lm1 lgnp lprc unemp linv @ c
' 1st identification restrictions (p.12 (12)-(17))
coef(6) b
var1.append(svar) @e1 = c(1)*@e2 + b(1)*@u1
var1.append(svar) @e2 = c(2)*@e3 + c(3)*@e4 + c(4)*@e1 + b(2)*@u2
var1.append(svar) @e3 = c(5)*@e1 + c(6)*@e6 + b(4)*@u4
var1.append(svar) @e4 = c(7)*@e1 + c(8)*@e3 + c(9)*@e6 + b(5)*@u5
var1.append(svar) @e5 = c(10)*@e1 + c(11)*@e3 + c(12)*@e6 + c(13)*@e4 + b(6)*@u6
var1.append(svar) @e6 = b(3)*@u3
' estimate structural factorization
freeze(svar1) var1.svar(rtype=text,conv=1e-5)
show svar1
' replicate Chart 2 (p.13)
freeze(chart2) var1.impulse(32,m,imp=struct,se=a) lgnp linv lprc lm1 unemp tb3sec @ 1 2 4 5 6 3
show chart2
' plot only some of money supply responses
freeze(chart2a) var1.impulse(32,m,imp=struct,se=a) lgnp lprc tb3sec @ 1
'----------------------------------------------------------------------
' replicate Chart 3 (p.13)
'----------------------------------------------------------------------
' clear previous restrictions
var1.cleartext(svar)
' 2nd identification restrictions (p.14 (18)-(23))
var1.append(svar) @e1 = c(1)*@e2 + b(1)*@u1
var1.append(svar) @e2 = c(2)*@e3 + c(3)*@e6 + c(4)*@e4 + c(5)*@e1 + b(2)*@u2
var1.append(svar) @e3 = c(6)*@e1 + c(7)*@e6 + b(3)*@u3
var1.append(svar) @e4 = c(8)*@e1 + c(9)*@e3 + c(10)*@e2 + b(5)*@u5
var1.append(svar) @e5 = c(11)*@e1 + c(12)*@e3 + c(13)*@e4 + c(14)*@e6 + b(6)*@u6
var1.append(svar) @e6 = b(4)*@u4
' estimate structural factorization
' default starting value does not converge
rndseed 123456
freeze(svar2) var1.svar(rtype=text,conv=1e-5,f0=u)
show svar2
' replicate Chart 3 (p.13)
freeze(chart3) var1.impulse(32,m,imp=struct,se=a) lgnp linv lprc lm1 unemp tb3sec @ 1 2 3 5 6 4
show chart3
' plot only some of money supply responses
freeze(chart3a) var1.impulse(32,m,imp=struct,se=a) lgnp lprc tb3sec @ 1
The following program replicates the long-run restriction model by
Blanchard-Quah (1989). (The data are not the same and the results do not
exactly match.) The DY (output growth) and U (unemployment) series are already
demeaned following the procedure in Blanchard-Quah (1989). In order to plot the
accumulated response and the level response in one graph, we first store the
impulse responses in separate matrices, extract and place the relevant columns into a
third matrix and use the line view of that matrix.
' replicate impulse response graph
' from long-run restriction (Blanchard-Quah 1989)
' for version 4
' 11/12/99 h
' last checked 10/27/2000 h
' include subroutine
include sub_rmaxdiff.prg
' create workfile
workfile blanquah q 1948:1 1987:4
' fetch data from database (already demeaned)
fetch(d=data_svar) dy u
' change to percentage
dy = 100*dy
' estimate VAR (no constant)
var var1.ls 1 8 dy u @
' impose long-run restrictions
' @u1 = aggregate demand shock
' @u2 = aggregate supply shock
var1.cleartext(svar)
var1.append(svar) @lr1(@u1)=0
' estimate factorization
freeze(tab1) var1.svar(rtype=text,conv=1e-5)
' set response horizon
!hrz = 40
' replicate Figures 3-4 (p.662)
' accumulate to get level response of output
freeze(fig3) var1.impulse(!hrz,imp=struct,se=a,a,matbys=rsp_acc) dy @ 1
freeze(fig4) var1.impulse(!hrz,imp=struct,se=a,a) dy @ 2
freeze(fig34) fig3 fig4
show fig34
' replicate Figures 5-6
freeze(fig5) var1.impulse(!hrz,imp=struct,se=a,matbys=rsp_lvl) u @ 1
freeze(fig6) var1.impulse(!hrz,imp=struct,se=a) u @ 2
freeze(fig56) fig5 fig6
show fig56
' replicate Figure 1
matrix(!hrz,2) mtmp
' get accumulated output response
vector v = @columnextract(rsp_acc,1)
colplace(mtmp,v,1)
' get level unemployment response
v = @columnextract(rsp_lvl,2)
colplace(mtmp,v,2)
freeze(fig1) mtmp.line
fig1.draw(line,left) 0
fig1.elem(1) legend(Output response to demand)
fig1.elem(2) legend(Unemployment response to demand)
show fig1
' replicate Figure 2
' get accumulated output response
vector v = @columnextract(rsp_acc,3)
colplace(mtmp,v,1)
' get level unemployment response
v = @columnextract(rsp_lvl,4)
colplace(mtmp,v,2)
freeze(fig2) mtmp.line
fig2.draw(line,left) 0
fig2.elem(1) legend(Output response to supply)
fig2.elem(2) legend(Unemployment response to supply)
show fig2
Currently EViews does not provide standard error bounds for impulse
responses from structural impulses identified by long-run restrictions.
The following program illustrates how to obtain bootstrap impulse response
confidence bounds for the Blanchard-Quah (1989) model. The bootstrap is
performed as follows:
- Resample from the VAR residuals with replacement.
- Obtain bootstrap data using the estimated VAR coefficients and
the bootstrap residuals.
- Estimate the VAR and structural factorization using the bootstrap
data from step 2.
- Store the impulse responses using estimates from step 3.
- Repeat steps 1-4 many times.
Step 2 is implemented by specifying the bootstrap residuals as (level shift)
add factors in model solve. Note also that estimation of the structural
factorization at each bootstrap step may fail to converge; in such cases
the code only stores results if convergence was achieved. The code below
does 100 bootstrap replications; you should probably increase this number
to get more reliable estimates but be forewarned that this make take a
very long time to execute.
' bootstrap structural VAR impulse responses
' from long-run restriction (Blanchard-Quah 1989)
' for version 4
' 11/28/2000 h
include sub_bootirf.prg
' create workfile
workfile blanquah q 1948:1 1987:4
' fetch data from database (already demeaned)
fetch(d="i:\hiroyuki\example files\var\data_svar") dy u
' change to percentage
dy = 100*dy
' estimate VAR (no constant)
var var1.ls 1 8 dy u @
' impose long-run restrictions
' @u1 = aggregate demand shock
' @u2 = aggregate supply shock
var1.cleartext(svar)
var1.append(svar) @lr1(@u1)=0
' estimate factorization
var1.svar(rtype=text,conv=1e-5)
' set response horizon
!hrz = 40
' store level and accumulated responses
do var1.impulse(!hrz,imp=struct,matbys=rsp_lvl) dy u @ 1
do var1.impulse(!hrz,imp=struct,a,matbys=rsp_acc) dy u @ 1
' get residuals to bootstrap
var1.makeresid(n=gres) res_dy res_u
' make model to generate bootstrap data
var1.makemodel(mod1)
' assign add factors for bootstrap residuals
mod1.addassign(i) @all
'set random number generator
rndseed(type=mt) 1234567
'set number of bootstrap replications
!reps = 100
'declare storage matrices
for !i=1 to @columns(rsp_lvl)
matrix(!hrz,!reps) brsp_lvl{!i}
matrix(!hrz,!reps) brsp_acc{!i}
next
'temporary working vector
vector vtmp
'count number of cases that converged
scalar cnt = 0
'bootstrap loop
for !i=1 to !reps
'bootstrap residuals
smpl @all if res_dy<>na
gres.resample dy_a u_a
'generate bootstrap data
mod1.solve
'estimate VAR (no constant) with bootstrap data
smpl @all
var var2.ls 1 8 dy_0 u_0 @
'impose long-run restrictions
var2.cleartext(svar)
var2.append(svar) @lr1(@u1)=0
'estimate factorization
var2.svar(rtype=text,conv=1e-5,maxiter=100,nostop)
'store results only if converged
if (var2.@svarcvgtype=0) then
'increment counter
cnt = cnt + 1
'get level and accumulated bootstrap responses
do var2.impulse(!hrz,imp=struct,matbys=tmp_lvl) dy_0 u_0 @ 1
do var2.impulse(!hrz,imp=struct,a,matbys=tmp_acc) dy_0 u_0 @ 1
'and store in appropriate position
for !j=1 to @columns(rsp_lvl)
vtmp = @columnextract(tmp_lvl,!j)
colplace(brsp_lvl{!j}, vtmp, cnt)
vtmp = @columnextract(tmp_acc,!j)
colplace(brsp_acc{!j}, vtmp, cnt)
next
delete tmp_lvl tmp_acc
endif
next
'reshape results matrix if there were non-converged cases
if (cnt <> !reps) then
for !i=1 to @columns(rsp_lvl)
brsp_lvl{!i} = @subextract(brsp_lvl{!i}, 1, 1, !hrz, cnt)
brsp_acc{!i} = @subextract(brsp_acc{!i}, 1, 1, !hrz, cnt)
next
'display warning in table
table warning
%msg = @str(!reps-cnt) + " case(s) out of " + @str(!reps) + " replications failed to converge"
warning(1,1) = %msg
show warning
endif
'find bootstrap 66 percentiles for each horizon
matrix(!hrz,@columns(rsp_lvl)) lvl_upp
matrix(!hrz,@columns(rsp_lvl)) lvl_low
matrix(!hrz,@columns(rsp_lvl)) acc_upp
matrix(!hrz,@columns(rsp_lvl)) acc_low
rowvector rtmp
for !i=1 to !hrz
for !j=1 to @columns(rsp_lvl)
rtmp = @rowextract(brsp_lvl{!j},!i)
lvl_upp(!i,!j) = @quantile(rtmp,0.833)
lvl_low(!i,!j) = @quantile(rtmp,0.167)
rtmp = @rowextract(brsp_acc{!j},!i)
acc_upp(!i,!j) = @quantile(rtmp,0.833)
acc_low(!i,!j) = @quantile(rtmp,0.167)
next
next
'plot impulse responses with bootstrap error bounds
%gname = "fig_lvl"
call sub_bootirf(rsp_lvl, lvl_upp, lvl_low, %gname)
' add zero line to all graphs
{%gname}.draw(line,left,rgb(180,180,180)) 0
' add title to graph
{%gname}.addtext(0.1,-0.2) Bootstrap 66% Confidence Bounds (Levels)
' realign
{%gname}.align(2,1,0.5)
show {%gname}
%gname = "fig_acc"
call sub_bootirf(rsp_acc, acc_upp, acc_low, %gname)
' add zero line to all graphs
{%gname}.draw(line,left,rgb(180,180,180)) 0
' add title to graph
{%gname}.addtext(0.1,-0.2) Bootstrap 66% Confidence Bounds (Accumulated Responses)
' realign
{%gname}.align(2,1,0.5)
show {%gname}