From 378f405395440176dd85056d9779c31b0496ea29 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Wed, 19 Feb 2025 16:27:51 +0100
Subject: [PATCH] WIP: perfect foresight conditional forecast
[skip ci]
---
.../perfect_foresight_solver.m | 19 +++++-
.../perfect_foresight_solver_core.m | 25 +++++++-
...foresight_with_expectation_errors_solver.m | 6 +-
matlab/perfect-foresight-models/sim1.m | 36 +++++++++++-
...ked_jacobian_substitute_controlled_paths.m | 48 +++++++++++++++
scripts/dynare.el | 6 +-
.../ramst_pf_controlled_paths.mod | 58 +++++++++++++++++++
7 files changed, 185 insertions(+), 13 deletions(-)
create mode 100644 matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m
create mode 100644 tests/deterministic_simulations/ramst_pf_controlled_paths.mod
diff --git a/matlab/perfect-foresight-models/perfect_foresight_solver.m b/matlab/perfect-foresight-models/perfect_foresight_solver.m
index cc4f2d60a..f55b6e3e2 100644
--- a/matlab/perfect-foresight-models/perfect_foresight_solver.m
+++ b/matlab/perfect-foresight-models/perfect_foresight_solver.m
@@ -23,7 +23,7 @@ function [oo_, ts]=perfect_foresight_solver(M_, options_, oo_, no_error_if_learn
% SPECIAL REQUIREMENTS
% none
-% Copyright © 1996-2024 Dynare Team
+% Copyright © 1996-2025 Dynare Team
%
% This file is part of Dynare.
%
@@ -158,13 +158,19 @@ if completed_share == 1
if options_.simul.endval_steady
oo_.steady_state = steady_state;
end
- % NB: no need to modify oo_.exo_simul and oo_.exo_steady_state, since we simulated 100% of the shock
+ % NB: no need to modify oo_.exo_simul and oo_.exo_steady_state, since we simulated 100% of the shock, except if controlled paths
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ oo_.exo_simul = exo_simul;
+ end
if ~options_.noprint
fprintf('Perfect foresight solution found.\n\n')
end
oo_.deterministic_simulation.status = true;
elseif options_.simul.homotopy_linearization_fallback && completed_share > 0
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ error('Option homotopy_linearization_fallback not available with perfect_foresight_controlled_paths')
+ end
oo_.deterministic_simulation.sim1.endo_simul = endo_simul;
oo_.deterministic_simulation.sim1.exo_simul = exo_simul;
oo_.deterministic_simulation.sim1.steady_state = steady_state;
@@ -188,6 +194,9 @@ elseif options_.simul.homotopy_linearization_fallback && completed_share > 0
oo_.deterministic_simulation.status = true;
oo_.deterministic_simulation.homotopy_linearization = true;
elseif options_.simul.homotopy_marginal_linearization_fallback > 0 && completed_share > options_.simul.homotopy_marginal_linearization_fallback
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ error('Option homotopy_marginal_linearization_fallback not available with perfect_foresight_controlled_paths')
+ end
oo_.deterministic_simulation.sim1.endo_simul = endo_simul;
oo_.deterministic_simulation.sim1.exo_simul = exo_simul;
oo_.deterministic_simulation.sim1.steady_state = steady_state;
@@ -344,7 +353,7 @@ while step > options_.simul.homotopy_min_step_size
end
% Solve for the paths of the endogenous variables.
- [endo_simul, success, maxerror, solver_iter, per_block_status] = perfect_foresight_solver_core(endo_simul, exo_simul, steady_state, exo_steady_state, M_, options_);
+ [endo_simul, success, maxerror, solver_iter, per_block_status, exo_simul] = perfect_foresight_solver_core(endo_simul, exo_simul, steady_state, exo_steady_state, M_, options_);
else
success = false;
maxerror = NaN;
@@ -362,6 +371,10 @@ while step > options_.simul.homotopy_min_step_size
break
end
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ error('Homotopy not available with perfect_foresight_controlled_paths')
+ end
+
if iteration == 1 && ~options_.noprint
fprintf('\nEntering the homotopy method iterations...\n')
iter_summary_table = { sprintf('\nIter. \t | Share \t | Status \t | Max. residual\t | Duration (sec)\n'),
diff --git a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m
index 05d20af69..8d38504c3 100644
--- a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m
+++ b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m
@@ -1,4 +1,4 @@
-function [y, success, maxerror, iter, per_block_status] = perfect_foresight_solver_core(y, exo_simul, steady_state, exo_steady_state, M_, options_)
+function [y, success, maxerror, iter, per_block_status, exo_simul] = perfect_foresight_solver_core(y, exo_simul, steady_state, exo_steady_state, M_, options_)
% Core function calling solvers for perfect foresight model
%
@@ -17,7 +17,7 @@ function [y, success, maxerror, iter, per_block_status] = perfect_foresight_solv
% - iter [integer] Number of iterations of the underlying nonlinear solver (empty for non-iterative methods)
% - per_block_status [struct] In the case of block decomposition, provides per-block solver status information (empty if no block decomposition)
-% Copyright © 2015-2024 Dynare Team
+% Copyright © 2015-2025 Dynare Team
%
% This file is part of Dynare.
%
@@ -55,6 +55,25 @@ if options_.linear && (isequal(options_.stack_solve_algo, 0) || isequal(options_
options_.linear_approximation = true;
end
+if ~isempty(M_.perfect_foresight_controlled_paths)
+ assert(nargout >= 6); % Ensure modified exos are used
+ if options_.bytecode
+ error('perfect_foresight_controlled_paths is not compatible with bytecode option')
+ end
+ if options_.block
+ error('perfect_foresight_controlled_paths is not compatible with block option')
+ end
+ if options_.linear
+ error('perfect_foresight_controlled_paths is not compatible with linear option')
+ end
+ if ~ismember(options_.stack_solve_algo, [0 2 3])
+ error('perfect_foresight_controlled_paths is only available with stack_solve_algo equal to 0, 2 or 3')
+ end
+ if M_.maximum_endo_lead == 0 || M_.maximum_endo_lag == 0
+ error('perfect_foresight_controlled_paths is not compatible with purely backward, purely forward or static models')
+ end
+end
+
maxerror = [];
iter = [];
per_block_status = [];
@@ -103,7 +122,7 @@ else
end
[y, success, maxerror] = sim1_linear(y, exo_simul, steady_state, exo_steady_state, M_, options_);
else
- [y, success, maxerror, iter] = sim1(y, exo_simul, steady_state, M_, options_);
+ [y, success, maxerror, iter, exo_simul] = sim1(y, exo_simul, steady_state, M_, options_);
end
case {1 6}
if options_.linear_approximation
diff --git a/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m b/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m
index 22a07c276..c69029c9d 100644
--- a/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m
+++ b/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m
@@ -7,7 +7,7 @@ function [oo_, ts] = perfect_foresight_with_expectation_errors_solver(M_, option
% OUTPUTS
% oo_ [structure] storing the results
-% Copyright © 2021-2024 Dynare Team
+% Copyright © 2021-2025 Dynare Team
%
% This file is part of Dynare.
%
@@ -28,6 +28,10 @@ if options_.pfwee.constant_simulation_length && ~isempty(options_.simul.last_sim
error('Options constant_simulation_length and last_simulation_period cannot be used together')
end
+if ~isempty(M_.perfect_foresight_controlled_paths)
+ error('perfect_foresight_with_expectation_errors_solver not available with perfect_foresight_controlled_paths')
+end
+
[periods, first_simulation_period] = get_simulation_periods(options_);
% Retrieve initial paths built by pfwee_setup
diff --git a/matlab/perfect-foresight-models/sim1.m b/matlab/perfect-foresight-models/sim1.m
index 55f5f6f63..3e5be836e 100644
--- a/matlab/perfect-foresight-models/sim1.m
+++ b/matlab/perfect-foresight-models/sim1.m
@@ -1,4 +1,4 @@
-function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_)
+function [endogenousvariables, success, err, iter, exogenousvariables] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_)
% [endogenousvariables, success, err, iter] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_)
% Performs deterministic simulations with lead or lag of one period, using
% a basic Newton solver on sparse matrices.
@@ -15,8 +15,10 @@ function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, e
% - success [logical] Whether a solution was found
% - err [double] ∞-norm of the residual
% - iter [integer] Number of iterations
+% - exogenousvariables [double] T*M array, paths for the exogenous variables
+% (may be modified if perfect_foresight_controlled_paths present)
-% Copyright © 1996-2024 Dynare Team
+% Copyright © 1996-2025 Dynare Team
%
% This file is part of Dynare.
%
@@ -39,6 +41,19 @@ ny = M_.endo_nbr;
periods = get_simulation_periods(options_);
vperiods = periods*ones(1,options_.simul.maxit);
+if ~isempty(M_.perfect_foresight_controlled_paths)
+ for i=1:length(M_.perfect_foresight_controlled_paths)
+ endogenousvariables(M_.perfect_foresight_controlled_paths(i).exogenize_id, M_.perfect_foresight_controlled_paths(i).periods + M_.maximum_lag) = M_.perfect_foresight_controlled_paths(i).value;
+ end
+
+ if options_.debug
+ error('Debugging not available with perfect_foresight_controlled_paths')
+ end
+ if options_.endogenous_terminal_period
+ error('The endogenous_terminal_period option not available with perfect_foresight_controlled_paths')
+ end
+end
+
if M_.maximum_lag > 0
y0 = endogenousvariables(:, M_.maximum_lag);
else
@@ -70,6 +85,11 @@ for iter = 1:options_.simul.maxit
[res, A] = perfect_foresight_problem(y, y0, yT, exogenousvariables, M_.params, steadystate, periods, M_, options_);
% A is the stacked Jacobian with period x equations alongs the rows and
% periods times variables (in declaration order) along the columns
+
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ A = stacked_jacobian_substitute_controlled_paths(A, y, y0, yT, exogenousvariables, steadystate, periods, M_);
+ end
+
if options_.debug && iter==1
[row,col]=find(A);
row=setdiff(1:periods*ny,row);
@@ -136,6 +156,16 @@ for iter = 1:options_.simul.maxit
end
end
y = y + dy;
+
+ if ~isempty(M_.perfect_foresight_controlled_paths)
+ for i = 1:length(M_.perfect_foresight_controlled_paths)
+ p = M_.perfect_foresight_controlled_paths(i).periods;
+ endogenize_id = M_.perfect_foresight_controlled_paths(i).endogenize_id;
+ exogenize_id = M_.perfect_foresight_controlled_paths(i).exogenize_id;
+ y(exogenize_id+(p-1)*M_.endo_nbr) = M_.perfect_foresight_controlled_paths(i).value;
+ exogenousvariables(p,endogenize_id) = exogenousvariables(p,endogenize_id) + dy(exogenize_id+(p-1)*M_.endo_nbr);
+ end
+ end
end
endogenousvariables(:, M_.maximum_lag+(1:periods)) = reshape(y, ny, periods);
@@ -278,7 +308,7 @@ end
function display_critical_variables(dyy, M_, noprint)
-if noprint
+if noprint || ~isempty(M_.perfect_foresight_controlled_paths)
return
end
diff --git a/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m b/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m
new file mode 100644
index 000000000..9597c0190
--- /dev/null
+++ b/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m
@@ -0,0 +1,48 @@
+function A = stacked_jacobian_substitute_controlled_paths(A, y, y0, yT, exogenousvariables, steadystate, periods, M_)
+
+% Copyright © 2025 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 <https://www.gnu.org/licenses/>.
+
+% Reorganize M_.perfect_foresight_controlled_paths period by period
+flips_by_period = struct('exogenize_id', cell(periods, 1), 'endogenize_id', cell(periods, 1));
+for i=1:length(M_.perfect_foresight_controlled_paths)
+ for p = M_.perfect_foresight_controlled_paths(i).periods
+ flips_by_period(p).exogenize_id(end+1) = M_.perfect_foresight_controlled_paths(i).exogenize_id;
+ flips_by_period(p).endogenize_id(end+1) = M_.perfect_foresight_controlled_paths(i).endogenize_id;
+ end
+end
+
+dynamic_g1 = str2func([M_.fname '.sparse.dynamic_g1']);
+
+for p = 1:periods
+ nsubst = length(flips_by_period(p).exogenize_id);
+ if nsubst > 0
+ if p == 1
+ y3n = [ y0; y(1:2*M_.endo_nbr) ];
+ elseif p == periods
+ y3n = [ y(end-2*M_.endo_nbr+1:end); yT ];
+ else
+ y3n = y((p-2)*M_.endo_nbr+(1:3*M_.endo_nbr));
+ end
+ g1 = dynamic_g1(y3n, exogenousvariables(p,:), M_.params, steadystate, M_.dynamic_g1_sparse_rowval, M_.dynamic_g1_sparse_colval, M_.dynamic_g1_sparse_colptr);
+
+ subst = [ zeros((p-1)*M_.endo_nbr, nsubst);
+ g1(:,3*M_.endo_nbr+flips_by_period(p).endogenize_id);
+ zeros((periods-p)*M_.endo_nbr, nsubst) ];
+ A(:,(p-1)*M_.endo_nbr+flips_by_period(p).exogenize_id) = subst;
+ end
+end
diff --git a/scripts/dynare.el b/scripts/dynare.el
index e723bf6b7..50d6e1d25 100644
--- a/scripts/dynare.el
+++ b/scripts/dynare.el
@@ -1,7 +1,7 @@
;;; dynare.el --- major mode for editing Dynare mod files
;; Copyright © 2010 Yannick Kalantzis
-;; Copyright © 2019-2024 Dynare Team
+;; Copyright © 2019-2025 Dynare Team
;;
;; This program is free software; you can redistribute it and/or modify
;; it under the terms of the GNU General Public License as published by
@@ -88,7 +88,7 @@
'("stderr" "values" "periods" "scales" "restriction" "exclusion"
"upper_cholesky" "lower_cholesky" "equation" "bind" "relax" "error_bind"
"error_relax" "add" "multiply" "target" "auxname_target_nonstationary"
- "component" "growth" "auxname" "kind" "weights")
+ "component" "growth" "auxname" "kind" "weights" "exogenize" "endogenize")
"Dynare statements-like keywords.")
;; Those keywords that makes the lexer enter the DYNARE_BLOCK start condition
@@ -106,7 +106,7 @@
"conditional_forecast_paths" "svar_identification" "moment_calibration"
"irf_calibration" "ramsey_constraints" "generate_irfs" "matched_moments"
"occbin_constraints" "model_replace" "pac_target_info" "matched_irfs"
- "matched_irfs_weights" "verbatim")
+ "matched_irfs_weights" "perfect_foresight_controlled_paths" "verbatim")
"Dynare block keywords."))
;; Mathematical functions and operators used in model equations (see "hand_side" in Bison file)
diff --git a/tests/deterministic_simulations/ramst_pf_controlled_paths.mod b/tests/deterministic_simulations/ramst_pf_controlled_paths.mod
new file mode 100644
index 000000000..d39f80ce2
--- /dev/null
+++ b/tests/deterministic_simulations/ramst_pf_controlled_paths.mod
@@ -0,0 +1,58 @@
+var c k y;
+varexo x e_y;
+
+parameters alph gam delt bet aa;
+alph=0.5;
+gam=0.5;
+delt=0.02;
+bet=0.05;
+aa=0.5;
+
+
+model;
+c + k - aa*x*k(-1)^alph - (1-delt)*k(-1);
+c^(-gam) - (1+bet)^(-1)*(aa*alph*x(+1)*k^(alph-1) + 1 - delt)*c(+1)^(-gam);
+y = 0.98*y(-1)+e_y;
+end;
+
+initval;
+x = 1;
+k = ((delt+bet)/(1.0*aa*alph))^(1/(alph-1));
+c = aa*k^alph-delt*k;
+y = 0;
+end;
+
+steady;
+
+check;
+
+perfect_foresight_controlled_paths;
+ exogenize c;
+ periods 2 4;
+ values 1.6 1.7;
+ endogenize x;
+
+ exogenize k;
+ periods 7;
+ values 13;
+ endogenize x;
+
+ exogenize y;
+ periods 3:5;
+ values 0.5;
+ endogenize e_y;
+end;
+
+shocks;
+var x;
+periods 1;
+values 1.2;
+end;
+
+
+
+perfect_foresight_setup(periods=200);
+perfect_foresight_solver;
+
+rplot c;
+rplot k;
--
GitLab