diff --git a/matlab/check_posterior_sampler_options.m b/matlab/check_posterior_sampler_options.m
index 5ed6cb80e..c296cc94a 100644
--- a/matlab/check_posterior_sampler_options.m
+++ b/matlab/check_posterior_sampler_options.m
@@ -1,193 +1,405 @@
-function [posterior_sampler_options, options_] = check_posterior_sampler_options(posterior_sampler_options, options_)
-
-% function posterior_sampler_options = check_posterior_sampler_options(posterior_sampler_options, options_)
-% initialization of posterior samplers
-%
-% INPUTS
-% posterior_sampler_options: posterior sampler options
-% options_: structure storing the options
-
-% OUTPUTS
-% posterior_sampler_options: checked posterior sampler options
-%
-% SPECIAL REQUIREMENTS
-% none
-
-% Copyright (C) 2015 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 .
-
-posterior_sampler_options.posterior_sampling_method = options_.posterior_sampling_method;
-posterior_sampler_options.proposal_distribution = options_.proposal_distribution;
-bounds = posterior_sampler_options.bounds;
-invhess = posterior_sampler_options.invhess;
-
-% here are all samplers requiring a proposal distribution
-if ~strcmp(posterior_sampler_options.posterior_sampling_method,'slice')
- if ~options_.cova_compute
- error('I Cannot start the MCMC because the Hessian of the posterior kernel at the mode was not computed.')
- end
- if options_.load_mh_file && options_.use_mh_covariance_matrix,
- [junk, invhess] = compute_mh_covariance_matrix;
- posterior_sampler_options.invhess = invhess;
- end
- posterior_sampler_options.parallel_bar_refresh_rate=50;
- posterior_sampler_options.serial_bar_refresh_rate=3;
- posterior_sampler_options.parallel_bar_title='MH';
- posterior_sampler_options.serial_bar_title='Metropolis-Hastings';
- posterior_sampler_options.save_tmp_file=1;
-end
-
-
-% check specific options for slice sampler
-if strcmp(options_.posterior_sampling_method,'slice')
- % by default, slice sampler should trigger
- % mode_compute=0 and
- % mh_replic=100 (much smaller than the default mh_replic=20000 of RWMH)
- % moreover slice must be associated to:
- % options_.mh_posterior_mode_estimation = 0;
- % this is done below, but perhaps preprocessing should do this?
-
- posterior_sampler_options.parallel_bar_refresh_rate=1;
- posterior_sampler_options.serial_bar_refresh_rate=1;
- posterior_sampler_options.parallel_bar_title='SLICE';
- posterior_sampler_options.serial_bar_title='SLICE';
- posterior_sampler_options.save_tmp_file=1;
- posterior_sampler_options = set_default_option(posterior_sampler_options,'rotated',0);
- posterior_sampler_options = set_default_option(posterior_sampler_options,'slice_initialize_with_mode',0);
- posterior_sampler_options = set_default_option(posterior_sampler_options,'use_mh_covariance_matrix',0);
- posterior_sampler_options = set_default_option(posterior_sampler_options,'WR',[]);
- if ~isfield(posterior_sampler_options,'mode'),
- posterior_sampler_options.mode = [];
- else % multimodal case
- posterior_sampler_options.rotated = 1;
- end
- posterior_sampler_options = set_default_option(posterior_sampler_options,'mode_files',[]);
- posterior_sampler_options = set_default_option(posterior_sampler_options,'W1',0.8*(bounds.ub-bounds.lb));
-
- if ~isempty(options_.optim_opt)
- options_list = read_key_value_string(options_.optim_opt);
- for i=1:rows(options_list)
- switch options_list{i,1}
- case 'rotated'
- % triggers rotated slice iterations using a covariance
- % matrix from initial burn-in iterations
- % must be associated with:
- % or
- % default = 0
- posterior_sampler_options.rotated = options_list{i,2};
-
- case 'mode'
- % for multimodal posteriors, provide the list of modes as a
- % matrix, ordered by column, i.e. [x1 x2 x3] for three
- % modes x1 x2 x3
- % MR note: not sure this is possible with the
- % read_key_value_string ???
- % if this is not possible does to job in any case
- % This will automatically trigger
- % default = []
- tmp_mode = options_list{i,2};
- for j=1:size(tmp_mode,2),
- posterior_sampler_options.mode(j).m = tmp_mode(:,j);
- end
-
- case 'mode_files'
- % for multimodal posteriors provide a list of mode files,
- % one per mode. With this info, the code will automatically
- % set the option. The mode files need only to
- % contain the xparam1 variable.
- % This will automatically trigger
- % default = []
- posterior_sampler_options.mode_files = options_list{i,2};
-
- case 'slice_initialize_with_mode'
- % the default for slice is to set mode_compute = 0 in the
- % preprocessor and start the chain(s) from a random location in the prior.
- % This option first runs the optimizer and then starts the
- % chain from the mode. Associated with optios , it will
- % use invhess from the mode to perform rotated slice
- % iterations.
- % default = 0
- posterior_sampler_options.slice_initialize_with_mode = options_list{i,2};
-
- case 'initial_step_size'
- % sets the initial size of the interval in the STEPPING-OUT PROCEDURE
- % the initial_step_size must be a real number in [0, 1],
- % and it sets the size as a proportion of the prior bounds,
- % i.e. the size will be initial_step_size*(UB-LB)
- % slice sampler requires prior_truncation > 0!
- % default = 0.8
- posterior_sampler_options.W1 = options_list{i,2}*(bounds.ub-bounds.lb);
-
- case 'use_mh_covariance_matrix'
- % in association with indicates to use the
- % covariance matrix from previous iterations to define the
- % rotated slice
- % default = 0
- posterior_sampler_options.use_mh_covariance_matrix = options_list{i,2};
-
- otherwise
- warning(['slice_sampler: Unknown option (' options_list{i,1} ')!'])
- end
- end
- end
- if options_.load_mh_file,
- posterior_sampler_options.slice_initialize_with_mode = 0;
- else
- if ~posterior_sampler_options.slice_initialize_with_mode,
- posterior_sampler_options.invhess=[];
- end
- end
-
- if posterior_sampler_options.rotated,
- if isempty(posterior_sampler_options.mode_files) && isempty(posterior_sampler_options.mode), % rotated unimodal
-
- if ~options_.cova_compute && ~(options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix)
- skipline()
- disp('I cannot start rotated slice sampler because')
- disp('there is no previous MCMC to load ')
- disp('or the Hessian at the mode is not computed.')
- error('Rotated slice cannot start')
- end
- if isempty(invhess)
- error('oops! This error should not occur, please contact developers.')
- end
- if options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix,
- [junk, invhess] = compute_mh_covariance_matrix;
- posterior_sampler_options.invhess = invhess;
- end
- [V1 D]=eig(invhess);
- posterior_sampler_options.V1=V1;
- posterior_sampler_options.WR=sqrt(diag(D))*3;
- end
- end
-
- if ~isempty(posterior_sampler_options.mode_files), % multimodal case
- modes = posterior_sampler_options.mode_files; % these can be also mean files from previous parallel slice chains
- for j=1:length( modes ),
- load(modes{j}, 'xparam1')
- mode(j).m=xparam1;
- end
- posterior_sampler_options.mode = mode;
- posterior_sampler_options.rotated = 1;
- posterior_sampler_options.WR=[];
- end
- options_.mh_posterior_mode_estimation = 0;
-
-end
-
-
-
+function [posterior_sampler_options, options_] = check_posterior_sampler_options(posterior_sampler_options, options_, bounds)
+
+% function [posterior_sampler_options, options_] = check_posterior_sampler_options(posterior_sampler_options, options_, bounds)
+% initialization of posterior samplers
+%
+% INPUTS
+% posterior_sampler_options: posterior sampler options
+% options_: structure storing the options
+
+% OUTPUTS
+% posterior_sampler_options: checked posterior sampler options
+%
+% SPECIAL REQUIREMENTS
+% none
+
+% Copyright (C) 2015 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 .
+
+
+init=0;
+if isempty(posterior_sampler_options),
+ init=1;
+end
+
+if init,
+ % set default options and user defined options
+ posterior_sampler_options.posterior_sampling_method = options_.posterior_sampler_options.posterior_sampling_method;
+ posterior_sampler_options.bounds = bounds;
+
+ switch posterior_sampler_options.posterior_sampling_method
+
+ case 'random_walk_metropolis_hastings'
+ posterior_sampler_options.parallel_bar_refresh_rate=50;
+ posterior_sampler_options.serial_bar_refresh_rate=3;
+ posterior_sampler_options.parallel_bar_title='RWMH';
+ posterior_sampler_options.serial_bar_title='RW Metropolis-Hastings';
+ posterior_sampler_options.save_tmp_file=1;
+
+ % default options
+ posterior_sampler_options = add_fields_(posterior_sampler_options,options_.posterior_sampler_options.rwmh);
+
+ % user defined options
+ if ~isempty(options_.posterior_sampler_options.sampling_opt)
+ options_list = read_key_value_string(options_.posterior_sampler_options.sampling_opt);
+ for i=1:rows(options_list)
+ switch options_list{i,1}
+
+ case 'proposal_distribution'
+ if ~(strcmpi(options_list{i,2}, 'rand_multivariate_student') || ...
+ strcmpi(options_list{i,2}, 'rand_multivariate_normal'))
+ error(['initial_estimation_checks:: the proposal_distribution option to estimation takes either ' ...
+ 'rand_multivariate_student or rand_multivariate_normal as options']);
+ else
+ posterior_sampler_options.proposal_distribution=options_list{i,2};
+ end
+
+
+ case 'student_degrees_of_freedom'
+ if options_list{i,2} <= 0
+ error('initial_estimation_checks:: the student_degrees_of_freedom takes a positive integer argument');
+ else
+ posterior_sampler_options.student_degrees_of_freedom=options_list{i,2};
+ end
+
+ case 'use_mh_covariance_matrix'
+ % indicates to use the covariance matrix from previous iterations to
+ % define the covariance of the proposal distribution
+ % default = 0
+ posterior_sampler_options.use_mh_covariance_matrix = options_list{i,2};
+ options_.use_mh_covariance_matrix = options_list{i,2};
+
+ otherwise
+ warning(['imh_sampler: Unknown option (' options_list{i,1} ')!'])
+ end
+ end
+ end
+
+ case 'tailored_random_block_metropolis_hastings'
+ posterior_sampler_options.parallel_bar_refresh_rate=5;
+ posterior_sampler_options.serial_bar_refresh_rate=1;
+ posterior_sampler_options.parallel_bar_title='TaRB-MH';
+ posterior_sampler_options.serial_bar_title='TaRB Metropolis-Hastings';
+ posterior_sampler_options.save_tmp_file=1;
+
+ % default options
+ posterior_sampler_options = add_fields_(posterior_sampler_options,options_.posterior_sampler_options.tarb);
+
+ % user defined options
+ if ~isempty(options_.posterior_sampler_options.sampling_opt)
+ options_list = read_key_value_string(options_.posterior_sampler_options.sampling_opt);
+ for i=1:rows(options_list)
+
+ switch options_list{i,1}
+
+ case 'proposal_distribution'
+ if ~(strcmpi(options_list{i,2}, 'rand_multivariate_student') || ...
+ strcmpi(options_list{i,2}, 'rand_multivariate_normal'))
+ error(['initial_estimation_checks:: the proposal_distribution option to estimation takes either ' ...
+ 'rand_multivariate_student or rand_multivariate_normal as options']);
+ else
+ posterior_sampler_options.proposal_distribution=options_list{i,2};
+ end
+
+
+ case 'student_degrees_of_freedom'
+ if options_list{i,2} <= 0
+ error('initial_estimation_checks:: the student_degrees_of_freedom takes a positive integer argument');
+ else
+ posterior_sampler_options.student_degrees_of_freedom=options_list{i,2};
+ end
+
+ case 'mode_compute'
+ posterior_sampler_options.mode_compute=options_list{i,2};
+
+ case 'new_block_probability'
+ if options_list{i,2}<0 || options_list{i,2}>1
+ error('check_posterior_sampler_options:: The tarb new_block_probability must be between 0 and 1!')
+ else
+ posterior_sampler_options.new_block_probability=options_list{i,2};
+ end
+
+ otherwise
+
+ end
+
+ end
+
+ end
+
+ case 'independent_metropolis_hastings'
+ posterior_sampler_options.parallel_bar_refresh_rate=50;
+ posterior_sampler_options.serial_bar_refresh_rate=3;
+ posterior_sampler_options.parallel_bar_title='IMH';
+ posterior_sampler_options.serial_bar_title='Ind. Metropolis-Hastings';
+ posterior_sampler_options.save_tmp_file=1;
+
+ % default options
+ posterior_sampler_options = add_fields_(posterior_sampler_options,options_.posterior_sampler_options.imh);
+
+ % user defined options
+ if ~isempty(options_.posterior_sampler_options.sampling_opt)
+ options_list = read_key_value_string(options_.posterior_sampler_options.sampling_opt);
+ for i=1:rows(options_list)
+ switch options_list{i,1}
+
+ case 'proposal_distribution'
+ if ~(strcmpi(options_list{i,2}, 'rand_multivariate_student') || ...
+ strcmpi(options_list{i,2}, 'rand_multivariate_normal'))
+ error(['initial_estimation_checks:: the proposal_distribution option to estimation takes either ' ...
+ 'rand_multivariate_student or rand_multivariate_normal as options']);
+ else
+ posterior_sampler_options.proposal_distribution=options_list{i,2};
+ end
+
+
+ case 'student_degrees_of_freedom'
+ if options_list{i,2} <= 0
+ error('initial_estimation_checks:: the student_degrees_of_freedom takes a positive integer argument');
+ else
+ posterior_sampler_options.student_degrees_of_freedom=options_list{i,2};
+ end
+
+ case 'use_mh_covariance_matrix'
+ % indicates to use the covariance matrix from previous iterations to
+ % define the covariance of the proposal distribution
+ % default = 0
+ posterior_sampler_options.use_mh_covariance_matrix = options_list{i,2};
+ options_.use_mh_covariance_matrix = options_list{i,2};
+
+ otherwise
+ warning(['imh_sampler: Unknown option (' options_list{i,1} ')!'])
+ end
+ end
+ end
+
+
+ case 'slice'
+ posterior_sampler_options.parallel_bar_refresh_rate=1;
+ posterior_sampler_options.serial_bar_refresh_rate=1;
+ posterior_sampler_options.parallel_bar_title='SLICE';
+ posterior_sampler_options.serial_bar_title='SLICE';
+ posterior_sampler_options.save_tmp_file=1;
+
+ % default options
+ posterior_sampler_options = add_fields_(posterior_sampler_options,options_.posterior_sampler_options.slice);
+
+ % user defined options
+ if ~isempty(options_.posterior_sampler_options.sampling_opt)
+ options_list = read_key_value_string(options_.posterior_sampler_options.sampling_opt);
+ for i=1:rows(options_list)
+ switch options_list{i,1}
+ case 'rotated'
+ % triggers rotated slice iterations using a covariance
+ % matrix from initial burn-in iterations
+ % must be associated with:
+ % or
+ % default = 0
+ posterior_sampler_options.rotated = options_list{i,2};
+
+ case 'mode'
+ % for multimodal posteriors, provide the list of modes as a
+ % matrix, ordered by column, i.e. [x1 x2 x3] for three
+ % modes x1 x2 x3
+ % MR note: not sure this is possible with the
+ % read_key_value_string ???
+ % if this is not possible does to job in any case
+ % This will automatically trigger
+ % default = []
+ tmp_mode = options_list{i,2};
+ for j=1:size(tmp_mode,2),
+ posterior_sampler_options.mode(j).m = tmp_mode(:,j);
+ end
+
+ case 'mode_files'
+ % for multimodal posteriors provide a list of mode files,
+ % one per mode. With this info, the code will automatically
+ % set the option. The mode files need only to
+ % contain the xparam1 variable.
+ % This will automatically trigger
+ % default = []
+ posterior_sampler_options.mode_files = options_list{i,2};
+
+ case 'slice_initialize_with_mode'
+ % the default for slice is to set mode_compute = 0 in the
+ % preprocessor and start the chain(s) from a random location in the prior.
+ % This option first runs the optimizer and then starts the
+ % chain from the mode. Associated with optios , it will
+ % use invhess from the mode to perform rotated slice
+ % iterations.
+ % default = 0
+ posterior_sampler_options.slice_initialize_with_mode = options_list{i,2};
+
+ case 'initial_step_size'
+ % sets the initial size of the interval in the STEPPING-OUT PROCEDURE
+ % the initial_step_size must be a real number in [0, 1],
+ % and it sets the size as a proportion of the prior bounds,
+ % i.e. the size will be initial_step_size*(UB-LB)
+ % slice sampler requires prior_truncation > 0!
+ % default = 0.8
+ if options_list{i,2}<=0 || options_list{i,2}>=1
+ error('check_posterior_sampler_options:: slice initial_step_size must be between 0 and 1')
+ else
+ posterior_sampler_options.initial_step_size=options_list{i,2};
+ end
+ case 'use_mh_covariance_matrix'
+ % in association with indicates to use the
+ % covariance matrix from previous iterations to define the
+ % rotated slice
+ % default = 0
+ posterior_sampler_options.use_mh_covariance_matrix = options_list{i,2};
+ options_.use_mh_covariance_matrix = options_list{i,2};
+
+ otherwise
+ warning(['slice_sampler: Unknown option (' options_list{i,1} ')!'])
+ end
+ end
+ end
+
+ % slice posterior sampler does not require mode or hessian to run
+ % needs to be set to 1 to skip parts in dynare_estimation_1.m
+ % requiring posterior maximization/calibrated smoother before MCMC
+ options_.mh_posterior_mode_estimation=1;
+
+ if ~ posterior_sampler_options.slice_initialize_with_mode
+ % by default, slice sampler should trigger
+ % mode_compute=0 and
+ % mh_replic=100 (much smaller than the default mh_replic=20000 of RWMH)
+ options_.mode_compute = 0;
+ options_.cova_compute = 0;
+ else
+ if (isequal(options_.mode_compute,0) && isempty(options_.mode_file) )
+ skipline()
+ disp('check_posterior_sampler_options:: You have specified the option "slice_initialize_with_mode"')
+ disp('check_posterior_sampler_options:: to initialize the slice sampler using mode information')
+ disp('check_posterior_sampler_options:: but no mode file nor posterior maximization is selected,')
+ error('check_posterior_sampler_options:: The option "slice_initialize_with_mode" is inconsistent with mode_compute=0 or empty mode_file.')
+ else
+ options_.mh_posterior_mode_estimation=0;
+ end
+ end
+
+ if any(isinf(bounds.lb)) || any(isinf(bounds.ub)),
+ skipline()
+ disp('some priors are unbounded and prior_trunc is set to zero')
+ error('The option "slice" is inconsistent with prior_trunc=0.')
+ end
+
+ % moreover slice must be associated to:
+ % options_.mh_posterior_mode_estimation = 0;
+ % this is done below, but perhaps preprocessing should do this?
+
+ if ~isfield(posterior_sampler_options,'mode'),
+ posterior_sampler_options.mode = [];
+ else % multimodal case
+ posterior_sampler_options.rotated = 1;
+ end
+ % posterior_sampler_options = set_default_option(posterior_sampler_options,'mode_files',[]);
+
+
+ posterior_sampler_options.W1=posterior_sampler_options.initial_step_size*(bounds.ub-bounds.lb);
+ if options_.load_mh_file,
+ posterior_sampler_options.slice_initialize_with_mode = 0;
+ else
+ if ~posterior_sampler_options.slice_initialize_with_mode,
+ posterior_sampler_options.invhess=[];
+ end
+ end
+
+ if ~isempty(posterior_sampler_options.mode_files), % multimodal case
+ modes = posterior_sampler_options.mode_files; % these can be also mean files from previous parallel slice chains
+ for j=1:length( modes ),
+ load(modes{j}, 'xparam1')
+ mode(j).m=xparam1;
+ end
+ posterior_sampler_options.mode = mode;
+ posterior_sampler_options.rotated = 1;
+ posterior_sampler_options.WR=[];
+ end
+
+ otherwise
+
+ end
+
+ return
+end
+
+% here are all samplers requiring a proposal distribution
+if ~strcmp(posterior_sampler_options.posterior_sampling_method,'slice')
+ if ~options_.cova_compute && ~(options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix)
+ skipline()
+ disp('check_posterior_sampler_options:: I cannot start the MCMC because the Hessian of the posterior kernel at the mode was not computed')
+ disp('check_posterior_sampler_options:: or there is no previous MCMC to load ')
+ error('check_posterior_sampler_options:: MCMC cannot start')
+ end
+end
+
+if options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix,
+ [junk, invhess] = compute_mh_covariance_matrix;
+ posterior_sampler_options.invhess = invhess;
+end
+
+
+
+% check specific options for slice sampler
+if strcmp(posterior_sampler_options.posterior_sampling_method,'slice')
+ invhess = posterior_sampler_options.invhess;
+
+ if posterior_sampler_options.rotated,
+ if isempty(posterior_sampler_options.mode_files) && isempty(posterior_sampler_options.mode), % rotated unimodal
+
+ if ~options_.cova_compute && ~(options_.load_mh_file && posterior_sampler_options.use_mh_covariance_matrix)
+ skipline()
+ disp('check_posterior_sampler_options:: I cannot start rotated slice sampler because')
+ disp('check_posterior_sampler_options:: there is no previous MCMC to load ')
+ disp('check_posterior_sampler_options:: or the Hessian at the mode is not computed.')
+ error('check_posterior_sampler_options:: Rotated slice cannot start')
+ end
+ if isempty(invhess)
+ error('check_posterior_sampler_options:: This error should not occur, please contact developers.')
+ end
+ % % % if options_.load_mh_file && options_.use_mh_covariance_matrix,
+ % % % [junk, invhess] = compute_mh_covariance_matrix;
+ % % % posterior_sampler_options.invhess = invhess;
+ % % % end
+ [V1 D]=eig(invhess);
+ posterior_sampler_options.V1=V1;
+ posterior_sampler_options.WR=sqrt(diag(D))*3;
+ end
+ end
+
+ % needs to be re-set to zero otherwise posterior analysis is filtered
+ % out in dynare_estimation_1.m
+ options_.mh_posterior_mode_estimation = 0;
+
+
+else
+
+
+end
+
+return
+
+function posterior_sampler_options = add_fields_(posterior_sampler_options, sampler_options)
+
+fnam = fieldnames(sampler_options);
+for j=1:length(fnam)
+ posterior_sampler_options.(fnam{j}) = sampler_options.(fnam{j});
+end
+
+
+
diff --git a/matlab/dynare_estimation_1.m b/matlab/dynare_estimation_1.m
index 61f8d57cc..14ab8e34d 100644
--- a/matlab/dynare_estimation_1.m
+++ b/matlab/dynare_estimation_1.m
@@ -420,15 +420,16 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
if options_.mh_replic
ana_deriv_old = options_.analytic_derivation;
options_.analytic_derivation = 0;
- posterior_sampler_options = options_.posterior_sampler_options;
- posterior_sampler_options.bounds = bounds;
+ posterior_sampler_options = options_.posterior_sampler_options.current_options;
posterior_sampler_options.invhess = invhess;
[posterior_sampler_options, options_] = check_posterior_sampler_options(posterior_sampler_options, options_);
- if strcmpi(options_.posterior_sampling_method,'adaptive_metropolis_hastings'), % keep old form only for this ...
+ % store current options in global
+ options_.posterior_sampler_options.current_options = posterior_sampler_options;
+ if strcmpi(options_.posterior_sampler_options.posterior_sampling_method,'adaptive_metropolis_hastings'), % keep old form only for this ...
invhess = posterior_sampler_options.invhess;
feval(options_.posterior_sampling_method,objective_function,options_.proposal_distribution,xparam1,invhess,bounds,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,oo_);
else
- posterior_sampler(objective_function,options_.proposal_distribution,xparam1,posterior_sampler_options,bounds,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,oo_);
+ posterior_sampler(objective_function,posterior_sampler_options.proposal_distribution,xparam1,posterior_sampler_options,bounds,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,oo_);
end
options_.analytic_derivation = ana_deriv_old;
end
diff --git a/matlab/dynare_estimation_init.m b/matlab/dynare_estimation_init.m
index c1b445709..9de2de7e6 100644
--- a/matlab/dynare_estimation_init.m
+++ b/matlab/dynare_estimation_init.m
@@ -591,28 +591,7 @@ else
end
-% slice posterior sampler does not require mode or hessian to run
-if strcmp(options_.posterior_sampling_method,'slice')
- options_.mh_posterior_mode_estimation=1;
- options_.posterior_sampler_options = set_default_option(options_.posterior_sampler_options,'slice_initialize_with_mode',0);
-
- if options_.posterior_sampler_options.slice_initialize_with_mode
- if (isequal(options_.mode_compute,0) && isempty(options_.mode_file) )
- skipline()
- disp('You have specified the option "slice_initialize_with_mode"')
- disp('to initialize the slice sampler using mode information')
- disp('but no mode file nor posterior maximization is selected,')
- error('The option "slice_initialize_with_mode" is inconsistent with mode_compute=0 or empty mode_file.')
- else
- options_.mh_posterior_mode_estimation=0;
- end
- end
-
- if any(isinf(bounds.lb)) || any(isinf(bounds.ub)),
- skipline()
- disp('some priors are unbounded and prior_trunc is set to zero')
- error('The option "slice" is inconsistent with prior_trunc=0.')
- end
-
+if options_.mh_replic
+ [current_options, options_] = check_posterior_sampler_options([], options_, bounds);
+ options_.posterior_sampler_options.current_options = current_options;
end
-
diff --git a/matlab/global_initialization.m b/matlab/global_initialization.m
index 88fa04bfa..6f8f5d1c7 100644
--- a/matlab/global_initialization.m
+++ b/matlab/global_initialization.m
@@ -435,11 +435,6 @@ options_.MCMC_jumping_covariance='hessian';
options_.use_calibration_initialization = 0;
options_.endo_vars_for_moment_computations_in_estimation=[];
-% Tailored Random Block Metropolis-Hastings
-options_.TaRB.use_TaRB = 0;
-options_.TaRB.mode_compute=4;
-options_.TaRB.new_block_probability=0.25; %probability that next parameter belongs to new block
-
% Run optimizer silently
options_.silent_optimizer=0;
@@ -460,15 +455,32 @@ options_.prior_trunc = 1e-10;
options_.smoother = 0;
options_.posterior_max_subsample_draws = 1200;
options_.sub_draws = [];
-options_.use_mh_covariance_matrix = 0;
+% options_.use_mh_covariance_matrix = 0;
options_.gradient_method = 2; %used by csminwel and newrat
options_.gradient_epsilon = 1e-6; %used by csminwel and newrat
-options_.posterior_sampler_options = []; %extended set of options for individual posterior samplers
+options_.posterior_sampler_options.sampling_opt = []; %extended set of options for individual posterior samplers
+% Random Walk Metropolis-Hastings
options_.posterior_sampler_options.posterior_sampling_method = 'random_walk_metropolis_hastings';
options_.posterior_sampler_options.rwmh.proposal_distribution = 'rand_multivariate_normal';
options_.posterior_sampler_options.rwmh.student_degrees_of_freedom = 3;
+options_.posterior_sampler_options.rwmh.use_mh_covariance_matrix=0;
+% Tailored Random Block Metropolis-Hastings
options_.posterior_sampler_options.tarb.proposal_distribution = 'rand_multivariate_normal';
options_.posterior_sampler_options.tarb.student_degrees_of_freedom = 3;
+options_.posterior_sampler_options.tarb.mode_compute=4;
+options_.posterior_sampler_options.tarb.new_block_probability=0.25; %probability that next parameter belongs to new block
+% Slice
+options_.posterior_sampler_options.slice.proposal_distribution = '';
+options_.posterior_sampler_options.slice.rotated=0;
+options_.posterior_sampler_options.slice.slice_initialize_with_mode=0;
+options_.posterior_sampler_options.slice.use_mh_covariance_matrix=0;
+options_.posterior_sampler_options.slice.WR=[];
+options_.posterior_sampler_options.slice.mode_files=[];
+options_.posterior_sampler_options.slice.initial_step_size=0.8;
+% Independent Metropolis-Hastings
+options_.posterior_sampler_options.imh.proposal_distribution = 'rand_multivariate_normal';
+options_.posterior_sampler_options.imh.use_mh_covariance_matrix=0;
+
options_.trace_plot_ma = 200;
options_.mh_autocorrelation_function_size = 30;
options_.plot_priors = 1;
diff --git a/matlab/initial_estimation_checks.m b/matlab/initial_estimation_checks.m
index fbab33579..84adbde0c 100644
--- a/matlab/initial_estimation_checks.m
+++ b/matlab/initial_estimation_checks.m
@@ -53,20 +53,6 @@ if maximum_number_non_missing_observations>length(find(diag(Model.Sigma_e)))+Est
error(['initial_estimation_checks:: Estimation can''t take place because too many shocks have been calibrated with a zero variance!'])
end
-if ~(strcmpi(DynareOptions.proposal_distribution, 'rand_multivariate_student') || ...
- strcmpi(DynareOptions.proposal_distribution, 'rand_multivariate_normal'))
- error(['initial_estimation_checks:: the proposal_distribution option to estimation takes either ' ...
- 'rand_multivariate_student or rand_multivariate_normal as options']);
-end
-
-if DynareOptions.student_degrees_of_freedom <= 0
- error('initial_estimation_checks:: the student_degrees_of_freedom takes a positive integer argument');
-end
-
-if DynareOptions.TaRB.use_TaRB && (DynareOptions.TaRB.new_block_probability<0 || DynareOptions.TaRB.new_block_probability>1)
- error(['initial_estimation_checks:: The tarb_new_block_probability must be between 0 and 1!'])
-end
-
if (any(BayesInfo.pshape >0 ) && DynareOptions.mh_replic) && DynareOptions.mh_nblck<1
error(['initial_estimation_checks:: Bayesian estimation cannot be conducted with mh_nblocks=0.'])
end
diff --git a/matlab/posterior_sampler.m b/matlab/posterior_sampler.m
index 57924ddb6..decea5d32 100644
--- a/matlab/posterior_sampler.m
+++ b/matlab/posterior_sampler.m
@@ -89,11 +89,6 @@ load_last_mh_history_file(MetropolisFolder, ModelName);
% on many cores). The mandatory variables for local/remote parallel
% computing are stored in the localVars struct.
-if options_.TaRB.use_TaRB
- options_.silent_optimizer=1; %locally set optimizer to silent mode
- sampler_options.posterior_sampling_method='tailored_random_block_metropolis_hastings';
-end
-
localVars = struct('TargetFun', TargetFun, ...
'ProposalFun', ProposalFun, ...
'xparam1', xparam1, ...
@@ -120,6 +115,9 @@ localVars = struct('TargetFun', TargetFun, ...
'oo_', oo_,...
'varargin',[]);
+if strcmp(sampler_options.posterior_sampling_method,'tailored_random_block_metropolis_hastings');
+ localVars.options_.silent_optimizer=1; %locally set optimizer to silent mode
+end
% User doesn't want to use parallel computing, or wants to compute a
% single chain compute Random walk Metropolis-Hastings algorithm sequentially.
@@ -145,11 +143,7 @@ else
end
% from where to get back results
% NamFileOutput(1,:) = {[M_.dname,'/metropolis/'],'*.*'};
- if options_.TaRB.use_TaRB
- [fout, nBlockPerCPU, totCPU] = masterParallel(options_.parallel, fblck, nblck,NamFileInput,'TaRB_metropolis_hastings_core', localVars, globalVars, options_.parallel_info);
- else
- [fout, nBlockPerCPU, totCPU] = masterParallel(options_.parallel, fblck, nblck,NamFileInput,'posterior_sampler_core', localVars, globalVars, options_.parallel_info);
- end
+ [fout, nBlockPerCPU, totCPU] = masterParallel(options_.parallel, fblck, nblck,NamFileInput,'posterior_sampler_core', localVars, globalVars, options_.parallel_info);
for j=1:totCPU,
offset = sum(nBlockPerCPU(1:j-1))+fblck-1;
record.LastLogPost(offset+1:sum(nBlockPerCPU(1:j)))=fout(j).record.LastLogPost(offset+1:sum(nBlockPerCPU(1:j)));
@@ -190,4 +184,4 @@ if max(record.FunctionEvalPerIteration)>1
disp([' Chain ' num2str(i) ': ' num2str(record.FunctionEvalPerIteration(i))])
end
end
-end
\ No newline at end of file
+end
diff --git a/matlab/posterior_sampler_initialization.m b/matlab/posterior_sampler_initialization.m
index 937683bfd..75df703c8 100644
--- a/matlab/posterior_sampler_initialization.m
+++ b/matlab/posterior_sampler_initialization.m
@@ -202,7 +202,7 @@ if ~options_.load_mh_file && ~options_.mh_recover
fprintf(['Estimation::mcmc: Write details about the MCMC... ']);
AnticipatedNumberOfFiles = ceil(nruns(1)/MAX_nruns);
AnticipatedNumberOfLinesInTheLastFile = nruns(1) - (AnticipatedNumberOfFiles-1)*MAX_nruns;
- record.Sampler = options_.posterior_sampling_method;
+ record.Sampler = options_.posterior_sampler_options.posterior_sampling_method;
record.Nblck = NumberOfBlocks;
record.MhDraws = zeros(1,3);
record.MhDraws(1,1) = nruns(1);
@@ -446,4 +446,4 @@ elseif options_.mh_recover
end
write_mh_history_file(MetropolisFolder, ModelName, record);
end
-end
\ No newline at end of file
+end
diff --git a/tests/estimation/fs2000_slice.mod b/tests/estimation/fs2000_slice.mod
index 664819b57..167ce5fda 100644
--- a/tests/estimation/fs2000_slice.mod
+++ b/tests/estimation/fs2000_slice.mod
@@ -1,76 +1,81 @@
-// See fs2000.mod in the examples/ directory for details on the model
-
-var m P c e W R k d n l gy_obs gp_obs y dA;
-varexo e_a e_m;
-
-parameters alp bet gam mst rho psi del;
-
-alp = 0.33;
-bet = 0.99;
-gam = 0.003;
-mst = 1.011;
-rho = 0.7;
-psi = 0.787;
-del = 0.02;
-
-model;
-dA = exp(gam+e_a);
-log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
--P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0;
-W = l/n;
--(psi/(1-psi))*(c*P/(1-n))+l/n = 0;
-R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W;
-1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0;
-c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1);
-P*c = m;
-m-1+d = l;
-e = exp(e_a);
-y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a));
-gy_obs = dA*y/y(-1);
-gp_obs = (P/P(-1))*m(-1)/dA;
-end;
-
-initval;
-k = 6;
-m = mst;
-P = 2.25;
-c = 0.45;
-e = 1;
-W = 4;
-R = 1.02;
-d = 0.85;
-n = 0.19;
-l = 0.86;
-y = 0.6;
-gy_obs = exp(gam);
-gp_obs = exp(-gam);
-dA = exp(gam);
-end;
-
-shocks;
-var e_a; stderr 0.014;
-var e_m; stderr 0.005;
-end;
-
-steady;
-
-check;
-
-estimated_params;
-alp, beta_pdf, 0.356, 0.02;
-bet, beta_pdf, 0.993, 0.002;
-gam, normal_pdf, 0.0085, 0.003;
-mst, normal_pdf, 1.0002, 0.007;
-rho, beta_pdf, 0.129, 0.223;
-psi, beta_pdf, 0.65, 0.05;
-del, beta_pdf, 0.01, 0.005;
-stderr e_a, inv_gamma_pdf, 0.035449, inf;
-stderr e_m, inv_gamma_pdf, 0.008862, inf;
-end;
-
-varobs gp_obs gy_obs;
-
-options_.posterior_sampling_method = 'slice';
-estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=50,mh_nblocks=2,mh_drop=0.2,mode_compute=0,cova_compute=0);
-// continue with rotated slice
-estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=100,mh_nblocks=2,mh_drop=0.5,load_mh_file,mode_compute=0,optim=('rotated',1,'use_mh_covariance_matrix',1));
+// See fs2000.mod in the examples/ directory for details on the model
+
+var m P c e W R k d n l gy_obs gp_obs y dA;
+varexo e_a e_m;
+
+parameters alp bet gam mst rho psi del;
+
+alp = 0.33;
+bet = 0.99;
+gam = 0.003;
+mst = 1.011;
+rho = 0.7;
+psi = 0.787;
+del = 0.02;
+
+model;
+dA = exp(gam+e_a);
+log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
+-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0;
+W = l/n;
+-(psi/(1-psi))*(c*P/(1-n))+l/n = 0;
+R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W;
+1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0;
+c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1);
+P*c = m;
+m-1+d = l;
+e = exp(e_a);
+y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a));
+gy_obs = dA*y/y(-1);
+gp_obs = (P/P(-1))*m(-1)/dA;
+end;
+
+initval;
+k = 6;
+m = mst;
+P = 2.25;
+c = 0.45;
+e = 1;
+W = 4;
+R = 1.02;
+d = 0.85;
+n = 0.19;
+l = 0.86;
+y = 0.6;
+gy_obs = exp(gam);
+gp_obs = exp(-gam);
+dA = exp(gam);
+end;
+
+shocks;
+var e_a; stderr 0.014;
+var e_m; stderr 0.005;
+end;
+
+steady;
+
+check;
+
+estimated_params;
+alp, beta_pdf, 0.356, 0.02;
+bet, beta_pdf, 0.993, 0.002;
+gam, normal_pdf, 0.0085, 0.003;
+mst, normal_pdf, 1.0002, 0.007;
+rho, beta_pdf, 0.129, 0.223;
+psi, beta_pdf, 0.65, 0.05;
+del, beta_pdf, 0.01, 0.005;
+stderr e_a, inv_gamma_pdf, 0.035449, inf;
+stderr e_m, inv_gamma_pdf, 0.008862, inf;
+end;
+
+varobs gp_obs gy_obs;
+
+//options_.posterior_sampling_method = 'slice';
+estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=50,mh_nblocks=2,mh_drop=0.2, //mode_compute=0,cova_compute=0,
+posterior_sampling_method='slice'
+);
+// continue with rotated slice
+estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=100,mh_nblocks=2,mh_drop=0.5,load_mh_file,//mode_compute=0,
+posterior_sampling_method='slice',
+posterior_sampler_options=('rotated',1,'use_mh_covariance_matrix',1)
+);