% TABLE8_line1.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;

% Exact monthly growth rates of U.S. real consumption
% Total real U.S. PCE, monthly 1970.1-2016.3, quantity index, SA
% Source: BEA
load rctotal.txt;

load transport.txt; load industrial.txt; load commercial.txt;
motorgas=(transport+industrial+commercial);
load gasimports.txt; importshare=gasimports./motorgas; mu=mean(importshare);
gasimportshare=deseasonal(importshare)+mu;

drctotal=(rctotal(2:end)-rctotal(1:end-1,:))./rctotal(1:end-1,:);

% Percent changes in real U.S. price for gasoline and other motor fuel
load pgasdefl.txt; % PCE deflator for gasoline and other motor fuel

% Modified
%load rameyfactor.txt; pgasdefl=pgasdefl.*rameyfactor;

load pcedefl.txt;  % Overall PCE deflator
rpgas=pgasdefl./pcedefl;
drpgas=(rpgas(2:end)-rpgas(1:end-1,:))./rpgas(1:end-1,:);

% 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;

% MER US data, 1973.1-2016.3
load expcrude.txt; load impcrude.txt; load prodcrude.txt;
ratio=(impcrude-expcrude)./(prodcrude+impcrude-expcrude); mu=mean(ratio); ratio=deseasonal(ratio)+mu;

% Purchasing power (PP) change 
time=1973+1/12:1/12:2016+3/12;
pploss=-drpgas.*share(2:end); pploss=pploss(36:end,1).*(gasimportshare+(1-gasimportshare).*ratio);

% Allow for different drift after 2008.12 in consumption growth
drctotal=[drctotal(1:467)-mean(drctotal(1:467)); drctotal(468:end)-mean(drctotal(468:end))];

% 1973.1-2016.3
drctotal=drctotal(36:end);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% VAR(6) for PP change and growth in real consumption, PP change
% predetermined, growth response to 1% shock (= 1std dev shock) in PP is 
% cumulated
y=[pploss drctotal]; [t,K]=size(y);  p=6;
[A,SIGMA,U,V]=olsvarc(y,p); h=21; 

% 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.7-2016.3
D1=sum(yhat1(528-36:end))*100; %Cumulative effect of PP shock
D2=sum(yhat2(528-36:end))*100; %Cumulative effect of residual shock
disp('2014')
disp(D1)

% 1986.1-1987.9
D1=sum(yhat1(151:171))*100; %Cumulative effect of PP shock
D2=sum(yhat2(151:171))*100; %Cumulative effect of residual shock
disp('1986')
disp(D1)

