Merge branch 'dynare-newrat'

Ref. !2193
kalman-mex
Sébastien Villemot 2023-10-11 09:33:33 -04:00
commit cf9ea686ef
No known key found for this signature in database
GPG Key ID: 2CECE9350ECEBE4A
40 changed files with 488 additions and 637 deletions

View File

@ -155,7 +155,7 @@ Dynare misc commands
|br|
.. matcomm:: prior [OPTIONS[, ...]];
.. matcomm:: prior [OPTIONS[ ...]];
Prints information about the prior distribution given the provided
options. If no options are provided, the command returns the list of

View File

@ -4498,6 +4498,11 @@ Computing the stochastic solution
Dont print moments of the endogenous variables (printing them
is the default).
.. option:: nomodelsummary
Dont print the model summary and the covariance of the exogenous shocks (printing them
is the default).
.. option:: nograph
Do not create graphs (which implies that they are not saved to

View File

@ -249,7 +249,7 @@ CheckPath('graphs',options_mom_.dirname);
options_mom_.mom.compute_derivs = false; % flag to compute derivs in objective function (might change for GMM with either analytic_standard_errors or analytic_jacobian (dependent on optimizer))
options_mom_.mom.vector_output = false; % specifies whether the objective function returns a vector
% decision rule
oo_.dr = set_state_space(oo_.dr,M_,options_mom_); % get state-space representation
oo_.dr = set_state_space(oo_.dr,M_); % get state-space representation
oo_.mom.obs_var = []; % create index of observed variables in DR order
for i = 1:options_mom_.obs_nbr
oo_.mom.obs_var = [oo_.mom.obs_var; find(strcmp(options_mom_.varobs{i}, M_.endo_names(oo_.dr.order_var)))];
@ -316,7 +316,7 @@ if strcmp(options_mom_.mom.mom_method,'GMM') || strcmp(options_mom_.mom.mom_meth
end
% check for calibrated covariances before updating parameters
estim_params_ = check_for_calibrated_covariances(xparam0,estim_params_,M_);
estim_params_ = check_for_calibrated_covariances(estim_params_,M_);
% checks on parameter calibration and initialization
xparam_calib = get_all_parameters(estim_params_,M_); % get calibrated parameters

View File

@ -105,7 +105,7 @@ if isfield(options_.osr,'maxit') || isfield(options_.osr,'tolf')
end
end
oo_.dr = set_state_space(oo_.dr,M_,options_);
oo_.dr = set_state_space(oo_.dr,M_);
par_0 = M_.params(i_params);
inv_order_var = oo_.dr.inv_order_var;

View File

@ -1,11 +1,11 @@
function [eigenvalues_,result,info] = check(M, options, oo)
function [eigenvalues_,result,info] = check(M_, options_, oo_)
%[eigenvalues_,result,info] = check(M_, options_, oo_)
% Checks determinacy conditions by computing the generalized eigenvalues.
%
% INPUTS
% - M [structure] Matlab's structure describing the model (M_).
% - options [structure] Matlab's structure describing the current options (options_).
% - oo [structure] Matlab's structure containing the results (oo_).
% - M_ [structure] Matlab's structure describing the model
% - options_ [structure] Matlab's structure describing the current options
% - oo_ [structure] Matlab's structure containing the results
%
% OUTPUTS
% - eigenvalues_ [double] vector, eigenvalues.
@ -30,40 +30,40 @@ function [eigenvalues_,result,info] = check(M, options, oo)
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
if ~options.initval_file && M.exo_nbr > 1
oo.exo_simul = ones(M.maximum_lead+M.maximum_lag+1,1)*oo.exo_steady_state';
if ~options_.initval_file && M_.exo_nbr > 1
oo_.exo_simul = ones(M_.maximum_lead+M_.maximum_lag+1,1)*oo_.exo_steady_state';
end
options.order = 1;
options_.order = 1;
if isempty(options.qz_criterium)
options.qz_criterium = 1+1e-6;
if isempty(options_.qz_criterium)
options_.qz_criterium = 1+1e-6;
end
oo.dr=set_state_space(oo.dr,M,options);
oo_.dr=set_state_space(oo_.dr,M_);
[dr,info] = resol(1,M,options,oo.dr ,oo.steady_state, oo.exo_steady_state, oo.exo_det_steady_state);
[dr,info] = resol(1,M_,options_,oo_.dr ,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
if info(1) ~= 0 && info(1) ~= 3 && info(1) ~= 4
print_info(info, 0, options);
print_info(info, 0, options_);
end
eigenvalues_ = dr.eigval;
[m_lambda,i]=sort(abs(eigenvalues_));
result = 0;
if (M.nsfwrd == dr.edim) && (dr.full_rank)
if (M_.nsfwrd == dr.edim) && (dr.full_rank)
result = 1;
end
if ~options.noprint
if ~options_.noprint
skipline()
disp('EIGENVALUES:')
disp(sprintf('%16s %16s %16s\n','Modulus','Real','Imaginary'))
fprintf('%16s %16s %16s\n','Modulus','Real','Imaginary');
z=[m_lambda real(eigenvalues_(i)) imag(eigenvalues_(i))]';
disp(sprintf('%16.4g %16.4g %16.4g\n',z))
disp(sprintf('\nThere are %d eigenvalue(s) larger than 1 in modulus ', dr.edim));
disp(sprintf('for %d forward-looking variable(s)', M.nsfwrd));
fprintf('%16.4g %16.4g %16.4g\n',z);
fprintf('\nThere are %d eigenvalue(s) larger than 1 in modulus ', dr.edim);
fprintf('for %d forward-looking variable(s)', M_.nsfwrd);
skipline()
if result
disp('The rank condition is verified.')

View File

@ -1,10 +1,9 @@
function estim_params=check_for_calibrated_covariances(xparam1,estim_params,M)
% function check_for_calibrated_covariances(xparam1,estim_params,M)
function estim_params=check_for_calibrated_covariances(estim_params,M_)
% function check_for_calibrated_covariances(estim_params,M)
% find calibrated covariances to consider during estimation
% Inputs
% -xparam1 [vector] parameters to be estimated
% -estim_params [structure] describing parameters to be estimated
% -M [structure] describing the model
% -M_ [structure] describing the model
%
% Outputs
% -estim_params [structure] describing parameters to be estimated
@ -12,7 +11,7 @@ function estim_params=check_for_calibrated_covariances(xparam1,estim_params,M)
% Notes: M is local to this function and not updated when calling
% set_all_parameters
% Copyright © 2013-2017 Dynare Team
% Copyright © 2013-2023 Dynare Team
%
% This file is part of Dynare.
%
@ -28,27 +27,108 @@ function estim_params=check_for_calibrated_covariances(xparam1,estim_params,M)
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
Sigma_e_calibrated=M.Sigma_e;
H_calibrated=M.H;
%check covariance for structural errors
covariance_pos=find(tril(Sigma_e_calibrated,-1)); %off-diagonal elements set by covariances before updating correlation matrix to reflect estimated covariances
covariance_pos_ME=find(tril(H_calibrated,-1)); %off-diagonal elements set by covariances before updating correlation matrix to reflect estimated covariances
%locally updated M
M = set_all_parameters(xparam1,estim_params,M);
correlation_pos=find(tril(M.Correlation_matrix,-1)); %off-diagonal elements set by correlations after accounting for estimation
calibrated_covariance_pos=covariance_pos(~ismember(covariance_pos,correlation_pos));
if any(calibrated_covariance_pos)
[rows, columns]=ind2sub(size(M.Sigma_e),calibrated_covariance_pos); %find linear indices of lower triangular covariance entries
estim_params.calibrated_covariances.position=[calibrated_covariance_pos;sub2ind(size(M.Sigma_e),columns,rows)]; %get linear entries of upper triangular parts
estim_params.calibrated_covariances.cov_value=Sigma_e_calibrated(estim_params.calibrated_covariances.position);
if isfield(estim_params,'calibrated_covariances')
estim_params = rmfield(estim_params,'calibrated_covariances'); %remove if already present
end
if isfield(estim_params,'calibrated_covariances_ME')
estim_params = rmfield(estim_params,'calibrated_covariances_ME'); %remove if already present
end
correlation_pos_ME=find(tril(M.Correlation_matrix_ME,-1)); %off-diagonal elements set by correlations after accounting for estimation
calibrated_covariance_pos_ME=covariance_pos_ME(~ismember(covariance_pos_ME,correlation_pos_ME));
if any(calibrated_covariance_pos_ME)
[rows, columns]=ind2sub(size(M.H),calibrated_covariance_pos_ME); %find linear indices of lower triangular covariance entries
estim_params.calibrated_covariances_ME.position=[calibrated_covariance_pos_ME;sub2ind(size(M.H),columns,rows)]; %get linear entries of upper triangular parts
estim_params.calibrated_covariances_ME.cov_value=H_calibrated(estim_params.calibrated_covariances_ME.position);
[rows_calibrated, columns_calibrated]=ind2sub(size(M_.Sigma_e),find(tril(M_.Sigma_e,-1))); %find linear indices of preset lower triangular covariance entries
if estim_params.ncx %delete preset entries actually estimated
for i=1:estim_params.ncx
shock_1 = estim_params.corrx(i,1);
shock_2 = estim_params.corrx(i,2);
estimated_corr_pos=find(rows_calibrated==shock_1 & columns_calibrated==shock_2);
if ~isempty(estimated_corr_pos)
rows_calibrated(estimated_corr_pos)=[];
columns_calibrated(estimated_corr_pos)=[];
end
estimated_corr_pos=find(rows_calibrated==shock_2 & columns_calibrated==shock_1);
if ~isempty(estimated_corr_pos)
rows_calibrated(estimated_corr_pos)=[];
columns_calibrated(estimated_corr_pos)=[];
end
end
if any(rows_calibrated)
estim_params.calibrated_covariances.position=[sub2ind(size(M_.Sigma_e),rows_calibrated,columns_calibrated);sub2ind(size(M_.Sigma_e),columns_calibrated,rows_calibrated)]; %get linear entries of upper triangular parts
estim_params.calibrated_covariances.cov_value=M_.Sigma_e(estim_params.calibrated_covariances.position);
end
end
[rows_calibrated, columns_calibrated]=ind2sub(size(M_.H),find(tril(M_.H,-1))); %find linear indices of preset lower triangular covariance entries
if estim_params.ncn %delete preset entries actually estimated
for i=1:estim_params.ncn
shock_1 = estim_params.corrn(i,1);
shock_2 = estim_params.corrn(i,2);
estimated_corr_pos=find(rows_calibrated==shock_1 & columns_calibrated==shock_2);
if ~isempty(estimated_corr_pos)
rows_calibrated(estimated_corr_pos)=[];
columns_calibrated(estimated_corr_pos)=[];
end
estimated_corr_pos=find(rows_calibrated==shock_2 & columns_calibrated==shock_1);
if ~isempty(estimated_corr_pos)
rows_calibrated(estimated_corr_pos)=[];
columns_calibrated(estimated_corr_pos)=[];
end
end
end
if any(rows_calibrated)
estim_params.calibrated_covariances_ME.position=[sub2ind(size(M_.H),rows_calibrated,columns_calibrated);sub2ind(size(M_.H),columns_calibrated,rows_calibrated)]; %get linear entries of upper triangular parts
estim_params.calibrated_covariances_ME.cov_value=M_.H(estim_params.calibrated_covariances_ME.position);
end
return % --*-- Unit tests --*--
%@test:1
M_.Sigma_e=[1 0; 0 1];
M_.H=[1 0; 0 1];
M_.Correlation_matrix= [1 -0.5; -0.5 1];
M_.Correlation_matrix_ME=[1 -0.5; -0.5 1];
estim_params.ncx=1;
estim_params.ncn=1;
estim_params.corrx=[2 1 NaN -1 1 3 0 0.2000 NaN NaN NaN];
estim_params.corrn=[2 1 NaN -1 1 3 0 0.2000 NaN NaN NaN];
estim_params=check_for_calibrated_covariances(estim_params,M_);
if isfield(estim_params,'calibrated_covariances_ME') || isfield(estim_params,'calibrated_covariances')
t(1)=false;
else
t(1)=true;
end
M_.Sigma_e=[1 -0.1; -0.1 1];
M_.H=[1 -0.1; -0.1 1];
M_.Correlation_matrix= [1 -0.5; -0.5 1];
M_.Correlation_matrix_ME=[1 0; 0 1];
estim_params.ncx=1;
estim_params.ncn=0;
estim_params.corrx=[2 1 NaN -1 1 3 0 0.2000 NaN NaN NaN];
estim_params.corrn=[];
estim_params=check_for_calibrated_covariances(estim_params,M_);
t(2)=isequal(estim_params.calibrated_covariances_ME.position,[2;3]);
t(3)=isequal(estim_params.calibrated_covariances_ME.cov_value,[-0.1;-0.1]);
M_.Sigma_e=[1 -0.1; -0.1 1];
M_.H=[1 -0.1; -0.1 1];
M_.Correlation_matrix= [1 -0.5; -0.5 1];
M_.Correlation_matrix_ME=[1 0; 0 1];
estim_params.ncx=1;
estim_params.ncn=1;
estim_params.corrx=[2 1 NaN -1 1 3 0 0.2000 NaN NaN NaN];
estim_params.corrn=[2 1 NaN -1 1 3 0 0.2000 NaN NaN NaN];
estim_params=check_for_calibrated_covariances(estim_params,M_);
if isfield(estim_params,'calibrated_covariances_ME') || isfield(estim_params,'calibrated_covariances')
t(4)=false;
else
t(4)=true;
end
T = all(t);
%@eof:1

View File

@ -1,6 +1,8 @@
function check_model(DynareModel)
function check_model(M_)
% check_model(M_)
% Performs various consistency checks on the model
% Copyright © 2005-2013 Dynare Team
% Copyright (C) 2005-2033 Dynare Team
%
% This file is part of Dynare.
%
@ -17,15 +19,21 @@ function check_model(DynareModel)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
if ~all(DynareModel.lead_lag_incidence(DynareModel.maximum_lag+1,:) > 0)
if ~all(M_.lead_lag_incidence(M_.maximum_lag+1,:) > 0)
warning('Problem in the model specification: some variables don''t appear as current. Check whether this is desired.');
end
if ~check_consistency_covariances(DynareModel.Sigma_e)
if any(diag(M_.Sigma_e)<0)
error('You have specified negative shock variances. That is not allowed.')
end
if ~check_consistency_covariances(M_.Sigma_e)
error('The specified covariances for the structural errors are not consistent with the variances as they imply a correlation larger than +-1')
end
if ~isequal(DynareModel.H,0)
if ~check_consistency_covariances(DynareModel.H)
if any(diag(M_.H)<0)
error('You have specified negative measurement error variances. That is not allowed.')
end
if ~isequal(M_.H,0)
if ~check_consistency_covariances(M_.H)
error('The specified covariances for the measurement errors are not consistent with the variances as they imply a correlation larger than +-1')
end
end

View File

@ -0,0 +1,10 @@
{
"_schemaVersion": "1.0.0",
"prior":
{
"inputs":
[
{"name":"options", "kind":"flag", "type":["char", "choices={'table', 'moments', 'moments(distribution)', 'optimize', 'simulate', 'plot'}"], "repeating":true}
]
}
}

View File

@ -135,7 +135,7 @@ if ismember('moments', varargin) % Prior simulations (2nd order moments).
check_model(Model);
% Compute state space representation of the model.
oo__ = oo_;
oo__.dr = set_state_space(oo__.dr, Model, options_);
oo__.dr = set_state_space(oo__.dr, Model);
% Solve model
[T,R,~,info,oo__.dr, Model.params] = dynare_resolve(Model , options_ , oo__.dr, oo__.steady_state, oo__.exo_steady_state, oo__.exo_det_steady_state,'restrict');
if ~info(1)

View File

@ -171,7 +171,11 @@ if options_.order==1
skipline();
end
if ~all(diag(M_.H)==0)
[observable_name_requested_vars, varlist_pos] = intersect(var_list_, options_.varobs, 'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[observable_name_requested_vars, varlist_pos] = intersect_stable(var_list_, options_.varobs);
else
[observable_name_requested_vars, varlist_pos] = intersect(var_list_, options_.varobs, 'stable');
end
if ~isempty(observable_name_requested_vars)
NumberOfObservedEndogenousVariables = length(observable_name_requested_vars);
temp = NaN(NumberOfObservedEndogenousVariables, NumberOfExogenousVariables+1);

View File

@ -1,4 +1,5 @@
function [ConditionalVarianceDecomposition, ConditionalVarianceDecomposition_ME]= conditional_variance_decomposition(M_,options_,dr, Steps, SubsetOfVariables)
% [ConditionalVarianceDecomposition, ConditionalVarianceDecomposition_ME]= conditional_variance_decomposition(M_,options_,dr, Steps, SubsetOfVariables)
% This function computes the conditional variance decomposition of a given state space model
% for a subset of endogenous variables.
%
@ -88,7 +89,11 @@ end
% Measurement error
if ~all(diag(M_.H)==0)
[observable_pos,index_subset,index_observables]=intersect(SubsetOfVariables,options_.varobs_id,'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[observable_pos,index_subset,index_observables]=intersect_stable(SubsetOfVariables,options_.varobs_id);
else
[observable_pos,index_subset,index_observables]=intersect(SubsetOfVariables,options_.varobs_id,'stable');
end
ME_Variance=diag(M_.H);
ConditionalVarianceDecomposition_ME = zeros(length(observable_pos),length(Steps),number_of_state_innovations+1);

View File

@ -1,5 +1,6 @@
function oo_ = ...
conditional_variance_decomposition_ME_mc_analysis(NumberOfSimulations, type, dname, fname, Steps, exonames, exo, var_list, endo, mh_conf_sig, oo_,options_)
%oo_ = conditional_variance_decomposition_ME_mc_analysis(NumberOfSimulations, type, dname, fname, Steps, exonames, exo, var_list, endo, mh_conf_sig, oo_,options_)
% This function analyses the (posterior or prior) distribution of the
% endogenous variables' conditional variance decomposition with measurement error.
%
@ -22,7 +23,7 @@ function oo_ = ...
% OUTPUTS
% oo_ [structure] Dynare structure where the results are saved.
% Copyright © 2017-2020 Dynare Team
% Copyright © 2017-2023 Dynare Team
%
% This file is part of Dynare.
%
@ -63,7 +64,12 @@ if isempty(exogenous_variable_index)
end
end
[observable_pos_requested_vars,index_subset,index_observables]=intersect(var_list,options_.varobs,'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[~,index_subset]=intersect_stable(var_list,options_.varobs);
else
[~,index_subset]=intersect(var_list,options_.varobs,'stable');
end
matrix_pos=strmatch(endo, var_list(index_subset),'exact');
name_1 = endo;
name_2 = exo;

View File

@ -164,6 +164,7 @@ options_.one_sided_hp_filter = 0;
options_.filtered_theoretical_moments_grid = 512;
options_.nodecomposition = false;
options_.nomoments = false;
options_.nomodelsummary = false;
options_.nocorr = false;
options_.periods = 0;
options_.noprint = false;

View File

@ -116,7 +116,7 @@ end
%write back solution to dr
dr.ys =ys;
dr=set_state_space(dr,M_,options_);
dr=set_state_space(dr,M_);
T=H(dr.order_var,dr.order_var);
dr.ghu=G(dr.order_var,:);
if M_.maximum_endo_lag

View File

@ -11,7 +11,7 @@ function oo_=disp_moments(y,var_list,M_,options_,oo_)
% OUTPUTS
% oo_ [structure] Dynare's results structure,
% Copyright © 2001-2021 Dynare Team
% Copyright © 2001-2023 Dynare Team
%
% This file is part of Dynare.
%
@ -50,7 +50,11 @@ y = y(ivar,options_.drop+1:end)';
ME_present=0;
if ~all(M_.H==0)
[observable_pos_requested_vars, index_subset, index_observables] = intersect(ivar, options_.varobs_id, 'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[observable_pos_requested_vars, index_subset, index_observables] = intersect_stable(ivar, options_.varobs_id);
else
[observable_pos_requested_vars, index_subset, index_observables] = intersect(ivar, options_.varobs_id, 'stable');
end
if ~isempty(observable_pos_requested_vars)
ME_present=1;
i_ME = setdiff([1:size(M_.H,1)],find(diag(M_.H) == 0)); % find ME with 0 variance

View File

@ -122,6 +122,11 @@ if isoctave || matlab_ver_less_than('8.4')
p{end+1} = '/missing/datetime';
end
% intersect with 'stable' flag is broken before Octave 8.4, bug #60347
if isoctave && octave_ver_less_than('8.4')
p{end+1} = '/missing/intersect_stable';
end
P = cellfun(@(c)[dynareroot(1:end-1) c], p, 'uni',false);
% Get mex files folder(s)

View File

@ -163,7 +163,7 @@ end
%check for calibrated covariances before updating parameters
if ~isempty(estim_params_) && ~(isfield(estim_params_,'nvx') && sum(estim_params_.nvx+estim_params_.nvn+estim_params_.ncx+estim_params_.ncn+estim_params_.np)==0)
estim_params_=check_for_calibrated_covariances(xparam1,estim_params_,M_);
estim_params_=check_for_calibrated_covariances(estim_params_,M_);
end
%%read out calibration that was set in mod-file and can be used for initialization
@ -268,7 +268,7 @@ else% Yes!
end
% Get informations about the variables of the model.
dr = set_state_space(oo_.dr,M_,options_);
dr = set_state_space(oo_.dr,M_);
oo_.dr = dr;
nstatic = M_.nstatic; % Number of static variables.
npred = M_.nspred; % Number of predetermined variables.

View File

@ -77,7 +77,7 @@ DynareOptions.stack_solve_algo = ep.stack_solve_algo;
dr = struct();
if ep.init
DynareOptions.order = 1;
DynareResults.dr=set_state_space(dr,DynareModel,DynareOptions);
DynareResults.dr=set_state_space(dr,DynareModel);
[DynareResults.dr,Info,DynareModel.params] = resol(0,DynareModel,DynareOptions,DynareResults.dr,DynareResults.steady_state, DynareResults.exo_steady_state, DynareResults.exo_det_steady_state);
end
@ -103,7 +103,7 @@ end
% hybrid correction
pfm.hybrid_order = ep.stochastic.hybrid_order;
if pfm.hybrid_order
DynareResults.dr = set_state_space(DynareResults.dr, DynareModel, DynareOptions);
DynareResults.dr = set_state_space(DynareResults.dr, DynareModel);
options = DynareOptions;
options.order = pfm.hybrid_order;
[pfm.dr, DynareModel.params] = resol(0, DynareModel, options, DynareResults.dr, DynareResults.steady_state, DynareResults.exo_steady_state, DynareResults.exo_det_steady_state);

View File

@ -0,0 +1,27 @@
{
"_schemaVersion": "1.0.0",
"dynare":
{
"inputs":
[
{"name":"fname", "kind":"required", "type":["file=*.mod,*.dyn"], "purpose":"Model file"},
{"name":"options", "kind":"flag", "type":["char", "choices={'debug', 'noclearall', 'onlyclearglobals', 'savemacro[=macro_file]', 'onlymacro', 'linemacro', 'notmpterms', 'nolog', 'warn_uninit', 'console', 'nograph', 'nointeractive', 'parallel[=cluster_name]', 'conffile=parallel_config_path_and_filename', 'parallel_slave_open_mode', 'parallel_test', 'parallel_use_psexec=true|false', '-D<variable>[=<value>]', '-I/path', 'nostrict', 'stochastic', 'fast', 'minimal_workspace', 'compute_xrefs', 'output=second|third', 'params_derivs_order=0|1|2', 'transform_unary_ops', 'exclude_eqs=<equation_tag_list_or_file>', 'include_eqs=<equation_tag_list_or_file>', 'json=parse|check|transform|compute', 'jsonstdout', 'onlyjson', 'jsonderivsimple', 'nopathchange', 'nopreprocessoroutput', 'mexext=<extension>', 'matlabroot=<path>', 'onlymodel', 'notime', 'use_dll', 'nocommutativity'}"], "repeating":true}
]
},
"internals":
{
"inputs":
[
{"name":"flag", "kind":"required", "type":["char", "choices={'--test'}"], "purpose":"Test routine"},
{"name":"fname", "kind":"required", "type":["file=*.m"], "purpose":"Routine name"}
]
},
"internals":
{
"inputs":
[
{"name":"flag", "kind":"required", "type":["char", "choices={'--display-mh-history','--load-mh-history'}"], "purpose":"Model file"},
{"name":"fname", "kind":"required", "type":["file=*.mod,*.dyn"], "purpose":"Model name"}
]
}
}

View File

@ -109,7 +109,7 @@ if strcmpi(flag,'--load-mh-history') || strcmpi(flag,'--display-mh-history')
assignin('caller','mcmc_informations',o);
else
oo = load_first_mh_history_file([dname filesep 'metropolis'],fname);
local = load([fname '_results'],'bayestopt_');
local = load([dname filesep 'Output' filesep fname '_results'],'bayestopt_');
names = local.bayestopt_.name; %evalin('base','bayestopt_.name');
str = ['MCMC set-up for ' fname ' mod file'];
ltr = length(str);

View File

@ -30,8 +30,8 @@ function [fval,info,exit_flag,fake_1,fake_2] = minus_logged_prior_density(xparam
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
fake_1 = 1;
fake_2 = 1;
fake_1 = [];
fake_2 = [];
exit_flag = 1;
info = zeros(4,1);

View File

@ -0,0 +1,91 @@
function [c, ia, ib] = intersect_stable(a, b)
% Crude implementation of intersect(…, 'stable'), which is missing before
% Octave 6.2 and buggy before 8.4 (#60347)
% Copyright (C) 2019-2023 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/>.
if nargin ~= 2
error('intersect_stable: needs exactly 2 input arguments');
end
if isnumeric (a) && isnumeric (b)
c = [];
elseif iscell (b)
c = {};
else
c = '';
end
ia = [];
ib = [];
if isempty (a) || isempty (b)
return
else
isrowvec = isrow (a) && isrow (b);
for i = 1:numel(a)
if iscellstr(c)
idx = strcmp(a(i), b);
else
idx = a(i) == b;
end
if any(idx) && ~ismember(a(i), c)
c = [c(:); a(i)];
if nargout > 1
ia = [ia, i];
ib = [ib, find(idx)];
end
end
end
%% Adjust output orientation for MATLAB compatibility
if isrowvec
c = c.';
end
end
end
%!test
%! a = [3 4 1 5];
%! b = [2 4 9 1 6];
%! [c,ia,ib]=intersect_stable(a,b);
%! assert(c, [4 1])
%! assert(ia, [2 3])
%! assert(ib, [2 4])
%! assert(a(ia), c)
%! assert(b(ib), c)
%!test
%! a = [3 4 1 5]';
%! b = [2 4 9 1 6]';
%! [c,ia,ib]=intersect_stable(a,b);
%! assert(c, [4 1]')
%! assert(ia, [2 3])
%! assert(ib, [2 4])
%! assert(a(ia), c)
%! assert(b(ib), c)
%!test
%! a = { 'defun', 'mapcar', 'let', 'eval-when'};
%! b = { 'setf', 'let', 'list', 'cdr', 'defun'};
%! [c,ia,ib]=intersect_stable(a,b);
%! assert(c, { 'defun', 'let' })
%! assert(ia, [1 3])
%! assert(ib, [5 2])
%! assert(a(ia), c)
%! assert(b(ib), c)

View File

@ -1,10 +1,10 @@
function model_diagnostics(M,options,oo)
% function model_diagnostics(M,options,oo)
function model_diagnostics(M_,options_,oo_)
% function model_diagnostics(M_,options_,oo_)
% computes various diagnostics on the model
% INPUTS
% M [matlab structure] Definition of the model.
% options [matlab structure] Global options.
% oo [matlab structure] Results
% M_ [matlab structure] Definition of the model.
% options_ [matlab structure] Global options.
% oo_ [matlab structure] Results
%
% OUTPUTS
% none
@ -33,22 +33,22 @@ function model_diagnostics(M,options,oo)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
endo_names = M.endo_names;
lead_lag_incidence = M.lead_lag_incidence;
maximum_endo_lag = M.maximum_endo_lag;
endo_names = M_.endo_names;
lead_lag_incidence = M_.lead_lag_incidence;
maximum_endo_lag = M_.maximum_endo_lag;
if options.ramsey_policy
if options_.ramsey_policy
%test whether specification matches
inst_nbr = size(options.instruments,1);
inst_nbr = size(options_.instruments,1);
if inst_nbr~=0
implied_inst_nbr = M.ramsey_orig_endo_nbr - M.ramsey_orig_eq_nbr;
implied_inst_nbr = M_.ramsey_orig_endo_nbr - M_.ramsey_orig_eq_nbr;
if inst_nbr>implied_inst_nbr
warning('You have specified more steady state instruments than there are omitted equations. While there are use cases for this setup, it is rather unusual. Check whether this is desired.')
elseif inst_nbr<implied_inst_nbr
warning('You have specified fewer steady state instruments than there are omitted equations. While there are use cases for this setup, it is rather unusual. Check whether this is desired.')
end
else
if options.steadystate_flag
if options_.steadystate_flag
warning('You have specified a steady state file, but not provided steady state instruments. In this case, you typically need to make sure to provide all steady state values, including the ones for the planner''s instrument(s).')
end
end
@ -57,20 +57,20 @@ end
problem_dummy=0;
%naming conflict in steady state file
if options.steadystate_flag == 1
if strmatch('ys',M.endo_names,'exact')
if options_.steadystate_flag == 1
if strmatch('ys',M_.endo_names,'exact')
disp(['MODEL_DIAGNOSTICS: using the name ys for an endogenous variable will typically conflict with the internal naming in user-defined steady state files.'])
problem_dummy=1;
end
if strmatch('ys',M.param_names,'exact')
if strmatch('ys',M_.param_names,'exact')
disp(['MODEL_DIAGNOSTICS: using the name ys for a parameter will typically conflict with the internal naming in user-defined steady state files.'])
problem_dummy=1;
end
if strmatch('M_',M.endo_names,'exact')
if strmatch('M_',M_.endo_names,'exact')
disp(['MODEL_DIAGNOSTICS: using the name M_ for an endogenous variable will typically conflict with the internal naming in user-defined steady state files.'])
problem_dummy=1;
end
if strmatch('M_',M.param_names,'exact')
if strmatch('M_',M_.param_names,'exact')
disp(['MODEL_DIAGNOSTICS: using the name M_ for a parameter will typically conflict with the internal naming in user-defined steady state files.'])
problem_dummy=1;
end
@ -94,24 +94,24 @@ end
%
info = 0;
if M.exo_nbr == 0
oo.exo_steady_state = [] ;
if M_.exo_nbr == 0
oo_.exo_steady_state = [] ;
end
info=test_for_deep_parameters_calibration(M);
info=test_for_deep_parameters_calibration(M_);
if info
problem_dummy=1;
end
% check if ys is steady state
options.debug=true; %locally set debug option to true
if options.logged_steady_state %if steady state was previously logged, undo this
oo.dr.ys=exp(oo.dr.ys);
oo.steady_state=exp(oo.steady_state);
options.logged_steady_state=0;
options_.debug=true; %locally set debug option to true
if options_.logged_steady_state %if steady state was previously logged, undo this
oo_.dr.ys=exp(oo_.dr.ys);
oo_.steady_state=exp(oo_.steady_state);
options_.logged_steady_state=0;
end
[dr.ys,M.params,check1]=evaluate_steady_state(oo.steady_state,[oo.exo_steady_state; oo.exo_det_steady_state],M,options,options.steadystate.nocheck);
[dr.ys,M_.params,check1]=evaluate_steady_state(oo_.steady_state,[oo_.exo_steady_state; oo_.exo_det_steady_state],M_,options_,options_.steadystate.nocheck);
% testing for problem
if check1(1)
@ -137,35 +137,35 @@ end
% singular Jacobian of static model
%
singularity_problem = 0;
if ~options.block
if ~options_.block
nb = 1;
else
nb = length(M.block_structure_stat.block);
nb = length(M_.block_structure_stat.block);
end
exo = [oo.exo_steady_state; oo.exo_det_steady_state];
exo = [oo_.exo_steady_state; oo_.exo_det_steady_state];
for b=1:nb
if options.bytecode
if options_.bytecode
if nb == 1
[res, jacob] = bytecode(dr.ys, exo, M.params, dr.ys, 1, exo, ...
[res, jacob] = bytecode(dr.ys, exo, M_.params, dr.ys, 1, exo, ...
'evaluate', 'static');
else
[res, jacob] = bytecode(dr.ys, exo, M.params, dr.ys, 1, exo, ...
[res, jacob] = bytecode(dr.ys, exo, M_.params, dr.ys, 1, exo, ...
'evaluate', 'static', 'block_decomposed', ['block=' ...
int2str(b)]);
end
else
if options.block
T = NaN(M.block_structure_stat.tmp_nbr, 1);
fh_static = str2func(sprintf('%s.sparse.block.static_%d', M.fname, b));
[~, ~,~, jacob] = fh_static(dr.ys, exo, M.params, M.block_structure_stat.block(b).g1_sparse_rowval, ...
M.block_structure_stat.block(b).g1_sparse_colval, ...
M.block_structure_stat.block(b).g1_sparse_colptr, T);
if options_.block
T = NaN(M_.block_structure_stat.tmp_nbr, 1);
fh_static = str2func(sprintf('%s.sparse.block.static_%d', M_.fname, b));
[~, ~,~, jacob] = fh_static(dr.ys, exo, M_.params, M_.block_structure_stat.block(b).g1_sparse_rowval, ...
M_.block_structure_stat.block(b).g1_sparse_colval, ...
M_.block_structure_stat.block(b).g1_sparse_colptr, T);
n_vars_jacob=size(jacob,2);
else
[res, T_order, T] = feval([M.fname '.sparse.static_resid'], dr.ys, exo, M.params);
jacob = feval([M.fname '.sparse.static_g1'], dr.ys, exo, M.params, M.static_g1_sparse_rowval, M.static_g1_sparse_colval, M.static_g1_sparse_colptr, T_order, T);
n_vars_jacob=M.endo_nbr;
[res, T_order, T] = feval([M_.fname '.sparse.static_resid'], dr.ys, exo, M_.params);
jacob = feval([M_.fname '.sparse.static_g1'], dr.ys, exo, M_.params, M_.static_g1_sparse_rowval, M_.static_g1_sparse_colval, M_.static_g1_sparse_colptr, T_order, T);
n_vars_jacob=M_.endo_nbr;
end
jacob=full(jacob);
end
@ -173,19 +173,19 @@ for b=1:nb
problem_dummy=1;
[infrow,infcol]=find(isinf(jacob) | isnan(jacob));
fprintf('\nMODEL_DIAGNOSTICS: The Jacobian of the static model contains Inf or NaN. The problem arises from: \n\n')
display_problematic_vars_Jacobian(infrow,infcol,M,dr.ys,'static','MODEL_DIAGNOSTICS: ')
display_problematic_vars_Jacobian(infrow,infcol,M_,dr.ys,'static','MODEL_DIAGNOSTICS: ')
end
if any(any(~isreal(jacob)))
problem_dummy=1;
[imagrow,imagcol]=find(abs(imag(jacob))>1e-15);
fprintf('\nMODEL_DIAGNOSTICS: The Jacobian of the static model contains imaginary parts. The problem arises from: \n\n')
display_problematic_vars_Jacobian(imagrow,imagcol,M,dr.ys,'static','MODEL_DIAGNOSTICS: ')
display_problematic_vars_Jacobian(imagrow,imagcol,M_,dr.ys,'static','MODEL_DIAGNOSTICS: ')
end
try
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options.jacobian_tolerance)
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options_.jacobian_tolerance)
rank_jacob = rank(jacob); %can sometimes fail
else
rank_jacob = rank(jacob,options.jacobian_tolerance); %can sometimes fail
rank_jacob = rank(jacob,options_.jacobian_tolerance); %can sometimes fail
end
catch
rank_jacob=size(jacob,1);
@ -197,10 +197,10 @@ for b=1:nb
'singular'])
disp(['MODEL_DIAGNOSTICS: there is ' num2str(n_vars_jacob-rank_jacob) ...
' collinear relationships between the variables and the equations'])
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options.jacobian_tolerance)
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options_.jacobian_tolerance)
ncol = null(jacob);
else
ncol = null(jacob,options.jacobian_tolerance); %can sometimes fail
ncol = null(jacob,options_.jacobian_tolerance); %can sometimes fail
end
n_rel = size(ncol,2);
for i = 1:n_rel
@ -214,16 +214,16 @@ for b=1:nb
break
end
end
if options.block && ~options.bytecode
fprintf('%s\n',endo_names{M.block_structure_stat.block(b).variable(k)})
if options_.block && ~options_.bytecode
fprintf('%s\n',endo_names{M_.block_structure_stat.block(b).variable(k)})
else
fprintf('%s\n',endo_names{k})
end
end
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options.jacobian_tolerance)
if (~isoctave && matlab_ver_less_than('9.12')) || isempty(options_.jacobian_tolerance)
neq = null(jacob'); %can sometimes fail
else
neq = null(jacob',options.jacobian_tolerance); %can sometimes fail
neq = null(jacob',options_.jacobian_tolerance); %can sometimes fail
end
n_rel = size(neq,2);
for i = 1:n_rel
@ -237,8 +237,8 @@ for b=1:nb
break
end
end
if options.block && ~options.bytecode
disp(M.block_structure_stat.block(b).equation(k))
if options_.block && ~options_.bytecode
disp(M_.block_structure_stat.block(b).equation(k))
else
disp(k')
end
@ -248,9 +248,9 @@ end
if singularity_problem
try
options_check=options;
options_check=options_;
options_check.noprint=1;
[eigenvalues_] = check(M, options_check, oo);
[eigenvalues_] = check(M_, options_check, oo_);
if any(abs(abs(eigenvalues_)-1)<1e-6)
fprintf('MODEL_DIAGNOSTICS: The singularity seems to be (partly) caused by the presence of a unit root\n')
fprintf('MODEL_DIAGNOSTICS: as the absolute value of one eigenvalue is in the range of +-1e-6 to 1.\n')
@ -265,32 +265,32 @@ if singularity_problem
end
%%check dynamic Jacobian
klen = M.maximum_lag + M.maximum_lead + 1;
exo_simul = [repmat(oo.exo_steady_state',klen,1) repmat(oo.exo_det_steady_state',klen,1)];
iyv = M.lead_lag_incidence';
klen = M_.maximum_lag + M_.maximum_lead + 1;
exo_simul = [repmat(oo_.exo_steady_state',klen,1) repmat(oo_.exo_det_steady_state',klen,1)];
iyv = M_.lead_lag_incidence';
iyv = iyv(:);
iyr0 = find(iyv) ;
it_ = M.maximum_lag + 1;
it_ = M_.maximum_lag + 1;
z = repmat(dr.ys,1,klen);
if options.order == 1
if (options.bytecode)
if options_.order == 1
if (options_.bytecode)
[~, loc_dr] = bytecode('dynamic','evaluate', z,exo_simul, ...
M.params, dr.ys, 1);
M_.params, dr.ys, 1);
jacobia_ = [loc_dr.g1 loc_dr.g1_x loc_dr.g1_xd];
else
[~,jacobia_] = feval([M.fname '.dynamic'],z(iyr0),exo_simul, ...
M.params, dr.ys, it_);
[~,jacobia_] = feval([M_.fname '.dynamic'],z(iyr0),exo_simul, ...
M_.params, dr.ys, it_);
end
elseif options.order >= 2
if (options.bytecode)
elseif options_.order >= 2
if (options_.bytecode)
[~, loc_dr] = bytecode('dynamic','evaluate', z,exo_simul, ...
M.params, dr.ys, 1);
M_.params, dr.ys, 1);
jacobia_ = [loc_dr.g1 loc_dr.g1_x];
else
[~,jacobia_,hessian1] = feval([M.fname '.dynamic'],z(iyr0),...
[~,jacobia_,hessian1] = feval([M_.fname '.dynamic'],z(iyr0),...
exo_simul, ...
M.params, dr.ys, it_);
M_.params, dr.ys, it_);
end
end
@ -298,14 +298,14 @@ if any(any(isinf(jacobia_) | isnan(jacobia_)))
problem_dummy=1;
[infrow,infcol]=find(isinf(jacobia_) | isnan(jacobia_));
fprintf('\nMODEL_DIAGNOSTICS: The Jacobian of the dynamic model contains Inf or NaN. The problem arises from: \n\n')
display_problematic_vars_Jacobian(infrow,infcol,M,dr.ys,'dynamic','MODEL_DIAGNOSTICS: ')
display_problematic_vars_Jacobian(infrow,infcol,M_,dr.ys,'dynamic','MODEL_DIAGNOSTICS: ')
end
if any(any(~isreal(jacobia_)))
[imagrow,imagcol]=find(abs(imag(jacobia_))>1e-15);
if ~isempty(imagrow)
problem_dummy=1;
fprintf('\nMODEL_DIAGNOSTICS: The Jacobian of the dynamic model contains imaginary parts. The problem arises from: \n\n')
display_problematic_vars_Jacobian(imagrow,imagcol,M,dr.ys,'dynamic','MODEL_DIAGNOSTICS: ')
display_problematic_vars_Jacobian(imagrow,imagcol,M_,dr.ys,'dynamic','MODEL_DIAGNOSTICS: ')
end
end
if exist('hessian1','var')

View File

@ -38,7 +38,11 @@ observable_pos_requested_vars=[];
ME_present=false;
if ~all(diag(M_.H)==0)
[observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(ivar,options_.varobs_id);
else
[observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable');
end
if ~isempty(observable_pos_requested_vars)
ME_present=true;
end

View File

@ -45,7 +45,11 @@ if M_.exo_nbr > 1
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, title, headers, labels, 100*oo_.gamma_y{options_.ar+2}(stationary_vars,:), lh, 8, 2);
if ME_present
[stationary_observables, pos_index_subset] = intersect(index_subset, stationary_vars, 'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
[stationary_observables, pos_index_subset] = intersect_stable(index_subset, stationary_vars);
else
[stationary_observables, pos_index_subset] = intersect(index_subset, stationary_vars, 'stable');
end
headers_ME = vertcat(headers, 'ME');
labels=get_labels_transformed_vars(M_.endo_names,ivar(stationary_observables),options_,false);
dyntable(options_, [title,' WITH MEASUREMENT ERROR'], headers_ME, labels, ...

View File

@ -82,11 +82,12 @@ if ischar(func0)
end
[fval0,exit_flag,gg,hh]=penalty_objective_function(x,func0,penalty,varargin{:});
fval=fval0;
if ~exit_flag
igg=NaN(nx);
disp_verbose('Bad initial parameter.',Verbose)
return
end
fval=fval0;
% initialize mr_gstep and mr_hessian

View File

@ -20,7 +20,7 @@ function optimize_prior(DynareOptions, ModelInfo, DynareResults, BayesInfo, Esti
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
% Initialize to the prior mean
DynareResults.dr = set_state_space(DynareResults.dr,ModelInfo,DynareOptions);
DynareResults.dr = set_state_space(DynareResults.dr,ModelInfo);
xparam1 = BayesInfo.p1;
% Pertubation of the initial condition.

View File

@ -157,7 +157,7 @@ if options_.debug
save([M_.dname filesep 'Output' filesep M_.fname '_debug.mat'],'jacobia_')
end
dr=set_state_space(dr,M_,options_);
dr=set_state_space(dr,M_);
kstate = dr.kstate;
nstatic = M_.nstatic;
nfwrd = M_.nfwrd;

View File

@ -45,7 +45,7 @@ end
if ~isfield(oo_,'dr') || ~isfield(oo_.dr,'ghx')
fprintf('computing the first order solution of the model as initial guess...');
dr = struct();
oo_.dr=set_state_space(dr,M_,options_);
oo_.dr=set_state_space(dr,M_);
options_.order = 1;
[oo_.dr,Info,M_.params] = resol(0,M_,options_,oo_.dr ,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
fprintf('done\n');

View File

@ -66,7 +66,11 @@ switch type
M_.exo_names,arg2,vartan,arg1,options_.mh_conf_sig,oo_,options_);
if ~all(diag(M_.H)==0)
if strmatch(arg1,options_.varobs,'exact')
[observable_name_requested_vars,index_subset,index_observables]=intersect(vartan,options_.varobs,'stable');
if isoctave && octave_ver_less_than('8.4') %Octave bug #60347
observable_name_requested_vars=intersect_stable(vartan,options_.varobs);
else
observable_name_requested_vars=intersect(vartan,options_.varobs,'stable');
end
oo_ = variance_decomposition_ME_mc_analysis(SampleSize,'posterior',M_.dname,M_.fname,...
[M_.exo_names;'ME'],arg2,observable_name_requested_vars,arg1,options_.mh_conf_sig,oo_,options_);
end

View File

@ -83,7 +83,7 @@ sampled_prior_covariance = zeros(NumberOfParameters,NumberOfParameters);
file_line_number = 0;
file_indx_number = 0;
oo_.dr=set_state_space(oo_.dr,M_,options_);
oo_.dr=set_state_space(oo_.dr,M_);
hh_fig = dyn_waitbar(0,'Please wait. Prior sampler...');
set(hh_fig,'Name','Prior sampler.');

View File

@ -55,8 +55,10 @@ 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;
if nvx || ncx
Sigma_e = M.Sigma_e;
Correlation_matrix = M.Correlation_matrix;
end
H = M.H;
Correlation_matrix_ME = M.Correlation_matrix_ME;
% setting shocks variance on the diagonal of Covariance matrix; used later
@ -93,13 +95,6 @@ if ncx
Correlation_matrix(k2,k1) = Correlation_matrix(k1,k2);
end
end
%build covariance matrix from correlation matrix and variances already on
%diagonal
Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e)));
%if calibrated covariances, set them now to their stored value
if isfield(estim_params,'calibrated_covariances')
Sigma_e(estim_params.calibrated_covariances.position)=estim_params.calibrated_covariances.cov_value;
end
% update offset
offset = nvx+nvn+ncx;
@ -113,13 +108,6 @@ if ncn
Correlation_matrix_ME(k2,k1) = Correlation_matrix_ME(k1,k2);
end
end
%build covariance matrix from correlation matrix and variances already on
%diagonal
H = diag(sqrt(diag(H)))*Correlation_matrix_ME*diag(sqrt(diag(H)));
%if calibrated covariances, set them now to their stored value
if isfield(estim_params,'calibrated_covariances_ME')
H(estim_params.calibrated_covariances_ME.position)=estim_params.calibrated_covariances_ME.cov_value;
end
% update offset
offset = nvx+ncx+nvn+ncn;
@ -131,10 +119,24 @@ end
% updating matrices in M
if nvx || ncx
%build covariance matrix from correlation matrix and variances already on
%diagonal
Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e)));
%if calibrated covariances, set them now to their stored value
if isfield(estim_params,'calibrated_covariances')
Sigma_e(estim_params.calibrated_covariances.position)=estim_params.calibrated_covariances.cov_value;
end
M.Sigma_e = Sigma_e;
M.Correlation_matrix=Correlation_matrix;
end
if nvn || ncn
%build covariance matrix from correlation matrix and variances already on
%diagonal
H = diag(sqrt(diag(H)))*Correlation_matrix_ME*diag(sqrt(diag(H)));
%if calibrated covariances, set them now to their stored value
if isfield(estim_params,'calibrated_covariances_ME')
H(estim_params.calibrated_covariances_ME.position)=estim_params.calibrated_covariances_ME.cov_value;
end
M.H = H;
M.Correlation_matrix_ME=Correlation_matrix_ME;
end

View File

@ -1,8 +1,9 @@
function dr=set_state_space(dr,DynareModel,DynareOptions)
function dr=set_state_space(dr,M_)
% dr=set_state_space(dr,M_)
% Write the state space representation of the reduced form solution.
%@info:
%! @deftypefn {Function File} {[@var{dr} =} set_state_space (@var{dr},@var{DynareModel},@var{DynareOptions})
%! @deftypefn {Function File} {[@var{dr} =} set_state_space (@var{dr},@var{M_})
%! @anchor{set_state_space}
%! @sp 1
%! Write the state space representation of the reduced form solution.
@ -12,10 +13,8 @@ function dr=set_state_space(dr,DynareModel,DynareOptions)
%! @table @ @var
%! @item dr
%! Matlab's structure describing decision and transition rules.
%! @item DynareModel
%! @item M_
%! Matlab's structure describing the model (initialized by dynare, see @ref{M_})
%! @item DynareOptions
%! Matlab's structure describing the current options (initialized by dynare, see @ref{options_}).
%! @end table
%! @sp 2
%! @strong{Outputs}
@ -51,10 +50,10 @@ function dr=set_state_space(dr,DynareModel,DynareOptions)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
max_lead = DynareModel.maximum_endo_lead;
max_lag = DynareModel.maximum_endo_lag;
endo_nbr = DynareModel.endo_nbr;
lead_lag_incidence = DynareModel.lead_lag_incidence;
max_lead = M_.maximum_endo_lead;
max_lag = M_.maximum_endo_lag;
endo_nbr = M_.endo_nbr;
lead_lag_incidence = M_.lead_lag_incidence;
klen = max_lag + max_lead + 1;
fwrd_var = find(lead_lag_incidence(max_lag+2:end,:))';

View File

@ -58,7 +58,7 @@ if M_.exo_nbr > 0
oo_.exo_simul= ones(max(options_.periods,1) + M_.maximum_lag + M_.maximum_lead,1) * oo_.exo_steady_state';
end
oo_.dr=set_state_space(oo_.dr,M_,options_);
oo_.dr=set_state_space(oo_.dr,M_);
if options_.logged_steady_state %if steady state was previously logged, undo this

View File

@ -103,7 +103,7 @@ end
check_model(M_);
oo_.dr=set_state_space(oo_.dr,M_,options_);
oo_.dr=set_state_space(oo_.dr,M_);
if PI_PCL_solver
[oo_.dr, info] = PCL_resol(oo_.steady_state,0);
@ -134,36 +134,38 @@ if info(1)
end
if ~options_.noprint
skipline()
disp('MODEL SUMMARY')
skipline()
disp([' Number of variables: ' int2str(M_.endo_nbr)])
disp([' Number of stochastic shocks: ' int2str(M_.exo_nbr)])
disp([' Number of state variables: ' int2str(M_.nspred)])
disp([' Number of jumpers: ' int2str(M_.nsfwrd)])
disp([' Number of static variables: ' int2str(M_.nstatic)])
my_title='MATRIX OF COVARIANCE OF EXOGENOUS SHOCKS';
labels = M_.exo_names;
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, my_title, headers, labels, M_.Sigma_e, lh, 10, 6);
if options_.TeX
labels = M_.exo_names_tex;
if ~options_.nomodelsummary
skipline()
disp('MODEL SUMMARY')
skipline()
disp([' Number of variables: ' int2str(M_.endo_nbr)])
disp([' Number of stochastic shocks: ' int2str(M_.exo_nbr)])
disp([' Number of state variables: ' int2str(M_.nspred)])
disp([' Number of jumpers: ' int2str(M_.nsfwrd)])
disp([' Number of static variables: ' int2str(M_.nstatic)])
my_title='MATRIX OF COVARIANCE OF EXOGENOUS SHOCKS';
labels = M_.exo_names;
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, my_title, 'covar_ex_shocks', headers, labels, M_.Sigma_e, lh, 10, 6);
end
if ~all(diag(M_.H)==0)
my_title='MATRIX OF COVARIANCE OF MEASUREMENT ERRORS';
labels = cellfun(@(x) horzcat('SE_', x), options_.varobs, 'UniformOutput', false);
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, my_title, headers, labels, M_.H, lh, 10, 6);
dyntable(options_, my_title, headers, labels, M_.Sigma_e, lh, 10, 6);
if options_.TeX
labels = M_.exo_names_tex;
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, my_title, 'covar_ME', headers, labels, M_.H, lh, 10, 6);
dyn_latex_table(M_, options_, my_title, 'covar_ex_shocks', headers, labels, M_.Sigma_e, lh, 10, 6);
end
if ~all(diag(M_.H)==0)
my_title='MATRIX OF COVARIANCE OF MEASUREMENT ERRORS';
labels = cellfun(@(x) horzcat('SE_', x), options_.varobs, 'UniformOutput', false);
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyntable(options_, my_title, headers, labels, M_.H, lh, 10, 6);
if options_.TeX
labels = M_.exo_names_tex;
headers = vertcat('Variables', labels);
lh = cellofchararraymaxlength(labels)+2;
dyn_latex_table(M_, options_, my_title, 'covar_ME', headers, labels, M_.H, lh, 10, 6);
end
end
end
if options_.partial_information

View File

@ -87,7 +87,7 @@ end
if options_.k_order_solver
orig_order = options_.order;
options_.order = local_order;
dr = set_state_space(dr,M_,options_);
dr = set_state_space(dr,M_);
[dr,info] = k_order_pert(dr,M_,options_);
options_.order = orig_order;
return

View File

@ -772,8 +772,7 @@ mod_and_m_tests = [
'pruning/AS_pruned_state_space_red_shock.mod' ],
'extra' : [ 'pruning/Andreasen_et_al_2018_Dynare44Pruning_v2.mat'] },
{ 'test' : [ 'pruning/fs2000_pruning.mod' ] },
{ 'test' : [ 'measurement_errors/fs2000_corr_me_ml_mcmc/fs2000_corr_ME.mod' ],
'extra' : [ 'measurement_errors/fs2000_corr_me_ml_mcmc/fsdat_simul.m'] },
{ 'test' : [ 'measurement_errors/fs2000_corr_me_ml_mcmc/fs2000_corr_ME.mod' ] },
{ 'test' : [ 'TeX/fs2000_corr_ME.mod' ],
'extra' : [ 'fs2000/fsdat_simul.m'] },
{ 'test' : [ 'estimation/MH_recover/fs2000_recover.mod',

@ -1 +1 @@
Subproject commit 084372a314a5f3081dc055ba83dd879947809576
Subproject commit 3c20b9862b87e894b00c70b10e7f5a425ba0372a

View File

@ -94,26 +94,32 @@ varobs gp_obs gy_obs;
shocks;
var e_a; stderr 0.014;
var e_m; stderr 0.005;
corr gy_obs,gp_obs = 0.5;
corr gy_obs,gp_obs = 0.25;
end;
steady;
stoch_simul(order=1,irf=0,nomoments,nofunctions,periods=1000);
verbatim;
temp=oo_.endo_simul([strmatch('gy_obs',M_.endo_names,'exact'),strmatch('gp_obs',M_.endo_names,'exact')],:);
exo_simul = randn(2,options_.periods)'*chol([0.005 0.25*sqrt(0.005*0.01); 0.25*sqrt(0.005*0.01) 0.01]);
temp=temp'+exo_simul;
gy_obs=temp(:,1);
gp_obs=temp(:,2);
save('fsdat_simul_corr.mat','gy_obs','gp_obs')
end;
estimated_params;
alp, 0.356;
gam, 0.0085;
del, 0.01;
stderr e_a, 0.035449;
stderr e_m, 0.008862;
corr e_m, e_a, 0;
stderr gp_obs, 1;
stderr gy_obs, 1;
corr gp_obs, gy_obs,0;
alp, 0.33, 0.2, 0.4;
gam, 0.0085, 0 , 0.1;
del, 0.01, 0, 0.1;
stderr e_a, 0.014, 0, Inf;
stderr e_m, 0.005, 0, Inf;
corr e_m, e_a, 0, 0, 1;
stderr gp_obs, 1, 0, Inf;
stderr gy_obs, 1, 0, Inf;
corr gp_obs, gy_obs,0, 0 , 1;
end;
warning('off','MATLAB:nearlySingularMatrix')
estimation(order=1,datafile=fsdat_simul,silent_optimizer,prior_trunc=0,mode_check,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20,tex) m P c e W R k d y gy_obs;
estimation(order=1,datafile=fsdat_simul_corr,silent_optimizer,prior_trunc=0,mode_check,smoother,filter_decomposition,forecast = 8,filtered_vars,filter_step_ahead=[1,3],irf=20,tex) m P c e W R k d y gy_obs;
@ -121,14 +127,14 @@ estimated_params(overwrite);
//alp, beta_pdf, 0.356, 0.02;
gam, normal_pdf, 0.0085, 0.003;
//del, beta_pdf, 0.01, 0.005;
stderr e_a, inv_gamma_pdf, 0.035449, inf;
stderr e_a, inv_gamma_pdf, 0.05449, inf;
stderr e_m, inv_gamma_pdf, 0.008862, inf;
corr e_m, e_a, normal_pdf, 0, 0.2;
stderr gp_obs, inv_gamma_pdf, 0.001, inf;
corr e_m, e_a, , -1, 1, normal_pdf, 0, 0.2;
stderr gp_obs, 0.1, inv_gamma_pdf, 0.001, inf;
//stderr gy_obs, inv_gamma_pdf, 0.001, inf;
//corr gp_obs, gy_obs,normal_pdf, 0, 0.2;
end;
estimation(mode_compute=5,silent_optimizer,order=1,datafile=fsdat_simul,mode_check,smoother,filter_decomposition,mh_replic=2000, mh_nblocks=1, mh_jscale=0.8,forecast = 8,bayesian_irf,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y;
estimation(mode_compute=5,silent_optimizer,order=1,datafile=fsdat_simul_corr,mode_check,smoother,filter_decomposition,mh_replic=2000, mh_nblocks=1, mh_jscale=0.8,forecast = 8,bayesian_irf,filtered_vars,filter_step_ahead=[1,3],irf=20) m P c e W R k d y;
shock_decomposition y W R;
//identification(advanced=1,max_dim_cova_group=3,prior_mc=250);

View File

@ -1,416 +0,0 @@
% Generated data, used by fs2000.mod
gy_obs =[
1.0030045
1.0002599
0.99104664
1.0321162
1.0223545
1.0043614
0.98626929
1.0092127
1.0357197
1.0150827
1.0051548
0.98465775
0.99132132
0.99904153
1.0044641
1.0179198
1.0113462
0.99409421
0.99904293
1.0448336
0.99932433
1.0057004
0.99619787
1.0267504
1.0077645
1.0058026
1.0025891
0.9939097
0.99604693
0.99908569
1.0151094
0.99348134
1.0039124
1.0145805
0.99800868
0.98578138
1.0065771
0.99843919
0.97979062
0.98413351
0.96468174
1.0273857
1.0225211
0.99958667
1.0111157
1.0099585
0.99480311
1.0079265
0.98924573
1.0070613
1.0075706
0.9937151
1.0224711
1.0018891
0.99051863
1.0042944
1.0184055
0.99419508
0.99756624
1.0015983
0.9845772
1.0004407
1.0116237
0.9861885
1.0073094
0.99273355
1.0013224
0.99777979
1.0301686
0.96809556
0.99917088
0.99949253
0.96590004
1.0083938
0.96662298
1.0221454
1.0069792
1.0343996
1.0066531
1.0072525
0.99743563
0.99723703
1.000372
0.99013917
1.0095223
0.98864268
0.98092242
0.98886488
1.0030341
1.01894
0.99155059
0.99533235
0.99734316
1.0047356
1.0082737
0.98425116
0.99949212
1.0055899
1.0065075
0.99385069
0.98867975
0.99804843
1.0184038
0.99301902
1.0177222
1.0051924
1.0187852
1.0098985
1.0097172
1.0145811
0.98721038
1.0361722
1.0105821
0.99469309
0.98626785
1.013871
0.99858924
0.99302637
1.0042186
0.99623745
0.98545708
1.0225435
1.0011861
1.0130321
0.97861347
1.0228193
0.99627435
1.0272779
1.0075172
1.0096762
1.0129306
0.99966549
1.0262882
1.0026914
1.0061475
1.009523
1.0036127
0.99762992
0.99092634
1.0058469
0.99887292
1.0060653
0.98673557
0.98895709
0.99111967
0.990118
0.99788054
0.97054709
1.0099157
1.0107431
0.99518695
1.0114048
0.99376019
1.0023369
0.98783327
1.0051727
1.0100462
0.98607387
1.0000064
0.99692442
1.012225
0.99574078
0.98642833
0.99008207
1.0197359
1.0112849
0.98711069
0.99402748
1.0242141
1.0135349
0.99842505
1.0130714
0.99887044
1.0059058
1.0185998
1.0073314
0.98687706
1.0084551
0.97698964
0.99482714
1.0015302
1.0105331
1.0261767
1.0232822
1.0084176
0.99785167
0.99619733
1.0055223
1.0076326
0.99205461
1.0030587
1.0137012
1.0145878
1.0190297
1.0000681
1.0153894
1.0140649
1.0007236
0.97961463
1.0125257
1.0169503
1.0197363
1.0221185
];
gp_obs =[
1.0079715
1.0115853
1.0167502
1.0068957
1.0138189
1.0258364
1.0243817
1.017373
1.0020171
1.0003742
1.0008974
1.0104804
1.0116393
1.0114294
0.99932124
0.99461459
1.0170349
1.0051446
1.020639
1.0051964
1.0093042
1.007068
1.01086
0.99590086
1.0014883
1.0117332
0.9990095
1.0108284
1.0103672
1.0036722
1.0005124
1.0190331
1.0130978
1.007842
1.0285436
1.0322054
1.0213403
1.0246486
1.0419306
1.0258867
1.0156316
0.99818589
0.9894107
1.0127584
1.0146882
1.0136529
1.0340107
1.0343652
1.02971
1.0077932
1.0198114
1.013971
1.0061083
1.0089573
1.0037926
1.0082071
0.99498155
0.99735772
0.98765026
1.006465
1.0196088
1.0053233
1.0119974
1.0188066
1.0029302
1.0183459
1.0034218
1.0158799
0.98824798
1.0274357
1.0168832
1.0180641
1.0294657
0.98864091
1.0358326
0.99889969
1.0178322
0.99813566
1.0073549
1.0215985
1.0084245
1.0080939
1.0157021
1.0075815
1.0032633
1.0117871
1.0209276
1.0077569
0.99680958
1.0120266
1.0017625
1.0138811
1.0198358
1.0059629
1.0115416
1.0319473
1.0167074
1.0116111
1.0048627
1.0217622
1.0125221
1.0142045
0.99792469
0.99823971
0.99561547
0.99850373
0.9898464
1.0030963
1.0051373
1.0004213
1.0144117
0.97185592
0.9959518
1.0073529
1.0051603
0.98642572
0.99433423
1.0112131
1.0007695
1.0176867
1.0134363
0.99926191
0.99879835
0.99878754
1.0331374
1.0077797
1.0127221
1.0047393
1.0074106
0.99784213
1.0056495
1.0057708
0.98817494
0.98742176
0.99930555
1.0000687
1.0129754
1.009529
1.0226731
1.0149534
1.0164295
1.0239469
1.0293458
1.026199
1.0197525
1.0126818
1.0054473
1.0254423
1.0069461
1.0153135
1.0337515
1.0178187
1.0240469
1.0079489
1.0186953
1.0008628
1.0113799
1.0140118
1.0168007
1.011441
0.98422774
0.98909729
1.0157859
1.0151586
0.99756232
0.99497777
1.0102841
1.0221659
0.9937759
0.99877193
1.0079433
0.99667692
1.0095959
1.0128804
1.0156949
1.0111951
1.0228887
1.0122083
1.0190197
1.0074927
1.0268096
0.99689352
0.98948474
1.0024938
1.0105543
1.014116
1.0141217
1.0056504
1.0101026
1.0105069
0.99619053
1.0059439
0.99449473
0.99482458
1.0037702
1.0068087
0.99575975
1.0030815
1.0334014
0.99879386
0.99625634
1.0171195
0.99233844
];