Commit 0f3b68ee authored by Johannes Pfeifer 's avatar Johannes Pfeifer

Fix get_posterior_parameters

Get rid of setting variables from the base workspace in the function, making the function read-only (as the name suggests)
parent 0f84dadb
function [xparams, logpost] = GetOneDraw(type) function [xparams, logpost] = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_)
% function [xparams, logpost] = GetOneDraw(type) % function [xparams, logpost] = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_)
% draws one parameter vector and its posterior from MCMC or the prior % draws one parameter vector and its posterior from MCMC or the prior
% %
% INPUTS % INPUTS
% type: [string] 'posterior': draw from MCMC draws % type: [string] 'posterior': draw from MCMC draws
% 'prior': draw from prior % 'prior': draw from prior
% M_ [structure] Definition of the model
% estim_params_ [structure] characterizing parameters to be estimated
% oo_ [structure] Storage of results
% options_ [structure] Options
% bayestopt_ [structure] describing the priors
% %
% OUTPUTS % OUTPUTS
% xparams: vector of estimated parameters (drawn from posterior or prior distribution) % xparams: vector of estimated parameters (drawn from posterior or prior distribution)
...@@ -36,6 +41,6 @@ switch type ...@@ -36,6 +41,6 @@ switch type
case 'prior' case 'prior'
xparams = prior_draw(); xparams = prior_draw();
if nargout>1 if nargout>1
logpost = evaluate_posterior_kernel(xparams'); logpost = evaluate_posterior_kernel(xparams',M_,estim_params_,oo_,options_,bayestopt_);
end end
end end
\ No newline at end of file
...@@ -174,7 +174,7 @@ localVars.type=type; ...@@ -174,7 +174,7 @@ localVars.type=type;
if strcmpi(type,'posterior') if strcmpi(type,'posterior')
while b<B while b<B
b = b + 1; b = b + 1;
x(b,:) = GetOneDraw(type); x(b,:) = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
end end
end end
......
...@@ -141,7 +141,7 @@ while fpar<B ...@@ -141,7 +141,7 @@ while fpar<B
irun = irun+1; irun = irun+1;
irun2 = irun2+1; irun2 = irun2+1;
if strcmpi(type,'prior') if strcmpi(type,'prior')
deep = GetOneDraw(type); deep = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
else else
deep = x(fpar,:); deep = x(fpar,:);
end end
......
...@@ -519,7 +519,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ... ...@@ -519,7 +519,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
end end
prior_posterior_statistics('posterior',dataset_,dataset_info); prior_posterior_statistics('posterior',dataset_,dataset_info);
end end
xparam1 = get_posterior_parameters('mean'); xparam1 = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
M_ = set_all_parameters(xparam1,estim_params_,M_); M_ = set_all_parameters(xparam1,estim_params_,M_);
end end
end end
......
...@@ -271,15 +271,15 @@ if iload <=0 ...@@ -271,15 +271,15 @@ if iload <=0
case 'posterior_mode' case 'posterior_mode'
parameters_TeX = 'Posterior mode'; parameters_TeX = 'Posterior mode';
disp('Testing posterior mode') disp('Testing posterior mode')
params(1,:) = get_posterior_parameters('mode'); params(1,:) = get_posterior_parameters('mode',M_,estim_params_,oo_,options_);
case 'posterior_mean' case 'posterior_mean'
parameters_TeX = 'Posterior mean'; parameters_TeX = 'Posterior mean';
disp('Testing posterior mean') disp('Testing posterior mean')
params(1,:) = get_posterior_parameters('mean'); params(1,:) = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
case 'posterior_median' case 'posterior_median'
parameters_TeX = 'Posterior median'; parameters_TeX = 'Posterior median';
disp('Testing posterior median') disp('Testing posterior median')
params(1,:) = get_posterior_parameters('median'); params(1,:) = get_posterior_parameters('median',M_,estim_params_,oo_,options_);
case 'prior_mode' case 'prior_mode'
parameters_TeX = 'Prior mode'; parameters_TeX = 'Prior mode';
disp('Testing prior mode') disp('Testing prior mode')
......
function [llik,parameters] = evaluate_likelihood(parameters) function [llik,parameters] = evaluate_likelihood(parameters,M_,estim_params_,oo_,options_,bayestopt_)
% Evaluate the logged likelihood at parameters. % Evaluate the logged likelihood at parameters.
% %
% INPUTS % INPUTS
% o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for % o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for
% the (estimated) parameters of the model. % the (estimated) parameters of the model.
% % o M_ [structure] Definition of the model
% o estim_params_ [structure] characterizing parameters to be estimated
% o oo_ [structure] Storage of results
% o options_ [structure] Options
% o bayestopt_ [structure] describing the priors
% %
% OUTPUTS % OUTPUTS
% o ldens [double] value of the sample logged density at parameters. % o ldens [double] value of the sample logged density at parameters.
...@@ -35,8 +39,6 @@ function [llik,parameters] = evaluate_likelihood(parameters) ...@@ -35,8 +39,6 @@ function [llik,parameters] = evaluate_likelihood(parameters)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>. % along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global options_ M_ bayestopt_ oo_ estim_params_
persistent dataset dataset_info persistent dataset dataset_info
if nargin==0 if nargin==0
...@@ -46,11 +48,11 @@ end ...@@ -46,11 +48,11 @@ end
if ischar(parameters) if ischar(parameters)
switch parameters switch parameters
case 'posterior mode' case 'posterior mode'
parameters = get_posterior_parameters('mode'); parameters = get_posterior_parameters('mode',M_,estim_params_,oo_,options_);
case 'posterior mean' case 'posterior mean'
parameters = get_posterior_parameters('mean'); parameters = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
case 'posterior median' case 'posterior median'
parameters = get_posterior_parameters('median'); parameters = get_posterior_parameters('median',M_,estim_params_,oo_,options_);
case 'prior mode' case 'prior mode'
parameters = bayestopt_.p5(:); parameters = bayestopt_.p5(:);
case 'prior mean' case 'prior mean'
...@@ -72,5 +74,5 @@ end ...@@ -72,5 +74,5 @@ end
options_=select_qz_criterium_value(options_); options_=select_qz_criterium_value(options_);
llik = -dsge_likelihood(parameters,dataset,dataset_info,options_,M_,estim_params_,bayestopt_,prior_bounds(bayestopt_,options_.prior_trunc),oo_); llik = -dsge_likelihood(parameters,dataset,dataset_info,options_,M_,estim_params_,bayestopt_,prior_bounds(bayestopt_,options_.prior_trunc),oo_);
ldens = evaluate_prior(parameters); ldens = evaluate_prior(parameters,M_,estim_params_,oo_,options_,bayestopt_);
llik = llik - ldens; llik = llik - ldens;
\ No newline at end of file
function lpkern = evaluate_posterior_kernel(parameters,llik) function lpkern = evaluate_posterior_kernel(parameters,M_,estim_params_,oo_,options_,bayestopt_,llik)
% Evaluate the prior density at parameters. % Evaluate the evaluate_posterior_kernel at parameters.
% %
% INPUTS % INPUTS
% o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for % o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for
% the (estimated) parameters of the model. % the (estimated) parameters of the model.
% % o M_ [structure] Definition of the model
% o estim_params_ [structure] characterizing parameters to be estimated
% o oo_ [structure] Storage of results
% o options_ [structure] Options
% o bayestopt_ [structure] describing the priors
% o llik [double] value of the logged likelihood if it
% should not be computed
% %
% OUTPUTS % OUTPUTS
% o lpkern [double] value of the logged posterior kernel. % o lpkern [double] value of the logged posterior kernel.
...@@ -34,8 +40,8 @@ function lpkern = evaluate_posterior_kernel(parameters,llik) ...@@ -34,8 +40,8 @@ function lpkern = evaluate_posterior_kernel(parameters,llik)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>. % along with Dynare. If not, see <http://www.gnu.org/licenses/>.
[ldens,parameters] = evaluate_prior(parameters); [ldens,parameters] = evaluate_prior(parameters,M_,estim_params_,oo_,options_,bayestopt_);
if nargin==1 if nargin==6 %llik provided as an input
llik = evaluate_likelihood(parameters); llik = evaluate_likelihood(parameters,M_,estim_params_,oo_,options_,bayestopt_);
end end
lpkern = ldens+llik; lpkern = ldens+llik;
\ No newline at end of file
function [ldens,parameters] = evaluate_prior(parameters) function [ldens,parameters] = evaluate_prior(parameters,M_,estim_params_,oo_,options_,bayestopt_)
% Evaluate the prior density at parameters. % Evaluate the prior density at parameters.
% %
% INPUTS % INPUTS
% o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for % o parameters a string ('posterior mode','posterior mean','posterior median','prior mode','prior mean') or a vector of values for
% the (estimated) parameters of the model. % the (estimated) parameters of the model.
% o M_ [structure] Definition of the model
% o oo_ [structure] Storage of results
% o options_ [structure] Options
% o bayestopt_ [structure] describing the priors
% o estim_params_ [structure] characterizing parameters to be estimated
% %
% %
% OUTPUTS % OUTPUTS
...@@ -33,8 +38,6 @@ function [ldens,parameters] = evaluate_prior(parameters) ...@@ -33,8 +38,6 @@ function [ldens,parameters] = evaluate_prior(parameters)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>. % along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global bayestopt_
if nargin==0 if nargin==0
parameters = 'posterior mode'; parameters = 'posterior mode';
end end
...@@ -42,11 +45,11 @@ end ...@@ -42,11 +45,11 @@ end
if ischar(parameters) if ischar(parameters)
switch parameters switch parameters
case 'posterior mode' case 'posterior mode'
parameters = get_posterior_parameters('mode'); parameters = get_posterior_parameters('mode',M_,estim_params_,oo_,options_);
case 'posterior mean' case 'posterior mean'
parameters = get_posterior_parameters('mean'); parameters = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
case 'posterior median' case 'posterior median'
parameters = get_posterior_parameters('median'); parameters = get_posterior_parameters('median',M_,estim_params_,oo_,options_);
case 'prior mode' case 'prior mode'
parameters = bayestopt_.p5(:); parameters = bayestopt_.p5(:);
case 'prior mean' case 'prior mean'
......
...@@ -73,13 +73,13 @@ end ...@@ -73,13 +73,13 @@ end
if ischar(parameters) if ischar(parameters)
switch parameters switch parameters
case 'posterior_mode' case 'posterior_mode'
parameters = get_posterior_parameters('mode'); parameters = get_posterior_parameters('mode',M_,estim_params_,oo_,options_);
case 'posterior_mean' case 'posterior_mean'
parameters = get_posterior_parameters('mean'); parameters = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
case 'posterior_median' case 'posterior_median'
parameters = get_posterior_parameters('median'); parameters = get_posterior_parameters('median',M_,estim_params_,oo_,options_);
case 'mle_mode' case 'mle_mode'
parameters = get_posterior_parameters('mode','mle_'); parameters = get_posterior_parameters('mode',M_,estim_params_,oo_,options_,'mle_');
case 'prior_mode' case 'prior_mode'
parameters = bayestopt_.p5(:); parameters = bayestopt_.p5(:);
case 'prior_mean' case 'prior_mean'
......
...@@ -68,11 +68,11 @@ else ...@@ -68,11 +68,11 @@ else
end end
%get draws for later use %get draws for later use
first_draw=GetOneDraw(type); first_draw=GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
parameter_mat=NaN(n_draws,length(first_draw)); parameter_mat=NaN(n_draws,length(first_draw));
parameter_mat(1,:)=first_draw; parameter_mat(1,:)=first_draw;
for draw_iter=2:n_draws for draw_iter=2:n_draws
parameter_mat(draw_iter,:) = GetOneDraw(type); parameter_mat(draw_iter,:) = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
end end
% get output size % get output size
......
...@@ -6,7 +6,7 @@ function xparam1=get_all_parameters(estim_params_,M_) ...@@ -6,7 +6,7 @@ function xparam1=get_all_parameters(estim_params_,M_)
% parameter values % parameter values
% %
% INPUTS % INPUTS
% estim_params_: Dynare structure describing the estimated parameters. % estim_params_: Dynare structure describing the estimated parameters.
% M_: Dynare structure describing the model. % M_: Dynare structure describing the model.
% %
% OUTPUTS % OUTPUTS
......
function xparam = get_posterior_parameters(type,field1) function xparam = get_posterior_parameters(type,M_,estim_params_,oo_,options_,field1)
% function xparam = get_posterior_parameters(type) % function xparam = get_posterior_parameters(type,M_,estim_params_,oo_,options_,field1)
% Selects (estimated) parameters (posterior mode or posterior mean). % Selects (estimated) parameters (posterior mode or posterior mean).
% %
% INPUTS % INPUTS
% o type [char] = 'mode' or 'mean'. % o type [char] = 'mode' or 'mean'.
% o field_1 [char] optional field like 'mle_'. % o M_: [structure] Dynare structure describing the model.
% o estim_params_: [structure] Dynare structure describing the estimated parameters.
% o field_1 [char] optional field like 'mle_'.
% %
% OUTPUTS % OUTPUTS
% o xparam vector of estimated parameters % o xparam vector of estimated parameters
...@@ -30,9 +32,7 @@ function xparam = get_posterior_parameters(type,field1) ...@@ -30,9 +32,7 @@ function xparam = get_posterior_parameters(type,field1)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>. % along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global estim_params_ oo_ options_ M_ if nargin<6
if nargin<2
field1='posterior_'; field1='posterior_';
end end
nvx = estim_params_.nvx; nvx = estim_params_.nvx;
...@@ -48,7 +48,6 @@ for i=1:nvx ...@@ -48,7 +48,6 @@ for i=1:nvx
k1 = estim_params_.var_exo(i,1); k1 = estim_params_.var_exo(i,1);
name1 = deblank(M_.exo_names(k1,:)); name1 = deblank(M_.exo_names(k1,:));
xparam(m) = oo_.([field1 type]).shocks_std.(name1); xparam(m) = oo_.([field1 type]).shocks_std.(name1);
M_.Sigma_e(k1,k1) = xparam(m)^2;
m = m+1; m = m+1;
end end
...@@ -65,8 +64,6 @@ for i=1:ncx ...@@ -65,8 +64,6 @@ for i=1:ncx
name1 = deblank(M_.exo_names(k1,:)); name1 = deblank(M_.exo_names(k1,:));
name2 = deblank(M_.exo_names(k2,:)); name2 = deblank(M_.exo_names(k2,:));
xparam(m) = oo_.([field1 type]).shocks_corr.([name1 '_' name2]); xparam(m) = oo_.([field1 type]).shocks_corr.([name1 '_' name2]);
M_.Sigma_e(k1,k2) = xparam(m);
M_.Sigma_e(k2,k1) = xparam(m);
m = m+1; m = m+1;
end end
...@@ -84,10 +81,5 @@ FirstDeep = m; ...@@ -84,10 +81,5 @@ FirstDeep = m;
for i=1:np for i=1:np
name1 = deblank(M_.param_names(estim_params_.param_vals(i,1),:)); name1 = deblank(M_.param_names(estim_params_.param_vals(i,1),:));
xparam(m) = oo_.([field1 type]).parameters.(name1); xparam(m) = oo_.([field1 type]).parameters.(name1);
%assignin('base',name1,xparam(m));% Useless with version 4 (except maybe for users)
m = m+1; m = m+1;
end
if np
M_.params(estim_params_.param_vals(:,1)) = xparam(FirstDeep:end);
end end
\ No newline at end of file
...@@ -80,16 +80,16 @@ if estimated_model ...@@ -80,16 +80,16 @@ if estimated_model
if ischar(options_cond_fcst.parameter_set) if ischar(options_cond_fcst.parameter_set)
switch options_cond_fcst.parameter_set switch options_cond_fcst.parameter_set
case 'posterior_mode' case 'posterior_mode'
xparam = get_posterior_parameters('mode'); xparam = get_posterior_parameters('mode',M_,estim_params_,oo_,options_);
graph_title='Posterior Mode'; graph_title='Posterior Mode';
case 'posterior_mean' case 'posterior_mean'
xparam = get_posterior_parameters('mean'); xparam = get_posterior_parameters('mean',M_,estim_params_,oo_,options_);
graph_title='Posterior Mean'; graph_title='Posterior Mean';
case 'posterior_median' case 'posterior_median'
xparam = get_posterior_parameters('median'); xparam = get_posterior_parameters('median',M_,estim_params_,oo_,options_);
graph_title='Posterior Median'; graph_title='Posterior Median';
case 'mle_mode' case 'mle_mode'
xparam = get_posterior_parameters('mode','mle_'); xparam = get_posterior_parameters('mode',M_,estim_params_,oo_,options_,'mle_');
graph_title='ML Mode'; graph_title='ML Mode';
case 'prior_mode' case 'prior_mode'
xparam = bayestopt_.p5(:); xparam = bayestopt_.p5(:);
......
...@@ -226,7 +226,7 @@ if strcmpi(type,'posterior') ...@@ -226,7 +226,7 @@ if strcmpi(type,'posterior')
else else
logpost=NaN(B,1); logpost=NaN(B,1);
for b=1:B for b=1:B
[x(b,:), logpost(b)] = GetOneDraw(type); [x(b,:), logpost(b)] = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
end end
end end
localVars.logpost=logpost; localVars.logpost=logpost;
......
...@@ -190,14 +190,14 @@ end ...@@ -190,14 +190,14 @@ end
for b=fpar:B for b=fpar:B
if strcmpi(type,'prior') if strcmpi(type,'prior')
[deep, logpo] = GetOneDraw(type); [deep, logpo] = GetOneDraw(type,M_,estim_params_,oo_,options_,bayestopt_);
else else
deep = x(b,:); deep = x(b,:);
if strcmpi(type,'posterior') if strcmpi(type,'posterior')
logpo = logpost(b); logpo = logpost(b);
else else
logpo = evaluate_posterior_kernel(deep'); logpo = evaluate_posterior_kernel(deep',M_,estim_params_,oo_,options_,bayestopt_);
end end
end end
M_ = set_all_parameters(deep,estim_params_,M_); M_ = set_all_parameters(deep,estim_params_,M_);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment