From 6cb41e0252eb8524b0d809fe4eb3140bdbded286 Mon Sep 17 00:00:00 2001 From: stepan Date: Tue, 7 Apr 2009 22:08:39 +0000 Subject: [PATCH] * Added info=19 in resol.m (problem in the steady state file). * Bug fix in prior_sampler. * Print more informations (BK conditions, steadys state or mjdgges * problems...) when running get_prior_info. git-svn-id: https://www.dynare.org/svn/dynare/trunk@2567 ac1d8469-bf42-47a9-8791-bf33cf982152 --- matlab/get_prior_info.m | 12 ++++++++++-- matlab/prior_sampler.m | 27 ++++++++++++++++++++------- matlab/resol.m | 18 ++++++++++-------- 3 files changed, 40 insertions(+), 17 deletions(-) diff --git a/matlab/get_prior_info.m b/matlab/get_prior_info.m index a2df20ad3..daabe6be1 100644 --- a/matlab/get_prior_info.m +++ b/matlab/get_prior_info.m @@ -94,8 +94,16 @@ function get_prior_info(info) M_.dname = M_.fname; if info% Prior simulations. - results = prior_sampler(1,M_,bayestopt_,options_,oo_); - results.prior.mass + results = prior_sampler(0,M_,bayestopt_,options_,oo_); + 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 diff --git a/matlab/prior_sampler.m b/matlab/prior_sampler.m index 32b8d2f1a..fa1e8ce1c 100644 --- a/matlab/prior_sampler.m +++ b/matlab/prior_sampler.m @@ -42,17 +42,21 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) count_bk_singularity = 0; count_static_var_def = 0; count_no_steadystate = 0; + count_steadystate_file_exit = 0; count_dll_problem = 0; + count_complex_jacobian = 0; + count_complex_steadystate = 0; count_unknown_problem = 0; NumberOfSimulations = options_.prior_mc; NumberOfParameters = length(bayestopt_.p1); NumberOfEndogenousVariables = size(M_.endo_names,1); NumberOfElementsPerFile = ceil(options_.MaxNumberOfBytes/NumberOfParameters/NumberOfEndogenousVariables/8) ; + if NumberOfSimulations <= NumberOfElementsPerFile TableOfInformations = [ 1 , NumberOfSimulations , 1] ; else NumberOfFiles = fix(NumberOfSimulations/NumberOfElementsPerFile) ; - NumberOfElementsInTheLastFile = NumberOfSimulations - NumberOfElementsPerFile*NumberOfFiles ; + NumberOfElementsInTheLastFile = NumberOfSimulations - NumberOfElementsPerFile*(NumberOfFiles-1) ; if ~isint(NumberOfSimulations/NumberOfElementsPerFile) NumberOfFiles = NumberOfFiles + 1 ; end @@ -63,11 +67,11 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) TableOfInformations(1,3) = 1; TableOfInformations(2:end,3) = cumsum(TableOfInformations(2:end,2))+1; end + pdraws = cell(TableOfInformations(1,2),drsave+1) ; sampled_prior_expectation = zeros(NumberOfParameters,1); sampled_prior_covariance = zeros(NumberOfParameters,NumberOfParameters); - % Simulations. while iteration <= NumberOfSimulations loop_indx = loop_indx+1; @@ -100,10 +104,14 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) count_bk_singularity = count_bk_singularity + 1 ; case 20 count_no_steadystate = count_no_steadystate + 1 ; + case 19 + count_steadystate_file_exit = count_steadystate_file_exit + 1 ; + case 6 + count_complex_jacobian = count_complex_jacobian + 1 ; + case 21 + count_complex_steadystate = count_complex_steadystate + 1 ; otherwise - % To be checked... count_unknown_problem = count_unknown_problem + 1 ; - continue end end @@ -116,18 +124,23 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) results.bk.singularity_share = count_bk_singularity/loop_indx; results.dll.problem_share = count_dll_problem/loop_indx; results.ss.problem_share = count_no_steadystate/loop_indx; - results.garbage_share = results.bk.indeterminacy_share + ... + results.ss.complex_share = count_complex_steadystate/loop_indx; + results.ass.problem_share = count_steadystate_file_exit/loop_indx; + results.jacobian.problem_share = count_complex_jacobian/loop_indx; + results.garbage_share = ... + results.bk.indeterminacy_share + ... results.bk.unstability_share + ... results.bk.singularity_share + ... results.dll.problem_share + ... results.ss.problem_share + ... - (count_unknown_problem/loop_indx) ; + results.ass.problem_share + ... + results.jacobian.problem_share + ... + count_unknown_problem/loop_indx ; results.prior.mean = sampled_prior_expectation; results.prior.variance = sampled_prior_covariance; results.prior.mass = 1-results.garbage_share; - function [mu,sigma] = recursive_prior_moments(m0,s0,newobs,iter) % Recursive estimation of order one and two moments (expectation and % covariance matrix). newobs should be a row vector. I do not use the diff --git a/matlab/resol.m b/matlab/resol.m index 4fef9ca36..944e47bfb 100644 --- a/matlab/resol.m +++ b/matlab/resol.m @@ -15,6 +15,7 @@ function [dr,info]=resol(ys,check_flag) % info=4: Blanchard Kahn conditions are not satisfied:'...' indeterminacy % info=5: Blanchard Kahn conditions are not satisfied:'...' indeterminacy due to rank failure % info=6: The jacobian evaluated at the steady state is complex. +% info=19: The steadystate file did not compute the steady state (inconsistent deep parameters). % info=20: can't find steady state info(2) contains sum of sqare residuals % info=21: steady state is complex % info(2) contains sum of sqare of @@ -92,14 +93,15 @@ else end % testing for problem if check1 - info(1)= 20; - if options_.steadystate_flag - resid = check1 ; - else - resid = feval(fh,ys,oo_.exo_steady_state, M_.params); - end - info(2) = resid'*resid ; % penalty... - return + if options_.steadystate_flag + info(1)= 19; + resid = check1 ; + else + info(1)= 20; + resid = feval(fh,ys,oo_.exo_steady_state, M_.params); + end + info(2) = resid'*resid ; + return end if ~isreal(dr.ys)