stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Deal with cases where too few draws are available for posterior moments

time-shift
Johannes Pfeifer 2014-03-06 08:10:51 +01:00 committed by Stéphane Adjemian (Scylla)
parent e7727ba2d3
commit bf674fa59f
2 changed files with 33 additions and 18 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
matlab/GetPosteriorParametersStatistics.m
View File

@ -75,6 +75,14 @@ catch
[marginal,oo_] = marginal_density(M_, options_, estim_params_, oo_);
disp(sprintf('Log data density is %f.',oo_.MarginalDensity.ModifiedHarmonicMean))
end
num_draws=length(NumberOfDraws*options_.mh_nblck);
hpd_draws = round((1-options_.mh_conf_sig)*num_draws);
if hpd_draws<2
fprintf('posterior_moments: There are not enough draws computes to compute HPD Intervals. Skipping their computation.\n')
end
if num_draws<9
fprintf('posterior_moments: There are not enough draws computes to compute deciles. Skipping their computation.\n')
end
if np
type = 'parameters';
if TeX

43
matlab/posterior_moments.m
View File

@ -48,25 +48,32 @@ post_var = var(xx);
number_of_draws = length(xx);
hpd_draws = round((1-mh_conf_sig)*number_of_draws);
kk = zeros(hpd_draws,1);
jj = number_of_draws-hpd_draws;
for ii = 1:hpd_draws
kk(ii) = xx(jj)-xx(ii);
jj = jj + 1;
if hpd_draws>2
kk = zeros(hpd_draws,1);
jj = number_of_draws-hpd_draws;
for ii = 1:hpd_draws
kk(ii) = xx(jj)-xx(ii);
jj = jj + 1;
end
[kmin,idx] = min(kk);
hpd_interval = [xx(idx) xx(idx)+kmin];
else
hpd_interval=NaN(1,2);
end
if length(xx)>9
post_deciles = xx([round(0.1*number_of_draws) ...
round(0.2*number_of_draws) ...
round(0.3*number_of_draws) ...
round(0.4*number_of_draws) ...
round(0.5*number_of_draws) ...
round(0.6*number_of_draws) ...
round(0.7*number_of_draws) ...
round(0.8*number_of_draws) ...
round(0.9*number_of_draws)]);
else
post_deciles=NaN(9,1);
end
[kmin,idx] = min(kk);
hpd_interval = [xx(idx) xx(idx)+kmin];
post_deciles = xx([round(0.1*number_of_draws) ...
round(0.2*number_of_draws) ...
round(0.3*number_of_draws) ...
round(0.4*number_of_draws) ...
round(0.5*number_of_draws) ...
round(0.6*number_of_draws) ...
round(0.7*number_of_draws) ...
round(0.8*number_of_draws) ...
round(0.9*number_of_draws)]);
density = [];
if info
number_of_grid_points = 2^9; % 2^9 = 512 !... Must be a power of two.