% TABULARWELFAREGIANDECOMP.M
%
% Produces the welfare gain table of "Production Networks and Monetary
% Policy in Distorted Economies". Run from the root of the replication
% package (the directory that contains data_clean/ and functions/). The main
% text figures are produced by main.m and the appendix figures by appendix.m.
%
%   Output                                         Paper
%   --------------------------------------------   ----------------
%   LaTeX table echoed to the Command Window       Table 2
%
% No file is written: latexTable.m builds the table and displays it, so the
% printed block is pasted into the paper source.
%
% Data inputs, all under data_clean/:
%   FPA/FPA_sum.mat          FPA                       66 x 1
%   IO_table/IO_sum.mat      IO_matrix_sum             79 x 94 x 25
%   markup/markup_sum.mat    mu_sum                    66 x 19 x 3
%   productivity/MFP.mat     Sigmaz                    66 x 66
%
% Functions called, all under functions/:
%   optimalpolicy.m, solveoutputgap.m, welfarefunction.m,
%   expectedwelfaregain.m, latexTable.m, weightedcov.m
%
% solveoutputgap.m and welfarefunction.m read the calibration from the
% global workspace rather than through their argument lists, which is why
% the global declaration below is required before they are called.
%
% Markup measures are stored in the third dimension of mu_sum:
%   (:,:,1) user cost (UC), (:,:,2) accounting profit (AP),
%   (:,:,3) production function (PF). The table reports UC markups.
%
% The Monte Carlo draws 10,000 productivity shocks and evaluates three
% policies at each draw. The seed is fixed below, so the reported numbers
% are exactly reproducible.

clc;
clear;

Path_FPA          = [cd, '/data_clean', '/FPA', '/FPA_sum.mat'];
Path_IO           = [cd, '/data_clean', '/IO_table', '/IO_sum.mat'];
Path_markup       = [cd, '/data_clean', '/markup', '/markup_sum.mat'];
Path_productivity = [cd, '/data_clean', '/productivity', '/MFP.mat'];


load(Path_FPA);
load(Path_IO);
load(Path_markup);
load(Path_productivity);

addpath('functions');


global N mu theta sigma_c sigma varepsilon omega_tilde alpha_tilde beta_tilde gamma varphi

%% Baseline calibration
% Unlike main.m and appendix.m, which evaluate the model year by year, the
% table is computed at a single calibration obtained by averaging the
% input-output coefficients, labor shares, final expenditure shares and
% markups over 1997-2015. The loop below accumulates those averages into
% omega_tilde, alpha_tilde, beta_tilde and mu.

% preference parameters
gamma          = 1;
varphi         = 2;

% number of industries

N              = 66;

% transforming the data from monthly to quarterly and calculating the
% degree of price flexibility of each industry
theta          = 1 - exp(-3*FPA);

% Calibrate the IO Table

IO             = IO_matrix_sum;
Final          = IO(:,93,:);
IO(isnan(IO))  = 0;
IO(IO < 0)     = 0;

% Labor as single input
grosscost      = squeeze(sum(IO(1:66,1:66,:)))+squeeze(IO(75,1:66,:));
labor          = squeeze(IO(75,1:66,:))./grosscost;
IO(67:end,:,:) = [];
IO(:,67:end,:) = [];


Omega_sum      = zeros(66,66);
beta_sum       = zeros(66,1);

for year = 1997:2015

    temp         = IO(:,:,year-1996);
    temp         = temp';
    Omega_sum    = Omega_sum + bsxfun(@rdivide, temp, grosscost(:,year-1996));
    beta_tilde   = Final(1:66,year-1996);
    beta_tilde(beta_tilde<0) = 0;
    beta_tilde   = beta_tilde/sum(beta_tilde);
    beta_sum     = beta_sum + beta_tilde;
end

omega_tilde  = Omega_sum/19;
alpha_tilde  = mean(labor(:,1:19),2);
beta_tilde   = beta_sum/19;
mu           = mean(mu_sum(:,:,1),2);

varepsilon   = 6;
sigma        = 0.5 * ones(N,1);
sigma_c      = 1;


%% Drawing log-productivity shocks with covariance matrix Sigmaz

% number of iterations
iteration = 10000;

% initialize the random number generator to make the results reproducable
rng('default');
rng(1);

% drawing productivity shocks with variance-covariance matrix Sigmaz
dlogz = mvnrnd(zeros(N,1),Sigmaz,iteration);
dlogz = dlogz';


