diff --git a/matlab/+occbin/DSGE_smoother.m b/matlab/+occbin/DSGE_smoother.m
index 364b328862aaf3b86735f10fcbe5c26117419457..f820e27555f49e91537bb778ddc16dc991d0bb05 100644
--- a/matlab/+occbin/DSGE_smoother.m
+++ b/matlab/+occbin/DSGE_smoother.m
@@ -85,7 +85,7 @@ if options_.occbin.smoother.linear_smoother && nargin==12
oo_.occbin.linear_smoother.alphahat0=alphahat0;
oo_.occbin.linear_smoother.state_uncertainty0=state_uncertainty0;
- fprintf('\nOccbin: linear smoother done.\n')
+ disp_verbose('Occbin: linear smoother done.',options_.verbosity)
options_.occbin.smoother.status=true;
end
% if init_mode
@@ -123,22 +123,20 @@ occbin_options.opts_simul = opts_simul; % this builds the opts_simul options fie
occbin_options.opts_regime.binding_indicator = options_.occbin.smoother.init_binding_indicator;
occbin_options.opts_regime.regime_history=options_.occbin.smoother.init_regime_history;
-error_indicator=false;
options_.noprint = true;
-try
- %blanket try-catch should be replaced be proper error handling, see https://git.dynare.org/Dynare/dynare/-/merge_requests/2226#note_20318
- [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T0,R0,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,occbin_options);% T1=TT;
-catch ME
- error_indicator=true;
- disp(ME.message)
- for iter = 1:numel(ME.stack)
- ME.stack(iter)
- end
+[alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T0,R0,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0,~,error_indicator] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,occbin_options);% T1=TT;
+
+if error_indicator(1) || isempty(alphahat0)
+ if ~options_.occbin.smoother.linear_smoother || nargin~=12 %make sure linear smoother results are set before using them
+ options_.occbin.smoother.status=false;
+ [~,etahat,~,~,~,~,~,~,~,~,~,~,~,~,~,~,alphahat0] = ...
+ DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_);
+ options_.occbin.smoother.status=true;
+ else
+ etahat= oo_.occbin.linear_smoother.etahat;
+ alphahat0= oo_.occbin.linear_smoother.alphahat0;
end
-if error_indicator || isempty(alphahat0)
- etahat= oo_.occbin.linear_smoother.etahat;
- alphahat0= oo_.occbin.linear_smoother.alphahat0;
base_regime = struct();
if M_.occbin.constraint_nbr==1
base_regime.regime = 0;
@@ -173,7 +171,7 @@ occbin_options.first_period_occbin_update = inf;
opts_regime.binding_indicator=[];
regime_history0 = regime_history;
-fprintf('Occbin smoother iteration 1.\n')
+disp_verbose('Occbin smoother iteration 1.',options_.verbosity)
opts_simul.SHOCKS = [etahat(:,1:end)'; zeros(1,M_.exo_nbr)];
opts_simul.exo_pos = 1:M_.exo_nbr;
opts_simul.endo_init = alphahat0(oo_.dr.inv_order_var,1);
@@ -182,8 +180,7 @@ opts_simul.periods = size(opts_simul.SHOCKS,1);
options_.occbin.simul=opts_simul;
[~, out, ss] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
if out.error_flag
- fprintf('Occbin smoother:: simulation within smoother did not converge.\n')
- print_info(out.error_flag, options_.noprint, options_)
+ disp_verbose('Occbin smoother:: simulation within smoother did not converge.',options_.verbosity)
oo_.occbin.smoother.error_flag=321;
return;
end
@@ -226,7 +223,7 @@ end
while is_changed && maxiter>iter && ~is_periodic
iter=iter+1;
- fprintf('Occbin smoother iteration %u.\n', iter)
+ disp_verbose(sprintf('Occbin smoother iteration %u.', iter),options_.verbosity)
occbin_options.opts_regime.regime_history=regime_history;
[alphahat,etahat,epsilonhat,~,SteadyState,trend_coeff,~,T0,R0,P,~,decomp,Trend,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0]...
= DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,occbin_options,TT,RR,CC);
@@ -245,8 +242,7 @@ while is_changed && maxiter>iter && ~is_periodic
options_.occbin.simul=opts_simul;
[~, out, ss] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
if out.error_flag
- fprintf('Occbin smoother:: simulation within smoother did not converge.\n')
- print_info(out.error_flag, false, options_)
+ disp_verbose('Occbin smoother:: simulation within smoother did not converge.',options_.verbosity)
oo_.occbin.smoother.error_flag=321;
return;
end
@@ -300,13 +296,13 @@ while is_changed && maxiter>iter && ~is_periodic
eee(:,k) = eig(TT(:,:,k));
end
if options_.debug
- err_eig(iter-1) = max(max(abs(sort(eee)-sort(sto_eee))));
- err_alphahat(iter-1) = max(max(max(abs(alphahat-sto_alphahat))));
- err_etahat(iter-1) = max(max(max(abs(etahat-sto_etahat{iter-1}))));
- err_CC(iter-1) = max(max(max(abs(CC-sto_CC))));
- err_RR(iter-1) = max(max(max(abs(RR-sto_RR))));
- err_TT(iter-1) = max(max(max(abs(TT-sto_TT))));
- end
+ err_eig(iter-1) = max(max(abs(sort(eee)-sort(sto_eee))));
+ err_alphahat(iter-1) = max(max(max(abs(alphahat-sto_alphahat))));
+ err_etahat(iter-1) = max(max(max(abs(etahat-sto_etahat{iter-1}))));
+ err_CC(iter-1) = max(max(max(abs(CC-sto_CC))));
+ err_RR(iter-1) = max(max(max(abs(RR-sto_RR))));
+ err_TT(iter-1) = max(max(max(abs(TT-sto_TT))));
+ end
end
if occbin_smoother_debug || is_periodic
@@ -391,22 +387,22 @@ if occbin_smoother_debug
end
if (maxiter==iter && is_changed) || is_periodic
- disp('occbin.DSGE_smoother: smoother did not converge.')
- fprintf('occbin.DSGE_smoother: The algorithm did not reach a fixed point for the smoothed regimes.\n')
+ disp_verbose('occbin.DSGE_smoother: smoother did not converge.',options_.verbosity)
+ disp_verbose('occbin.DSGE_smoother: The algorithm did not reach a fixed point for the smoothed regimes.',options_.verbosity)
if is_periodic
oo_.occbin.smoother.error_flag=0;
- fprintf('occbin.DSGE_smoother: For the periods indicated above, regimes loops between the "regime_" and the "regime_new_" pattern displayed above.\n')
- fprintf('occbin.DSGE_smoother: We provide smoothed shocks consistent with "regime_" in oo_.\n')
+ disp_verbose('occbin.DSGE_smoother: For the periods indicated above, regimes loops between the "regime_" and the "regime_new_" pattern displayed above.',options_.verbosity)
+ disp_verbose('occbin.DSGE_smoother: We provide smoothed shocks consistent with "regime_" in oo_.',options_.verbosity)
else
- fprintf('occbin.DSGE_smoother: The respective fields in oo_ will be left empty.\n')
+ disp_verbose('occbin.DSGE_smoother: The respective fields in oo_ will be left empty.',options_.verbosity)
oo_.occbin.smoother=[];
oo_.occbin.smoother.error_flag=322;
end
else
- disp('occbin.DSGE_smoother: smoother converged.')
+ disp_verbose('occbin.DSGE_smoother: smoother converged.',options_.verbosity)
oo_.occbin.smoother.error_flag=0;
if occbin_smoother_fast && is_changed_start
- disp('occbin.DSGE_smoother: WARNING: fast algo is used, regime duration was not forced to converge')
+ disp_verbose('occbin.DSGE_smoother: WARNING: fast algo is used, regime duration was not forced to converge',options_.verbosity)
end
end
if (~is_changed || occbin_smoother_debug) && nargin==12
@@ -484,9 +480,9 @@ if (~is_changed || occbin_smoother_debug) && nargin==12
fprintf(fidTeX,'\\label{Fig:smoothedshocks_occbin:%s}\n',int2str(ifig));
fprintf(fidTeX,'\\end{figure}\n');
fprintf(fidTeX,' \n');
+ end
end
end
- end
end
if mod(j1,9)~=0 && j==M_.exo_nbr
diff --git a/matlab/estimation/dynare_estimation_init.m b/matlab/estimation/dynare_estimation_init.m
index fe246cbc26c4bf4e9e51a47d6032a2153e1325e9..ffa982c0c29a434a867a655288db6e00f0d217b1 100644
--- a/matlab/estimation/dynare_estimation_init.m
+++ b/matlab/estimation/dynare_estimation_init.m
@@ -120,13 +120,10 @@ end
options_=select_qz_criterium_value(options_);
% Set options related to filtered variables.
-if isequal(options_.filtered_vars,0) && ~isempty(options_.filter_step_ahead)
+if isequal(options_.filtered_vars,0) && ~isempty(options_.filter_step_ahead) %require filtered var because filter_step_ahead was set
options_.filtered_vars = 1;
end
-if ~isequal(options_.filtered_vars,0) && isempty(options_.filter_step_ahead)
- options_.filter_step_ahead = 1;
-end
-if ~isequal(options_.filtered_vars,0) && isequal(options_.filter_step_ahead,0)
+if ~isequal(options_.filtered_vars,0) && (isempty(options_.filter_step_ahead) || isequal(options_.filter_step_ahead,0)) %require filter_step_ahead because filtered_vars was requested
options_.filter_step_ahead = 1;
end
if ~isequal(options_.filter_step_ahead,0)
@@ -533,12 +530,17 @@ if (options_.occbin.likelihood.status && options_.occbin.likelihood.inversion_fi
error('IVF-filter: an observable is mapped to a zero variance shock.')
end
end
+ if ~isequal(M_.H,0)
+ error('IVF-filter: Measurement erros are not allowed with the inversion filter.')
+ end
end
if options_.occbin.smoother.status && options_.occbin.smoother.inversion_filter
if ~isempty(options_.nk)
- fprintf('dynare_estimation_init: the inversion filter does not support filter_step_ahead. Disabling the option.\n')
+ fprintf('dynare_estimation_init: the inversion filter does not support filter_step_ahead and filtered_variables. Disabling the option.\n')
options_.nk=[];
+ options_.filter_step_ahead=[];
+ options_.filtered_vars=0;
end
if options_.filter_covariance
fprintf('dynare_estimation_init: the inversion filter does not support filter_covariance. Disabling the option.\n')
diff --git a/matlab/estimation/prior_posterior_statistics.m b/matlab/estimation/prior_posterior_statistics.m
index 36349dc6b872f26b5fad5b60dd1e6b8b39b5dec5..a0fb6128cc73a382992d562163771867004122db 100644
--- a/matlab/estimation/prior_posterior_statistics.m
+++ b/matlab/estimation/prior_posterior_statistics.m
@@ -141,9 +141,9 @@ if options_.filter_covariance
filter_covariance=1;
end
-smoothed_state_uncertainty=0;
+smoothed_state_uncertainty=false;
if options_.smoothed_state_uncertainty
- smoothed_state_uncertainty=1;
+ smoothed_state_uncertainty=true;
end
% Store the variable mandatory for local/remote parallel computing.
@@ -327,6 +327,9 @@ if options_.smoother
varlist,'UpdatedVariables',DirectoryName, ...
'_update',dispString);
if smoothed_state_uncertainty
+ if isfield(oo_,'Smoother') && isfield(oo_.Smoother,'State_uncertainty')
+ oo_.Smoother=rmfield(oo_.Smoother,'State_uncertainty'); %needs to be removed as classical smoother field has a different format
+ end
oo_=pm3(M_,options_,oo_,endo_nbr,endo_nbr,ifil(13),B,'State Uncertainty',...
varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'StateUncertainty',DirectoryName,'_state_uncert',dispString);
diff --git a/matlab/estimation/prior_posterior_statistics_core.m b/matlab/estimation/prior_posterior_statistics_core.m
index 5cf8f364e5ff56e0b1de55f0842b981effb04afa..478dab2e60eda6405d0995d83c786f226f79e138 100644
--- a/matlab/estimation/prior_posterior_statistics_core.m
+++ b/matlab/estimation/prior_posterior_statistics_core.m
@@ -230,11 +230,28 @@ for b=fpar:B
opts_local.occbin.simul.waitbar=0;
opts_local.occbin.smoother.waitbar = 0;
opts_local.occbin.smoother.linear_smoother=false; % speed-up
- [alphahat,etahat,epsilonhat,alphatilde,SteadyState,trend_coeff,aK,~,~,P,~,~,trend_addition,state_uncertainty,oo_,bayestopt_] = ...
- occbin.DSGE_smoother(deep,gend,Y,data_index,missing_value,M_,oo_,opts_local,bayestopt_,estim_params_);
- if oo_.occbin.smoother.error_flag(1)
- message=get_error_message(oo_.occbin.smoother.error_flag,opts_local);
- fprintf('\nprior_posterior_statistics: One of the draws failed with the error:\n%s\n',message)
+ if options_.occbin.smoother.inversion_filter
+ dataset_.data=Y';
+ [~, info, ~, ~, ~, ~, ~, ~, ~, ~, oo_.dr, alphahat, etahat] = ...
+ occbin.IVF_posterior(deep,dataset_,[],options_,M_,estim_params_,bayestopt_,prior_bounds(bayestopt_,options_.prior_trunc),oo_.dr, oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
+ if info(1)
+ message=get_error_message(info,opts_local);
+ fprintf('\nprior_posterior_statistics: IVF smoother failed for one of the draws:\n%s\n',message)
+ else
+ alphatilde = alphahat*nan;
+ SteadyState=oo_.dr.ys;
+ trend_coeff = zeros(length(options_.varobs_id),1);
+ trend_addition=zeros(options_.number_of_observed_variables,gend);
+ end
+ %epsilonhat not available as no measurement error allowed
+ else
+ opts_local.verbosity=0;
+ [alphahat,etahat,epsilonhat,alphatilde,SteadyState,trend_coeff,aK,~,~,P,~,~,trend_addition,state_uncertainty,oo_,bayestopt_] = ...
+ occbin.DSGE_smoother(deep,gend,Y,data_index,missing_value,M_,oo_,opts_local,bayestopt_,estim_params_);
+ if oo_.occbin.smoother.error_flag(1)
+ message=get_error_message(oo_.occbin.smoother.error_flag,opts_local);
+ fprintf('\nprior_posterior_statistics: One of the draws failed with the error:\n%s\n',message)
+ end
end
else
[alphahat,etahat,epsilonhat,alphatilde,SteadyState,trend_coeff,aK,~,~,P,~,~,trend_addition,state_uncertainty,oo_,bayestopt_] = ...
diff --git a/matlab/kalman/DsgeSmoother.m b/matlab/kalman/DsgeSmoother.m
index 7a9ecd279a9bc0714c53f34a76daf981e491e8ea..7784ab7ae7c4e072246ca5feb2ef712f2d9e48df 100644
--- a/matlab/kalman/DsgeSmoother.m
+++ b/matlab/kalman/DsgeSmoother.m
@@ -1,5 +1,5 @@
-function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,decomp,trend_addition,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0,d] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,varargin)
-% [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,decomp,trend_addition,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0,d] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,varargin)
+function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,decomp,trend_addition,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0,d,info] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,varargin)
+% [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,decomp,trend_addition,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0,d,info] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,varargin)
% Estimation of the smoothed variables and innovations.
%
% INPUTS
@@ -35,7 +35,14 @@ function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,de
% about the smoothed state (decision-rule order)
% o oo_ [structure] storing the results
% o bayestopt_ [structure] describing the priors
-%
+% o alphahat0 [double] (m*1) matrix, smoothed endogenous variables
+% (a_{0}) for initial period from PKF
+% o state_uncertainty0 [double] (K,K) matrix storing the uncertainty about
+% the smoothed state for the initial
+% period from the PKF
+% o d [integer] number of diffuse periods
+% o info [1 by 4 double] error code and penalty
+
% Notes:
% m: number of endogenous variables (M_.endo_nbr)
% T: number of Time periods (options_.nobs)
@@ -74,9 +81,9 @@ function [alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T,R,P,PK,de
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
-alphahat = [];
-etahat = [];
-epsilonhat = [];
+alphahat = [];
+etahat = [];
+epsilonhat = [];
ahat = [];
SteadyState = [];
trend_coeff = [];
@@ -86,8 +93,11 @@ R = [];
P = [];
PK = [];
decomp = [];
-vobs = length(options_.varobs);
+vobs = length(options_.varobs);
smpl = size(Y,2);
+alphahat0 = [];
+state_uncertainty0 =[];
+d = 0;
if ~isempty(xparam1) %not calibrated model
M_ = set_all_parameters(xparam1,estim_params_,M_);
diff --git a/tests/occbin/filter/NKM.mod b/tests/occbin/filter/NKM.mod
index a54142f94ad7f1c928af1183ac1c8b716bd12ba1..96f0250fe4720a500ebd3f909bd91c53e86d5120 100644
--- a/tests/occbin/filter/NKM.mod
+++ b/tests/occbin/filter/NKM.mod
@@ -389,13 +389,22 @@ varobs yg inom pi;
else
disp('smoother redux successfully recovers full k-step ahead variables results!')
end
+
+ //run PKF with MCMC
+ options_.smoother_redux=false;
+ estimation(
+ datafile=dataobsfile2, mode_file=NKM_mh_mode_saved,
+ mode_compute=0, nobs=120, first_obs=1,
+ mh_replic=50, plot_priors=0, smoother,
+ consider_all_endogenous,heteroskedastic_filter,filter_step_ahead=[1:8],smoothed_state_uncertainty);
+
// use inversion filter (note that IF provides smoother together with likelihood)
occbin_setup(likelihood_inversion_filter,smoother_inversion_filter);
estimation(
datafile=dataobsfile2, mode_file=NKM_mh_mode_saved,
mode_compute=0, nobs=120, first_obs=1,
- mh_replic=0, plot_priors=0, smoother,
+ mh_replic=50, plot_priors=0, smoother,
consider_all_endogenous,heteroskedastic_filter,filter_step_ahead=[1:8],smoothed_state_uncertainty);
// show initial condition effect of IF
@@ -417,6 +426,12 @@ varobs yg inom pi;
legend('PKF','IF')
occbin_write_regimes(smoother);
+ estimation(
+ datafile=dataobsfile2, mode_file=NKM_mh_mode_saved,
+ mode_compute=0, nobs=120, first_obs=1,
+ mh_replic=50, plot_priors=0, smoother,
+ consider_all_endogenous,heteroskedastic_filter,filter_step_ahead=[1:8],smoothed_state_uncertainty);
+
write_latex_dynamic_model;
collect_latex_files;
[status, cmdout]=system(['pdflatex -halt-on-error -interaction=nonstopmode ' M_.fname '_TeX_binder.tex']);