From 73cdaa3091b109564c55bd538f8769461bf94fa9 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Charybdis=29?=
<stephane.adjemian@univ-lemans.fr>
Date: Fri, 27 Feb 2015 17:12:04 +0100
Subject: [PATCH] Replaced get_prior_info routine by prior routine.
This new routine gives the possibility to obtain a description of the
prior from the command line.
>> prior table
will display a table with the prior mean, prior mode, prior std, prior
bounds and prior HPD interval.
>> prior plot
will plot the prior densities.
>> prior optimize
will trigger the maximization of the prior density and display the
results (often, because of BK conditions or other issues, Dynare is even
not able to get the prior mode).
>> prior moments
will display the moments of the endogenous variables at the prior mode.
>> prior simulate
will run a Monte-Carlo, by sampling from the prior and discarding vector
of parameters such that the steady state does not exist or the BK
conditions are not met, and return an estimate of the effective prior
mass.
---
matlab/cli/prior.m | 129 ++++++++++++++++++++
matlab/dynare_config.m | 1 +
matlab/get_prior_info.m | 253 ----------------------------------------
3 files changed, 130 insertions(+), 253 deletions(-)
create mode 100644 matlab/cli/prior.m
delete mode 100644 matlab/get_prior_info.m
diff --git a/matlab/cli/prior.m b/matlab/cli/prior.m
new file mode 100644
index 0000000000..f490a570b3
--- /dev/null
+++ b/matlab/cli/prior.m
@@ -0,0 +1,129 @@
+function varargout = prior(varargin)
+
+% Computes various prior statistics and display them in the command window.
+%
+% INPUTS
+% 'table', 'moments', 'optimize', 'simulate', 'plot'
+%
+% OUTPUTS
+% none
+%
+% 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 <http://www.gnu.org/licenses/>.
+
+global options_ M_ estim_params_ bayestopt_ oo_ objective_function_penalty_base
+
+donesomething = false;
+
+% Temporarly change qz_criterium option value
+changed_qz_criterium_flag = 0;
+if isempty(options_.qz_criterium)
+ options_.qz_criterium = 1+1e-9;
+ changed_qz_criterium_flag = 1;
+end
+
+M_.dname = M_.fname;
+
+% Temporarly set options_.order equal to one
+order = options_.order;
+options_.order = 1;
+
+[xparam1,estim_params_,bayestopt_,lb,ub,M_] = set_prior(estim_params_,M_,options_);
+
+if ismember('plot', varargin)
+ plot_priors(bayestopt_,M_,estim_params_,options_)
+ donesomething = true;
+end
+
+if ismember('table', varargin)
+ print_table_prior(lb, ub, options_, M_, bayestopt_, estim_params_);
+ donesomething = true;
+end
+
+if ismember('simulate', varargin) % Prior simulations (BK).
+ results = prior_sampler(0,M_,bayestopt_,options_,oo_,estim_params_);
+ % Display prior mass info
+ skipline(2)
+ disp(['Prior mass = ' num2str(results.prior.mass)])
+ disp(['BK indeterminacy share = ' num2str(results.bk.indeterminacy_share)])
+ disp(['BK unstability share = ' num2str(results.bk.unstability_share)])
+ disp(['BK singularity share = ' num2str(results.bk.singularity_share)])
+ disp(['Complex jacobian share = ' num2str(results.jacobian.problem_share)])
+ disp(['mjdgges crash share = ' num2str(results.dll.problem_share)])
+ disp(['Steady state problem share = ' num2str(results.ss.problem_share)])
+ disp(['Complex steady state share = ' num2str(results.ss.complex_share)])
+ disp(['Analytical steady state problem share = ' num2str(results.ass.problem_share)])
+ skipline(2)
+ donesomething = true;
+end
+
+if ismember('optimize', varargin) % Prior optimization.
+ optimize_prior(options_, M_, oo_, bayestopt_, estim_params_);
+ donesomething = true;
+end
+
+if ismember('moments', varargin) % Prior simulations (2nd order moments).
+ % Set estimated parameters to the prior mode...
+ xparam1 = bayestopt_.p5;
+ % ... Except for uniform priors!
+ k = find(isnan(xparam1));
+ xparam1(k) = bayestopt_.p5(k);
+ % Update vector of parameters and covariance matrices
+ M_ = set_all_parameters(xparam1,estim_params_,M_);
+ % Check model.
+ check_model(M_);
+ % Compute state space representation of the model.
+ oo_.dr=set_state_space(oo_.dr, M_, options_);
+ % Solve model
+ [dr,info, M_ ,options_ , oo_] = resol(0, M_ , options_ ,oo_);
+ % Compute and display second order moments
+ disp_th_moments(oo_.dr,[]);
+ skipline(2)
+ donesomething = true;
+end
+
+if changed_qz_criterium_flag
+ options_.qz_criterium = [];
+end
+
+options_.order = order;
+
+if ~donesomething
+ error('prior: Unexpected arguments!')
+end
+
+function format_string = build_format_string(PriorStandardDeviation,LowerBound,UpperBound)
+format_string = ['%s & %s & %6.4f &'];
+if ~isnumeric(PriorStandardDeviation)
+ format_string = [ format_string , ' %s &'];
+else
+ format_string = [ format_string , ' %6.4f &'];
+end
+if ~isnumeric(LowerBound)
+ format_string = [ format_string , ' %s &'];
+else
+ format_string = [ format_string , ' %6.4f &'];
+end
+if ~isnumeric(UpperBound)
+ format_string = [ format_string , ' %s &'];
+else
+ format_string = [ format_string , ' %6.4f &'];
+end
+format_string = [ format_string , ' %6.4f & %6.4f \\\\ \n'];
\ No newline at end of file
diff --git a/matlab/dynare_config.m b/matlab/dynare_config.m
index 24543ddcde..160cd62220 100644
--- a/matlab/dynare_config.m
+++ b/matlab/dynare_config.m
@@ -58,6 +58,7 @@ addpath([dynareroot '/parallel/'])
addpath([dynareroot '/particle/'])
addpath([dynareroot '/gsa/'])
addpath([dynareroot '/ep/'])
+addpath([dynareroot '/cli/'])
addpath([dynareroot '/lmmcp/'])
addpath([dynareroot '/optimization/'])
addpath([dynareroot '/modules/dates/src/'])
diff --git a/matlab/get_prior_info.m b/matlab/get_prior_info.m
deleted file mode 100644
index 6a82442bf2..0000000000
--- a/matlab/get_prior_info.m
+++ /dev/null
@@ -1,253 +0,0 @@
-function results = get_prior_info(info,plt_flag)
-% Computes various prior statistics.
-%
-% INPUTS
-% info [integer] scalar specifying what has to be done.
-%
-% OUTPUTS
-% none
-%
-% SPECIAL REQUIREMENTS
-% none
-
-% Copyright (C) 2009-2012 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 <http://www.gnu.org/licenses/>.
-global options_ M_ estim_params_ oo_ objective_function_penalty_base
-
-if ~nargin
- info = 0;
- plt_flag = 0;
-end
-
-if nargin==1
- plt_flag = 1;
-end
-
-% Initialize returned variable.
-results = [];
-
-changed_qz_criterium_flag = 0;
-if isempty(options_.qz_criterium)
- options_.qz_criterium = 1+1e-9;
- changed_qz_criterium_flag = 1;
-end
-
-M_.dname = M_.fname;
-
-% Temporarly set options_.order equal to one
-order = options_.order;
-options_.order = 1;
-
-[xparam1,estim_params_,bayestopt_,lb,ub,M_] = set_prior(estim_params_,M_,options_);
-
-if plt_flag
- plot_priors(bayestopt_,M_,estim_params_,options_)
-end
-
-PriorNames = { 'Beta' , 'Gamma' , 'Gaussian' , 'Inverted Gamma' , 'Uniform' , 'Inverted Gamma -- 2' };
-
-if size(M_.param_names,1)==size(M_.param_names_tex,1)% All the parameters have a TeX name.
- fidTeX = fopen('priors_data.tex','w+');
- fprintf(fidTeX,'%% TeX-table generated by get_prior_info (Dynare).\n');
- fprintf(fidTeX,'%% Prior Information\n');
- fprintf(fidTeX,['%% ' datestr(now,0)]);
- fprintf(fidTeX,' \n');
- fprintf(fidTeX,' \n');
- fprintf(fidTeX,'\\begin{center}\n');
- fprintf(fidTeX,'\\begin{longtable}{l|ccccccc} \n');
- fprintf(fidTeX,'\\caption{Prior information (parameters)}\\\\\n ');
- fprintf(fidTeX,'\\label{Table:Prior}\\\\\n');
- fprintf(fidTeX,'\\hline\\hline \\\\ \n');
- fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Lower Bound & Upper Bound & LB Untrunc. 80\\%% HPDI & UB Untrunc. 80\\%% HPDI \\\\ \n');
- fprintf(fidTeX,'\\hline \\endfirsthead \n');
- fprintf(fidTeX,'\\caption{(continued)}\\\\\n ');
- fprintf(fidTeX,'\\hline\\hline \\\\ \n');
- fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Lower Bound & Upper Bound & LB Untrunc. 80\\%% HPDI & UB Untrunc. 80\\%% HPDI \\\\ \n');
- fprintf(fidTeX,'\\hline \\endhead \n');
- fprintf(fidTeX,'\\hline \\multicolumn{8}{r}{(Continued on next page)} \\\\ \\hline \\endfoot \n');
- fprintf(fidTeX,'\\hline \\hline \\endlastfoot \n');
-
- % Column 1: a string for the name of the prior distribution.
- % Column 2: the prior mean.
- % Column 3: the prior standard deviation.
- % Column 4: the lower bound of the prior density support.
- % Column 5: the upper bound of the prior density support.
- % Column 6: the lower bound of the interval containing 80% of the prior mass.
- % Column 7: the upper bound of the interval containing 80% of the prior mass.
- prior_trunc_backup = options_.prior_trunc ;
- options_.prior_trunc = (1-options_.prior_interval)/2 ;
- PriorIntervals = prior_bounds(bayestopt_,options_) ;
- options_.prior_trunc = prior_trunc_backup ;
- for i=1:size(bayestopt_.name,1)
- [tmp,TexName] = get_the_name(i,1,M_,estim_params_,options_);
- PriorShape = PriorNames{ bayestopt_.pshape(i) };
- PriorMean = bayestopt_.p1(i);
- PriorStandardDeviation = bayestopt_.p2(i);
- switch bayestopt_.pshape(i)
- case { 1 , 5 }
- LowerBound = bayestopt_.p3(i);
- UpperBound = bayestopt_.p4(i);
- if ~isinf(lb(i))
- LowerBound=max(LowerBound,lb(i));
- end
- if ~isinf(ub(i))
- UpperBound=min(UpperBound,ub(i));
- end
- case { 2 , 4 , 6 }
- LowerBound = bayestopt_.p3(i);
- if ~isinf(lb(i))
- LowerBound=max(LowerBound,lb(i));
- end
- if ~isinf(ub(i))
- UpperBound=ub(i);
- else
- UpperBound = '$\infty$';
- end
- case 3
- if isinf(bayestopt_.p3(i)) && isinf(lb(i))
- LowerBound = '$-\infty$';
- else
- LowerBound = bayestopt_.p3(i);
- if ~isinf(lb(i))
- LowerBound=max(LowerBound,lb(i));
- end
- end
- if isinf(bayestopt_.p4(i)) && isinf(ub(i))
- UpperBound = '$\infty$';
- else
- UpperBound = bayestopt_.p4(i);
- if ~isinf(ub(i))
- UpperBound=min(UpperBound,ub(i));
- end
- end
- otherwise
- error('get_prior_info:: Dynare bug!')
- end
- format_string = build_format_string(PriorStandardDeviation,LowerBound,UpperBound);
- fprintf(fidTeX,format_string, ...
- TexName, ...
- PriorShape, ...
- PriorMean, ...
- PriorStandardDeviation, ...
- LowerBound, ...
- UpperBound, ...
- PriorIntervals.lb(i), ...
- PriorIntervals.ub(i) );
- end
- fprintf(fidTeX,'\\end{longtable}\n ');
- fprintf(fidTeX,'\\end{center}\n');
- fprintf(fidTeX,'%% End of TeX file.\n');
- fclose(fidTeX);
-end
-
-M_.dname = M_.fname;
-
-if info==1% Prior simulations (BK).
- results = prior_sampler(0,M_,bayestopt_,options_,oo_,estim_params_);
- % Display prior mass info
- disp(['Prior mass = ' num2str(results.prior.mass)])
- disp(['BK indeterminacy share = ' num2str(results.bk.indeterminacy_share)])
- disp(['BK unstability share = ' num2str(results.bk.unstability_share)])
- disp(['BK singularity share = ' num2str(results.bk.singularity_share)])
- disp(['Complex jacobian share = ' num2str(results.jacobian.problem_share)])
- disp(['mjdgges crash share = ' num2str(results.dll.problem_share)])
- disp(['Steady state problem share = ' num2str(results.ss.problem_share)])
- disp(['Complex steady state share = ' num2str(results.ss.complex_share)])
- disp(['Analytical steady state problem share = ' num2str(results.ass.problem_share)])
-end
-
-if info==2% Prior optimization.
- % Initialize to the prior mode if possible
- oo_.dr=set_state_space(oo_.dr,M_,options_);
- k = find(~isnan(bayestopt_.p5));
- xparam1(k) = bayestopt_.p5(k);
- % Pertubation of the initial condition.
- look_for_admissible_initial_condition = 1;
- scale = 1.0;
- iter = 0;
- while look_for_admissible_initial_condition
- xinit = xparam1+scale*randn(size(xparam1));
- if all(xinit(:)>bayestopt_.p3) && all(xinit(:)<bayestopt_.p4)
- M_ = set_all_parameters(xinit,estim_params_,M_);
- [dr,INFO,M_,options_,oo_] = resol(0,M_,options_,oo_);
- if ~INFO(1)
- look_for_admissible_initial_condition = 0;
- end
- else
- if iter == 2000;
- scale = scale/1.1;
- iter = 0;
- else
- iter = iter+1;
- end
- end
- end
- objective_function_penalty_base = minus_logged_prior_density(xinit, bayestopt_.pshape, ...
- bayestopt_.p6, ...
- bayestopt_.p7, ...
- bayestopt_.p3, ...
- bayestopt_.p4,options_,M_,estim_params_,oo_);
- % Maximization
- [xparams,lpd,hessian] = ...
- maximize_prior_density(xinit, bayestopt_.pshape, ...
- bayestopt_.p6, ...
- bayestopt_.p7, ...
- bayestopt_.p3, ...
- bayestopt_.p4,options_,M_,bayestopt_,estim_params_,oo_);
- % Display the results.
- skipline(2)
- disp('------------------')
- disp('PRIOR OPTIMIZATION')
- disp('------------------')
- skipline()
- for i = 1:length(xparams)
- disp(['deep parameter ' int2str(i) ': ' get_the_name(i,0,M_,estim_params_,options_) '.'])
- disp([' Initial condition ....... ' num2str(xinit(i)) '.'])
- disp([' Prior mode .............. ' num2str(bayestopt_.p5(i)) '.'])
- disp([' Optimized prior mode .... ' num2str(xparams(i)) '.'])
- skipline()
- end
-end
-
-if info==3% Prior simulations (2nd order moments).
- oo_ = compute_moments_varendo('prior',options_,M_,oo_);
-end
-
-if changed_qz_criterium_flag
- options_.qz_criterium = [];
-end
-
-options_.order = order;
-
-function format_string = build_format_string(PriorStandardDeviation,LowerBound,UpperBound)
-format_string = ['%s & %s & %6.4f &'];
-if ~isnumeric(PriorStandardDeviation)
- format_string = [ format_string , ' %s &'];
-else
- format_string = [ format_string , ' %6.4f &'];
-end
-if ~isnumeric(LowerBound)
- format_string = [ format_string , ' %s &'];
-else
- format_string = [ format_string , ' %6.4f &'];
-end
-if ~isnumeric(UpperBound)
- format_string = [ format_string , ' %s &'];
-else
- format_string = [ format_string , ' %6.4f &'];
-end
-format_string = [ format_string , ' %6.4f & %6.4f \\\\ \n'];
\ No newline at end of file
--
GitLab