1function [post_mean, post_median, post_var, hpd_interval, post_deciles, density] = posterior_moments(xx,info,mh_conf_sig,kernel_options) 2% Computes posterior mean, median, variance, HPD interval, deciles, and density from posterior draws. 3% 4% INPUTS 5% xx [double] Vector of posterior draws (or prior draws ;-) 6% info [integer] If equal to one the posterior density is estimated. 7% mh_config_sig [double] Scalar between 0 and 1 specifying the size of the HPD interval. 8% 9% 10% OUTPUTS 11% post_mean [double] Scalar, posterior mean. 12% post_median [double] Scalar, posterior median. 13% post_var [double] Scalar, posterior variance. 14% hpd_interval [double] Vector (1*2), Highest Probability Density interval 15% post_deciles [double] Vector (9*1), deciles of the posterior distribution. 16% density [double] Matrix (n*2), non parametric estimate of the posterior density. First and second 17% columns are respectively abscissa and ordinate coordinates. 18% 19% SPECIAL REQUIREMENTS 20% Other matlab routines distributed with Dynare: mh_optimal_bandwidth.m 21% kernel_density_estimate.m. 22% 23 24% Copyright (C) 2005-2017 Dynare Team 25% 26% This file is part of Dynare. 27% 28% Dynare is free software: you can redistribute it and/or modify 29% it under the terms of the GNU General Public License as published by 30% the Free Software Foundation, either version 3 of the License, or 31% (at your option) any later version. 32% 33% Dynare is distributed in the hope that it will be useful, 34% but WITHOUT ANY WARRANTY; without even the implied warranty of 35% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 36% GNU General Public License for more details. 37% 38% You should have received a copy of the GNU General Public License 39% along with Dynare. If not, see <http://www.gnu.org/licenses/>. 40 41if nargin<4 42 number_of_grid_points = 2^9; % 2^9 = 512 !... Must be a power of two. 43 bandwidth = 0; % Rule of thumb optimal bandwidth parameter. 44 kernel_function = 'gaussian'; % Gaussian kernel for Fast Fourrier Transform approximaton. 45else 46 number_of_grid_points = kernel_options.gridpoints; 47 bandwidth = kernel_options.bandwidth; 48 kernel_function = kernel_options.kernel_function; 49end 50 51xx = xx(:); 52xx = sort(xx); 53 54 55post_mean = mean(xx); 56post_median = median(xx); 57post_var = var(xx); 58 59number_of_draws = length(xx); 60hpd_draws = round((1-mh_conf_sig)*number_of_draws); 61 62if hpd_draws>2 63 kk = zeros(hpd_draws,1); 64 jj = number_of_draws-hpd_draws; 65 for ii = 1:hpd_draws 66 kk(ii) = xx(jj)-xx(ii); 67 jj = jj + 1; 68 end 69 [kmin,idx] = min(kk); 70 hpd_interval = [xx(idx) xx(idx)+kmin]; 71else 72 hpd_interval=NaN(1,2); 73end 74if length(xx)>9 75 post_deciles = xx([round(0.1*number_of_draws) ... 76 round(0.2*number_of_draws) ... 77 round(0.3*number_of_draws) ... 78 round(0.4*number_of_draws) ... 79 round(0.5*number_of_draws) ... 80 round(0.6*number_of_draws) ... 81 round(0.7*number_of_draws) ... 82 round(0.8*number_of_draws) ... 83 round(0.9*number_of_draws)]); 84else 85 post_deciles=NaN(9,1); 86end 87density = []; 88if info 89 if post_var > 1e-12 90 optimal_bandwidth = mh_optimal_bandwidth(xx,number_of_draws,bandwidth,kernel_function); 91 [density(:,1),density(:,2)] = kernel_density_estimate(xx,number_of_grid_points,... 92 number_of_draws,optimal_bandwidth,kernel_function); 93 else 94 density = NaN(number_of_grid_points,2); 95 end 96end