Replication for SVAR

[Jay Jeong]

Dear all,

I am studying VAR, Orthogonalized VAR, Structural VAR.
I successfully replicated VAR and OVAR, but I failed to do for SVAR. Could anyone suggest how to replicate it?

Question1 = I am not sure which one should I use for SVAR: oirf or sirf?
Question2 = My Replication does not match the Stata built-in results. How can I fix it?

clear all
import delimited "https://raw.githubusercontent.com/jayjeo/public/master/warehouse/SVARpractice.csv", varnames(1) clear


******** oirf
tsset t
mat A=(1,.\0,1)
mat B=(.,0\0,.)
svar y z, lags(1) aeq(A) beq(B)
irf create order1, set(var2.irf) replace step(10)
irf graph oirf, xlabel(0(2)10) irf(order1) noci yline(0,lcolor(black)) byopts(yrescale)
irf drop order1


******** sirf
tsset t
mat A=(1,.\0,1)
mat B=(.,0\0,.)
svar y z, lags(1) aeq(A) beq(B)
irf create order1, set(var2.irf) replace step(10)
irf graph sirf, xlabel(0(2)10) irf(order1) noci yline(0,lcolor(black)) byopts(yrescale)
irf drop order1


******** My Replication
tsset t
sureg (y L.y L.z)(z L.y L.z)
gen ybly=_b[y:L.y]
gen yblz=_b[y:L.z]
gen zbly=_b[z:L.y]
gen zblz=_b[z:L.z]
putmata ybly yblz zbly zblz
mata A1=ybly[1],yblz[1]\zbly[1],zblz[1]
mata
    B=1,-0.3864334590141\0,1
    Gamma=B*A1
    pi_svar=I(2)
    pi11_svar=pi_svar[1,1]
    pi12_svar=pi_svar[1,2]
    pi21_svar=pi_svar[2,1]
    pi22_svar=pi_svar[2,2]
end
forvalues i=2(1)10000 {
    mata: pi_svar=Gamma*pi_svar
    mata: pi11_svar=pi11_svar\pi_svar[1,1]
    mata: pi12_svar=pi12_svar\pi_svar[1,2]
    mata: pi21_svar=pi21_svar\pi_svar[2,1]
    mata: pi22_svar=pi22_svar\pi_svar[2,2]
}
getmata pi11_svar pi12_svar pi21_svar pi22_svar
preserve
    tsset t
    keep if t<=10
    tsline pi11_svar, name(pi11_svar)
    tsline pi12_svar, name(pi12_svar)
    tsline pi21_svar, name(pi21_svar)
    tsline pi22_svar, name(pi22_svar)
    graph combine pi11_svar pi21_svar pi12_svar pi22_svar
restore

[Jay Jeong]

I made a lot of progress on replication the program with STATA, but it still does not match exactly. I also posted an updated question on Stackoverflow: https://stackoverflow.com/questions/...ation-for-svar

My Replication does not exactly match the R built-in results. How can I fix it?


Below is the built-in R code programmed by Danne, Christian. VARsignR: Estimating VARs Using Sign Restrictions in R, 2015.

 ###### Install required packages
install.packages("minqa")  
install.packages("HI")  
install.packages("mvnfast")
install.packages("lubridate")  
install.packages("VARsignR")    

###### Import data
rm(list = ls())
set.seed(12345)
library(VARsignR)  
SVARdata <- read.csv("https://raw.githubusercontent.com/jayjeo/public/master/LaborShortage/SVARdata_seasonadjusted2_practice.csv")
SVARdata <- ts (SVARdata, frequency = 12, start = c(2012, 1))  

###### set sign restrictions
constr <- c(-1,+2)    

###### Uhlig’s (2005) Rejection Method
model <- uhlig.reject(Y=SVARdata, nlags=1, draws=200, subdraws=100, nkeep=1000, KMIN=1, KMAX=6, constrained=constr, constant=TRUE, steps=120)
irfs <- model$IRFS  
vl <- c("y","z")
irfplot(irfdraws=irfs, type="median", labels=vl, save=FALSE, bands=NULL, grid=TRUE, bw=TRUE)

Below is a STATA code that I wrote.

 *** SVAR Sign-restrictions.  
clear all  

capture program drop subprogram
program subprogram      
      mata mata clear    
     putmata ybly yblz zbly zblz    
     mata A1=ybly[1],yblz[1]\zbly[1],zblz[1]    
     mata sigma7=7600000, 7.76655  \7.76655, .007221    
     mata L=cholesky(sigma7)    
     mata A=rnormal(1,1,0,1),rnormal(1,1,0,1)\rnormal(1,1,0,1),rnormal(1,1,0,1)      
     mata Q=R=.    
     mata qrd(A,Q,R)    
     mata pi_svar=L*Q    
     mata pi_svar1=pi_svar[1,1]    
     mata pi_svar2=pi_svar[2,1]    

     forvalues i=2(1)$obs {        
          mata: pi_svar=A1*pi_svar            
          mata: pi_svar1=pi_svar1\pi_svar[1,1]        
          mata: pi_svar2=pi_svar2\pi_svar[2,1]    
     }     

     getmata pi_svar1 pi_svar2      

     qui tsset t    
     qui keep if t<=120    
     qui keep t pi_svar1 pi_svar2    
     qui rename (pi_svar1 pi_svar2)(p1 p2)    
     qui gen d=0    
     qui replace d=1 if p1[1]<0&p2[1]>0&p1[2]<0&p2[2]>0&p1[3]<0&p2[3]>0&p1[4]<0&p2[4]>0&p1[5]<0&p2[5]>0&p1[6]<0&p2[6]>0
end  

import delimited "https://raw.githubusercontent.com/jayjeo/public/master/LaborShortage/SVARdata_seasonadjusted2_practice.csv", varnames(1) clear  
global obs=120
gen t=_n
tsset t
sureg (y L.y L.z)(z L.y L.z)
gen ybly=_b[y:L.y]
gen yblz=_b[y:L.z]
gen zbly=_b[z:L.y]
gen zblz=_b[z:L.z]  

reg y L.y L.z
predict ey, residual
reg z L.y L.z
predict ez, residual
corr ey ez, covariance  
*mata sigma7=7600000, 7.76655  \7.76655, .007221  

forvalues i=1(1)100 {      
      preserve                  
          qui subprogram        
          qui save pp`i', replace          
          scalar process=`i'        
          di process    
     restore  
}  

use pp1, clear
forvalues i=2(1)100 {      
      append using pp`i'
}
drop if d==0
collapse (median) p1 p2, by(t)  

tsset t
tsline p1, name(p1) yline(0)
tsline p2, name(p2) yline(0)
graph combine p1 p2, rows(2)

Below figure is the result using the built-in R program by Danne, Christian. "VARsignR"

Below figure is the result using the STATA code that I wrote (which does not match exactly).