diff --git a/matlab/perfect-foresight-models/det_cond_forecast.m b/matlab/perfect-foresight-models/det_cond_forecast.m
index 4652c6c5829df7fdaf4cfdca218456665ad3afd3..4a43627e5d34e11473e18be22de383cf97e68d3e 100644
--- a/matlab/perfect-foresight-models/det_cond_forecast.m
+++ b/matlab/perfect-foresight-models/det_cond_forecast.m
@@ -1,18 +1,16 @@
-function data_set = det_cond_forecast(varargin)
+function data_set = det_cond_forecast(plan, dset, range)
% Computes conditional forecasts using the extended path method.
%
% INPUTS
% o plan [structure] A structure describing the different shocks and the endogenous varibales, the date of the shocks and the path of the shock.
% The plan structure is created by the functions init_plan, basic_plan and flip_plan
-% o [dataset] [dseries] A dserie containing the initial values of the shocks and the endogenous variables (usually the dseries generated by the smoother).
-% o [starting_date] [dates] The first date of the forecast.
-%
+% o dset [dseries] The initial values of the shocks and the endogenous variables (usually the dseries generated by the smoother).
+% o range [dates] The dates of the forecast.
%
% OUTPUTS
% dataset [dseries] Returns a dseries containing the forecasted endgenous variables and shocks
-%
-% Copyright © 2013-2023 Dynare Team
+% Copyright © 2013-2024 Dynare Team
%
% This file is part of Dynare.
%
@@ -56,206 +54,200 @@ surprise = 0;
pf = 0;
is_shock = [];
is_constraint = [];
-if length(varargin) > 3
- % regular way to call
- constrained_paths = varargin{1};
- max_periods_simulation = size(constrained_paths, 2);
- constrained_vars = varargin{2};
- options_cond_fcst = varargin{3};
- constrained_perfect_foresight = varargin{4};
- constraint_index = cell(max_periods_simulation,1);
- nvars = length(constrained_vars);
- for i = 1:max_periods_simulation
- constraint_index{i} = 1:nvars;
- end
- direct_mode = 0;
- shocks_present = 0;
- controlled_varexo = options_cond_fcst.controlled_varexo;
- nvarexo = size(controlled_varexo, 1);
- options_cond_fcst.controlled_varexo = zeros(nvarexo,1);
- exo_names = M_.exo_names;
- for i = 1:nvarexo
- j = find(strcmp(controlled_varexo(i,:), exo_names));
- if ~isempty(j)
- options_cond_fcst.controlled_varexo(i) = j;
- else
- error(['Unknown exogenous variable ' controlled_varexo(i,:)]);
- end
- end
-else
- % alternative way to call: plan, dset, dates_of_frcst
- plan = varargin{1};
- if length(varargin) >= 2
- dset = varargin{2};
- if ~isa(dset,'dseries')
- error('det_cond_forecast: The second argument should be a dseries');
- end
- if length(varargin) >= 3
- range = varargin{3};
- if ~isa(range,'dates')
- error('det_cond_forecast: The third argument should be a dates object');
- end
- %if (range(range.ndat) > dset.time(dset.nobs) )
- if (range(range.ndat) > dset.dates(dset.nobs)+1 )
- s1 = strings(dset.dates(dset.nobs));
- s2 = strings(range(range.ndat));
- error(['det_cond_forecast: The dseries ' inputname(2) ' finishes at time ' s1{1} ' before the last period of forecast ' s2{1}]);
- end
+if ~isa(dset,'dseries')
+ error('det_cond_forecast: The second argument should be a dseries');
+end
- sym_dset = dset(dates(-range(1)):dates(range(range.ndat)));
- periods = options_.periods + M_.maximum_lag + M_.maximum_lead;
- total_periods = periods + range.ndat;
- if isfield(oo_, 'exo_simul')
- if size(oo_.exo_simul, 1) ~= total_periods
- oo_.exo_simul = repmat(oo_.exo_steady_state',total_periods,1);
- end
- else
- oo_.exo_simul = repmat(oo_.exo_steady_state',total_periods,1);
- end
+if ~isa(range,'dates')
+ error('det_cond_forecast: The third argument should be a dates object');
+end
- oo_.endo_simul = repmat(oo_.steady_state, 1, total_periods);
+if (range(range.ndat) > dset.dates(dset.nobs)+1 )
+ s1 = strings(dset.dates(dset.nobs));
+ s2 = strings(range(range.ndat));
+ error(['det_cond_forecast: The dseries ' inputname(2) ' finishes at time ' s1{1} ' before the last period of forecast ' s2{1}]);
+end
- for i = 1:sym_dset.vobs
- iy = find(strcmp(strtrim(sym_dset.name{i}), strtrim(plan.endo_names)));
- if ~isempty(iy)
- oo_.endo_simul(iy,1:sym_dset.nobs) = sym_dset.data(:, i);
- initial_conditions(iy) = sym_dset.data(1, i);
- else
- ix = find(strcmp(strtrim(sym_dset.name{i}), strtrim(plan.exo_names)));
- if ~isempty(ix)
- oo_.exo_simul(1, ix) = sym_dset.data(1, i)';
- else
- %warning(['The variable ' sym_dset.name{i} ' in the dataset ' inputname(2) ' is not a endogenous neither an exogenous variable!!']);
- end
- end
- end
- for i = 1:length(M_.aux_vars)
- if M_.aux_vars(i).type == 1 %lag variable
- iy = find(strcmp(M_.endo_names{M_.aux_vars(i).orig_index}, sym_dset.name));
- if ~isempty(iy)
- oo_.endo_simul(M_.aux_vars(i).endo_index, 1:sym_dset.nobs) = dset(dates(range(1) + (M_.aux_vars(i).orig_lead_lag - 1))).data(:,iy);
- initial_conditions(M_.aux_vars(i).endo_index) = dset(dates(range(1) + (M_.aux_vars(i).orig_lead_lag - 1))).data(:,iy);
- else
- warning(['det_cond_forecast: The variable auxiliary ' M_.endo_names{M_.aux_vars(i).endo_index} ' associated with the variable ' M_.endo_names{M_.aux_vars(i).orig_index} ' does not appear in the dataset']);
- end
- else
- oo_.endo_simul(M_.aux_vars(i).endo_index, 1:sym_dset.nobs) = repmat(oo_.steady_state(M_.aux_vars(i).endo_index), 1, range.ndat + 1);
- end
- end
- %Compute the initial path using the the steady-state
- % steady-state
- %for jj = 2 : (options_.periods + 2)
- for jj = 2 : (range.ndat + 2)
- oo_.endo_simul(:, jj) = oo_.steady_state;
- end
- missings = isnan(oo_.endo_simul(:,1));
- if any(missings)
- for jj = 1:M_.endo_nbr
- if missings(jj)
- oo_.endo_simul(jj,1) = oo_.steady_state(jj,1);
- end
- end
- end
+sym_dset = dset(dates(-range(1)):dates(range(range.ndat)));
+periods = options_.periods + M_.maximum_lag + M_.maximum_lead;
+total_periods = periods + range.ndat;
+if isfield(oo_, 'exo_simul')
+ if size(oo_.exo_simul, 1) ~= total_periods
+ oo_.exo_simul = repmat(oo_.exo_steady_state',total_periods,1);
+ end
+else
+ oo_.exo_simul = repmat(oo_.exo_steady_state',total_periods,1);
+end
- if options_.bytecode
- save_options_dynatol_f = options_.dynatol.f;
- options_.dynatol.f = 1e-7;
- if options_.block
- [endo, exo] = bytecode('extended_path', plan, 'block_decomposed', M_, options_, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
- else
- [endo, exo] = bytecode('extended_path', plan, M_, options_, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
- end
- options_.dynatol.f = save_options_dynatol_f;
-
- oo_.endo_simul = endo;
- oo_.exo_simul = exo;
-
- endo = endo';
- endo_l = size(endo(1+M_.maximum_lag:end,:),1);
- jrng = dates(plan.date(1)):dates(plan.date(1)+endo_l);
- data_set = dseries(nan(endo_l, dset.vobs), plan.date(1), dset.name);
- for i = 1:length(dset.name)
- pos = find(strcmp(dset.name{i},plan.endo_names));
- if ~isempty(pos)
- data_set.(dset.name{i}) = dseries(endo(1+M_.maximum_lag:end,pos), plan.date(1), dset.name{i});
- else
- pos = find(strcmp(dset.name{i},plan.exo_names));
- if ~isempty(pos)
- data_set{dset.name{i}} = dseries(exo(1+M_.maximum_lag:end,pos), plan.date(1),dset.name{i});
- end
- end
- end
- data_set = [dset(dset.dates(1):(plan.date(1)-1)) ; data_set];
- for i=1:M_.exo_nbr
- pos = find(strcmp(M_.exo_names{i}, dset.name));
- if isempty(pos)
- data_set{M_.exo_names{i}} = dseries(exo(1+M_.maximum_lag:end,i), plan.date(1), M_.exo_names{i});
- else
- data_set{M_.exo_names{i}}(plan.date(1):plan.date(1)+ (size(exo, 1) - M_.maximum_lag)) = exo(1+M_.maximum_lag:end,i);
- end
- end
- data_set = merge(dset(dset.dates(1):(plan.date(1)-1)), data_set);
- return
- union_names = union(data_set.name, dset.name);
- dif = setdiff(union_names, data_set.name);
- data_set_nobs = data_set.nobs;
- for i = 1:length(dif)
- data_set{dif{i}} = dseries(nan(data_set_nobs,1),plan.date(1), dif(i), dif(i));
- end
- dif = setdiff(union_names, dset.name);
- dset_nobs = dset.nobs;
- for i = 1:length(dif)
- dset{dif{i}} = dseries(nan(dset_nobs,1),dset.dates(1), dif(i), dif(i));
- end
- data_set = [dset(dset.dates(1):(plan.date(1)-1)) ; data_set];
- return
- end
+oo_.endo_simul = repmat(oo_.steady_state, 1, total_periods);
+
+for i = 1:sym_dset.vobs
+ iy = find(strcmp(strtrim(sym_dset.name{i}), strtrim(plan.endo_names)));
+ if ~isempty(iy)
+ oo_.endo_simul(iy,1:sym_dset.nobs) = sym_dset.data(:, i);
+ initial_conditions(iy) = sym_dset.data(1, i);
+ else
+ ix = find(strcmp(strtrim(sym_dset.name{i}), strtrim(plan.exo_names)));
+ if ~isempty(ix)
+ oo_.exo_simul(1, ix) = sym_dset.data(1, i)';
else
- error('det_cond_forecast: unaccounted case, please contact the developers');
+ %%warning(['The variable ' sym_dset.name{i} ' in the dataset ' inputname(2) ' is not a endogenous neither an exogenous variable!!']);
+ end
+ end
+end
+for i = 1:length(M_.aux_vars)
+ if M_.aux_vars(i).type == 1 %lag variable
+ iy = find(strcmp(M_.endo_names{M_.aux_vars(i).orig_index}, sym_dset.name));
+ if ~isempty(iy)
+ oo_.endo_simul(M_.aux_vars(i).endo_index, 1:sym_dset.nobs) = dset(dates(range(1) + (M_.aux_vars(i).orig_lead_lag - 1))).data(:,iy);
+ initial_conditions(M_.aux_vars(i).endo_index) = dset(dates(range(1) + (M_.aux_vars(i).orig_lead_lag - 1))).data(:,iy);
+ else
+ warning(['det_cond_forecast: The variable auxiliary ' M_.endo_names{M_.aux_vars(i).endo_index} ' associated with the variable ' M_.endo_names{M_.aux_vars(i).orig_index} ' does not appear in the dataset']);
end
-
else
- oo_.exo_simul = repmat(oo_.exo_steady_state',options_.periods+2,1);
- oo_.endo_simul = repmat(oo_.steady_state, 1, options_.periods+2);
+ oo_.endo_simul(M_.aux_vars(i).endo_index, 1:sym_dset.nobs) = repmat(oo_.steady_state(M_.aux_vars(i).endo_index), 1, range.ndat + 1);
+ end
+end
+%% Compute the initial path using the the steady-state
+for jj = 2 : (range.ndat + 2)
+ oo_.endo_simul(:, jj) = oo_.steady_state;
+end
+missings = isnan(oo_.endo_simul(:,1));
+if any(missings)
+ for jj = 1:M_.endo_nbr
+ if missings(jj)
+ oo_.endo_simul(jj,1) = oo_.steady_state(jj,1);
+ end
end
+end
- direct_mode = 1;
- constrained_paths = plan.constrained_paths_;
- constrained_vars = plan.constrained_vars_;
- options_cond_fcst = plan.options_cond_fcst_;
- constrained_perfect_foresight = plan.constrained_perfect_foresight_;
- constrained_periods = plan.constrained_date_;
- if ~isempty(plan.shock_paths_)
- shock_paths = plan.shock_paths_;
- shock_vars = plan.shock_vars_;
- shock_perfect_foresight = plan.shock_perfect_foresight_;
- shock_periods = plan.shock_date_;
- shocks_present = 1;
+if options_.bytecode
+ save_options_dynatol_f = options_.dynatol.f;
+ options_.dynatol.f = 1e-7;
+ if options_.block
+ [endo, exo] = bytecode('extended_path', plan, 'block_decomposed', M_, options_, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
else
- shocks_present = 0;
+ [endo, exo] = bytecode('extended_path', plan, M_, options_, oo_.endo_simul, oo_.exo_simul, M_.params, oo_.steady_state, options_.periods);
+ end
+ options_.dynatol.f = save_options_dynatol_f;
+
+ oo_.endo_simul = endo;
+ oo_.exo_simul = exo;
+
+ endo = endo';
+ endo_l = size(endo(1+M_.maximum_lag:end,:),1);
+ jrng = dates(plan.date(1)):dates(plan.date(1)+endo_l);
+ data_set = dseries(nan(endo_l, dset.vobs), plan.date(1), dset.name);
+ for i = 1:length(dset.name)
+ pos = find(strcmp(dset.name{i},plan.endo_names));
+ if ~isempty(pos)
+ data_set.(dset.name{i}) = dseries(endo(1+M_.maximum_lag:end,pos), plan.date(1), dset.name{i});
+ else
+ pos = find(strcmp(dset.name{i},plan.exo_names));
+ if ~isempty(pos)
+ data_set{dset.name{i}} = dseries(exo(1+M_.maximum_lag:end,pos), plan.date(1),dset.name{i});
+ end
+ end
end
+ data_set = [dset(dset.dates(1):(plan.date(1)-1)) ; data_set];
+ for i=1:M_.exo_nbr
+ pos = find(strcmp(M_.exo_names{i}, dset.name));
+ if isempty(pos)
+ data_set{M_.exo_names{i}} = dseries(exo(1+M_.maximum_lag:end,i), plan.date(1), M_.exo_names{i});
+ else
+ data_set{M_.exo_names{i}}(plan.date(1):plan.date(1)+ (size(exo, 1) - M_.maximum_lag)) = exo(1+M_.maximum_lag:end,i);
+ end
+ end
+ data_set = merge(dset(dset.dates(1):(plan.date(1)-1)), data_set);
+ return
+ union_names = union(data_set.name, dset.name);
+ dif = setdiff(union_names, data_set.name);
+ data_set_nobs = data_set.nobs;
+ for i = 1:length(dif)
+ data_set{dif{i}} = dseries(nan(data_set_nobs,1),plan.date(1), dif(i), dif(i));
+ end
+ dif = setdiff(union_names, dset.name);
+ dset_nobs = dset.nobs;
+ for i = 1:length(dif)
+ dset{dif{i}} = dseries(nan(dset_nobs,1),dset.dates(1), dif(i), dif(i));
+ end
+ data_set = [dset(dset.dates(1):(plan.date(1)-1)) ; data_set];
+ return
+end
- total_periods = plan.date;
+constrained_paths = plan.constrained_paths_;
+constrained_vars = plan.constrained_vars_;
+options_cond_fcst = plan.options_cond_fcst_;
+constrained_perfect_foresight = plan.constrained_perfect_foresight_;
+constrained_periods = plan.constrained_date_;
+if ~isempty(plan.shock_paths_)
+ shock_paths = plan.shock_paths_;
+ shock_vars = plan.shock_vars_;
+ shock_perfect_foresight = plan.shock_perfect_foresight_;
+ shock_periods = plan.shock_date_;
+ shocks_present = 1;
+else
+ shocks_present = 0;
end
+total_periods = plan.date;
+
if ~isfield(options_cond_fcst,'periods') || isempty(options_cond_fcst.periods)
options_cond_fcst.periods = 100;
end
options_.periods = 10;
-if direct_mode == 1
- n_periods = length(constrained_periods);
- is_constraint = zeros(length(total_periods), n_periods);
- constrained_paths_cell = constrained_paths;
- clear constrained_paths;
- constrained_paths = zeros(n_periods, length(total_periods));
- max_periods_simulation = 0;
+n_periods = length(constrained_periods);
+is_constraint = zeros(length(total_periods), n_periods);
+constrained_paths_cell = constrained_paths;
+clear constrained_paths;
+constrained_paths = zeros(n_periods, length(total_periods));
+max_periods_simulation = 0;
+for i = 1:n_periods
+ period_i = constrained_periods{i};
+ %period_i
+ tp = total_periods(1);
+ if size(period_i) > 1
+ init_periods = period_i(1);
+ tp_end = period_i(end);
+ else
+ init_periods = period_i;
+ tp_end = period_i;
+ end
+ tp0 = tp;
+ while tp < init_periods
+ tp = tp + 1;
+ end
+ j = 0;
+ while tp <= tp_end
+ is_constraint(tp - tp0 + 1, i) = 1;
+ constrained_paths(i, tp - tp0 + 1) = constrained_paths_cell{i}(j + 1);
+ tp = tp + 1;
+ j = j + 1;
+ end
+ if tp - tp0 > max_periods_simulation
+ max_periods_simulation = tp - tp0;
+ end
+end
+n_nnz = length(sum(is_constraint,2));
+if n_nnz > 0
+ constraint_index = cell(n_nnz,1);
+ for i= 1:n_nnz
+ constraint_index{i} = find(is_constraint(i,:));
+ end
+end
+if shocks_present
+ n_periods = length(shock_periods);
+ shock_paths_cell = shock_paths;
+ clear shock_paths;
+ shock_paths = zeros(n_periods, length(total_periods));
+ is_shock = zeros(length(total_periods), n_periods);
for i = 1:n_periods
- period_i = constrained_periods{i};
- %period_i
+ period_i = shock_periods{i};
+ %period_i
tp = total_periods(1);
if size(period_i) > 1
init_periods = period_i(1);
@@ -270,8 +262,8 @@ if direct_mode == 1
end
j = 0;
while tp <= tp_end
- is_constraint(tp - tp0 + 1, i) = 1;
- constrained_paths(i, tp - tp0 + 1) = constrained_paths_cell{i}(j + 1);
+ is_shock(tp - tp0 + 1, i) = 1;
+ shock_paths(i, tp - tp0 + 1) = shock_paths_cell{i}(j + 1);
tp = tp + 1;
j = j + 1;
end
@@ -279,59 +271,16 @@ if direct_mode == 1
max_periods_simulation = tp - tp0;
end
end
- n_nnz = length(sum(is_constraint,2));
+ n_nnz = length(sum(is_shock,2));
if n_nnz > 0
- constraint_index = cell(n_nnz,1);
+ shock_index = cell(n_nnz, 1);
for i= 1:n_nnz
- constraint_index{i} = find(is_constraint(i,:));
- end
- end
- if shocks_present
- n_periods = length(shock_periods);
- shock_paths_cell = shock_paths;
- clear shock_paths;
- shock_paths = zeros(n_periods, length(total_periods));
- is_shock = zeros(length(total_periods), n_periods);
- for i = 1:n_periods
- period_i = shock_periods{i};
- %period_i
- tp = total_periods(1);
- if size(period_i) > 1
- init_periods = period_i(1);
- tp_end = period_i(end);
- else
- init_periods = period_i;
- tp_end = period_i;
- end
- tp0 = tp;
- while tp < init_periods
- tp = tp + 1;
- end
- j = 0;
- while tp <= tp_end
- is_shock(tp - tp0 + 1, i) = 1;
- shock_paths(i, tp - tp0 + 1) = shock_paths_cell{i}(j + 1);
- tp = tp + 1;
- j = j + 1;
- end
- if tp - tp0 > max_periods_simulation
- max_periods_simulation = tp - tp0;
- end
- end
- n_nnz = length(sum(is_shock,2));
- if n_nnz > 0
- shock_index = cell(n_nnz, 1);
- for i= 1:n_nnz
- shock_index{i} = find(is_shock(i,:));
- end
+ shock_index{i} = find(is_shock(i,:));
end
end
-else
- is_constraint = ones(size(constrained_paths));
end
-
maximum_lag = M_.maximum_lag;
ys = oo_.steady_state;
@@ -350,24 +299,7 @@ else
error('det_cond_forecast: to run a deterministic conditional forecast you have to specify the controlled exogenous variables using the option controlled_varexo in forecast command');
end
-% if n_endo_constrained == 0
-% options_.ep.use_bytecode = options_.bytecode;
-% data_set = extended_path(initial_conditions, max_periods_simulation);
-% end
-
-if length(varargin) >= 1
- controlled_varexo = options_cond_fcst.controlled_varexo;
-else
- exo_names = M_.exo_names;
- controlled_varexo = zeros(1,n_control_exo);
- for i = 1:nx
- for j=1:n_control_exo
- if strcmp(exo_names{i}, deblank(options_cond_fcst.controlled_varexo(j,:)))
- controlled_varexo(j) = i;
- end
- end
- end
-end
+controlled_varexo = options_cond_fcst.controlled_varexo;
%todo check if zero => error message
@@ -572,7 +504,7 @@ if pf && ~surprise
else
for t = 1:constrained_periods
- if direct_mode && ~isempty(is_constraint)
+ if ~isempty(is_constraint)
[pos_constrained_pf, ~] = find(constrained_perfect_foresight .* is_constraint(t, :)');
indx_endo_solve_pf = constrained_vars(pos_constrained_pf);
if isempty(indx_endo_solve_pf)
@@ -591,7 +523,7 @@ else
end
end
- if direct_mode && ~isempty(is_shock)
+ if ~isempty(is_shock)
[pos_shock_pf, ~] = find(shock_perfect_foresight .* is_shock(t, :)');
indx_endo_solve_pf = shock_vars(pos_shock_pf);
if isempty(indx_endo_solve_pf)
@@ -628,17 +560,15 @@ else
time_index_constraint = maximum_lag + 1:maximum_lag + constrained_periods - t + 1;
- if direct_mode
- nb_shocks = length(plan.shock_vars_);
- if nb_shocks > 0
- shock_index_t = shock_index{t};
- shock_vars_t = shock_vars(shock_index_t);
- for shock_indx = shock_index_t
- if shock_perfect_foresight(shock_indx)
- oo_.exo_simul(time_index_constraint,shock_vars_t(shock_indx)) = shock_paths(shock_indx,t:constrained_periods);
- else
- oo_.exo_simul(maximum_lag + 1,shock_vars_t(shock_indx)) = shock_paths(shock_indx,t);
- end
+ nb_shocks = length(plan.shock_vars_);
+ if nb_shocks > 0
+ shock_index_t = shock_index{t};
+ shock_vars_t = shock_vars(shock_index_t);
+ for shock_indx = shock_index_t
+ if shock_perfect_foresight(shock_indx)
+ oo_.exo_simul(time_index_constraint,shock_vars_t(shock_indx)) = shock_paths(shock_indx,t:constrained_periods);
+ else
+ oo_.exo_simul(maximum_lag + 1,shock_vars_t(shock_indx)) = shock_paths(shock_indx,t);
end
end
end
@@ -913,6 +843,4 @@ options_.initval_file = save_options_initval_file;
options_.verbosity = verbosity;
oo_.endo_simul = endo';
oo_.exo_simul = exo;
-if direct_mode
- data_set = dseries([endo(2:end,1:M_.orig_endo_nbr) exo(2:end,:)], total_periods(1), {plan.endo_names{:} plan.exo_names{:}}, {plan.endo_names{:} plan.exo_names{:}});
-end
+data_set = dseries([endo(2:end,1:M_.orig_endo_nbr) exo(2:end,:)], total_periods(1), {plan.endo_names{:} plan.exo_names{:}}, {plan.endo_names{:} plan.exo_names{:}});