% TABLE8_line4_2014.m
% in: Baumeister, C. and L. Kilian, "Lower Oil Prices and the U.S. 
%     Economy: Is This Time Different?" Brookings Papers on Economic
%     Activity, Fall 2016

clear;
clc;

% 1970Q2-2016Q1
% Nonresidential investment excluding all oil investment and railroad 
% equipment investment, quarterly
load nriexoil.txt;

% 1970Q3-2016Q1, 183 obs, quarterly growth rate
drinvq=(nriexoil(2:end)-nriexoil(1:end-1,:))./nriexoil(1:end-1,:);

% 1970.1-2016.3, monthly, level
% Percent changes in real U.S. price for gasoline and other motor fuel
load pgasdefl.txt; % PCE deflator for gasoline and other motor fuel

load pcedefl.txt;  % Overall PCE deflator
rpgas=pgasdefl./pcedefl;
for i=1:183
    drpgasq(i,1)=(rpgas(9+(i-1)*3,1)-rpgas(6+(i-1)*3,1))./rpgas(6+(i-1)*3,1);
end;

% Share of gasoline and motor fuel in total PCE expenditures
load pcetotal.txt; % Total nominal consumer expenditures
load pcegas.txt; % Gasoline and motor fuel consumer expenditures
share=pcegas./pcetotal;
for i=1:183
    shareq(i,1)=mean([share(7+(i-1)*3,1) share(8+(i-1)*3,1) share(9+(i-1)*3,1)]);
end;

% Purchasing power (PP) change analogous to Edelstein and Kilian (JME 2009):
% defined as: - drpgas x (share of gasline in total PCE)
time=1970+3/4:1/4:2016+1/4;
pploss=-drpgasq.*shareq;

% Allow for different drift after 2008.12 in investment growth
drinvq=[drinvq(1:154)-mean(drinvq(1:154)); drinvq(155:end)-mean(drinvq(155:end))];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% VAR(4) for PP change and growth in real consumption
y=[pploss drinvq]; [t,K]=size(y); p=4; 
[A,SIGMA,U,V]=olsvarc(y,p); h=7; 

% Restricted VAR(6) for PP change and growth in real consumption with PP
% change as measure of exogenous shock
A=[zeros(1,p*K); A(2:K,:)]; A=[V [A; eye(K*(p-1),K*(p-1)) zeros(K*(p-1),K)]];
y=y';
Y=y(:,p:t);	
for i=1:p-1
 	Y=[Y; y(:,p-i:t-i)];
end;
X=[ones(1,t-p); Y(:,1:t-p)]; Y=Y(:,2:t-p+1); Uhat=Y-A*X; SIGMA=Uhat*Uhat'/(t-p);	
A=A(:,2:K*p+1);

% Compute historical decomposition from structural MA representation
[IRF]=irfvar(A,SIGMA,p,t-p,K); 

% Compute structural shocks Ehat from reduced form shocks Uhat
Ehat=inv(chol(SIGMA(1:K,1:K))')*Uhat(1:K,:);

% Cross-multiply the weights for the effect of a given shock on the real
% oil price (given by the relevant row of IRF) with the structural shock
% in question
yhat1=zeros(t-p,1); yhat2=zeros(t-p,1); 
for i=1:t-p
	yhat1(i,:)=dot(IRF(2,1:i),Ehat(1,i:-1:1));
    yhat2(i,:)=dot(IRF(4,1:i),Ehat(2,i:-1:1));
end;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Bar charts for cumulatuve effect 
% 2014.Q3-2016.Q1
D1=sum(yhat1(end-6:end))*100; %Cumulative effect of PP shock
D2=sum(yhat2(end-6:end))*100; %Cumulative effect of residual shock
disp('2014')
disp(D1)

% 1986.Q1-1987Q3
D1=sum(yhat1(59:65))*100; %Cumulative effect of PP shock
D2=sum(yhat2(59:65))*100; %Cumulative effect of residual shock
disp('1986')
disp(D1)