%% Monte Carlo simulation
% Three policies are compared, each described by a vector of sectoral
% inflation weights and an inflation bias:
%   optimal policy         phi_star, pi_star   (Proposition 5)
%   output-gap targeting   phi_og,   0
%   CPI targeting          phi_cpi,  0         weights are the final
%                                              expenditure shares beta_tilde
% For every productivity draw, solveoutputgap.m returns the output gap the
% policy implies and welfarefunction.m evaluates welfare at that gap. The
% seven columns of w_cpi, w_og and w_star are, in order: total welfare, the
% first-order aggregate wedge and supply-side terms, and the second-order
% output gap, within-sector price dispersion, cross-sector price dispersion
% and allocative efficiency terms. w_dc is preallocated but not used: the
% divine coincidence policy is not evaluated here.

[phi_star,pi_star,~,~,phi_og,~,~,~,phi_dc, phi_dc_eff]  = optimalpolicy(mu, theta, sigma_c, sigma, varepsilon, omega_tilde, alpha_tilde, beta_tilde, gamma, varphi);

phi_cpi      = beta_tilde';

iteration = size(dlogz,2);


w_cpi        = zeros(iteration,7);
w_og         = zeros(iteration,7);
w_star       = zeros(iteration,7);
w_dc         = zeros(iteration,7);

for i=1:iteration
    y_tilde_cpi          = solveoutputgap(phi_cpi, 0, dlogz(:,i));
    y_tilde_og           = solveoutputgap(phi_og, 0, dlogz(:,i));
    y_tilde_star         = solveoutputgap(phi_star, pi_star, dlogz(:,i));


    [w_cpi(i,1),w_cpi(i,2),w_cpi(i,3),w_cpi(i,4),w_cpi(i,5),w_cpi(i,6),w_cpi(i,7)]...
        = welfarefunction (y_tilde_cpi, dlogz(:,i));

    [w_og(i,1),w_og(i,2),w_og(i,3),w_og(i,4),w_og(i,5),w_og(i,6),w_og(i,7)]...
        = welfarefunction (y_tilde_og, dlogz(:,i));

    [w_star(i,1),w_star(i,2),w_star(i,3),w_star(i,4),w_star(i,5),w_star(i,6),w_star(i,7)]...
        = welfarefunction (y_tilde_star, dlogz(:,i));

end

%% cosine_similarity to optimal policy
% How closely the inflation index of each simple rule is aligned with the
% optimal one, reported in the last row of the table.

CS_og             = dot(phi_star, phi_og) / (norm(phi_star) * norm(phi_og));
CS_cpi            = dot(phi_star, phi_cpi) / (norm(phi_star) * norm(phi_cpi));


%% Calculating the welfare gain for each policy

EWG_cpi  = mean(w_cpi);
EWG_og   = mean(w_og);
EWG_star = mean(w_star);

%% Expected welfare gain measured as a fraction of steady-state consumption
% Converts the expected welfare terms into consumption equivalents, in
% percent, applying the transformation componentwise.

EWG_GDP_cpi   = expectedwelfaregain(EWG_cpi,gamma);
EWG_GDP_og    = expectedwelfaregain(EWG_og,gamma);
EWG_GDP_star  = expectedwelfaregain(EWG_star,gamma);

input.data    = [EWG_GDP_star; EWG_GDP_og; EWG_GDP_cpi];
input.data    = input.data';
input.data    = [input.data; 1, CS_og, CS_cpi];



%% Generates a LaTeX table
% Columns are the three policies, rows are the seven welfare components in
% the order returned by welfarefunction.m followed by the cosine similarity
% row appended to input.data above.

% Optional fields (if not set default values will be used):

% Set the position of the table in the LaTex document using h, t, p, b, H or !
input.tablePositioning = 'h';

% Set column labels (use empty string for no label):
input.tableColLabels = {'Optimal Policy', 'Output-Gap Targeting', 'CPI Targeting'};
% Set row labels (use empty string for no label):
input.tableRowLabels = {'Welfare Gain (Percent GDP)','\hline \quad First order welfare gain\\ \quad \quad Aggregate Wedge',...
    '\quad \quad Supply-side effect','\hline \quad Second order welfare gain\\ \quad \quad Output gap volatility',...
    '\quad \quad Within-sector price dispersion','\quad \quad Cross-sector price dispersion',...
    '\quad \quad Variation in allocative efficiency','\hline Cosine similarity to optimal policy'};

% Switch transposing/pivoting your table:
input.transposeTable = 0;
input.dataFormat = {'%.3f',3};
input.dataNanString = '-';
input.tableColumnAlignment = 'c';
input.tableBorders = 1;
input.booktabs = 1;

% LaTex table:
input.tableCaption = 'Welfare Gain under Various Policies';
input.tableLabel = 'welfare-gain-decomp';
input.makeCompleteLatexDocument = 1;
latex = latexTable(input);





















