% FIGURE1A.M 

clear;

% Load posterior draws
load posterior_arrw0; z=za; clear za; 

% Shells for impulse responses of admissible structural models
IRF11=zeros(size(z,1),H+1); IRF12=zeros(size(z,1),H+1); IRF13=zeros(size(z,1),H+1); 
IRF21=zeros(size(z,1),H+1); IRF22=zeros(size(z,1),H+1); IRF23=zeros(size(z,1),H+1); 
IRF31=zeros(size(z,1),H+1); IRF32=zeros(size(z,1),H+1); IRF33=zeros(size(z,1),H+1);

for j=1:size(z,1)
        
    A      = (reshape(z(j,1:(1+n*p)*n),1+n*p,n))';   % RF level parameters
    B0inv  = reshape(z(j,(1+n*p)*n+n*n+1:(1+n*p)*n+2*n*n),n,n);    
    irf=irfvar([A(:,2:p*n+1); eye((p-1)*n,(p-1)*n) zeros((p-1)*n,n)],B0inv,p,H);	                    
    IRF11(j,:)=irf(1,:);   IRF21(j,:)=irf(2,:);   IRF31(j,:)=irf(3,:); 
    IRF12(j,:)=irf(4,:);   IRF22(j,:)=irf(5,:);   IRF32(j,:)=irf(6,:);  
    IRF13(j,:)=irf(7,:);   IRF23(j,:)=irf(8,:);   IRF33(j,:)=irf(9,:);

end

% Stack posterior draws for impulse responses
IRF=[IRF11 IRF12 IRF13 IRF21 IRF22 IRF23 IRF31 IRF32 IRF33]; 

clear A B0inv irf IRF11 IRF12 IRF13 IRF21 IRF22 IRF23 IRF31 IRF32 IRF33;

% Bayes estimate of vector of inpulse responses for joint inference and 68% joint credible set
[irfbayes,jcs,index]=jointL1(IRF); clear IRF;    

IRF12bayes=irfbayes(1,1*(H+1)+1:2*(H+1));
IRF22bayes=irfbayes(1,4*(H+1)+1:5*(H+1));
IRF32bayes=irfbayes(1,7*(H+1)+1:8*(H+1));

IRF12jcs=jcs(:,1*(H+1)+1:2*(H+1));
IRF22jcs=jcs(:,4*(H+1)+1:5*(H+1));
IRF32jcs=jcs(:,7*(H+1)+1:8*(H+1));

clear IRF jcs irfbayes index;

% Rescale responses of first variable to facilitate interpretation of IRF plot
IRF12jcs=IRF12jcs*100; IRF12bayes=IRF12bayes*100;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Plot structural impulse responses 
subplot(3,3,1); plot(0:1:H,IRF12jcs,'r-',0:1:H,IRF12bayes,'k-',0:1:H,zeros(1,H+1),'k-','linewidth',3); axis([0 H  -9 9])
ylabel('Real gas price','fontsize',15)
title('Nominal gasoline price shock','fontsize',15)
subplot(3,3,2); plot(0:1:H,IRF22jcs,'r-',0:1:H,IRF22bayes,'k-',0:1:H,zeros(1,H+1),'k-','linewidth',3); axis([0 H -0.1 0.25])
ylabel('Headline inflation','fontsize',15)
title('Nominal gasoline price shock','fontsize',15)
subplot(3,3,3); plot(0:1:H,IRF32jcs,'r-',0:1:H,IRF32bayes,'k-',0:1:H,zeros(1,H+1),'k-','linewidth',3); axis([0 H -0.1 0.5])
ylabel('Inflation expectation','fontsize',15)
title('Nominal gasoline price shock','fontsize',15)


