Merge branch 'discretion' into 'master'

evaluate_planner_objective.m: fix output for linear-quadratic problems solved at second order

See merge request Dynare/dynare!2056
mr#2067
Sébastien Villemot 2022-07-26 13:03:45 +00:00
commit b93d04646e
1 changed files with 264 additions and 268 deletions

View File

@ -59,7 +59,7 @@ function planner_objective_value = evaluate_planner_objective(M_,options_,oo_)
% 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.
% Copyright © 2007-2021 Dynare Team
% Copyright © 2007-2022 Dynare Team
%
% This file is part of Dynare.
%
@ -90,274 +90,271 @@ if beta>=1
fprintf('evaluate_planner_objective: the planner discount factor is not strictly smaller than 1. Unconditional welfare will not be finite.\n')
end
if options_.ramsey_policy
if oo_.gui.ran_perfect_foresight
T = size(oo_.endo_simul,2);
[U_term] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,T-M_.maximum_lead),oo_.exo_simul(T-M_.maximum_lead,:), M_.params);
EW = U_term/(1-beta);
W = EW;
for t=T-M_.maximum_lead:-1:1+M_.maximum_lag
[U] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,t),oo_.exo_simul(t,:), M_.params);
W = U + beta*W;
end
planner_objective_value = struct('conditional', W, 'unconditional', EW);
else
planner_objective_value = struct('conditional', struct('zero_initial_multiplier', 0., 'steady_initial_multiplier', 0.), 'unconditional', 0.);
if isempty(oo_.dr) || ~isfield(oo_.dr,'ys')
error('evaluate_planner_objective requires decision rules to have previously been computed (e.g. by stoch_simul)')
else
ys = oo_.dr.ys;
end
if options_.order == 1 || M_.hessian_eq_zero
[U,Uy] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
Gy = dr.ghx(nstatic+(1:nspred),:);
Gu = dr.ghu(nstatic+(1:nspred),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
%% Unconditional welfare
EW = U/(1-beta);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state (i) with Lagrange multipliers set to their steady-state value (ii) with Lagrange multipliers set to 0
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u;
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u;
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
elseif options_.order == 2 && ~M_.hessian_eq_zero
[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);
Uyygugy = A_times_B_kronecker_C(Uyy,gu,gy);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
if options_.pruning
fprintf('evaluate_planner_objective: pruning option is not supported and will be ignored\n')
end
oo_=disp_th_moments(dr,var_list,M_,options_,oo_);
if ~old_noprint
options_.noprint = 0;
end
if any(isnan(oo_.mean)) || any(any(isnan(oo_.var)))
fprintf('evaluate_planner_objective: encountered NaN moments in the endogenous variables often associated\n')
fprintf('evaluate_planner_objective: with either non-stationary variables or singularity due e.g. including\n')
fprintf('evaluate_planner_objective: the planner objective function (or additive parts of it) in the model.\n')
fprintf('evaluate_planner_objective: I will replace the NaN with a large number, but tread carefully,\n')
fprintf('evaluate_planner_objective: check your model, and watch out for strange results.\n')
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)));
Eyhatyhat = oo_.var(:);
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);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state (i) with Lagrange multipliers set to their steady-state value (ii) with Lagrange multipliers set to 0
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
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);
Wyygugy = A_times_B_kronecker_C(Wyy,Gu,Gy);
Wuu = Uyygugu + Uy*guu + beta*(Wyygugu + Wy*Guu);
Wss = (Uy*gss + beta*(Wy*Gss + Wuu*M_.Sigma_e(:)))/(1-beta);
Wyu = Uyygugy + Uy*gyu + beta*(Wyygugy + Wy*Gyu);
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
Wyu_yu_L_SS = A_times_B_kronecker_C(Wyu,yhat_L_SS,u);
Wyy_yy_L_SS = A_times_B_kronecker_C(Wyy,yhat_L_SS,yhat_L_SS);
Wuu_uu_L_SS = A_times_B_kronecker_C(Wuu,u,u);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u+Wyu_yu_L_SS+0.5*(Wss+Wyy_yy_L_SS+Wuu_uu_L_SS);
Wyu_yu_L_0 = A_times_B_kronecker_C(Wyu,yhat_L_0,u);
Wyy_yy_L_0 = A_times_B_kronecker_C(Wyy,yhat_L_0,yhat_L_0);
Wuu_uu_L_0 = A_times_B_kronecker_C(Wuu,u,u);
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u+Wyu_yu_L_0+0.5*(Wss+Wyy_yy_L_0+Wuu_uu_L_0);
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
else
% Computes the welfare decision rule
[W] = k_order_welfare(dr,M_,options_);
% Appends the welfare decision rule to the endogenous variables decision
% rule
for i=0:options_.order
dr.(['g_' num2str(i)]) = [dr.(['g_' num2str(i)]); W.(['W_' num2str(i)])];
end
% Amends the steady-state vector accordingly
[U] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
ysteady = [ys(oo_.dr.order_var); U/(1-beta)];
% Generates the sequence of shocks to compute unconditional welfare
i_exo_var = setdiff([1:M_.exo_nbr],find(diag(M_.Sigma_e) == 0));
nxs = length(i_exo_var);
chol_S = chol(M_.Sigma_e(i_exo_var,i_exo_var));
exo_simul = zeros(M_.exo_nbr,options_.ramsey.periods);
if ~isempty(M_.Sigma_e)
exo_simul(i_exo_var,:) = chol_S*randn(nxs,options_.ramsey.periods);
end
yhat_start = zeros(M_.endo_nbr+1,1);
[moment] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, options_.ramsey.drop, yhat_start, exo_simul, ysteady, dr);
% Stores the result for unconditional welfare
planner_objective_value.unconditional = moment(end);
% Conditional welfare
% Gets initial values
[yhat_L_SS,yhat_L_0, u] = get_initial_state(ys,M_,dr,oo_);
% Conditional welfare (i) with Lagrange multipliers set to their
% steady-state values
yhat_start(M_.nstatic+1:M_.nstatic+M_.npred+M_.nboth) = yhat_L_SS;
[~,sim] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, 0, yhat_start, u, ysteady, dr);
planner_objective_value.conditional.steady_initial_multiplier = sim(end,1);
% Conditional welfare (ii) with Lagrange multipliers set to 0
yhat_start(M_.nstatic+1:M_.nstatic+M_.npred+M_.nboth) = yhat_L_0;
[~,sim] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, 0, yhat_start, u, ysteady, dr);
planner_objective_value.conditional.zero_initial_multiplier = sim(end,1);
end
if options_.ramsey_policy && oo_.gui.ran_perfect_foresight
T = size(oo_.endo_simul,2);
[U_term] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,T-M_.maximum_lead),oo_.exo_simul(T-M_.maximum_lead,:), M_.params);
EW = U_term/(1-beta);
W = EW;
for t=T-M_.maximum_lead:-1:1+M_.maximum_lag
[U] = feval([M_.fname '.objective.static'],oo_.endo_simul(:,t),oo_.exo_simul(t,:), M_.params);
W = U + beta*W;
end
elseif options_.discretionary_policy
ys = oo_.dr.ys;
planner_objective_value = struct('conditional', W, 'unconditional', EW);
else
planner_objective_value = struct('conditional', struct('zero_initial_multiplier', 0., 'steady_initial_multiplier', 0.), 'unconditional', 0.);
[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),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
Uyy = full(Uyy);
Uyygygy = A_times_B_kronecker_C(Uyy,gy,gy);
Uyygugu = A_times_B_kronecker_C(Uyy,gu,gu);
Uyygugy = A_times_B_kronecker_C(Uyy,gu,gy);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
if isempty(oo_.dr) || ~isfield(oo_.dr,'ys')
error('evaluate_planner_objective requires decision rules to have previously been computed (e.g. by stoch_simul or discretionary_policy)')
else
ys = oo_.dr.ys;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
oo_=disp_th_moments(dr,var_list,M_,options_,oo_);
if ~old_noprint
options_.noprint = 0;
if options_.order == 1 && ~options_.discretionary_policy
[U,Uy] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
Gy = dr.ghx(nstatic+(1:nspred),:);
Gu = dr.ghu(nstatic+(1:nspred),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
%% Unconditional welfare
EW = U/(1-beta);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state (i) with Lagrange multipliers set to their steady-state value (ii) with Lagrange multipliers set to 0
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u;
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u;
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
elseif options_.order == 2 && ~M_.hessian_eq_zero %full second order approximation
[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);
Uyygugy = A_times_B_kronecker_C(Uyy,gu,gy);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
if options_.pruning
fprintf('evaluate_planner_objective: pruning option is not supported and will be ignored\n')
end
oo_=disp_th_moments(dr,var_list,M_,options_,oo_);
if ~old_noprint
options_.noprint = 0;
end
if any(isnan(oo_.mean)) || any(any(isnan(oo_.var)))
fprintf('evaluate_planner_objective: encountered NaN moments in the endogenous variables often associated\n')
fprintf('evaluate_planner_objective: with either non-stationary variables or singularity due e.g. including\n')
fprintf('evaluate_planner_objective: the planner objective function (or additive parts of it) in the model.\n')
fprintf('evaluate_planner_objective: I will replace the NaN with a large number, but tread carefully,\n')
fprintf('evaluate_planner_objective: check your model, and watch out for strange results.\n')
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)));
Eyhatyhat = oo_.var(:);
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);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state (i) with Lagrange multipliers set to their steady-state value (ii) with Lagrange multipliers set to 0
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
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);
Wyygugy = A_times_B_kronecker_C(Wyy,Gu,Gy);
Wuu = Uyygugu + Uy*guu + beta*(Wyygugu + Wy*Guu);
Wss = (Uy*gss + beta*(Wy*Gss + Wuu*M_.Sigma_e(:)))/(1-beta);
Wyu = Uyygugy + Uy*gyu + beta*(Wyygugy + Wy*Gyu);
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
Wyu_yu_L_SS = A_times_B_kronecker_C(Wyu,yhat_L_SS,u);
Wyy_yy_L_SS = A_times_B_kronecker_C(Wyy,yhat_L_SS,yhat_L_SS);
Wuu_uu_L_SS = A_times_B_kronecker_C(Wuu,u,u);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u+Wyu_yu_L_SS+0.5*(Wss+Wyy_yy_L_SS+Wuu_uu_L_SS);
Wyu_yu_L_0 = A_times_B_kronecker_C(Wyu,yhat_L_0,u);
Wyy_yy_L_0 = A_times_B_kronecker_C(Wyy,yhat_L_0,yhat_L_0);
Wuu_uu_L_0 = A_times_B_kronecker_C(Wuu,u,u);
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u+Wyu_yu_L_0+0.5*(Wss+Wyy_yy_L_0+Wuu_uu_L_0);
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
elseif (options_.order == 2 && M_.hessian_eq_zero) || options_.discretionary_policy %linear quadratic problem
[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),:);
gy(dr.order_var,:) = dr.ghx;
gu(dr.order_var,:) = dr.ghu;
Uyy = full(Uyy);
Uyygygy = A_times_B_kronecker_C(Uyy,gy,gy);
Uyygugu = A_times_B_kronecker_C(Uyy,gu,gu);
Uyygugy = A_times_B_kronecker_C(Uyy,gu,gy);
%% Unconditional welfare
old_noprint = options_.noprint;
if ~old_noprint
options_.noprint = 1;
end
var_list = M_.endo_names(dr.order_var(nstatic+(1:nspred)));
oo_=disp_th_moments(dr,var_list,M_,options_,oo_);
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)));
Eyhatyhat = oo_.var(:);
Euu = M_.Sigma_e(:);
EU = U + Uy*gy*Eyhat + 0.5*(Uyygygy*Eyhatyhat + Uyygugu*Euu);
EW = EU/(1-beta);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
%solve Lyapunuv equation Wyy=gy'*Uyy*gy+beta*Gy'Wyy*Gy
Wyy = reshape(lyapunov_symm(sqrt(beta)*Gy',reshape(Uyygygy,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);
Wyygugy = A_times_B_kronecker_C(Wyy,Gu,Gy);
Wuu = Uyygugu + beta*Wyygugu;
Wss = beta*Wuu*M_.Sigma_e(:)/(1-beta);
Wyu = Uyygugy + beta*Wyygugy;
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
Wyu_yu_L_SS = A_times_B_kronecker_C(Wyu,yhat_L_SS,u);
Wyy_yy_L_SS = A_times_B_kronecker_C(Wyy,yhat_L_SS,yhat_L_SS);
Wuu_uu_L_SS = A_times_B_kronecker_C(Wuu,u,u);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u+Wyu_yu_L_SS+0.5*(Wss+Wyy_yy_L_SS+Wuu_uu_L_SS);
Wyu_yu_L_0 = A_times_B_kronecker_C(Wyu,yhat_L_0,u);
Wyy_yy_L_0 = A_times_B_kronecker_C(Wyy,yhat_L_0,yhat_L_0);
Wuu_uu_L_0 = A_times_B_kronecker_C(Wuu,u,u);
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u+Wyu_yu_L_0+0.5*(Wss+Wyy_yy_L_0+Wuu_uu_L_0);
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
elseif options_.order > 2 || ~options_.discretionary_policy
% Computes the welfare decision rule
[W] = k_order_welfare(dr,M_,options_);
% Appends the welfare decision rule to the endogenous variables decision
% rule
for i=0:options_.order
dr.(['g_' num2str(i)]) = [dr.(['g_' num2str(i)]); W.(['W_' num2str(i)])];
end
% Amends the steady-state vector accordingly
[U] = feval([M_.fname '.objective.static'],ys,zeros(1,exo_nbr), M_.params);
ysteady = [ys(oo_.dr.order_var); U/(1-beta)];
% Generates the sequence of shocks to compute unconditional welfare
i_exo_var = setdiff([1:M_.exo_nbr],find(diag(M_.Sigma_e) == 0));
nxs = length(i_exo_var);
chol_S = chol(M_.Sigma_e(i_exo_var,i_exo_var));
exo_simul = zeros(M_.exo_nbr,options_.ramsey.periods);
if ~isempty(M_.Sigma_e)
exo_simul(i_exo_var,:) = chol_S*randn(nxs,options_.ramsey.periods);
end
yhat_start = zeros(M_.endo_nbr+1,1);
[moment] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, options_.ramsey.drop, yhat_start, exo_simul, ysteady, dr);
% Stores the result for unconditional welfare
planner_objective_value.unconditional = moment(end);
% Conditional welfare
% Gets initial values
[yhat_L_SS,yhat_L_0, u] = get_initial_state(ys,M_,dr,oo_);
% Conditional welfare (i) with Lagrange multipliers set to their
% steady-state values
yhat_start(M_.nstatic+1:M_.nstatic+M_.npred+M_.nboth) = yhat_L_SS;
[~,sim] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, 0, yhat_start, u, ysteady, dr);
planner_objective_value.conditional.steady_initial_multiplier = sim(end,1);
% Conditional welfare (ii) with Lagrange multipliers set to 0
yhat_start(M_.nstatic+1:M_.nstatic+M_.npred+M_.nboth) = yhat_L_0;
[~,sim] = k_order_mean(options_.order, M_.nstatic, M_.npred, M_.nboth, M_.nfwrd+1, M_.exo_nbr, 1, 0, yhat_start, u, ysteady, dr);
planner_objective_value.conditional.zero_initial_multiplier = sim(end,1);
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)));
Eyhatyhat = oo_.var(:);
Euu = M_.Sigma_e(:);
EU = U + Uy*gy*Eyhat + 0.5*(Uyygygy*Eyhatyhat + Uyygugu*Euu);
EW = EU/(1-beta);
planner_objective_value.unconditional = EW;
%% Conditional welfare starting from the non-stochastic steady-state
Wbar = U/(1-beta);
Wy = Uy*gy/(eye(nspred)-beta*Gy);
Wu = Uy*gu + beta*Wy*Gu;
%solve Lyapunuv equation Wyy=gy'*Uyy*gy+beta*Gy'Wyy*Gy
Wyy = reshape(lyapunov_symm(sqrt(beta)*Gy',reshape(Uyygygy,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);
Wyygugy = A_times_B_kronecker_C(Wyy,Gu,Gy);
Wuu = Uyygugu + beta*Wyygugu;
Wss = beta*Wuu*M_.Sigma_e(:)/(1-beta);
Wyu = Uyygugy + beta*Wyygugy;
[yhat_L_SS,yhat_L_0, u]=get_initial_state(ys,M_,dr,oo_);
Wyu_yu_L_SS = A_times_B_kronecker_C(Wyu,yhat_L_SS,u);
Wyy_yy_L_SS = A_times_B_kronecker_C(Wyy,yhat_L_SS,yhat_L_SS);
Wuu_uu_L_SS = A_times_B_kronecker_C(Wuu,u,u);
W_L_SS = Wbar+Wy*yhat_L_SS+Wu*u+Wyu_yu_L_SS+0.5*(Wss+Wyy_yy_L_SS+Wuu_uu_L_SS);
Wyu_yu_L_0 = A_times_B_kronecker_C(Wyu,yhat_L_0,u);
Wyy_yy_L_0 = A_times_B_kronecker_C(Wyy,yhat_L_0,yhat_L_0);
Wuu_uu_L_0 = A_times_B_kronecker_C(Wuu,u,u);
W_L_0 = Wbar+Wy*yhat_L_0+Wu*u+Wyu_yu_L_0+0.5*(Wss+Wyy_yy_L_0+Wuu_uu_L_0);
planner_objective_value.conditional.steady_initial_multiplier = W_L_SS;
planner_objective_value.conditional.zero_initial_multiplier = W_L_0;
end
if ~options_.noprint
if options_.ramsey_policy
if oo_.gui.ran_perfect_foresight
fprintf('\nSimulated value of unconditional welfare: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nSimulated value of conditional welfare: %10.8f\n', planner_objective_value.conditional)
fprintf('\nSimulated value of unconditional welfare: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nSimulated value of conditional welfare: %10.8f\n', planner_objective_value.conditional)
else
fprintf('\nApproximated value of unconditional welfare: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nApproximated value of conditional welfare:\n')
fprintf(' - with initial Lagrange multipliers set to 0: %10.8f\n', planner_objective_value.conditional.zero_initial_multiplier)
fprintf(' - with initial Lagrange multipliers set to steady state: %10.8f\n\n', planner_objective_value.conditional.steady_initial_multiplier)
fprintf('\nApproximated value of unconditional welfare: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nApproximated value of conditional welfare:\n')
fprintf(' - with initial Lagrange multipliers set to 0: %10.8f\n', planner_objective_value.conditional.zero_initial_multiplier)
fprintf(' - with initial Lagrange multipliers set to steady state: %10.8f\n\n', planner_objective_value.conditional.steady_initial_multiplier)
end
elseif options_.discretionary_policy
fprintf('\nApproximated value of unconditional welfare with discretionary policy: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nApproximated value of conditional welfare with discretionary policy:\n')
fprintf(' - with initial Lagrange multipliers set to 0: %10.8f\n', planner_objective_value.conditional.zero_initial_multiplier)
fprintf(' - with initial Lagrange multipliers set to steady state: %10.8f\n\n', planner_objective_value.conditional.steady_initial_multiplier)
fprintf('\nApproximated value of unconditional welfare with discretionary policy: %10.8f\n', planner_objective_value.unconditional)
fprintf('\nApproximated value of conditional welfare with discretionary policy:\n')
fprintf(' - with initial Lagrange multipliers set to 0: %10.8f\n', planner_objective_value.conditional.zero_initial_multiplier)
fprintf(' - with initial Lagrange multipliers set to steady state: %10.8f\n\n', planner_objective_value.conditional.steady_initial_multiplier)
end
end
@ -368,8 +365,8 @@ yhat_L_SS = ys;
% initialize Lagrange multipliers to 0 in yhat_L_0
yhat_L_0 = zeros(M_.endo_nbr,1);
if ~isempty(M_.aux_vars)
mult_indicator=([M_.aux_vars(:).type]==6);
mult_indices=[M_.aux_vars(mult_indicator).endo_index];
mult_indicator=([M_.aux_vars(:).type]==6);
mult_indices=[M_.aux_vars(mult_indicator).endo_index];
else
mult_indices=[];
end
@ -396,7 +393,7 @@ if ~isempty(M_.det_shocks)
end
if any(periods>1)
fprintf(['\nevaluate_planner_objective: Shock values for periods not contained in the initial information set (t=1) have been provided.\n' ...
'evaluate_planner_objective: Note that they will be ignored.\n'])
'evaluate_planner_objective: Note that they will be ignored.\n'])
end
shock_indices=find(periods==1);
if any([M_.det_shocks(shock_indices).multiplicative])
@ -406,4 +403,3 @@ if ~isempty(M_.det_shocks)
else
u = oo_.exo_simul(1,:)'; %first value of simulation series (set by simult.m if periods>0), 1 otherwise
end