Commit b03697b3 authored by Johannes Pfeifer's avatar Johannes Pfeifer
Browse files

Add possibility to initialize parameters from calibrated model

parent e1b92dc5
function check_prior_bounds(xparam1,bounds,M_,estim_params_,options_,bayestopt_)
% function check_prior_bounds(xparam1,bounds,M_,estim_params_,options_)
% checks the parameter vector of violations of the prior bounds
% Inputs:
% -xparam1 [double] vector of parameters to be estimated (initial values)
% -bounds [vector] vector containing the lower and upper
% bounds
% -M_ [structure] characterizing the model.
% -estim_params_ [structure] characterizing parameters to be estimated
% -options_ [structure] characterizing the options
% -bayestopt_ [structure] characterizing priors
% Copyright (C) 2013 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/>.
outside_bound_pars=find(xparam1 < bounds(:,1) | xparam1 > bounds(:,2));
if ~isempty(outside_bound_pars)
for ii=1:length(outside_bound_pars)
outside_bound_par_names{ii,1}=get_the_name(outside_bound_pars(ii),0,M_,estim_params_,options_);
end
disp_string=[outside_bound_par_names{1,:}];
for ii=2:size(outside_bound_par_names,1)
disp_string=[disp_string,', ',outside_bound_par_names{ii,:}];
end
error(['Initial value(s) of ', disp_string ,' are outside parameter bounds. Potentially, you should set prior_trunc=0. If you used the mode_file-option, check whether your mode-file is consistent with the priors.'])
end
inadmissible_inverse_gamma_values=find(bayestopt_.pshape==4 & xparam1 == 0);
if ~isempty(inadmissible_inverse_gamma_values)
for ii=1:length(inadmissible_inverse_gamma_values)
inadmissible_inverse_gamma_par_names{ii,1}=get_the_name(inadmissible_inverse_gamma_values(ii),0,M_,estim_params_,options_);
end
disp_string=[inadmissible_inverse_gamma_par_names{1,:}];
for ii=2:size(inadmissible_inverse_gamma_par_names,1)
disp_string=[disp_string,', ',inadmissible_inverse_gamma_par_names{ii,:}];
end
error(['Initial value(s) of ', disp_string ,' is zero. This is not allowed when using an inverse gamma prior.\n'])
end
\ No newline at end of file
function [xparam1,estim_params_,xparam1_explicitly_initialized,xparam1_properly_calibrated]=do_parameter_initialization(estim_params_,xparam1_calib,xparam1_NaN_set_to_prior_mean)
% function [xparam1,estim_params_]=get_initialized_parameters(estim_params_,xparam1_calib)
% gets explicitly initialized variables and properly calibrated parameters
%
% INPUTS
% o estim_params_ [structure] characterizing parameters to be estimated.
% o xparam1_calib [double] vector of parameters to be estimated, with parameters
% initialized from calibration using get_all_parameters
%
% o xparam1_NaN_set_to_prior_mean [double] vector of parameters to be estimated, with parameters
% initialized using dynare_estimation_init; not explicitly initialized
% parameters are at prior mean
% OUTPUTS
% o xparam1 [double] vector of initialized parameters; uses the hierarchy: 1) explicitly initialized parameters,
% 2) calibrated parameters, 3) prior mean
% o estim_params_ [structure] characterizing parameters to be estimated; it is
% updated here to reflect calibrated parameters
% o xparam1_explicitly_initialized [double] vector of parameters to be estimated that
% were explicitly initialized
% o xparam1_properly_calibrated [double] vector of parameters to be estimated that
% were properly calibrated
%
% SPECIAL REQUIREMENTS
% None
% Copyright (C) 2013 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/>.
nvx = size(estim_params_.var_exo,1);
nvn = size(estim_params_.var_endo,1);
ncx = size(estim_params_.corrx,1);
ncn = size(estim_params_.corrn,1);
np = size(estim_params_.param_vals,1);
estim_params_.nvx = nvx; %exogenous shock variances
estim_params_.nvn = nvn; %endogenous variances, i.e. measurement error
estim_params_.ncx = ncx; %exogenous shock correlations
estim_params_.ncn = ncn; % correlation between endogenous variables, i.e. measurement error.
estim_params_.np = np; % other parameters of the model
xparam1_explicitly_initialized = NaN(nvx+nvn+ncx+ncn+np,1);
xparam1_properly_calibrated = NaN(nvx+nvn+ncx+ncn+np,1);
offset=0;
if nvx
initialized_par_index=find(~isnan(estim_params_.var_exo(:,2)));
calibrated_par_index=find(isnan(estim_params_.var_exo(:,2)) & ~isnan(xparam1_calib(offset+1:offset+nvx,1)));
uninitialized_par_index=find(isnan(estim_params_.var_exo(:,2)) & isnan(xparam1_calib(offset+1:offset+nvx,1)));
xparam1_explicitly_initialized(offset+initialized_par_index,1) = estim_params_.var_exo(initialized_par_index,2);
%update estim_params_ with calibrated starting values
estim_params_.var_exo(calibrated_par_index,2)=xparam1_calib(offset+calibrated_par_index,1);
%find parameters that are calibrated and do not violate inverse gamma prior
xparam1_properly_calibrated(offset+calibrated_par_index,1) = xparam1_calib(offset+calibrated_par_index,1);
inv_gamma_violation=find(estim_params_.var_exo(calibrated_par_index,2)==0 & estim_params_.var_exo(calibrated_par_index,4)==0);
if inv_gamma_violation
estim_params_.var_exo(calibrated_par_index(inv_gamma_violation),2)=NaN;
xparam1_properly_calibrated(offset+calibrated_par_index(inv_gamma_violation),1)=NaN;
fprintf('PARAMETER INITIALIZATION: Some standard deviations of shocks of the calibrated model are 0 and\n')
fprintf('PARAMETER INITIALIZATION: violate the inverse gamma prior. They will instead be initialized with the prior mean.\n')
end
if uninitialized_par_index
fprintf('PARAMETER INITIALIZATION: Warning, some estimated standard deviations of shocks are not\n')
fprintf('PARAMETER INITIALIZATION: initialized. They will be initialized with the prior mean.\n')
end
end
offset=offset+nvx;
if nvn
initialized_par_index=find(~isnan(estim_params_.var_endo(:,2)));
calibrated_par_index=find(isnan(estim_params_.var_endo(:,2)) & ~isnan(xparam1_calib(offset+1:offset+nvn,1)));
uninitialized_par_index=find(isnan(estim_params_.var_endo(:,2)) & isnan(xparam1_calib(offset+1:offset+nvn,1)));
xparam1_explicitly_initialized(offset+initialized_par_index,1) = estim_params_.var_endo(initialized_par_index,2);
estim_params_.var_endo(calibrated_par_index,2)=xparam1_calib(offset+calibrated_par_index,1);
%find parameters that are calibrated and do not violate inverse gamma prior
xparam1_properly_calibrated(offset+calibrated_par_index,1) = xparam1_calib(offset+calibrated_par_index,1);
inv_gamma_violation=find(estim_params_.var_endo(calibrated_par_index,2)==0 & estim_params_.var_endo(calibrated_par_index,5)==4);
if inv_gamma_violation
estim_params_.var_endo(calibrated_par_index(inv_gamma_violation),2)=NaN;
xparam1_properly_calibrated(offset+calibrated_par_index(inv_gamma_violation),1)=NaN;
fprintf('PARAMETER INITIALIZATION: Some measurement errors of the calibrated model are 0 and violate the\n')
fprintf('PARAMETER INITIALIZATION: inverse gamma prior. They will instead be initialized with the prior mean.\n')
end
if uninitialized_par_index
fprintf('PARAMETER INITIALIZATION: Warning, some measurement errors are not initialized. They will be initialized\n')
fprintf('PARAMETER INITIALIZATION: with the prior mean.\n')
end
end
offset=offset+nvn;
if ncx
initialized_par_index=find(~isnan(estim_params_.corrx(:,3)));
calibrated_par_index=find(isnan(estim_params_.corrx(:,3)) & ~isnan(xparam1_calib(offset+1:offset+ncx,1)));
uninitialized_par_index=find(isnan(estim_params_.corrx(:,3)) & isnan(xparam1_calib(offset+1:offset+ncx,1)));
xparam1_explicitly_initialized(offset+initialized_par_index,1) = estim_params_.corrx(initialized_par_index,3);
estim_params_.corrx(calibrated_par_index,3)=xparam1_calib(offset+calibrated_par_index,1);
xparam1_properly_calibrated(offset+calibrated_par_index,1) = xparam1_calib(offset+calibrated_par_index,1);
if uninitialized_par_index
fprintf('PARAMETER INITIALIZATION: Warning, some correlations between structural shocks are not initialized.\n')
fprintf('PARAMETER INITIALIZATION: They will be initialized with the prior mean.\n')
end
end
offset=offset+ncx;
if ncn
initialized_par_index=find(~isnan(estim_params_.corrn(:,3)));
calibrated_par_index=find(isnan(estim_params_.corrn(:,3)) & ~isnan(xparam1_calib(offset+1:offset+ncn,1)));
uninitialized_par_index=find(isnan(estim_params_.corrn(:,3)) & isnan(xparam1_calib(offset+1:offset+ncn,1)));
xparam1_explicitly_initialized(offset+initialized_par_index,1) = estim_params_.corrn(initialized_par_index,3);
estim_params_.corrn(calibrated_par_index,3)=xparam1_calib(offset+calibrated_par_index,1);
xparam1_properly_calibrated(offset+calibrated_par_index,1) = xparam1_calib(offset+calibrated_par_index,1);
if uninitialized_par_index
fprintf('PARAMETER INITIALIZATION: Warning, some correlations between measurement errors are not initialized.\n')
fprintf('PARAMETER INITIALIZATION: They will be initialized with the prior mean.\n')
end
end
offset=offset+ncn;
if np
initialized_par_index=find(~isnan(estim_params_.param_vals(:,2)));
calibrated_par_index=find(isnan(estim_params_.param_vals(:,2)) & ~isnan(xparam1_calib(offset+1:offset+np,1)));
uninitialized_par_index=find(isnan(estim_params_.param_vals(:,2)) & isnan(xparam1_calib(offset+1:offset+np,1)));
xparam1_explicitly_initialized(offset+initialized_par_index,1) = estim_params_.param_vals(initialized_par_index,2);
estim_params_.param_vals(calibrated_par_index,2)=xparam1_calib(offset+calibrated_par_index,1);
xparam1_properly_calibrated(offset+calibrated_par_index,1) = xparam1_calib(offset+calibrated_par_index,1);
if uninitialized_par_index
fprintf('PARAMETER INITIALIZATION: Warning, some deep parameters are not initialized. They will be\n')
fprintf('PARAMETER INITIALIZATION: initialized with the prior mean.\n')
end
end
xparam1=xparam1_explicitly_initialized;
xparam1(isnan(xparam1))=xparam1_properly_calibrated(isnan(xparam1)); %set not explicitly initialized parameters that do not obviously violate prior distribution to calibrated parameter values
xparam1(isnan(xparam1))=xparam1_NaN_set_to_prior_mean(isnan(xparam1)); %set not yet initialized parameters to prior mean coming from dynare_estimation_init
......@@ -80,6 +80,23 @@ end
[dataset_,xparam1, hh, M_, options_, oo_, estim_params_,bayestopt_] = dynare_estimation_init(var_list_, dname, [], M_, options_, oo_, estim_params_, bayestopt_);
%%read out calibration that was set in mod-file and can be used for initialization
xparam1_calib=get_all_parameters(estim_params_,M_); %get calibrated parameters
if ~any(isnan(xparam1_calib)) %all estimated parameters are calibrated
full_calibration_detected=1;
else
full_calibration_detected=0;
end
if options_.use_calibration_initialization %set calibration as starting values
[xparam1,estim_params_]=do_parameter_initialization(estim_params_,xparam1_calib,xparam1); %get explicitly initialized parameters that have precedence to calibrated values
try
check_prior_bounds(xparam1,[bayestopt_.lb bayestopt_.ub],M_,estim_params_,options_,bayestopt_); %check whether calibration satisfies prior bounds
catch prior_bound_check_error
fprintf('Cannot use parameter values from calibration as they violate the prior bounds.')
rethrow(prior_bound_check_error);
end
end
% Set sigma_e_is_diagonal flag (needed if the shocks block is not declared in the mod file).
M_.sigma_e_is_diagonal = 1;
if estim_params_.ncx || any(nnz(tril(M_.Sigma_e,-1)))
......@@ -174,7 +191,22 @@ if options_.dsge_var
end
oo_ = initial_estimation_checks(objective_function,xparam1,dataset_,M_,estim_params_,options_,bayestopt_,oo_);
%% perform initial estimation checks;
try
oo_ = initial_estimation_checks(objective_function,xparam1,dataset_,M_,estim_params_,options_,bayestopt_,oo_);
catch initial_estimation_checks_fail % if check fails, provide info on using calibration if present
if full_calibration_detected %calibrated model present and no explicit starting values
skipline(1);
fprintf('ESTIMATION_CHECKS: There was an error in computing the likelihood for initial parameter values.\n')
fprintf('ESTIMATION_CHECKS: You should try using the calibrated version of the model as starting values. To do\n')
fprintf('ESTIMATION_CHECKS: this, add the following estimated_params_init-block immediately before the estimation\n')
fprintf('ESTIMATION_CHECKS: command (and after the estimated_params-block so that it does not get overwritten):\n');
skipline(1);
generate_estimated_params_init_block(xparam1_calib,estim_params_,M_,options_)
skipline(2);
end
rethrow(initial_estimation_checks_fail);
end
if isequal(options_.mode_compute,0) && isempty(options_.mode_file) && options_.mh_posterior_mode_estimation==0
if options_.smoother == 1
......
......@@ -242,28 +242,7 @@ if ~isempty(estim_params_)
bounds(:,2) = ub;
end
% Test if initial values of the estimated parameters are all between the prior lower and upper bounds.
outside_bound_pars=find(xparam1 < bounds(:,1) | xparam1 > bounds(:,2));
if ~isempty(outside_bound_pars)
for ii=1:length(outside_bound_pars)
outside_bound_par_names{ii,1}=get_the_name(outside_bound_pars(ii),0,M_,estim_params_,options_);
end
disp_string=[outside_bound_par_names{1,:}];
for ii=2:size(outside_bound_par_names,1)
disp_string=[disp_string,', ',outside_bound_par_names{ii,:}];
end
error(['Initial value(s) of ', disp_string ,' are outside parameter bounds. Potentially, you should set prior_trunc=0. If you used the mode_file-option, check whether your mode-file is consistent with the priors.'])
end
inadmissible_inverse_gamma_values=find(bayestopt_.pshape==4 & xparam1 == 0);
if ~isempty(inadmissible_inverse_gamma_values)
for ii=1:length(inadmissible_inverse_gamma_values)
inadmissible_inverse_gamma_par_names{ii,1}=get_the_name(inadmissible_inverse_gamma_values(ii),0,M_,estim_params_,options_);
end
disp_string=[inadmissible_inverse_gamma_par_names{1,:}];
for ii=2:size(inadmissible_inverse_gamma_par_names,1)
disp_string=[disp_string,', ',inadmissible_inverse_gamma_par_names{ii,:}];
end
error(['Initial value(s) of ', disp_string ,' is zero. This is not allowed when using an inverse gamma prior.\n'])
end
check_prior_bounds(xparam1,bounds,M_,estim_params_,options_,bayestopt_)
lb = bounds(:,1);
ub = bounds(:,2);
......
function generate_estimated_params_init_block(xparam1,estim_params_,M_,options_)
% function generate_estimated_params_init_block(xparam1,estim_params_,M_,options_)
% writes parameter values from xparam_calib derived from get_all_parameters into
% an estimated_params_init-block that can be used to start the estimation.
%
% INPUTS
% xparam1: Parameter vector from which to create block
% estim_params_: Dynare structure describing the estimated parameters.
% M_: Dynare structure describing the model.
% options_: Dynare options structure
%
% OUTPUTS
% none
%
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2013 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/>.
nvx = estim_params_.nvx;
ncx = estim_params_.ncx;
nvn = estim_params_.nvn;
ncn = estim_params_.ncn;
np = estim_params_.np;
% stderrs of the exogenous shocks
fprintf('\nestimated_params_init;\n')
if nvx
for ii=1:nvx
vname = deblank(M_.exo_names(estim_params_.var_exo(ii,1),:));
fprintf('stderr %s, %f;\n', vname,xparam1(ii));
end
end
% update offset
offset = nvx;
% setting measument error variance
if nvn
for ii=1:nvn
vname = deblank(options_.varobs(estim_params_.nvn_observable_correspondence(ii,1),:));
fprintf('stderr %s, %f;\n', vname,xparam1(offset+ii));
end
end
% update offset
offset = nvx+nvn;
% correlations among shocks (ncx)
if ncx
corrx = estim_params_.corrx;
for ii=1:ncx
k1 = corrx(ii,1);
k2 = corrx(ii,2);
vname1 = deblank(M_.exo_names(k1,:));
vname2 = deblank(M_.exo_names(k2,:));
fprintf('corr %s, %s, %f;\n', vname1,vname2,xparam1(offset+ii));
end
end
% update offset
offset = nvx+nvn+ncx;
if ncn
for ii=1:ncn
vname1 = deblank(options_.varobs(estim_params_.corrn_observable_correspondence(ii,1),:));
vname2 = deblank(options_.varobs(estim_params_.corrn_observable_correspondence(ii,2),:));
fprintf('corr %s, %s, %f;\n', vname1,vname2,xparam1(offset+ii));
end
end
% update offset
offset = nvx+ncx+nvn+ncn;
% structural parameters
if np
for ii=1:np
jj1 = estim_params_.param_vals(ii,1);
vname = deblank(M_.param_names(jj1,:));
fprintf('%s, %f;\n',vname,xparam1(offset+ii));
end
end
fprintf('end;\n')
function xparam1=get_all_parameters(estim_params_,M_)
% function xparam1=get_parameters
% gets parameters values from M_.params into xparam1 (inverse mapping to set_all_parameters)
% This is called if a model was calibrated before estimation to back out
% parameter values
%
% INPUTS
% estim_params_: Dynare structure describing the estimated parameters.
% M_: Dynare structure describing the model.
%
% OUTPUTS
% xparam1: N*1 double vector of parameters from calibrated model that are to be estimated
%
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2013 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/>.
nvx = estim_params_.nvx;
ncx = estim_params_.ncx;
nvn = estim_params_.nvn;
ncn = estim_params_.ncn;
np = estim_params_.np;
Sigma_e = M_.Sigma_e;
Correlation_matrix = M_.Correlation_matrix;
H = M_.H;
Correlation_matrix_ME = M_.Correlation_matrix_ME;
xparam1=NaN(nvx+ncx+nvn+ncn+np,1);
% stderrs of the exogenous shocks
if nvx
var_exo = estim_params_.var_exo;
for i=1:nvx
k = var_exo(i,1);
xparam1(i)=sqrt(Sigma_e(k,k));
end
end
% update offset
offset = nvx;
% setting measument error variance
if nvn
for i=1:nvn
k = estim_params_.nvn_observable_correspondence(i,1);
xparam1(offset+i)=sqrt(H(k,k));
end
end
% update offset
offset = nvx+nvn;
% correlations among shocks (ncx)
if ncx
corrx = estim_params_.corrx;
for i=1:ncx
k1 = corrx(i,1);
k2 = corrx(i,2);
xparam1(i+offset)=Correlation_matrix(k1,k2);
end
end
% update offset
offset = nvx+nvn+ncx;
if ncn
corrn_observable_correspondence = estim_params_.corrn_observable_correspondence;
for i=1:ncn
k1 = corrn_observable_correspondence(i,1);
k2 = corrn_observable_correspondence(i,2);
xparam1(i+offset)=Correlation_matrix_ME(k1,k2);
end
end
% update offset
offset = nvx+ncx+nvn+ncn;
% structural parameters
if np
xparam1(offset+1:end)=M_.params(estim_params_.param_vals(:,1));
end
\ No newline at end of file
......@@ -381,6 +381,7 @@ options_.mh_nblck = 2;
options_.mh_recover = 0;
options_.mh_replic = 20000;
options_.recursive_estimation_restart = 0;
options_.use_calibration_initialization = 0;
options_.mode_compute = 4;
options_.mode_file = '';
......
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