Skip to content
GitLab
Commit 9db12658 authored by Normann Rion's avatar Normann Rion Committed by NormannR
Browse files

Asssesing welfare in perfect-foresight frameworks 

Ref. #1680
parent 9066d31d
Hide whitespace changes
Inline Side-by-side
matlab/evaluate_planner_objective.m
View file @ 9db12658
......@@ -49,6 +49,8 @@ function planner_objective_value = evaluate_planner_objective(M_,options_,oo_)
% E(W) = (1 - beta)^{-1} ( Ubar + 0.5 ( U_xx h_y^2 E(yhat^2) + U_xx h_u^2 E(u^2) )
% As for the conditional welfare, the second-order formula above is still valid, but the derivatives of W no longer contain any second-order derivatives of the policy and transition functions h and g.
% In the deterministic case, resorting to approximations for welfare is no longer required as it is possible to simulate the model given initial conditions for pre-determined variables and terminal conditions for forward-looking variables, whether these initial and terminal conditions are explicitly or implicitly specified. Assuming that the number of simulated periods is high enough for the new steady-state to be reached, the new unconditional welfare is thus the last period's welfare. As for the conditional welfare, it can be derived using backward recursions on the equation W = U + beta*W(+1) starting from the final unconditional steady-state welfare.
% INPUTS
% M_: (structure) model description
% options_: (structure) options
......@@ -81,87 +83,101 @@ nstatic = M_.nstatic;
nspred = M_.nspred;
beta = get_optimal_policy_discount_factor(M_.params, M_.param_names);
ys = oo_.dr.ys;
planner_objective_value = zeros(2,1);
if options_.ramsey_policy
if options_.order == 1
[U] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
planner_objective_value(1) = U/(1-beta);
planner_objective_value(2) = U/(1-beta);
elseif options_.order == 2
[U,Uy,Uyy] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
Gy = dr.ghx(nstatic+(1:nspred),:);
Gu = dr.ghu(nstatic+(1:nspred),:);
Gyy = dr.ghxx(nstatic+(1:nspred),:);
Gyu = dr.ghxu(nstatic+(1:nspred),:);
Guu = dr.ghuu(nstatic+(1:nspred),:);
Gss = dr.ghs2(nstatic+(1:nspred),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
gyy(dr.order_var,:) = dr.ghxx;
gyu(dr.order_var,:) = dr.ghxu;
guu(dr.order_var,:) = dr.ghuu;
gss(dr.order_var,:) = dr.ghs2;
Uyy = full(Uyy);
Uyygygy = A_times_B_kronecker_C(Uyy,gy,gy);
Uyygugu = A_times_B_kronecker_C(Uyy,gu,gu);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
[info, oo_, options_] = stoch_simul(M_, options_, oo_, var_list); %get decision rules and moments
if ~old_noprint
options_.noprint = 0;
if oo_.gui.ran_perfect_foresight
T = size(oo_.endo_simul,2);
[U_term] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,T),oo_.exo_simul(T,:), M_.params);
EW = U_term/(1-beta);
W = EW;
for t=T:-1:2
[U] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,t),oo_.exo_simul(t,:), M_.params);
W = U + beta*W;
end
oo_.mean(isnan(oo_.mean)) = options_.huge_number;
oo_.var(isnan(oo_.var)) = options_.huge_number;
Ey = oo_.mean;
Eyhat = Ey - ys(dr.order_var(nstatic+(1:nspred)));
var_corr = Eyhat*Eyhat';
Eyhatyhat = oo_.var(:) + var_corr(:);
Euu = M_.Sigma_e(:);
EU = U + Uy*gy*Eyhat + 0.5*((Uyygygy + Uy*gyy)*Eyhatyhat + (Uyygugu + Uy*guu)*Euu + Uy*gss);
EW = EU/(1-beta);
%% Conditional welfare starting from the non-stochastic steady-state
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
if isempty(options_.qz_criterium)
options_.qz_criterium = 1+1e-6;
end
%solve Lyapunuv equation Wyy=gy'*Uyy*gy+Uy*gyy+beta*Wy*Gyy+beta*Gy'Wyy*Gy
Wyy = reshape(lyapunov_symm(sqrt(beta)*Gy',reshape(Uyygygy + Uy*gyy + beta*Wy*Gyy,nspred,nspred),options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold, 3, options_.debug),1,nspred*nspred);
Wyygugu = A_times_B_kronecker_C(Wyy,Gu,Gu);
Wuu = Uyygugu + Uy*guu + beta*(Wyygugu + Wy*Guu);
Wss = (Uy*gss + beta*(Wy*Gss + Wuu*M_.Sigma_e(:)))/(1-beta);
W = Wbar + 0.5*Wss;
planner_objective_value(1) = EW;
planner_objective_value(2) = W;
else
%Order k code will go here!
fprintf('\nevaluate_planner_objective: order>2 not yet supported\n')
planner_objective_value(1) = NaN;
planner_objective_value(2) = NaN;
return
ys = oo_.dr.ys;
if options_.order == 1
[U] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
planner_objective_value(1) = U/(1-beta);
planner_objective_value(2) = U/(1-beta);
elseif options_.order == 2
[U,Uy,Uyy] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
Gy = dr.ghx(nstatic+(1:nspred),:);
Gu = dr.ghu(nstatic+(1:nspred),:);
Gyy = dr.ghxx(nstatic+(1:nspred),:);
Gyu = dr.ghxu(nstatic+(1:nspred),:);
Guu = dr.ghuu(nstatic+(1:nspred),:);
Gss = dr.ghs2(nstatic+(1:nspred),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
gyy(dr.order_var,:) = dr.ghxx;
gyu(dr.order_var,:) = dr.ghxu;
guu(dr.order_var,:) = dr.ghuu;
gss(dr.order_var,:) = dr.ghs2;
Uyy = full(Uyy);
Uyygygy = A_times_B_kronecker_C(Uyy,gy,gy);
Uyygugu = A_times_B_kronecker_C(Uyy,gu,gu);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
[info, oo_, options_] = stoch_simul(M_, options_, oo_, var_list); %get decision rules and moments
if ~old_noprint
options_.noprint = 0;
end
oo_.mean(isnan(oo_.mean)) = options_.huge_number;
oo_.var(isnan(oo_.var)) = options_.huge_number;
Ey = oo_.mean;
Eyhat = Ey - ys(dr.order_var(nstatic+(1:nspred)));
var_corr = Eyhat*Eyhat';
Eyhatyhat = oo_.var(:) + var_corr(:);
Euu = M_.Sigma_e(:);
EU = U + Uy*gy*Eyhat + 0.5*((Uyygygy + Uy*gyy)*Eyhatyhat + (Uyygugu + Uy*guu)*Euu + Uy*gss);
EW = EU/(1-beta);
%% Conditional welfare starting from the non-stochastic steady-state
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
if isempty(options_.qz_criterium)
options_.qz_criterium = 1+1e-6;
end
%solve Lyapunuv equation Wyy=gy'*Uyy*gy+Uy*gyy+beta*Wy*Gyy+beta*Gy'Wyy*Gy
Wyy = reshape(lyapunov_symm(sqrt(beta)*Gy',reshape(Uyygygy + Uy*gyy + beta*Wy*Gyy,nspred,nspred),options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold, 3, options_.debug),1,nspred*nspred);
Wyygugu = A_times_B_kronecker_C(Wyy,Gu,Gu);
Wuu = Uyygugu + Uy*guu + beta*(Wyygugu + Wy*Guu);
Wss = (Uy*gss + beta*(Wy*Gss + Wuu*M_.Sigma_e(:)))/(1-beta);
W = Wbar + 0.5*Wss;
planner_objective_value(1) = EW;
planner_objective_value(2) = W;
else
%Order k code will go here!
fprintf('\nevaluate_planner_objective: order>2 not yet supported\n')
planner_objective_value(1) = NaN;
planner_objective_value(2) = NaN;
return
end
end
elseif options_.discretionary_policy
ys = oo_.dr.ys;
[U,Uy,Uyy] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
Gy = dr.ghx(nstatic+(1:nspred),:);
......@@ -221,11 +237,15 @@ elseif options_.discretionary_policy
end
if ~options_.noprint
if options_.ramsey_policy
fprintf('\nApproximated value of unconditional welfare: %10.8f\n', planner_objective_value(1))
fprintf('\nApproximated value of conditional welfare: %10.8f\n', planner_objective_value(2))
if oo_.gui.ran_perfect_foresight
fprintf('\nSimulated value of unconditional welfare: %10.8f\n', planner_objective_value(1))
fprintf('\nSimulated value of conditional welfare: %10.8f\n', planner_objective_value(2))
else
fprintf('\nApproximated value of unconditional welfare: %10.8f\n', planner_objective_value(1))
fprintf('\nApproximated value of conditional welfare: %10.8f\n', planner_objective_value(2))
end
elseif options_.discretionary_policy
fprintf('\nApproximated value of unconditional welfare with discretionary policy: %10.8f\n', planner_objective_value(1))
fprintf('\nApproximated value of conditional welfare with discretionary policy: %10.8f\n', planner_objective_value(2))
end
end
tests/Makefile.am
View file @ 9db12658
......@@ -580,6 +580,9 @@ MFILES = histval_initval_file/ramst_initval_file_data.m
optimal_policy/neo_growth_ramsey.m.trs: optimal_policy/neo_growth.m.trs
optimal_policy/neo_growth_ramsey.o.trs: optimal_policy/neo_growth.o.trs
optimal_policy/neo_growth_ramsey_foresight.m.trs: optimal_policy/neo_growth_foresight.m.trs
optimal_policy/neo_growth_ramsey_foresight.o.trs: optimal_policy/neo_growth_foresight.o.trs
example1_use_dll.m.trs: example1.m.trs
example1_use_dll.o.trs: example1.o.trs
......@@ -1185,6 +1188,8 @@ EXTRA_DIST = \
observation_trends_and_prefiltering/Trend_model_calib_no_prefilter_common.inc \
observation_trends_and_prefiltering/Trend_load_data_common.inc \
observation_trends_and_prefiltering/Trend_no_prefilter_conditional_forecast.inc \
optimal_policy/neo_growth_common.inc \
optimal_policy/neo_growth_ramsey_common.inc \
optimal_policy/Ramsey/oo_ramsey_policy_initval.mat \
optimizers/optimizer_function_wrapper.m \
optimizers/fs2000.common.inc \
......
tests/optimal_policy/neo_growth.mod
View file @ 9db12658
......@@ -21,34 +21,7 @@
* It is called by neo_growth_ramsey.mod to compare by-hand calculations of unconditional
* and conditional welfares and the output of the evaluate_planner_objective function.
*/
var U k z c W;
varexo e;
parameters beta gamma alpha delta rho s;
beta = 0.987;
gamma = 1;
delta = 0.012;
alpha = 0.4;
rho = 0.95;
s = 0.007;
model;
c^(-gamma)=beta*c(+1)^(-gamma)*(alpha*exp(z(+1))*k^(alpha-1)+1-delta);
W = U + beta*W(+1);
k=exp(z)*k(-1)^(alpha)-c+(1-delta)*k(-1);
z=rho*z(-1)+s*e;
U=ln(c);
end;
steady_state_model;
k = ((1/beta-(1-delta))/alpha)^(1/(alpha-1));
c = k^alpha-delta*k;
z = 0;
U = ln(c);
W = U/(1-beta);
end;
@#include "neo_growth_common.inc"
shocks;
var e;
......@@ -57,4 +30,5 @@ end;
steady;
resid;
stoch_simul(order=2, irf=0);
tests/optimal_policy/neo_growth_common.inc 0 → 100644
View file @ 9db12658
var U k z c W;
varexo e;
parameters beta gamma alpha delta rho s;
beta = 0.987;
gamma = 1;
delta = 0.012;
alpha = 0.4;
rho = 0.95;
s = 0.007;
model;
c^(-gamma)=beta*c(+1)^(-gamma)*(alpha*exp(z(+1))*k^(alpha-1)+1-delta);
W = U + beta*W(+1);
k=exp(z)*k(-1)^(alpha)-c+(1-delta)*k(-1);
z=rho*z(-1)+s*e;
U=ln(c);
end;
steady_state_model;
k = ((1/beta-(1-delta))/alpha)^(1/(alpha-1));
c = k^alpha-delta*k;
z = 0;
U = ln(c);
W = U/(1-beta);
end;
tests/optimal_policy/neo_growth_foresight.mod 0 → 100644
View file @ 9db12658
/*
* Copyright (C) 2021 Dynare Team
*
* This file is part of Dynare.
*
* Dynare is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* Dynare is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with Dynare. If not, see <http://www.gnu.org/licenses/>.
*/
/*
* This file simulates the neo-classical growth model in a perfect foresight framework.
* It is called by neo_growth_ramsey_foresight.mod to compare by-hand calculations of unconditional
* and conditional welfares and the output of the evaluate_planner_objective function.
*/
@#include "neo_growth_common.inc"
initval;
z = 0;
end;
steady;
shocks;
var e;
periods 1;
values 1;
end;
resid;
perfect_foresight_setup(periods=200);
perfect_foresight_solver;
\ No newline at end of file
tests/optimal_policy/neo_growth_ramsey.mod
View file @ 9db12658
......@@ -22,43 +22,13 @@
* and compares them to a by-hand assessment stemming from the results the model neo_growth.mod incur.
*/
var k z c;
varexo e;
parameters beta gamma alpha delta rho s;
beta = 0.987;
gamma = 1;
delta = 0.012;
alpha = 0.4;
rho = 0.95;
s = 0.007;
model;
k=exp(z)*k(-1)^(alpha)-c+(1-delta)*k(-1);
z=rho*z(-1)+s*e;
end;
steady_state_model;
z = 0;
end;
@#include "neo_growth_ramsey_common.inc"
shocks;
var e;
stderr 1;
end;
planner_objective ln(c);
ramsey_model(instruments=(k,c), planner_discount=beta);
initval;
k = ((1/beta-(1-delta))/alpha)^(1/(alpha-1));
c = k^alpha-delta*k;
end;
steady;
resid;
stoch_simul(order=2, irf=0);
planner_objective_value = evaluate_planner_objective(M_, options_, oo_);
......
tests/optimal_policy/neo_growth_ramsey_common.inc 0 → 100644
View file @ 9db12658
var k z c;
varexo e;
parameters beta gamma alpha delta rho s;
beta = 0.987;
gamma = 1;
delta = 0.012;
alpha = 0.4;
rho = 0.95;
s = 0.007;
model;
k=exp(z)*k(-1)^(alpha)-c+(1-delta)*k(-1);
z=rho*z(-1)+s*e;
end;
steady_state_model;
z = 0;
end;
planner_objective ln(c);
ramsey_model(instruments=(k,c), planner_discount=beta);
initval;
k = ((1/beta-(1-delta))/alpha)^(1/(alpha-1));
c = k^alpha-delta*k;
end;
steady;
resid;
tests/optimal_policy/neo_growth_ramsey_foresight.mod 0 → 100644
View file @ 9db12658
/*
* Copyright (C) 2021 Dynare Team
*
* This file is part of Dynare.
*
* Dynare is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* Dynare is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with Dynare. If not, see <http://www.gnu.org/licenses/>.
*/
/*
* This file simulates a perfect-foresight version of the neo-classical growth model.
* It assesses the conditional and unconditional welfares computed by the evaluate_planner_objective function
* and compares them to a by-hand assessment stemming from the results of the model neo_growth_foresight.mod
*/
@#include "neo_growth_ramsey_common.inc"
shocks;
var e;
periods 1;
values 1;
end;
perfect_foresight_setup(periods=200);
perfect_foresight_solver;
planner_objective_value = evaluate_planner_objective(M_, options_, oo_);
if ~exist('neo_growth_foresight_results.mat','file');
error('neo_growth_foresight must be run first');
end;
oo1 = load('neo_growth_foresight_results','oo_');
M1 = load('neo_growth_foresight_results','M_');
options1 = load('neo_growth_foresight_results','options_');
cond_W_hand = oo1.oo_.endo_simul(strmatch('W',M1.M_.endo_names,'exact'),2);
unc_W_hand = oo1.oo_.endo_simul(strmatch('W',M1.M_.endo_names,'exact'),end);
if abs((unc_W_hand - planner_objective_value(1))/unc_W_hand) > 1e-6;
error('Inaccurate unconditional welfare assessment');
end;
if abs((cond_W_hand - planner_objective_value(2))/cond_W_hand) > 1e-6;
error('Inaccurate conditional welfare assessment');
end;
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment