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

Merge branch 'MoM_testsuite' into 'master' 

Method of Moments: Updates to testsuite

See merge request Dynare/dynare!1799
time-shift
Sébastien Villemot 2021-01-07 19:52:13 +00:00
parent be16f29254 b3e3501a6d
commit 4434edae0b
25 changed files with 2199 additions and 498 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

109
matlab/method_of_moments/method_of_moments.m
View File

@ -63,7 +63,7 @@ function [oo_, options_mom_, M_] = method_of_moments(bayestopt_, options_, oo_,
% o set_all_parameters.m % o set_all_parameters.m
% o test_for_deep_parameters_calibration.m % o test_for_deep_parameters_calibration.m
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -92,6 +92,7 @@ function [oo_, options_mom_, M_] = method_of_moments(bayestopt_, options_, oo_,
% - [ ] SMM with extended path % - [ ] SMM with extended path
% - [ ] deal with measurement errors (once @wmutschl has implemented this in identification toolbox) % - [ ] deal with measurement errors (once @wmutschl has implemented this in identification toolbox)
% - [ ] improve check for duplicate moments by using the cellfun and unique functions % - [ ] improve check for duplicate moments by using the cellfun and unique functions
% - [ ] dirname option to save output to different directory not yet implemented
% ------------------------------------------------------------------------- % -------------------------------------------------------------------------
% Step 0: Check if required structures and options exist % Step 0: Check if required structures and options exist
% ------------------------------------------------------------------------- % -------------------------------------------------------------------------
@ -133,9 +134,9 @@ if strcmp(options_mom_.mom.mom_method,'GMM') || strcmp(options_mom_.mom.mom_meth
options_mom_.mom = set_default_option(options_mom_.mom,'bartlett_kernel_lag',20); % bandwith in optimal weighting matrix options_mom_.mom = set_default_option(options_mom_.mom,'bartlett_kernel_lag',20); % bandwith in optimal weighting matrix
options_mom_.mom = set_default_option(options_mom_.mom,'penalized_estimator',false); % include deviation from prior mean as additional moment restriction and use prior precision as weight options_mom_.mom = set_default_option(options_mom_.mom,'penalized_estimator',false); % include deviation from prior mean as additional moment restriction and use prior precision as weight
options_mom_.mom = set_default_option(options_mom_.mom,'verbose',false); % display and store intermediate estimation results options_mom_.mom = set_default_option(options_mom_.mom,'verbose',false); % display and store intermediate estimation results
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix',{'DIAGONAL'; 'DIAGONAL'}); % weighting matrix in moments distance objective function at each iteration of estimation; cell of strings with options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix',{'DIAGONAL'; 'DIAGONAL'}); % weighting matrix in moments distance objective function at each iteration of estimation;
% possible values are 'OPTIMAL', 'IDENTITY_MATRIX' ,'DIAGONAL' or a filename. Size of cell determines stages in iterated estimation. % possible values are 'OPTIMAL', 'IDENTITY_MATRIX' ,'DIAGONAL' or a filename. Size of cell determines stages in iterated estimation.
options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix_scaling_factor',1); % scaling of weighting matrix options_mom_.mom = set_default_option(options_mom_.mom,'weighting_matrix_scaling_factor',1); % scaling of weighting matrix in objective function
options_mom_.mom = set_default_option(options_mom_.mom,'se_tolx',1e-5); % step size for numerical computation of standard errors options_mom_.mom = set_default_option(options_mom_.mom,'se_tolx',1e-5); % step size for numerical computation of standard errors
options_mom_ = set_default_option(options_mom_,'order',1); % order of Taylor approximation in perturbation options_mom_ = set_default_option(options_mom_,'order',1); % order of Taylor approximation in perturbation
options_mom_ = set_default_option(options_mom_,'pruning',false); % use pruned state space system at higher-order options_mom_ = set_default_option(options_mom_,'pruning',false); % use pruned state space system at higher-order
@ -169,6 +170,7 @@ options_mom_.mom.compute_derivs = false;% flag to compute derivs in objective fu
% General options that can be set by the user in the mod file, otherwise default values are provided % General options that can be set by the user in the mod file, otherwise default values are provided
options_mom_ = set_default_option(options_mom_,'dirname',M_.dname); % specify directory in which to store estimation output [not yet working]
options_mom_ = set_default_option(options_mom_,'graph_format','eps'); % specify the file format(s) for graphs saved to disk options_mom_ = set_default_option(options_mom_,'graph_format','eps'); % specify the file format(s) for graphs saved to disk
options_mom_ = set_default_option(options_mom_,'nodisplay',false); % do not display the graphs, but still save them to disk options_mom_ = set_default_option(options_mom_,'nodisplay',false); % do not display the graphs, but still save them to disk
options_mom_ = set_default_option(options_mom_,'nograph',false); % do not create graphs (which implies that they are not saved to the disk nor displayed) options_mom_ = set_default_option(options_mom_,'nograph',false); % do not create graphs (which implies that they are not saved to the disk nor displayed)
@ -200,7 +202,7 @@ options_mom_ = set_default_option(options_mom_,'optim_opt',[]);
options_mom_ = set_default_option(options_mom_,'silent_optimizer',false); % run minimization of moments distance silently without displaying results or saving files in between options_mom_ = set_default_option(options_mom_,'silent_optimizer',false); % run minimization of moments distance silently without displaying results or saving files in between
% Check plot options that can be set by the user in the mod file, otherwise default values are provided % Check plot options that can be set by the user in the mod file, otherwise default values are provided
options_mom_.mode_check.nolik = false; % we don't do likelihood (also this initializes mode_check substructure) options_mom_.mode_check.nolik = false; % we don't do likelihood (also this initializes mode_check substructure)
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'status',false); % plot the target function for values around the computed minimum for each estimated parameter in turn. This is helpful to diagnose problems with the optimizer. options_mom_.mode_check = set_default_option(options_mom_.mode_check,'status',false); % plot the target function for values around the computed minimum for each estimated parameter in turn. This is helpful to diagnose problems with the optimizer.
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'neighbourhood_size',.5); % width of the window around the computed minimum to be displayed on the diagnostic plots. This width is expressed in percentage deviation. The Inf value is allowed, and will trigger a plot over the entire domain options_mom_.mode_check = set_default_option(options_mom_.mode_check,'neighbourhood_size',.5); % width of the window around the computed minimum to be displayed on the diagnostic plots. This width is expressed in percentage deviation. The Inf value is allowed, and will trigger a plot over the entire domain
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'symmetric_plots',true); % ensure that the check plots are symmetric around the minimum. A value of 0 allows to have asymmetric plots, which can be useful if the minimum is close to a domain boundary, or in conjunction with neighbourhood_size = Inf when the domain is not the entire real line options_mom_.mode_check = set_default_option(options_mom_.mode_check,'symmetric_plots',true); % ensure that the check plots are symmetric around the minimum. A value of 0 allows to have asymmetric plots, which can be useful if the minimum is close to a domain boundary, or in conjunction with neighbourhood_size = Inf when the domain is not the entire real line
options_mom_.mode_check = set_default_option(options_mom_.mode_check,'number_of_points',20); % number of points around the minimum where the target function is evaluated (for each parameter) options_mom_.mode_check = set_default_option(options_mom_.mode_check,'number_of_points',20); % number of points around the minimum where the target function is evaluated (for each parameter)
@ -738,11 +740,6 @@ catch last_error% if check fails, provide info on using calibration if present
rethrow(last_error); rethrow(last_error);
end end
if options_mom_.mode_compute == 0 %We only report value of moments distance at initial value of the parameters
fprintf('No minimization of moments distance due to ''mode_compute=0''\n')
return
end
% ------------------------------------------------------------------------- % -------------------------------------------------------------------------
% Step 7a: Method of moments estimation: print some info % Step 7a: Method of moments estimation: print some info
% ------------------------------------------------------------------------- % -------------------------------------------------------------------------
@ -760,27 +757,56 @@ end
if options_mom_.mom.penalized_estimator if options_mom_.mom.penalized_estimator
fprintf('\n - penalized estimation using deviation from prior mean and weighted with prior precision'); fprintf('\n - penalized estimation using deviation from prior mean and weighted with prior precision');
end end
if options_mom_.mode_compute == 1; fprintf('\n - optimizer (mode_compute=1): fmincon'); optimizer_vec=[options_mom_.mode_compute;num2cell(options_mom_.additional_optimizer_steps)]; % at each stage one can possibly use different optimizers sequentially
elseif options_mom_.mode_compute == 2; fprintf('\n - optimizer (mode_compute=2): continuous simulated annealing'); for i = 1:length(optimizer_vec)
elseif options_mom_.mode_compute == 3; fprintf('\n - optimizer (mode_compute=3): fminunc'); if i == 1
elseif options_mom_.mode_compute == 4; fprintf('\n - optimizer (mode_compute=4): csminwel'); str = '- optimizer (mode_compute';
elseif options_mom_.mode_compute == 5; fprintf('\n - optimizer (mode_compute=5): newrat'); else
elseif options_mom_.mode_compute == 6; fprintf('\n - optimizer (mode_compute=6): gmhmaxlik'); str = ' (additional_optimizer_steps';
elseif options_mom_.mode_compute == 7; fprintf('\n - optimizer (mode_compute=7): fminsearch'); end
elseif options_mom_.mode_compute == 8; fprintf('\n - optimizer (mode_compute=8): Dynare Nelder-Mead simplex'); switch optimizer_vec{i}
elseif options_mom_.mode_compute == 9; fprintf('\n - optimizer (mode_compute=9): CMA-ES'); case 0
elseif options_mom_.mode_compute == 10; fprintf('\n - optimizer (mode_compute=10): simpsa'); fprintf('\n %s=0): no minimization',str);
elseif options_mom_.mode_compute == 11; fprintf('\n - optimizer (mode_compute=11): online_auxiliary_filter'); case 1
elseif options_mom_.mode_compute == 12; fprintf('\n - optimizer (mode_compute=12): particleswarm'); fprintf('\n %s=1): fmincon',str);
elseif options_mom_.mode_compute == 101; fprintf('\n - optimizer (mode_compute=101): SolveOpt'); case 2
elseif options_mom_.mode_compute == 102; fprintf('\n - optimizer (mode_compute=102): simulannealbnd'); fprintf('\n %s=2): continuous simulated annealing',str);
elseif options_mom_.mode_compute == 13; fprintf('\n - optimizer (mode_compute=13): lsqnonlin'); case 3
elseif ischar(minimizer_algorithm); fprintf(['\n - user-defined optimizer: ' minimizer_algorithm]); fprintf('\n %s=3): fminunc',str);
else case 4
error('method_of_moments: Unknown optimizer, please contact the developers ') fprintf('\n %s=4): csminwel',str);
end case 5
if options_mom_.silent_optimizer fprintf('\n %s=5): newrat',str);
fprintf(' (silent)'); case 6
fprintf('\n %s=6): gmhmaxlik',str);
case 7
fprintf('\n %s=7): fminsearch',str);
case 8
fprintf('\n %s=8): Dynare Nelder-Mead simplex',str);
case 9
fprintf('\n %s=9): CMA-ES',str);
case 10
fprintf('\n %s=10): simpsa',str);
case 11
fprintf('\n %s=11): online_auxiliary_filter',str);
case 12
fprintf('\n %s=12): particleswarm',str);
case 101
fprintf('\n %s=101): SolveOpt',str);
case 102
fprintf('\n %s=102): simulannealbnd',str);
case 13
fprintf('\n %s=13): lsqnonlin',str);
otherwise
if ischar(optimizer_vec{i})
fprintf('\n %s=%s): user-defined',str,optimizer_vec{i});
else
error('method_of_moments: Unknown optimizer, please contact the developers ')
end
end
if options_mom_.silent_optimizer
fprintf(' (silent)');
end
end end
fprintf('\n - perturbation order: %d', options_mom_.order) fprintf('\n - perturbation order: %d', options_mom_.order)
if options_mom_.order > 1 && options_mom_.pruning if options_mom_.order > 1 && options_mom_.pruning
@ -802,8 +828,6 @@ if size(options_mom_.mom.weighting_matrix,1)>1 && ~(any(strcmpi('diagonal',optio
fprintf('\nYou did not specify the use of an optimal or diagonal weighting matrix. There is no point in running an iterated method of moments.\n') fprintf('\nYou did not specify the use of an optimal or diagonal weighting matrix. There is no point in running an iterated method of moments.\n')
end end
optimizer_vec=[options_mom_.mode_compute,options_mom_.additional_optimizer_steps]; % at each stage one can possibly use different optimizers sequentially
for stage_iter=1:size(options_mom_.mom.weighting_matrix,1) for stage_iter=1:size(options_mom_.mom.weighting_matrix,1)
fprintf('Estimation stage %u\n',stage_iter); fprintf('Estimation stage %u\n',stage_iter);
Woptflag = false; Woptflag = false;
@ -849,15 +873,20 @@ for stage_iter=1:size(options_mom_.mom.weighting_matrix,1)
end end
for optim_iter= 1:length(optimizer_vec) for optim_iter= 1:length(optimizer_vec)
if optimizer_vec(optim_iter)==13 if optimizer_vec{optim_iter}==0
options_mom_.vector_output = true; xparam1=xparam0; %no minimization, evaluate objective at current values
fval = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_);
else else
options_mom_.vector_output = false; if optimizer_vec{optim_iter}==13
end options_mom_.vector_output = true;
[xparam1, fval, exitflag] = dynare_minimize_objective(objective_function, xparam0, optimizer_vec(optim_iter), options_mom_, [Bounds.lb Bounds.ub], bayestopt_laplace.name, bayestopt_laplace, [],... else
Bounds, oo_, estim_params_, M_, options_mom_); options_mom_.vector_output = false;
if options_mom_.vector_output end
fval = fval'*fval; [xparam1, fval, exitflag] = dynare_minimize_objective(objective_function, xparam0, optimizer_vec{optim_iter}, options_mom_, [Bounds.lb Bounds.ub], bayestopt_laplace.name, bayestopt_laplace, [],...
Bounds, oo_, estim_params_, M_, options_mom_);
if options_mom_.vector_output
fval = fval'*fval;
end
end end
fprintf('\nStage %d Iteration %d: value of minimized moment distance objective function: %12.10f.\n',stage_iter,optim_iter,fval) fprintf('\nStage %d Iteration %d: value of minimized moment distance objective function: %12.10f.\n',stage_iter,optim_iter,fval)
if options_mom_.mom.verbose if options_mom_.mom.verbose

2
matlab/method_of_moments/method_of_moments_check_plot.m
View File

@ -2,7 +2,7 @@ function method_of_moments_check_plot(fun,xparam,SE_vec,options_,M_,estim_params
% Checks the estimated local minimum of the moment's distance objective % Checks the estimated local minimum of the moment's distance objective
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %

2
matlab/method_of_moments/method_of_moments_data_moments.m
View File

@ -16,7 +16,7 @@ function [dataMoments, m_data] = method_of_moments_data_moments(data, oo_, match
% o method_of_moments.m % o method_of_moments.m
% o method_of_moments_objective_function.m % o method_of_moments_objective_function.m
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %

4
matlab/method_of_moments/method_of_moments_objective_function.m
View File

@ -31,7 +31,7 @@ function [fval, info, exit_flag, junk1, junk2, oo_, M_, options_mom_] = method_o
% o resol % o resol
% o set_all_parameters % o set_all_parameters
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -232,7 +232,7 @@ elseif strcmp(options_mom_.mom.mom_method,'SMM')
i_ME = setdiff([1:size(M_.H,1)],find(diag(M_.H) == 0)); % find ME with 0 variance i_ME = setdiff([1:size(M_.H,1)],find(diag(M_.H) == 0)); % find ME with 0 variance
chol_S = chol(M_.H(i_ME,i_ME)); %decompose rest chol_S = chol(M_.H(i_ME,i_ME)); %decompose rest
shock_mat=zeros(size(options_mom_.mom.ME_shock_series)); %initialize shock_mat=zeros(size(options_mom_.mom.ME_shock_series)); %initialize
shock_mat(:,i_ME)=options_mom_.mom.ME_shock_series(:,i_exo_var)*chol_S; shock_mat(:,i_ME)=options_mom_.mom.ME_shock_series(:,i_ME)*chol_S;
y_sim = y_sim+shock_mat; y_sim = y_sim+shock_mat;
end end

2
matlab/method_of_moments/method_of_moments_optimal_weighting_matrix.m
View File

@ -19,7 +19,7 @@ function W_opt = method_of_moments_optimal_weighting_matrix(m_data, moments, q_l
% This function calls: % This function calls:
% o CorrMatrix (embedded) % o CorrMatrix (embedded)
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %

2
matlab/method_of_moments/method_of_moments_standard_errors.m
View File

@ -29,7 +29,7 @@ function [SE_values, Asympt_Var] = method_of_moments_standard_errors(xparam, obj
% o SMM_objective_function.m % o SMM_objective_function.m
% o method_of_moments_optimal_weighting_matrix % o method_of_moments_optimal_weighting_matrix
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %

6
tests/.gitignore vendored
View File

@ -50,8 +50,10 @@ wsOct
!/ep/mean_preserving_spread.m !/ep/mean_preserving_spread.m
!/ep/rbcii_steady_state.m !/ep/rbcii_steady_state.m
!/estimation/fsdat_simul.m !/estimation/fsdat_simul.m
!/estimation/method_of_moments/RBC_MoM_steady_helper.m !/estimation/method_of_moments/RBC/RBC_MoM_steady_helper.m
!/estimation/method_of_moments/RBC_Andreasen_Data_2.mat !/estimation/method_of_moments/RBC/RBC_Andreasen_Data_2.mat
!/estimation/method_of_moments/AFVRR/AFVRR_data.mat
!/estimation/method_of_moments/AFVRR/AFVRR_steady_helper.m
!/expectations/expectation_ss_old_steadystate.m !/expectations/expectation_ss_old_steadystate.m
!/external_function/extFunDeriv.m !/external_function/extFunDeriv.m
!/external_function/extFunNoDerivs.m !/external_function/extFunNoDerivs.m

22
tests/Makefile.am
View File

@ -50,10 +50,14 @@ MODFILES = \
estimation/MH_recover/fs2000_recover_3.mod \ estimation/MH_recover/fs2000_recover_3.mod \
estimation/t_proposal/fs2000_student.mod \ estimation/t_proposal/fs2000_student.mod \
estimation/tune_mh_jscale/fs2000.mod \ estimation/tune_mh_jscale/fs2000.mod \
estimation/method_of_moments/AnScho_MoM.mod \ estimation/method_of_moments/AnScho/AnScho_MoM.mod \
estimation/method_of_moments/RBC_MoM_Andreasen.mod \ estimation/method_of_moments/RBC/RBC_MoM_Andreasen.mod \
estimation/method_of_moments/RBC_MoM_SMM_ME.mod \ estimation/method_of_moments/RBC/RBC_MoM_SMM_ME.mod \
estimation/method_of_moments/RBC_MoM_prefilter.mod \ estimation/method_of_moments/RBC/RBC_MoM_prefilter.mod \
estimation/method_of_moments/RBC/RBC_MoM_optimizer.mod \
estimation/method_of_moments/AFVRR/AFVRR_M0.mod \
estimation/method_of_moments/AFVRR/AFVRR_MFB.mod \
estimation/method_of_moments/AFVRR/AFVRR_MFB_RRA.mod \
moments/example1_var_decomp.mod \ moments/example1_var_decomp.mod \
moments/example1_bp_test.mod \ moments/example1_bp_test.mod \
moments/test_AR1_spectral_density.mod \ moments/test_AR1_spectral_density.mod \
@ -835,6 +839,10 @@ particle: m/particle o/particle
m/particle: $(patsubst %.mod, %.m.trs, $(PARTICLEFILES)) m/particle: $(patsubst %.mod, %.m.trs, $(PARTICLEFILES))
o/particle: $(patsubst %.mod, %.o.trs, $(PARTICLEFILES)) o/particle: $(patsubst %.mod, %.o.trs, $(PARTICLEFILES))
method_of_moments: m/method_of_moments o/method_of_moments
m/method_of_moments: $(patsubst %.mod, %.m.trs, $(filter estimation/method_of_moments/%.mod, $(MODFILES)))
o/method_of_moments: $(patsubst %.mod, %.o.trs, $(filter estimation/method_of_moments/%.mod, $(MODFILES)))
# Matlab TRS Files # Matlab TRS Files
M_TRS_FILES = $(patsubst %.mod, %.m.trs, $(MODFILES)) M_TRS_FILES = $(patsubst %.mod, %.m.trs, $(MODFILES))
M_TRS_FILES += run_block_byte_tests_matlab.m.trs \ M_TRS_FILES += run_block_byte_tests_matlab.m.trs \
@ -984,8 +992,10 @@ EXTRA_DIST = \
lmmcp/sw-common-header.inc \ lmmcp/sw-common-header.inc \
lmmcp/sw-common-footer.inc \ lmmcp/sw-common-footer.inc \
estimation/tune_mh_jscale/fs2000.inc \ estimation/tune_mh_jscale/fs2000.inc \
estimation/method_of_moments/RBC_MoM_common.inc \ estimation/method_of_moments/RBC/RBC_MoM_common.inc \
estimation/method_of_moments/RBC_MoM_steady_helper.m \ estimation/method_of_moments/RBC/RBC_MoM_steady_helper.m \
estimation/method_of_moments/AFVRR/AFVRR_common.inc \
estimation/method_of_moments/AFVRR/AFVRR_steady_helper.m \
histval_initval_file_unit_tests.m \ histval_initval_file_unit_tests.m \
histval_initval_file/my_assert.m \ histval_initval_file/my_assert.m \
histval_initval_file/ramst_data.xls \ histval_initval_file/ramst_data.xls \

299
tests/estimation/method_of_moments/AFVRR/AFVRR_M0.mod Normal file
View File

@ -0,0 +1,299 @@
% DSGE model based on replication files of
% Andreasen, Fernandez-Villaverde, Rubio-Ramirez (2018), The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications, Review of Economic Studies, 85, p. 1-49
% Adapted for Dynare by Willi Mutschler (@wmutschl, willi@mutschler.eu), Jan 2021
% =========================================================================
% Copyright (C) 2021 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/>.
% =========================================================================
% This is the benchmark model with no feedback M_0
% Original code RunGMM_standardModel_RRA.m by Martin M. Andreasen, Jan 2016
@#include "AFVRR_common.inc"
%--------------------------------------------------------------------------
% Parameter calibration taken from RunGMM_standardModel_RRA.m
%--------------------------------------------------------------------------
% fixed parameters
INHABIT = 1;
PHI1 = 4;
PHI4 = 1;
KAPAone = 0;
DELTA = 0.025;
THETA = 0.36;
ETA = 6;
CHI = 0;
CONSxhr40 = 0;
BETTAxhr = 0;
BETTAxhr40= 0;
RHOD = 0;
GAMA = 0.9999;
CONSxhr20 = 0;
% estimated parameters
BETTA = 0.999544966118000;
B = 0.668859504661000;
H = 0.342483445196000;
PHI2 = 0.997924964981000;
RRA = 662.7953149595370;
KAPAtwo = 5.516226495551000;
ALFA = 0.809462321180000;
RHOR = 0.643873352513000;
BETTAPAI = 1.270087844103000;
BETTAY = 0.031812764291000;
MYYPS = 1.001189151180000;
MYZ = 1.005286347928000;
RHOA = 0.743239127127000;
RHOG = 0.793929380230000;
PAI = 1.012163659169000;
GoY = 0.206594858866000;
STDA = 0.016586292524000;
STDG = 0.041220613851000;
STDD = 0.013534473123000;
% endogenous parameters set via steady state, no need to initialize
%PHIzero = ;
%AA = ;
%PHI3 = ;
%negVf = ;
model_diagnostics;
% Model diagnostics show that some parameters are endogenously determined
% via the steady state, so we run steady to calibrate all parameters
steady;
model_diagnostics;
% Now all parameters are determined
resid;
check;
%--------------------------------------------------------------------------
% Shock distribution
%--------------------------------------------------------------------------
shocks;
var eps_a = STDA^2;
var eps_d = STDD^2;
var eps_g = STDG^2;
end;
%--------------------------------------------------------------------------
% Estimated Params block - these parameters will be estimated, we
% initialize at calibrated values
%--------------------------------------------------------------------------
estimated_params;
BETTA;
B;
H;
PHI2;
RRA;
KAPAtwo;
ALFA;
RHOR;
BETTAPAI;
BETTAY;
MYYPS;
MYZ;
RHOA;
RHOG;
PAI;
GoY;
stderr eps_a;
stderr eps_g;
stderr eps_d;
end;
estimated_params_init(use_calibration);
end;
%--------------------------------------------------------------------------
% Compare whether toolbox yields equivalent moments at second order
%--------------------------------------------------------------------------
% Note that we compare results for orderApp=1|2 and not for orderApp=3, because
% there is a small error in the replication files of the original article in the
% computation of the covariance matrix of the extended innovations vector.
% The authors have been contacted, fixed it, and report that the results
% change only slightly at orderApp=3 to what they report in the paper. At
% orderApp=2 all is correct and so the following part tests whether we get
% the same model moments at the calibrated parameters (we do not optimize).
% We compare it to the replication file RunGMM_standardModel_RRA.m with the
% following settings: orderApp=1|2, seOn=0, q_lag=10, weighting=1;
% scaled=0; optimizer=0; estimator=1; momentSet=2;
%
% Output of the replication files for orderApp=1
AndreasenEtAl.Q1 = 23893.072;
AndreasenEtAl.moments1 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023764' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.028517' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.048361' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.073945' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.073945' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.577' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.042861' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.0011816' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0016052' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.00090947' }
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.0016016' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.0017076' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.0013997' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0055317' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'0.00050106' }
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0018178' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0020186' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0064471' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0030519' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0042181' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.0039217' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.0019975' }
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0061403' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.0058343' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'-0.00089501'}
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0056883' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'-0.00041184'}
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.016255' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4919' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018384' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00065543' }
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.0033626' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0029033' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.006112' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.005683' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'3.3307e-16' }
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4912' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.0018378' }
];
% Output of the replication files for orderApp=2
AndreasenEtAl.Q2 = 65.8269;
AndreasenEtAl.moments2 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023764' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.028517' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.034882' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.056542' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.070145' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0.020825' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.5748' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.04335' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.001205' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0016067' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.00059406'}
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.0011949' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.0016104' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.0020245' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0060254' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'8.3563e-05'}
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0013176' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0019042' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0064261' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0020735' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0027621' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.0029257' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.0012165'}
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0040235' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.0044702' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'0.00030542'}
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0052718' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'0.0010045' }
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.018416' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4853' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018806' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00067309'}
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.0033293' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0019223' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.0039949' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.0052659' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'0.0004337' }
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4846' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.00188' }
];
@#for orderApp in 1:2
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 0 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
% Check results
fprintf('****************************************************************\n')
fprintf('Compare Results for perturbation order @{orderApp}\n')
fprintf('****************************************************************\n')
dev_Q = AndreasenEtAl.Q@{orderApp} - oo_.mom.Q;
dev_datamoments = str2double(AndreasenEtAl.moments@{orderApp}(:,5)) - oo_.mom.data_moments;
dev_modelmoments = str2double(AndreasenEtAl.moments@{orderApp}(:,6)) - oo_.mom.model_moments;
table([AndreasenEtAl.Q@{orderApp} ; str2double(AndreasenEtAl.moments@{orderApp}(:,5)) ; str2double(AndreasenEtAl.moments@{orderApp}(:,6))],...
[oo_.mom.Q ; oo_.mom.data_moments ; oo_.mom.model_moments ],...
[dev_Q ; dev_datamoments ; dev_modelmoments ],...
'VariableNames', {'Andreasen et al', 'Dynare', 'dev'},...
'RowNames', ['Q'; strcat('Data_', M_.matched_moments(:,4)); strcat('Model_', M_.matched_moments(:,4))])
if norm(dev_modelmoments)> 1e-4
error('Something wrong in the computation of moments at order @{orderApp}')
end
@#endfor
%--------------------------------------------------------------------------
% Replicate estimation at orderApp=3
%--------------------------------------------------------------------------
@#ifdef DoEstimation
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = 3 % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL', 'OPTIMAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 13 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, additional_optimizer_steps = [13]
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
@#endif

300
tests/estimation/method_of_moments/AFVRR/AFVRR_MFB.mod Normal file
View File

@ -0,0 +1,300 @@
% DSGE model based on replication files of
% Andreasen, Fernandez-Villaverde, Rubio-Ramirez (2018), The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications, Review of Economic Studies, 85, p. 1-49
% Adapted for Dynare by Willi Mutschler (@wmutschl, willi@mutschler.eu), Jan 2021
% =========================================================================
% Copyright (C) 2021 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/>.
% =========================================================================
% This is the model with Feedback M_FB
% Original code RunGMM_Feedback_estim_RRA.m by Martin M. Andreasen, Jan 2016
@#include "AFVRR_common.inc"
%--------------------------------------------------------------------------
% Parameter calibration taken from RunGMM_Feedback_estim_RRA.m
%--------------------------------------------------------------------------
% fixed parameters
INHABIT = 1;
PHI1 = 4;
PHI4 = 1;
KAPAone = 0;
DELTA = 0.025;
THETA = 0.36;
ETA = 6;
CHI = 0;
BETTAxhr = 0;
BETTAxhr40= 0;
RHOD = 0;
GAMA = 0.9999;
CONSxhr20 = 0;
% estimated parameters
BETTA = 0.997007023687000;
B = 0.692501768577000;
H = 0.339214495653000;
PHI2 = 0.688555040951000;
RRA = 24.346514272871001;
KAPAtwo = 10.018421876923000;
ALFA = 0.792507553312000;
RHOR = 0.849194030384000;
BETTAPAI = 2.060579322980000;
BETTAY = 0.220573712342000;
MYYPS = 1.001016690133000;
MYZ = 1.005356313981000;
RHOA = 0.784141391843000;
RHOG = 0.816924540497000;
PAI = 1.011924196487000;
CONSxhr40 = 0.878774662208000;
GoY = 0.207110300602000;
STDA = 0.013024450606000;
STDG = 0.051049871928000;
STDD = 0.008877423780000;
% endogenous parameters set via steady state, no need to initialize
%PHIzero = ;
%AA = ;
%PHI3 = ;
%negVf = ;
model_diagnostics;
% Model diagnostics show that some parameters are endogenously determined
% via the steady state, so we run steady to calibrate all parameters
steady;
model_diagnostics;
% Now all parameters are determined
resid;
check;
%--------------------------------------------------------------------------
% Shock distribution
%--------------------------------------------------------------------------
shocks;
var eps_a = STDA^2;
var eps_d = STDD^2;
var eps_g = STDG^2;
end;
%--------------------------------------------------------------------------
% Estimated Params block - these parameters will be estimated, we
% initialize at calibrated values
%--------------------------------------------------------------------------
estimated_params;
BETTA;
B;
H;
PHI2;
RRA;
KAPAtwo;
ALFA;
RHOR;
BETTAPAI;
BETTAY;
MYYPS;
MYZ;
RHOA;
RHOG;
PAI;
CONSxhr40;
GoY;
stderr eps_a;
stderr eps_g;
stderr eps_d;
end;
estimated_params_init(use_calibration);
end;
%--------------------------------------------------------------------------
% Compare whether toolbox yields equivalent moments at second order
%--------------------------------------------------------------------------
% Note that we compare results for orderApp=1|2 and not for orderApp=3, because
% there is a small error in the replication files of the original article in the
% computation of the covariance matrix of the extended innovations vector.
% The authors have been contacted, fixed it, and report that the results
% change only slightly at orderApp=3 to what they report in the paper. At
% orderApp=2 all is correct and so the following part tests whether we get
% the same model moments at the calibrated parameters (we do not optimize).
% We compare it to the replication file RunGMM_Feedback_estim_RRA.m with the
% following settings: orderApp=1|2, seOn=0, q_lag=10, weighting=1;
% scaled=0; optimizer=0; estimator=1; momentSet=2;
%
% Output of the replication files for orderApp=1
AndreasenEtAl.Q1 = 201778.9697;
AndreasenEtAl.moments1 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023654' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.027719' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.047415' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.083059' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.083059' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.5745' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.043245' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.0012253' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0015117' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.00080078' }
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.00182' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.001913' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.0016326' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0040112' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'0.00060604' }
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0021426' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0022348' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0039852' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0030058' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0044951' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.0042225' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.0021222' }
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0074776' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.0071906' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'-0.0006736' }
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0070599' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'-0.00036735'}
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.014516' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4866' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018713' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00076856' }
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.002163' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0028078' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.0074583' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.0070551' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'7.2164e-16' }
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4856' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.0018708' }
];
% Output of the replication files for orderApp=2
AndreasenEtAl.Q2 = 59.3323;
AndreasenEtAl.moments2 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023654' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.027719' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.034565' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.056419' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.07087' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0.01517' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.5743' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.043352' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.0012464' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0015247' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.0004867' }
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.0011867' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.0016146' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.0021395' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0043272' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'0.00021752'}
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0013919' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0018899' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0037854' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0021043' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0026571' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.0028566' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.0016279'}
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0039136' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.0044118' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'0.00016791'}
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0052851' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'0.00062143'}
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.018126' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4863' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018806' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00078586'}
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.0021519' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0019046' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.0038939' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.0052792' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'0.00023012'}
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4852' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.0018801' }
];
@#for orderApp in 1:2
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 0 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
% Check results
fprintf('****************************************************************\n')
fprintf('Compare Results for perturbation order @{orderApp}\n')
fprintf('****************************************************************\n')
dev_Q = AndreasenEtAl.Q@{orderApp} - oo_.mom.Q;
dev_datamoments = str2double(AndreasenEtAl.moments@{orderApp}(:,5)) - oo_.mom.data_moments;
dev_modelmoments = str2double(AndreasenEtAl.moments@{orderApp}(:,6)) - oo_.mom.model_moments;
table([AndreasenEtAl.Q@{orderApp} ; str2double(AndreasenEtAl.moments@{orderApp}(:,5)) ; str2double(AndreasenEtAl.moments@{orderApp}(:,6))],...
[oo_.mom.Q ; oo_.mom.data_moments ; oo_.mom.model_moments ],...
[dev_Q ; dev_datamoments ; dev_modelmoments ],...
'VariableNames', {'Andreasen et al', 'Dynare', 'dev'},...
'RowNames', ['Q'; strcat('Data_', M_.matched_moments(:,4)); strcat('Model_', M_.matched_moments(:,4))])
if norm(dev_modelmoments)> 1e-4
warning('Something wrong in the computation of moments at order @{orderApp}')
end
@#endfor
%--------------------------------------------------------------------------
% Replicate estimation at orderApp=3
%--------------------------------------------------------------------------
@#ifdef DoEstimation
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = 3 % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL', 'Optimal'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 13 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, additional_optimizer_steps = [13]
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
@#endif

299
tests/estimation/method_of_moments/AFVRR/AFVRR_MFB_RRA.mod Normal file
View File

@ -0,0 +1,299 @@
% DSGE model based on replication files of
% Andreasen, Fernandez-Villaverde, Rubio-Ramirez (2018), The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications, Review of Economic Studies, 85, p. 1-49
% Adapted for Dynare by Willi Mutschler (@wmutschl, willi@mutschler.eu), Jan 2021
% =========================================================================
% Copyright (C) 2021 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/>.
% =========================================================================
% This is the model with feedback and calibrated RRA
% Original code RunGMM_Feedback_estim_RRA_5.m by Martin M. Andreasen, Jan 2016
@#include "AFVRR_common.inc"
%--------------------------------------------------------------------------
% Parameter calibration taken from RunGMM_Feedback_estim_RRA_5.m
%--------------------------------------------------------------------------
% fixed parameters
INHABIT = 1;
PHI1 = 4;
PHI4 = 1;
KAPAone = 0;
DELTA = 0.025;
THETA = 0.36;
ETA = 6;
CHI = 0;
BETTAxhr = 0;
BETTAxhr40= 0;
RHOD = 0;
GAMA = 0.9999;
CONSxhr20 = 0;
RRA = 5;
% estimated parameters
BETTA = 0.996850651147000;
B = 0.684201133923000;
H = 0.338754441432000;
PHI2 = 0.738293581320000;
KAPAtwo = 11.664785970704999;
ALFA = 0.831836572237000;
RHOR = 0.772754520116000;
BETTAPAI = 3.020381242896000;
BETTAY = 0.288367683973000;
MYYPS = 1.000911709188000;
MYZ = 1.005433723022000;
RHOA = 0.749465413198000;
RHOG = 0.847225569814000;
PAI = 1.010428794858000;
CONSxhr40 = 0.992863217133000;
GoY = 0.207099399789000;
STDA = 0.015621059978000;
STDG = 0.047539390956000;
STDD = 0.008623441943000;
% endogenous parameters set via steady state, no need to initialize
%PHIzero = ;
%AA = ;
%PHI3 = ;
%negVf = ;
model_diagnostics;
% Model diagnostics show that some parameters are endogenously determined
% via the steady state, so we run steady to calibrate all parameters
steady;
model_diagnostics;
% Now all parameters are determined
resid;
check;
%--------------------------------------------------------------------------
% Shock distribution
%--------------------------------------------------------------------------
shocks;
var eps_a = STDA^2;
var eps_d = STDD^2;
var eps_g = STDG^2;
end;
%--------------------------------------------------------------------------
% Estimated Params block - these parameters will be estimated, we
% initialize at calibrated values
%--------------------------------------------------------------------------
estimated_params;
BETTA;
B;
H;
PHI2;
KAPAtwo;
ALFA;
RHOR;
BETTAPAI;
BETTAY;
MYYPS;
MYZ;
RHOA;
RHOG;
PAI;
CONSxhr40;
GoY;
stderr eps_a;
stderr eps_g;
stderr eps_d;
end;
estimated_params_init(use_calibration);
end;
%--------------------------------------------------------------------------
% Compare whether toolbox yields equivalent moments at second order
%--------------------------------------------------------------------------
% Note that we compare results for orderApp=1|2 and not for orderApp=3, because
% there is a small error in the replication files of the original article in the
% computation of the covariance matrix of the extended innovations vector.
% The authors have been contacted, fixed it, and report that the results
% change only slightly at orderApp=3 to what they report in the paper. At
% orderApp=2 all is correct and so the following part tests whether we get
% the same model moments at the calibrated parameters (we do not optimize).
% We compare it to the replication file RunGMM_Feedback_estim_RRA.m with the
% following settings: orderApp=1|2, seOn=1, q_lag=10, weighting=1+1;
% scaled=0; optimizer=0; estimator=1; momentSet=2;
%
% Output of the replication files for orderApp=1
AndreasenEtAl.Q1 = 60275.3715;
AndreasenEtAl.moments1 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023726' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.027372' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.041499' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.077843' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.077843' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.5746' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.043299' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.0012763' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0017759' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.00077354' }
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.0016538' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.0017949' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.0017847' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0053424' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'0.00064897' }
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0019533' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0020602' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0064856' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0020922' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0036375' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.0034139' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.0011665' }
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0066074' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.0062959' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'-0.00075499'}
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0061801' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'-0.00030456'}
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.012048' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4872' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018759' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00080528' }
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.0017036' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0020185' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.0065788' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.0061762' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'-4.5519e-15'}
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4863' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.0018755' }
];
% Output of the replication files for orderApp=2
AndreasenEtAl.Q2 = 140.8954;
AndreasenEtAl.moments2 =[ % note that we reshuffeled to be compatible with our matched moments block
{[ 1]} {'Ex' } {'Gr_C '} {' ' } {'0.024388' } {'0.023726' }
{[ 2]} {'Ex' } {'Gr_I '} {' ' } {'0.031046' } {'0.027372' }
{[ 3]} {'Ex' } {'Infl ' } {' ' } {'0.03757' } {'0.034618' }
{[ 4]} {'Ex' } {'r1 ' } {' ' } {'0.056048' } {'0.056437' }
{[ 5]} {'Ex' } {'r40 ' } {' ' } {'0.069929' } {'0.07051' }
{[ 6]} {'Ex' } {'xhr40 '} {' ' } {'0.017237' } {'0.014242' }
{[ 7]} {'Ex' } {'GoY '} {' ' } {'-1.5745' } {'-1.574' }
{[ 8]} {'Ex' } {'hours '} {' ' } {'-0.043353' } {'-0.043351' }
{[ 9]} {'Exx' } {'Gr_C '} {'Gr_C '} {'0.0013159' } {'0.0012917' }
{[17]} {'Exx' } {'Gr_C '} {'Gr_I '} {'0.0021789' } {'0.0017862' }
{[18]} {'Exx' } {'Gr_C '} {'Infl ' } {'0.00067495' } {'0.00061078' }
{[19]} {'Exx' } {'Gr_C '} {'r1 ' } {'0.0011655' } {'0.0011494' }
{[20]} {'Exx' } {'Gr_C '} {'r40 ' } {'0.0015906' } {'0.0016149' }
{[21]} {'Exx' } {'Gr_C '} {'xhr40 '} {'0.0020911' } {'0.002203' }
{[10]} {'Exx' } {'Gr_I '} {'Gr_I '} {'0.0089104' } {'0.0054317' }
{[22]} {'Exx' } {'Gr_I '} {'Infl ' } {'0.00063139' } {'0.00045278' }
{[23]} {'Exx' } {'Gr_I '} {'r1 ' } {'0.0011031' } {'0.0013672' }
{[24]} {'Exx' } {'Gr_I '} {'r40 ' } {'0.0018445' } {'0.0018557' }
{[25]} {'Exx' } {'Gr_I '} {'xhr40 '} {'0.00095556' } {'0.0067742' }
{[11]} {'Exx' } {'Infl ' } {'Infl ' } {'0.0020268' } {'0.0016583' }
{[26]} {'Exx' } {'Infl ' } {'r1 ' } {'0.0025263' } {'0.0024521' }
{[27]} {'Exx' } {'Infl ' } {'r40 ' } {'0.0029126' } {'0.002705' }
{[28]} {'Exx' } {'Infl ' } {'xhr40 '} {'-0.00077101'} {'-0.00065007'}
{[12]} {'Exx' } {'r1 ' } {'r1 ' } {'0.0038708' } {'0.0038274' }
{[29]} {'Exx' } {'r1 ' } {'r40 ' } {'0.0044773' } {'0.004297' }
{[30]} {'Exx' } {'r1 ' } {'xhr40 '} {'-0.00048202'} {'6.3243e-05' }
{[13]} {'Exx' } {'r40 ' } {'r40 ' } {'0.0054664' } {'0.0051686' }
{[31]} {'Exx' } {'r40 ' } {'xhr40 '} {'0.00053864' } {'0.00066645' }
{[14]} {'Exx' } {'xhr40 '} {'xhr40 '} {'0.053097' } {'0.013543' }
{[15]} {'Exx' } {'GoY '} {'GoY '} {'2.4863' } {'2.4858' }
{[16]} {'Exx' } {'hours '} {'hours '} {'0.0018799' } {'0.0018804' }
{[32]} {'Exx1'} {'Gr_C '} {'Gr_C '} {'0.00077917' } {'0.00081772' }
{[33]} {'Exx1'} {'Gr_I '} {'Gr_I '} {'0.0050104' } {'0.0017106' }
{[34]} {'Exx1'} {'Infl ' } {'Infl ' } {'0.0019503' } {'0.0015835' }
{[35]} {'Exx1'} {'r1 ' } {'r1 ' } {'0.0038509' } {'0.0037985' }
{[36]} {'Exx1'} {'r40 ' } {'r40 ' } {'0.0054699' } {'0.0051642' }
{[37]} {'Exx1'} {'xhr40 '} {'xhr40 '} {'-0.00098295'} {'0.00020285' }
{[38]} {'Exx1'} {'GoY '} {'GoY '} {'2.4868' } {'2.4848' }
{[39]} {'Exx1'} {'hours '} {'hours '} {'0.0018799' } {'0.0018799' }
];
@#for orderApp in 1:2
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 0 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
% Check results
fprintf('****************************************************************\n')
fprintf('Compare Results for perturbation order @{orderApp}\n')
fprintf('****************************************************************\n')
dev_Q = AndreasenEtAl.Q@{orderApp} - oo_.mom.Q;
dev_datamoments = str2double(AndreasenEtAl.moments@{orderApp}(:,5)) - oo_.mom.data_moments;
dev_modelmoments = str2double(AndreasenEtAl.moments@{orderApp}(:,6)) - oo_.mom.model_moments;
table([AndreasenEtAl.Q@{orderApp} ; str2double(AndreasenEtAl.moments@{orderApp}(:,5)) ; str2double(AndreasenEtAl.moments@{orderApp}(:,6))],...
[oo_.mom.Q ; oo_.mom.data_moments ; oo_.mom.model_moments ],...
[dev_Q ; dev_datamoments ; dev_modelmoments ],...
'VariableNames', {'Andreasen et al', 'Dynare', 'dev'},...
'RowNames', ['Q'; strcat('Data_', M_.matched_moments(:,4)); strcat('Model_', M_.matched_moments(:,4))])
if norm(dev_modelmoments)> 1e-4
warning('Something wrong in the computation of moments at order @{orderApp}')
end
@#endfor
%--------------------------------------------------------------------------
% Replicate estimation at orderApp=3
%--------------------------------------------------------------------------
@#ifdef DoEstimation
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'AFVRR_data.mat' % name of filename with data
, bartlett_kernel_lag = 10 % bandwith in optimal weighting matrix
, order = 3 % order of Taylor approximation in perturbation
, pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL', 'Optimal'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
% , TeX % print TeX tables and graphics
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 13 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, additional_optimizer_steps = [13]
, optim = ('TolFun', 1e-6
,'TolX', 1e-6
,'MaxIter', 3000
,'MaxFunEvals', 1D6
,'UseParallel' , 1
%,'Jacobian' , 'on'
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
%, analytic_standard_errors
, se_tolx=1e-10
);
@#endif

540
tests/estimation/method_of_moments/AFVRR/AFVRR_common.inc Normal file
View File

@ -0,0 +1,540 @@
% DSGE model based on replication files of
% Andreasen, Fernandez-Villaverde, Rubio-Ramirez (2018), The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications, Review of Economic Studies, 85, p. 1-49
% Original code by Martin M. Andreasen, Jan 2016
% Adapted for Dynare by Willi Mutschler (@wmutschl, willi@mutschler.eu), Jan 2021
% =========================================================================
% Copyright (C) 2021 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/>.
% =========================================================================
%--------------------------------------------------------------------------
% Variable declaration
%--------------------------------------------------------------------------
var
ln_k
ln_s
ln_a
ln_g
ln_d
ln_c
ln_r
ln_pai
ln_h
ln_q
ln_evf
ln_iv
ln_x2
ln_la
ln_goy
ln_Esdf
xhr20
xhr40
Exhr
@#for i in 1:40
ln_p@{i}
@#endfor
Obs_Gr_C
Obs_Gr_I
Obs_Infl
Obs_r1
Obs_r40
Obs_xhr40
Obs_GoY
Obs_hours
;
predetermined_variables ln_k ln_s;
varobs Obs_Gr_C Obs_Gr_I Obs_Infl Obs_r1 Obs_r40 Obs_xhr40 Obs_GoY Obs_hours;
%--------------------------------------------------------------------------
% Exogenous shocks
%--------------------------------------------------------------------------
varexo
eps_a
eps_d
eps_g
;
%--------------------------------------------------------------------------
% Parameter declaration
%--------------------------------------------------------------------------
parameters
BETTA
B
INHABIT
H
PHI1
PHI2
RRA
PHI4
KAPAone
KAPAtwo
DELTA
THETA
ETA
ALFA
CHI
RHOR
BETTAPAI
BETTAY
MYYPS
MYZ
RHOA
%STDA
RHOG
%STDG
RHOD
%STDD
CONSxhr40
BETTAxhr
BETTAxhr40
CONSxhr20
PAI
GAMA
GoY
%auxiliary
PHIzero
AA
PHI3
negVf
;
%--------------------------------------------------------------------------
% Model equations
%--------------------------------------------------------------------------
% Based on DSGE_model_NegVf_yieldCurve.m and DSGE_model_PosVf_yieldCurve.m
% Note that we include an auxiliary parameter negVf to distinguish whether
% the steady state value function is positive (negVf=0) or negative (negVf=1).
% This parameter is endogenously determined in the steady_state_model block.
model;
%--------------------------------------------------------------------------
% Auxiliary expressions
%--------------------------------------------------------------------------
% do exp transform such that variables are logged variables
@#for var in [ "k", "s", "c", "r", "a", "g", "d", "pai", "h", "q", "evf", "iv", "x2", "la", "goy", "Esdf" ]
#@{var}_ba1 = exp(ln_@{var}(-1));
#@{var}_cu = exp(ln_@{var});
#@{var}_cup = exp(ln_@{var}(+1));
@#endfor
@#for i in 1:40
#p@{i}_cu = exp(ln_p@{i});
#p@{i}_cup = exp(ln_p@{i}(+1));
@#endfor
% these variables are not transformed
#xhr20_cu = xhr20;
#xhr20_cup = xhr20(+1);
#xhr40_cu = xhr40;
#xhr40_cup = xhr40(+1);
#Exhr_cu = Exhr;
#Exhr_cup = Exhr(+1);
% auxiliary steady state variables
#K = exp(steady_state(ln_k));
#IV = exp(steady_state(ln_iv));
#C = exp(steady_state(ln_c));
#Y = (C + IV)/(1-GoY);
#R = exp(steady_state(ln_r));
#G = Y-C-IV;
#removeMeanXhr = 1;
% The atemporal relations if possible
% No stochastic trend in investment specific shocks
#myyps_cu = MYYPS;
#myyps_cup = MYYPS;
% No stochastic trend in non-stationary technology shocks
#myz_cu = MYZ;
#myz_cup = MYZ;
% Defining myzstar
#MYZSTAR = MYYPS^(THETA/(1-THETA))*MYZ;
#myzstar_cu = myyps_cu ^(THETA/(1-THETA))*myz_cu;
#myzstar_cup= myyps_cup^(THETA/(1-THETA))*myz_cup;
% The expression for the value function (only valid for deterministic trends!)
% Note that we make use of auxiliary parameter negVf to switch signs
#mvf_cup = -negVf*(d_cup/(1-PHI2)*((c_cup-B*c_cu*MYZSTAR^-1)^(1-PHI2)-1) + d_cup*PHIzero/(1-PHI1)*(1-h_cup)^(1-PHI1) - negVf* BETTA*MYZSTAR^((1-PHI4)*(1-PHI2))*AA*evf_cup^(1/(1-PHI3)));
% The growth rate in lambda
#myla_cup = (la_cup/la_cu)*(AA*evf_cu^(1/(1-PHI3))/mvf_cup)^PHI3*myzstar_cup^(-PHI2*(1-PHI4)-PHI4);
% The relation between the optimal price for the firms and the pris and inflation
%ptil_cu = ((1-ALFA*(pai_ba1^CHI/pai_cu )^(1-ETA))/(1-ALFA))^(1/(1-ETA));
%ptil_cup = ((1-ALFA*(pai_cu ^CHI/pai_cup)^(1-ETA))/(1-ALFA))^(1/(1-ETA));
#ptil_cu = ((1-ALFA*(1/pai_cu )^(1-ETA))/(1-ALFA))^(1/(1-ETA));
#ptil_cup = ((1-ALFA*(1/pai_cup)^(1-ETA))/(1-ALFA))^(1/(1-ETA));
% From the households' FOC for labor
#w_cu = d_cu*PHIzero*(1-h_cu )^(-PHI1)/la_cu;
#w_cup = d_cu*PHIzero*(1-h_cup)^(-PHI1)/la_cup;
% Shouldn't w_cup include d_cup? Let's stick to the original (wrong) code in the replication files as results don't change dramatically... [@wmutschl]
% The firms' FOC for labor
#mc_cu = w_cu /((1-THETA)*a_cu *myyps_cu ^(-THETA/(1-THETA))*myz_cu ^-THETA *k_cu ^THETA*h_cu ^(-THETA));
#mc_cup = w_cup/((1-THETA)*a_cup*myyps_cup^(-THETA/(1-THETA))*myz_cup^-THETA *k_cup^THETA*h_cup^(-THETA));
% The firms' FOC for capital
#rk_cu = mc_cu *THETA* a_cu *myyps_cu *myz_cu ^(1-THETA)*k_cu ^(THETA-1)*h_cu ^(1-THETA);
#rk_cup = mc_cup*THETA* a_cup*myyps_cup*myz_cup^(1-THETA)*k_cup^(THETA-1)*h_cup^(1-THETA);
% The income identity
#y_cu = c_cu + iv_cu + g_cu;
%--------------------------------------------------------------------------
% Actual model equations
%--------------------------------------------------------------------------
[name='Expected value of the value function']
0 = -evf_cu + (mvf_cup/AA)^(1-PHI3);
[name='Households FOC for capital']
0 = -q_cu+BETTA*myla_cup/myyps_cup*(rk_cup+q_cup*(1-DELTA) -q_cup*KAPAtwo/2*(iv_cup/k_cup*myyps_cup*myzstar_cup - IV/K*MYYPS*MYZSTAR)^2 +q_cup*KAPAtwo*(iv_cup/k_cup*myyps_cup*myzstar_cup - IV/K*MYYPS*MYZSTAR)*iv_cup/k_cup*myyps_cup*myzstar_cup);
[name='Households FOC for investments']
0 = -1+q_cu*(1-KAPAone/2*(iv_cu/IV-1)^2-iv_cu/IV*KAPAone*(iv_cu/IV-1)-KAPAtwo*(iv_cu/k_cu*myyps_cu*myzstar_cu - IV/K*MYYPS*MYZSTAR));
[name='Euler equation for consumption']
0 = -1+BETTA*r_cu*exp(CONSxhr40*xhr40_cu + CONSxhr20*xhr20_cu)*myla_cup/pai_cup;
[name='Households FOC for consumption']
0 = -la_cu + d_cu*(c_cu -B*c_ba1*myzstar_cu^-1)^(-PHI2) -INHABIT*B*BETTA*d_cup*(AA*evf_cu^(1/(1-PHI3))/mvf_cup)^PHI3*(c_cup -B*c_cu*myzstar_cup^-1)^(-PHI2)*myzstar_cup^(-PHI2*(1-PHI4)-PHI4);
[name='Nonlinear pricing, relation for x1 = (ETA-1)/ETA*x2']
0= -(ETA-1)/ETA*x2_cu+y_cu*mc_cu*ptil_cu^(-ETA-1) +ALFA*BETTA*myla_cup*(ptil_cu/ptil_cup)^(-ETA-1)*(1/pai_cup)^(-ETA)*(ETA-1)/ETA*x2_cup*myzstar_cup;
[name='Nonlinear pricing, relation for x2']
0=-x2_cu+y_cu*ptil_cu^-ETA +ALFA*BETTA*myla_cup*(ptil_cu/ptil_cup)^(-ETA)*(1/pai_cup)^(1-ETA)*x2_cup*myzstar_cup;
[name='Nonlinear pricing, relation for s']
0= -s_cup+(1-ALFA)*ptil_cu^(-ETA)+ALFA*(pai_cu/1)^ETA*s_cu;
[name='Interest rate rule']
0 = -log(r_cu/R)+RHOR*log(r_ba1/R)+(1-RHOR)*(BETTAPAI*log(pai_cu/PAI)+BETTAY*log(y_cu/Y) + BETTAxhr*(BETTAxhr40*xhr40_cu - removeMeanXhr*Exhr_cu));
[name='Production function']
0 = -y_cu*s_cup + a_cu *(k_cu *myyps_cu ^(-1/(1-THETA))*myz_cu ^-1)^THETA*h_cu ^(1-THETA);
[name='Relation for physical capital stock']
0= -k_cup + (1-DELTA)*k_cu*(myyps_cu*myzstar_cu)^-1 + iv_cu - iv_cu*KAPAone/2*(iv_cu/IV-1)^2 - k_cu*(myyps_cu*myzstar_cu)^-1*KAPAtwo/2*(iv_cu/k_cu*myyps_cu*myzstar_cu - IV/K*MYYPS*MYZSTAR)^2;
[name='Goverment spending over output']
0=-goy_cu + g_cu/y_cu;
[name='The yield curve: p1']
0= -p1_cu + 1/r_cu;
@#for i in 2:40
[name='The yield curve: p@{i}']
0= -p@{i}_cu + BETTA*myla_cup/pai_cup*p@{i-1}_cup;
@#endfor
[name='Stochastic discount factor']
0= -Esdf_cu+ BETTA*myla_cup/pai_cup;
[name='Expected 5 year excess holding period return']
0= -xhr20_cu+ log(p19_cup) - log(p20_cu) - log(r_cu);
[name='Expected 10 year excess holding period return']
0= -xhr40_cu+ log(p39_cup) - log(p40_cu) - log(r_cu);
[name='Mean of expected excess holding period return in Taylor rule']
0= -Exhr_cu + (1-GAMA)*(BETTAxhr40*xhr40_cu) + GAMA*Exhr_cup;
[name='Exogenous process for productivity']
0 = -log(a_cu)+RHOA*log(a_ba1) + eps_a;
[name='Exogenous process for government spending']
0 = -log(g_cu/G)+RHOG*log(g_ba1/G) + eps_g;
[name='Exogenous process for discount factor shifter']
0 = -log(d_cu)+RHOD*log(d_ba1) + eps_d;
[name='Observable annualized consumption growth']
Obs_Gr_C = 4*( ln_c -ln_c(-1) + log(MYZSTAR));
[name='Observable annualized investment growth']
Obs_Gr_I = 4*( ln_iv - ln_iv(-1) + log(MYZSTAR)+log(MYYPS));
[name='Observable annualized inflation']
Obs_Infl = 4*( ln_pai);
[name='Observable annualized one-quarter nominal yield']
Obs_r1 = 4*( ln_r);
[name='Observable annualized 10-year nominal yield']
Obs_r40 = 4*( -1/40*ln_p40);
[name='Observable annualized 10-year ex post excess holding period return']
Obs_xhr40 = 4*( ln_p39 - ln_p40(-1) - ln_r(-1) );
[name='Observable annualized log ratio of government spending to GDP']
Obs_GoY = 4*( 1/4*ln_goy);
[name='Observable annualized log of hours']
Obs_hours = 4*( 1/100*ln_h);
end;
%--------------------------------------------------------------------------
% Steady State Computations
%--------------------------------------------------------------------------
% Based on DSGE_model_yieldCurve_ss.m, getPHI3.m, ObjectGMM.m
% Note that we include an auxiliary parameter negVf to distinguish whether
% the steady state value function is positive (negVf=0) or negative (negVf=1).
% This parameter is endogenously determined in the steady_state_model block.
steady_state_model;
% The growth rate in the firms' fixed costs
MYZSTARBAR = MYYPS^(THETA/(1-THETA))*MYZ;
% The growth rate for lampda
MYLABAR = MYZSTARBAR^(-PHI2*(1-PHI4)-PHI4);
% The relative optimal price for firms
PTILBAR = ((1-ALFA*PAI^((CHI-1)*(1-ETA)))/(1-ALFA))^(1/(1-ETA));
% The state variable s for distortions between output and produktion
SBAR = ((1-ALFA)*PTILBAR^(-ETA))/(1-ALFA*PAI^((1-CHI)*ETA));
% The 1-period interest rate
RBAR = PAI/(BETTA*MYLABAR);
% The market price of capital
QBAR = 1;
% The real price of renting capital
RKBAR = QBAR*(MYYPS/(BETTA*MYLABAR)-(1-DELTA));
% The marginal costs in the firms
MCBAR = (1-ALFA*BETTA*MYLABAR*PAI^((1-CHI)*ETA)*MYZSTARBAR)*(ETA-1)/ETA*PTILBAR/(1-ALFA*BETTA*MYLABAR*PAI^((CHI-1)*(1-ETA))*MYZSTARBAR);
% The capital stock
KBAR = H*(RKBAR/(MCBAR*THETA*MYYPS*MYZ^(1-THETA)))^(1/(THETA-1));
% The wage level
WBAR = MCBAR*(1-THETA)*MYYPS^(-THETA/(1-THETA))*MYZ^-THETA*(KBAR/H)^THETA;
% The level of investment
IVBAR = KBAR - (1-DELTA)*KBAR*MYYPS^(-1/(1-THETA))*MYZ^-1;
% The consumption level
CBAR = ((1-GoY)*(KBAR*MYYPS^(-1/(1-THETA))*MYZ^-1)^THETA*H^(1-THETA))/SBAR-IVBAR;
% The output level
YBAR = (CBAR + IVBAR)/(1-GoY);
% The value of lambda
LABAR = (CBAR-B*CBAR*MYZSTARBAR^-1)^-PHI2 - INHABIT*B*BETTA*(CBAR-B*CBAR*MYZSTARBAR^-1)^-PHI2*MYZSTARBAR^(-PHI2*(1-PHI4)-PHI4);
% The value of PHIzero
PHIzero = LABAR*WBAR*(1-H)^PHI1;
% The level of the value function
VFBAR = 1/(1-BETTA*MYZSTARBAR^((1-PHI4)*(1-PHI2)))*(1/(1-PHI2)*((CBAR-B*CBAR*MYZSTARBAR^-1)^(1-PHI2)-1)+PHIzero/(1-PHI1)*(1-H)^(1-PHI1));
UBAR = 1/(1-PHI2)*((CBAR-B*CBAR*MYZSTARBAR^-1)^(1-PHI2)-1)+PHIzero/(1-PHI1)*(1-H)^(1-PHI1);
[AA, EVFBAR, PHI3, negVf, info]= AFVRR_steady_helper(VFBAR,RBAR,IVBAR,CBAR,KBAR,LABAR,QBAR,YBAR, BETTA,B,PAI,H,PHIzero,PHI1,PHI2,THETA,MYYPS,MYZ,INHABIT,RRA,CONSxhr40);
% The value of X2
X2BAR = YBAR*PTILBAR^(-ETA)/(1-BETTA*ALFA*MYLABAR*PAI^((CHI-1)*(1-ETA))*MYZSTARBAR);
% Government spending
GBAR = GoY*YBAR;
%**************************************************************************
% map into model variables
ln_k = log(KBAR);
ln_s = log(SBAR);
ln_c_ba1 = log(CBAR);
ln_r_ba1 = log(RBAR);
ln_a = log(1);
ln_g = log(GBAR);
ln_d = log(1);
ln_c = log(CBAR);
ln_r = log(RBAR);
ln_pai = log(PAI);
ln_h = log(H);
ln_q = log(QBAR);
ln_evf = log(EVFBAR);
ln_iv = log(IVBAR);
ln_x2 = log(X2BAR);
ln_la = log(LABAR);
ln_goy = log(GoY);
ln_Esdf = log(1/RBAR);
xhr20 = 0;
xhr40 = 0;
Exhr = 0;
% The yield curve
ln_p1 = log((1/RBAR)^1);
ln_p2 = log((1/RBAR)^2);
ln_p3 = log((1/RBAR)^3);
ln_p4 = log((1/RBAR)^4);
ln_p5 = log((1/RBAR)^5);
ln_p6 = log((1/RBAR)^6);
ln_p7 = log((1/RBAR)^7);
ln_p8 = log((1/RBAR)^8);
ln_p9 = log((1/RBAR)^9);
ln_p10 = log((1/RBAR)^10);
ln_p11 = log((1/RBAR)^11);
ln_p12 = log((1/RBAR)^12);
ln_p13 = log((1/RBAR)^13);
ln_p14 = log((1/RBAR)^14);
ln_p15 = log((1/RBAR)^15);
ln_p16 = log((1/RBAR)^16);
ln_p17 = log((1/RBAR)^17);
ln_p18 = log((1/RBAR)^18);
ln_p19 = log((1/RBAR)^19);
ln_p20 = log((1/RBAR)^20);
ln_p21 = log((1/RBAR)^21);
ln_p22 = log((1/RBAR)^22);
ln_p23 = log((1/RBAR)^23);
ln_p24 = log((1/RBAR)^24);
ln_p25 = log((1/RBAR)^25);
ln_p26 = log((1/RBAR)^26);
ln_p27 = log((1/RBAR)^27);
ln_p28 = log((1/RBAR)^28);
ln_p29 = log((1/RBAR)^29);
ln_p30 = log((1/RBAR)^30);
ln_p31 = log((1/RBAR)^31);
ln_p32 = log((1/RBAR)^32);
ln_p33 = log((1/RBAR)^33);
ln_p34 = log((1/RBAR)^34);
ln_p35 = log((1/RBAR)^35);
ln_p36 = log((1/RBAR)^36);
ln_p37 = log((1/RBAR)^37);
ln_p38 = log((1/RBAR)^38);
ln_p39 = log((1/RBAR)^39);
ln_p40 = log((1/RBAR)^40);
Obs_Gr_C = 4*( log(MYZSTARBAR) );
Obs_Gr_I = 4*( log(MYZSTARBAR)+log(MYYPS) );
Obs_Infl = 4*( ln_pai );
Obs_r1 = 4*( ln_r );
Obs_r40 = 4*( -1/40*ln_p40 );
Obs_xhr40 = 4*( xhr40 );
Obs_GoY = 4*( 1/4*ln_goy );
Obs_hours = 4*( 1/100*ln_h );
end;
%--------------------------------------------------------------------------
% Declare moments to use in estimation
%--------------------------------------------------------------------------
% These are the moments used in the paper; corresponds to momentSet=2 in the replication files
matched_moments;
%mean
Obs_Gr_C;
Obs_Gr_I;
Obs_Infl;
Obs_r1;
Obs_r40;
Obs_xhr40;
Obs_GoY;
Obs_hours;
% all variances
Obs_Gr_C*Obs_Gr_C;
Obs_Gr_I*Obs_Gr_I;
Obs_Infl*Obs_Infl;
Obs_r1*Obs_r1;
Obs_r40*Obs_r40;
Obs_xhr40*Obs_xhr40;
Obs_GoY*Obs_GoY;
Obs_hours*Obs_hours;
% covariance excluding GoY and hours
Obs_Gr_C*Obs_Gr_I;
Obs_Gr_C*Obs_Infl;
Obs_Gr_C*Obs_r1;
Obs_Gr_C*Obs_r40;
Obs_Gr_C*Obs_xhr40;
%Obs_Gr_C*Obs_GoY;
%Obs_Gr_C*Obs_hours;
Obs_Gr_I*Obs_Infl;
Obs_Gr_I*Obs_r1;
Obs_Gr_I*Obs_r40;
Obs_Gr_I*Obs_xhr40;
%Obs_Gr_I*Obs_GoY;
%Obs_Gr_I*Obs_hours;
Obs_Infl*Obs_r1;
Obs_Infl*Obs_r40;
Obs_Infl*Obs_xhr40;
%Obs_Infl*Obs_GoY;
%Obs_Infl*Obs_hours;
Obs_r1*Obs_r40;
Obs_r1*Obs_xhr40;
%Obs_r1*Obs_GoY;
%Obs_r1*Obs_hours;
Obs_r40*Obs_xhr40;
%Obs_r40*Obs_GoY;
%Obs_r40*Obs_hours;
%Obs_xhr40*Obs_GoY;
%Obs_xhr40*Obs_hours;
%Obs_GoY*Obs_hours;
%first autocovariance
Obs_Gr_C*Obs_Gr_C(-1);
Obs_Gr_I*Obs_Gr_I(-1);
Obs_Infl*Obs_Infl(-1);
Obs_r1*Obs_r1(-1);
Obs_r40*Obs_r40(-1);
Obs_xhr40*Obs_xhr40(-1);
Obs_GoY*Obs_GoY(-1);
Obs_hours*Obs_hours(-1);
end;
%--------------------------------------------------------------------------
% Create Data
%--------------------------------------------------------------------------
@#ifdef CreateData
verbatim;
% From 1961Q3 to 2007Q4
DataUS = xlsread('Data_PruningPaper_v5.xlsx','Data_used','E3:M188');
% ANNUALIZED (except for hours and GoY)
% 1 2 3 4 5 6 7 8 9
% Lables: Date Gr_C Gr_I GoY hours Infl_C r1 r40 xhr40
%label_data = {'Gr_C ', 'Gr_I ','Infl ', 'r1 ', 'r40 ', 'xhr40 ','GoY ', 'hours '};
%DataUS = [DataUS(:,2:3) DataUS(:,6:8) DataUS(:,9) log(DataUS(:,4)) 4*log(DataUS(:,5))/100];
Obs_Gr_C = DataUS(:,2);
Obs_Gr_I = DataUS(:,3);
Obs_Infl = DataUS(:,6);
Obs_r1 = DataUS(:,7);
Obs_r40 = DataUS(:,8);
Obs_xhr40 = DataUS(:,9);
Obs_GoY = log(DataUS(:,4));
Obs_hours = 4*log(DataUS(:,5))/100;
save('AFVRR_data.mat','Obs_Gr_C','Obs_Gr_I','Obs_Infl','Obs_r1','Obs_r40','Obs_xhr40','Obs_GoY','Obs_hours');
pause(1);
end;
@#endif

BIN
tests/estimation/method_of_moments/AFVRR/AFVRR_data.mat Normal file
View File

Binary file not shown.

80
tests/estimation/method_of_moments/AFVRR/AFVRR_steady_helper.m Normal file
View File

@ -0,0 +1,80 @@
% DSGE model based on replication files of
% Andreasen, Fernandez-Villaverde, Rubio-Ramirez (2018), The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications, Review of Economic Studies, 85, p. 1-49
% Adapted for Dynare by Willi Mutschler (@wmutschl, willi@mutschler.eu), Jan 2021
% =========================================================================
% Copyright (C) 2021 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/>.
% =========================================================================
% This is a helper function to compute steady state values and endogenous parameters
% Based on DSGE_model_yieldCurve_ss.m, getPHI3.m, ObjectGMM.m
function [AA, EVFBAR, PHI3, negVf, info]= AFVRR_steady_helper(VFBAR,RBAR,IVBAR,CBAR,KBAR,LABAR,QBAR,YBAR, BETTA,B,PAI,H,PHIzero,PHI1,PHI2,THETA,MYYPS,MYZ,INHABIT,RRA,CONSxhr40)
% We get nice values of EVF by setting AA app. equal to VF.
% The value of the expected value function raised to the power 1-PHI3
% Also we check bounds on other variables
% % Adding PHI3 to params. Note that PHI3 only affects the value function in
% % steady state, hence the value we assign to PHI3 is irrelevant
% PHI3 = -100;
info=0;
AA = NaN;
EVFBAR = NaN;
PHI3 = NaN;
negVf = NaN;
MYZSTAR = MYYPS^(THETA/(1-THETA))*MYZ;
% The wage level
WBAR = PHIzero*(1-H)^(-PHI1)/LABAR;
RRAc = RRA;
if INHABIT == 1
PHI3 = (RRAc - PHI2/((1-B*MYZSTAR^-1)/(1-BETTA*B)+PHI2/PHI1*WBAR*(1-H)/CBAR))/((1-PHI2)/((1-B*MYZSTAR^-1)/(1-BETTA*B)-(CBAR-B*CBAR*MYZSTAR^-1)^PHI2/((1-BETTA*B)*CBAR)+WBAR*(1-H)/CBAR*(1-PHI2)/(1-PHI1)));
else
PHI3 = (RRAc - PHI2/(1-B*MYZSTAR^-1+PHI2/PHI1*WBAR*(1-H)/CBAR))/((1-PHI2)/(1-B*MYZSTAR^-1-(CBAR-B*CBAR*MYZSTAR^-1)^PHI2/((1-BETTA*B)*CBAR)+WBAR*(1-H)/CBAR*(1-PHI2)/(1-PHI1)));
end
if abs(PHI3) > 30000
disp('abs of PHI3 exceeds 30000')
info=1;
return
end
if CONSxhr40 > 1
info=1;
return
end
if VFBAR < 0
AA = -VFBAR;
EVFBAR = (-VFBAR/AA)^(1-PHI3);
negVf = 1;
else
AA = VFBAR;
EVFBAR = (VFBAR/AA)^(1-PHI3);
negVf = -1;
disp('Positive Value Function');
end
if RBAR < 1 || IVBAR < 0 || CBAR < 0 || KBAR < 0 || PAI < 1 || H < 0 || H > 1 || QBAR < 0 || YBAR < 0
info = 1;
end
end

95
tests/estimation/method_of_moments/AnScho_MoM.mod → tests/estimation/method_of_moments/AnScho/AnScho_MoM.mod
View File

@ -1,8 +1,8 @@
% DSGE model used in replication files of % DSGE model used in replication files of
% An, Sungbae and Schorfheide, Frank, (2007), Bayesian Analysis of DSGE Models, Econometric Reviews, 26, issue 2-4, p. 113-172. % An, Sungbae and Schorfheide, Frank, (2007), Bayesian Analysis of DSGE Models, Econometric Reviews, 26, issue 2-4, p. 113-172.
% Adapted by Willi Mutschler (@wmutschl, willi@mutschler.eu) % Adapted by Willi Mutschler (@wmutschl, willi@mutschler.eu)
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -203,28 +203,33 @@ end
@#for mommethod in ["GMM", "SMM"] @#for mommethod in ["GMM", "SMM"]
method_of_moments( method_of_moments(
% Necessery options % Necessery options
mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM
, datafile = 'AnScho_MoM_data_@{orderApp}.mat' % name of filename with data , datafile = 'AnScho_MoM_data_@{orderApp}.mat' % name of filename with data
% Options for both GMM and SMM % Options for both GMM and SMM
% , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix % , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation , order = @{orderApp} % order of Taylor approximation in perturbation
% , penalized_estimator % use penalized optimization % , penalized_estimator % include deviation from prior mean as additional moment restriction and use prior precision as weight
, pruning % use pruned state space system at higher-order , pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results % , verbose % display and store intermediate estimation results
, weighting_matrix = ['optimal'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename , weighting_matrix = ['optimal'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename. Size of cell determines stages in iterated estimation, e.g. two state with ['DIAGONAL','OPTIMAL']
, additional_optimizer_steps = [4] % vector of numbers for the iterations in the 2-step feasible method of moments %, weighting_matrix_scaling_factor=1 % scaling of weighting matrix in objective function
% , prefilter=0 % demean each data series by its empirical mean and use centered moments , se_tolx=1e-6 % step size for numerical computation of standard errors
%
% Options for SMM % Options for SMM
% , burnin=500 % number of periods dropped at beginning of simulation
% , bounded_shock_support % trim shocks in simulation to +- 2 stdev % , bounded_shock_support % trim shocks in simulation to +- 2 stdev
% , drop = 500 % number of periods dropped at beginning of simulation
% , seed = 24051986 % seed used in simulations % , seed = 24051986 % seed used in simulations
% , simulation_multiple = 5 % multiple of the data length used for simulation % , simulation_multiple = 5 % multiple of the data length used for simulation
%
% General options % Options for GMM
%, dirname = 'MM' % directory in which to store estimation output @#if mommethod == "GMM"
, analytic_standard_errors % compute standard errors using analytical derivatives
@#endif
% General options
% , dirname = 'MM' % directory in which to store estimation output
% , graph_format = EPS % specify the file format(s) for graphs saved to disk % , graph_format = EPS % specify the file format(s) for graphs saved to disk
% , nodisplay % do not display the graphs, but still save them to disk % , nodisplay % do not display the graphs, but still save them to disk
% , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed) % , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
@ -232,44 +237,50 @@ end
% , plot_priors = 1 % control plotting of priors % , plot_priors = 1 % control plotting of priors
% , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters % , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters
% , TeX % print TeX tables and graphics % , TeX % print TeX tables and graphics
%
% Data and model options % Data and model options
%, first_obs = 501 % number of first observation % , first_obs = 501 % number of first observation
% , logdata % if loglinear is set, this option is necessary if the user provides data already in logs, otherwise the log transformation will be applied twice (this may result in complex data) % , logdata % if data is already in logs
% , loglinear % computes a log-linear approximation of the model instead of a linear approximation , nobs = 250 % number of observations
, nobs = 250 % number of observations % , prefilter=0 % demean each data series by its empirical mean and use centered moments
% , xls_sheet = willi % name of sheet with data in Excel
% , xls_sheet = data % name/number of sheet with data in Excel
% , xls_range = B2:D200 % range of data in Excel sheet % , xls_range = B2:D200 % range of data in Excel sheet
%
% Optimization options that can be set by the user in the mod file, otherwise default values are provided % Optimization options that can be set by the user in the mod file, otherwise default values are provided
% , analytic_derivation % uses analytic derivatives to compute standard errors for GMM % , huge_number=1e7 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons , mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance
, mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer , additional_optimizer_steps = [1] % vector of additional mode-finders run after mode_compute
%, optim = ('TolFun', 1e-5 % optim: a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute, some exemplary common options:
% ,'TolX', 1e-6 , optim = ('TolFun' , 1e-6 % termination tolerance on the function value, a positive scalar
% ) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute ,'TolX' , 1e-6 % termination tolerance on x, a positive scalar
, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between ,'MaxIter' , 3000 % maximum number of iterations allowed, a positive integer
% , tolf = 1e-5 % convergence criterion on function value for numerical differentiation ,'MaxFunEvals' , 1D6 % maximum number of function evaluations allowed, a positive integer
% , tolx = 1e-6 % convergence criterion on funciton input for numerical differentiation % ,'UseParallel' , 1 % when true (and supported by optimizer) solver estimates gradients in parallel (using Matlab/Octave's parallel toolbox)
% % ,'Jacobian' , 'off' % when 'off' gradient-based solvers approximate Jacobian using finite differences; for GMM we can also pass the analytical Jacobian to gradient-based solvers by setting this 'on'
% % Numerical algorithms options )
, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
% Numerical algorithms options
% , aim_solver % Use AIM algorithm to compute perturbation approximation % , aim_solver % Use AIM algorithm to compute perturbation approximation
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION % , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION
% , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm % , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm % , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the cycle reduction algorithm
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER % , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER
% , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver % , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver
% , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver % , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver
% , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm % , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm
% , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT % , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT
% , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver % , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver
% , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems [IS THIS CORRET @wmutschl] % , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems
% , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition % , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
@#if mommethod == "GMM" % , schur_vec_tol=1e-11 % tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix
, analytic_standard_errors % , mode_check % plot the target function for values around the computed minimum for each estimated parameter in turn
@#endif % , mode_check_neighbourhood_size = 5 % width of the window (expressed in percentage deviation) around the computed minimum to be displayed on the diagnostic plots
% , mode_check_symmetric_plots=1 % ensure that the check plots are symmetric around the minimum
% , mode_check_number_of_points = 20 % number of points around the minimum where the target function is evaluated (for each parameter)
); );
@#endfor @#endfor

0
tests/estimation/method_of_moments/RBC_Andreasen_Data_2.mat → tests/estimation/method_of_moments/RBC/RBC_Andreasen_Data_2.mat
View File

230
tests/estimation/method_of_moments/RBC/RBC_MoM_Andreasen.mod Normal file
View File

@ -0,0 +1,230 @@
% Tests SMM and GMM routines
%
% Copyright (C) 2020-2021 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/>.
% =========================================================================
% Define testscenario
@#define orderApp = 2
@#define estimParams = 1
% Note that we will set the numerical optimization tolerance levels very large to speed up the testsuite
@#define optimizer = 13
@#include "RBC_MoM_common.inc"
shocks;
var u_a; stderr 0.0072;
end;
varobs c iv n;
@#if estimParams == 0
estimated_params;
DELTA, 0.025;
BETTA, 0.984;
B, 0.5;
ETAc, 2;
ALFA, 0.667;
RHOA, 0.979;
stderr u_a, 0.0072;
end;
@#endif
@#if estimParams == 1
estimated_params;
DELTA, , 0, 1;
BETTA, , 0, 1;
B, , 0, 1;
ETAc, , 0, 10;
ALFA, , 0, 1;
RHOA, , 0, 1;
stderr u_a, , 0, 1;
end;
@#endif
@#if estimParams == 2
estimated_params;
DELTA, 0.025, 0, 1, normal_pdf, 0.02, 0.5;
BETTA, 0.98, 0, 1, beta_pdf, 0.90, 0.25;
B, 0.45, 0, 1, normal_pdf, 0.40, 0.5;
%ETAl, 1, 0, 10, normal_pdf, 0.25, 0.0.1;
ETAc, 1.8, 0, 10, normal_pdf, 1.80, 0.5;
ALFA, 0.65, 0, 1, normal_pdf, 0.60, 0.5;
RHOA, 0.95, 0, 1, normal_pdf, 0.90, 0.5;
stderr u_a, 0.01, 0, 1, normal_pdf, 0.01, 0.5;
%THETA, 3.48, 0, 10, normal_pdf, 0.25, 0.0.1;
end;
@#endif
% Simulate data
%stoch_simul(order=@{orderApp},pruning,nodisplay,nomoments,periods=500);
%save('RBC_MoM_data_@{orderApp}.mat', options_.varobs{:} );
%pause(1);
estimated_params_init(use_calibration);
end;
%--------------------------------------------------------------------------
% Method of Moments Estimation
%--------------------------------------------------------------------------
matched_moments;
c;
n;
iv;
c*c;
c*iv;
iv*n;
iv*iv;
n*c;
n*n;
c*c(-1);
n*n(-1);
iv*iv(-1);
c*c(-3);
n*n(-3);
iv*iv(-3);
c*c(-5);
n*n(-5);
iv*iv(-5);
end;
% get indices in declaration order
ic = strmatch('c', M_.endo_names,'exact');
iiv = strmatch('iv', M_.endo_names,'exact');
in = strmatch('n', M_.endo_names,'exact');
% first entry: number of variable in declaration order
% second entry: lag
% third entry: power
matched_moments_ = {
[ic ] [0 ], [1 ];
[in ] [0 ], [1 ];
[iiv ] [0 ], [1 ];
[ic ic ] [0 0], [1 1];
[ic iiv] [0 0], [1 1];
%[ic in ] [0 0], [1 1];
%[iiv ic ] [0 0], [1 1];
[in iiv] [0 0], [1 1];
[iiv iiv] [0 0], [1 1];
[ic in] [0 0], [1 1];
%[in iiv] [0 0], [1 1];
[in in ] [0 0], [1 1];
[ic ic ] [0 -1], [1 1];
[in in ] [0 -1], [1 1];
[iiv iiv] [0 -1], [1 1];
[ic ic ] [0 -3], [1 1];
[in in ] [0 -3], [1 1];
[iiv iiv] [0 -3], [1 1];
[ic ic ] [0 -5], [1 1];
[in in ] [0 -5], [1 1];
[iiv iiv] [0 -5], [1 1];
};
if ~isequal(M_.matched_moments,matched_moments_)
error('Translation to matched_moments-block failed')
end
method_of_moments(
% Necessery options
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'RBC_Andreasen_Data_2.mat' % name of filename with data
% Options for both GMM and SMM
% , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation
% , penalized_estimator % include deviation from prior mean as additional moment restriction and use prior precision as weight
% , pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL','OPTIMAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename. Size of cell determines stages in iterated estimation, e.g. two state with ['DIAGONAL','OPTIMAL']
% , weighting_matrix_scaling_factor=1 % scaling of weighting matrix in objective function
, se_tolx=1e-6 % step size for numerical computation of standard errors
% Options for SMM
% , burnin=500 % number of periods dropped at beginning of simulation
% , bounded_shock_support % trim shocks in simulation to +- 2 stdev
% , seed = 24051986 % seed used in simulations
% , simulation_multiple = 5 % multiple of the data length used for simulation
% Options for GMM
% , analytic_standard_errors % compute standard errors using analytical derivatives
% General options
% , dirname = 'MM' % directory in which to store estimation output
% , graph_format = EPS % specify the file format(s) for graphs saved to disk
% , nodisplay % do not display the graphs, but still save them to disk
% , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
% , noprint % do not print stuff to console
% , plot_priors = 1 % control plotting of priors
% , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters
, TeX % print TeX tables and graphics
% Data and model options
% , first_obs = 501 % number of first observation
% , logdata % if data is already in logs
% , nobs = 250 % number of observations
% , prefilter=0 % demean each data series by its empirical mean and use centered moments
% , xls_sheet = data % name/number of sheet with data in Excel
% , xls_range = B2:D200 % range of data in Excel sheet
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
% , huge_number=1e7 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 3 % specifies the optimizer for minimization of moments distance
, additional_optimizer_steps = [13] % vector of additional mode-finders run after mode_compute
% optim: a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute, some exemplary common options:
, optim = ('TolFun' , 1D-6 % termination tolerance on the function value, a positive scalar
,'TolX' , 1e-6 % termination tolerance on x, a positive scalar
% ,'MaxIter' , 3000 % maximum number of iterations allowed, a positive integer
% ,'MaxFunEvals' , 1D6 % maximum number of function evaluations allowed, a positive integer
% ,'UseParallel' , 1 % when true (and supported by optimizer) solver estimates gradients in parallel (using Matlab/Octave's parallel toolbox)
% ,'Jacobian' , 'off' % when 'off' gradient-based solvers approximate Jacobian using finite differences; for GMM we can also pass the analytical Jacobian to gradient-based solvers by setting this 'on'
)
% , silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
% Numerical algorithms options
% , aim_solver % Use AIM algorithm to compute perturbation approximation
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION
% , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm
% , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER
% , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver
% , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver
% , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm
% , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT
% , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver
% , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems
% , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
% , schur_vec_tol=1e-11 % tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix
, mode_check % plot the target function for values around the computed minimum for each estimated parameter in turn
% , mode_check_neighbourhood_size = 5 % width of the window (expressed in percentage deviation) around the computed minimum to be displayed on the diagnostic plots
% , mode_check_symmetric_plots=1 % ensure that the check plots are symmetric around the minimum
% , mode_check_number_of_points = 20 % number of points around the minimum where the target function is evaluated (for each parameter)
);

89
tests/estimation/method_of_moments/RBC_MoM_SMM_ME.mod → tests/estimation/method_of_moments/RBC/RBC_MoM_SMM_ME.mod
View File

@ -1,3 +1,5 @@
% =========================================================================
% Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -137,30 +139,31 @@ end
@#for mommethod in ["SMM"] @#for mommethod in ["SMM"]
method_of_moments( method_of_moments(
% Necessery options % Necessery options
mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM
, datafile = 'RBC_MoM_data_@{orderApp}.mat' % name of filename with data , datafile = 'RBC_MoM_data_@{orderApp}.mat' % name of filename with data
% Options for both GMM and SMM % Options for both GMM and SMM
% , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix % , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation , order = @{orderApp} % order of Taylor approximation in perturbation
% , penalized_estimator % use penalized optimization % , penalized_estimator % include deviation from prior mean as additional moment restriction and use prior precision as weight
, pruning % use pruned state space system at higher-order , pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results % , verbose % display and store intermediate estimation results
, weighting_matrix = ['identity_matrix'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename , weighting_matrix = ['identity_matrix'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename. Size of cell determines stages in iterated estimation, e.g. two state with ['DIAGONAL','OPTIMAL']
, weighting_matrix_scaling_factor = 10 , weighting_matrix_scaling_factor=10 % scaling of weighting matrix in objective function
, burnin=250 % , se_tolx=1e-6 % step size for numerical computation of standard errors
%, additional_optimizer_steps = [4] % vector of additional mode-finders run after mode_compute
% , prefilter=0 % demean each data series by its empirical mean and use centered moments % Options for SMM
% , burnin=250 % number of periods dropped at beginning of simulation
% Options for SMM
% , bounded_shock_support % trim shocks in simulation to +- 2 stdev % , bounded_shock_support % trim shocks in simulation to +- 2 stdev
% , drop = 500 % number of periods dropped at beginning of simulation
% , seed = 24051986 % seed used in simulations % , seed = 24051986 % seed used in simulations
% , simulation_multiple = 5 % multiple of the data length used for simulation % , simulation_multiple = 5 % multiple of the data length used for simulation
%
% General options % Options for GMM
%, dirname = 'MM' % directory in which to store estimation output % , analytic_standard_errors % compute standard errors using analytical derivatives
% General options
% , dirname = 'MM' % directory in which to store estimation output
% , graph_format = EPS % specify the file format(s) for graphs saved to disk % , graph_format = EPS % specify the file format(s) for graphs saved to disk
% , nodisplay % do not display the graphs, but still save them to disk % , nodisplay % do not display the graphs, but still save them to disk
% , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed) % , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
@ -168,41 +171,49 @@ end
% , plot_priors = 1 % control plotting of priors % , plot_priors = 1 % control plotting of priors
% , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters % , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters
% , TeX % print TeX tables and graphics % , TeX % print TeX tables and graphics
%
% Data and model options % Data and model options
%, first_obs = 501 % number of first observation % , first_obs = 501 % number of first observation
% , logdata % if loglinear is set, this option is necessary if the user provides data already in logs, otherwise the log transformation will be applied twice (this may result in complex data) % , logdata % if data is already in logs
% , loglinear % computes a log-linear approximation of the model instead of a linear approximation % , nobs = 250 % number of observations
%, nobs = 500 % number of observations % , prefilter=0 % demean each data series by its empirical mean and use centered moments
% , xls_sheet = willi % name of sheet with data in Excel % , xls_sheet = data % name/number of sheet with data in Excel
% , xls_range = B2:D200 % range of data in Excel sheet % , xls_range = B2:D200 % range of data in Excel sheet
%
% Optimization options that can be set by the user in the mod file, otherwise default values are provided % Optimization options that can be set by the user in the mod file, otherwise default values are provided
% , analytic_derivation % uses analytic derivatives to compute standard errors for GMM % , huge_number=1e7 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons , mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance
, mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer %, additional_optimizer_steps = [1 2 3 4] % vector of additional mode-finders run after mode_compute
%, optim = ('TolFun', 1e-3 % optim: a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute, some exemplary common options:
% ,'TolX', 1e-5 % , optim = ('TolFun' , 1e-6 % termination tolerance on the function value, a positive scalar
% ) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute % ,'TolX' , 1e-6 % termination tolerance on x, a positive scalar
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between % ,'MaxIter' , 3000 % maximum number of iterations allowed, a positive integer
% , tolf = 1e-5 % convergence criterion on function value for numerical differentiation % ,'MaxFunEvals' , 1D6 % maximum number of function evaluations allowed, a positive integer
% , tolx = 1e-6 % convergence criterion on funciton input for numerical differentiation % ,'UseParallel' , 1 % when true (and supported by optimizer) solver estimates gradients in parallel (using Matlab/Octave's parallel toolbox)
% % ,'Jacobian' , 'off' % when 'off' gradient-based solvers approximate Jacobian using finite differences; for GMM we can also pass the analytical Jacobian to gradient-based solvers by setting this 'on'
% % Numerical algorithms options % )
% , silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
% Numerical algorithms options
% , aim_solver % Use AIM algorithm to compute perturbation approximation % , aim_solver % Use AIM algorithm to compute perturbation approximation
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION % , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION
% , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm % , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm % , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the cycle reduction algorithm
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER % , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER
% , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver % , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver
% , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver % , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver
% , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm % , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm
% , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT % , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT
% , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver % , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver
% , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems [IS THIS CORRET @wmutschl] % , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems
% , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition % , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
% , schur_vec_tol=1e-11 % tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix
% , mode_check % plot the target function for values around the computed minimum for each estimated parameter in turn
% , mode_check_neighbourhood_size = 5 % width of the window (expressed in percentage deviation) around the computed minimum to be displayed on the diagnostic plots
% , mode_check_symmetric_plots=1 % ensure that the check plots are symmetric around the minimum
% , mode_check_number_of_points = 20 % number of points around the minimum where the target function is evaluated (for each parameter)
); );
@#endfor @#endfor

18
tests/estimation/method_of_moments/RBC_MoM_common.inc → tests/estimation/method_of_moments/RBC/RBC_MoM_common.inc
View File

@ -2,7 +2,23 @@
% Andreasen, Fernández-Villaverde, Rubio-Ramírez (2018): "The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications", Review of Economic Studies, 85(1):1-49. % Andreasen, Fernández-Villaverde, Rubio-Ramírez (2018): "The Pruned State-Space System for Non-Linear DSGE Models: Theory and Empirical Applications", Review of Economic Studies, 85(1):1-49.
% Adapted by Willi Mutschler (@wmutschl, willi@mutschler.eu) % Adapted by Willi Mutschler (@wmutschl, willi@mutschler.eu)
% ========================================================================= % =========================================================================
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 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/>.
% =========================================================================
var k $K$ var k $K$
c $C$ c $C$

146
tests/estimation/method_of_moments/RBC/RBC_MoM_optimizer.mod Normal file
View File

@ -0,0 +1,146 @@
% Test optimizers
%
% Copyright (C) 2020-2021 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/>.
% =========================================================================
% TO DO
% [ ] fix optimizers 11 and 12;
% note that 12 and 102 require GADS_Toolbox which is not available on servers, but need to be tested locally
% Define testscenario
@#define orderApp = 2
% Note that we will set the numerical optimization tolerance levels very large to speed up the testsuite
@#include "RBC_MoM_common.inc"
shocks;
var u_a; stderr 0.0072;
end;
varobs c iv n;
%--------------------------------------------------------------------------
% Method of Moments Estimation
%--------------------------------------------------------------------------
matched_moments;
c;
n;
iv;
c*c;
c*iv;
iv*n;
iv*iv;
n*c;
n*n;
c*c(-1);
n*n(-1);
iv*iv(-1);
end;
% reduce options to speed up testsuite
options_.newrat.maxiter = 10;
options_.newrat.tolerance.f = 1e-2;
options_.newrat.tolerance.f_analytic = 1e-2;
options_.mh_jscale = 0.6;
options_.gmhmaxlik.iterations=1;
options_.gmhmaxlik.number=2000;
options_.gmhmaxlik.nclimb=2000;
options_.gmhmaxlik.nscale=2000;
options_.gmhmaxlik.target=0.5;
options_.solveopt.MaxIter=300;
options_.solveopt.LBGradientStep=1e-3;
options_.solveopt.TolFun = 1e-3;
options_.solveopt.TolX = 1e-3;
options_.solveopt.TolXConstraint=1e-3;
@#for estimParams in [0, 1, 2]
clear estim_params_;
@#if estimParams == 0
estimated_params;
%DELTA, 0.025;
%BETTA, 0.984;
%B, 0.5;
%ETAc, 2;
ALFA, 0.667;
RHOA, 0.979;
stderr u_a, 0.0072;
end;
@#define OPTIMIZERS = [1, 2, 3, 4, 5, 7, 8, 9, 10, 13, 101]
@#endif
@#if estimParams == 1
estimated_params;
%DELTA, , 0, 1;
%BETTA, , 0, 1;
%B, , 0, 1;
%ETAc, , 0, 10;
ALFA, , 0, 1;
RHOA, , 0, 1;
stderr u_a, , 0, 1;
end;
@#define OPTIMIZERS = [1, 2, 3, 4, 7, 8, 9, 10, 13, 101]
@#endif
@#if estimParams == 2
estimated_params;
%DELTA, 0.025, 0, 1, normal_pdf, 0.02, 0.5;
%BETTA, 0.98, 0, 1, beta_pdf, 0.90, 0.25;
%B, 0.45, 0, 1, normal_pdf, 0.40, 0.5;
%ETAl, 1, 0, 10, normal_pdf, 0.25, 0.0.1;
%ETAc, 1.8, 0, 10, normal_pdf, 1.80, 0.5;
ALFA, 0.65, 0, 1, normal_pdf, 0.60, 0.5;
RHOA, 0.95, 0, 1, normal_pdf, 0.90, 0.5;
stderr u_a, 0.01, 0, 1, normal_pdf, 0.01, 0.5;
%THETA, 3.48, 0, 10, normal_pdf, 0.25, 0.0.1;
end;
@#define OPTIMIZERS = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 13, 101]
@#endif
estimated_params_init(use_calibration);
end;
@#for optimizer in OPTIMIZERS
method_of_moments(
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'RBC_Andreasen_Data_2.mat' % name of filename with data
, order = @{orderApp} % order of Taylor approximation in perturbation
, weighting_matrix = ['OPTIMAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename. Size of cell determines stages in iterated estimation, e.g. two state with ['DIAGONAL','OPTIMAL']
, nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
, mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance
@#if optimizer == 102
, optim = ('TolFun' , 1D-3 % termination tolerance on the function value, a positive scalar
,'MaxIter' , 300 % maximum number of iterations allowed, a positive integer
,'MaxFunEvals' , 1D3 % maximum number of function evaluations allowed, a positive integer
)
@#else
, optim = ('TolFun' , 1D-3 % termination tolerance on the function value, a positive scalar
,'TolX' , 1e-3 % termination tolerance on x, a positive scalar
,'MaxIter' , 300 % maximum number of iterations allowed, a positive integer
,'MaxFunEvals' , 1D3 % maximum number of function evaluations allowed, a positive integer
)
@#endif
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
);
@#endfor
@#endfor

89
tests/estimation/method_of_moments/RBC_MoM_prefilter.mod → tests/estimation/method_of_moments/RBC/RBC_MoM_prefilter.mod
View File

@ -1,6 +1,6 @@
% Tests SMM and GMM routines with prefilter, explicit initialization, and estimated_params_init(use_calibration); % Tests SMM and GMM routines with prefilter, explicit initialization, and estimated_params_init(use_calibration);
% %
% Copyright (C) 2020 Dynare Team % Copyright (C) 2020-2021 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -110,31 +110,31 @@ save('test_matrix.mat','weighting_matrix')
@#for mommethod in ["GMM", "SMM"] @#for mommethod in ["GMM", "SMM"]
method_of_moments( method_of_moments(
% Necessery options % Necessery options
mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM mom_method = @{mommethod} % method of moments method; possible values: GMM|SMM
, datafile = 'RBC_MoM_data_@{orderApp}.mat' % name of filename with data , datafile = 'RBC_MoM_data_@{orderApp}.mat' % name of filename with data
% Options for both GMM and SMM % Options for both GMM and SMM
% , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix % , bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation , order = @{orderApp} % order of Taylor approximation in perturbation
% , penalized_estimator % use penalized optimization % , penalized_estimator % include deviation from prior mean as additional moment restriction and use prior precision as weight
, pruning % use pruned state space system at higher-order , pruning % use pruned state space system at higher-order
% , verbose % display and store intermediate estimation results % , verbose % display and store intermediate estimation results
% , weighting_matrix = 'test_matrix.mat' % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename , weighting_matrix = ['test_matrix.mat','optimal'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename. Size of cell determines stages in iterated estimation, e.g. two state with ['DIAGONAL','OPTIMAL']
, weighting_matrix =['test_matrix.mat','optimal'] %, weighting_matrix_scaling_factor=1 % scaling of weighting matrix in objective function
%, weighting_matrix = optimal % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename , se_tolx=1e-5 % step size for numerical computation of standard errors
%, additional_optimizer_steps = [4] % vector of additional mode-finders run after mode_compute
, prefilter=1 % demean each data series by its empirical mean and use centered moments % Options for SMM
, se_tolx=1e-5 , burnin=500 % number of periods dropped at beginning of simulation
%
% Options for SMM
% , bounded_shock_support % trim shocks in simulation to +- 2 stdev % , bounded_shock_support % trim shocks in simulation to +- 2 stdev
, burnin = 500 % number of periods dropped at beginning of simulation
% , seed = 24051986 % seed used in simulations % , seed = 24051986 % seed used in simulations
% , simulation_multiple = 5 % multiple of the data length used for simulation % , simulation_multiple = 5 % multiple of the data length used for simulation
%
% General options % Options for GMM
%, dirname = 'MM' % directory in which to store estimation output % , analytic_standard_errors % compute standard errors using analytical derivatives
% General options
% , dirname = 'MM' % directory in which to store estimation output
% , graph_format = EPS % specify the file format(s) for graphs saved to disk % , graph_format = EPS % specify the file format(s) for graphs saved to disk
% , nodisplay % do not display the graphs, but still save them to disk % , nodisplay % do not display the graphs, but still save them to disk
% , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed) % , nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
@ -142,38 +142,49 @@ save('test_matrix.mat','weighting_matrix')
% , plot_priors = 1 % control plotting of priors % , plot_priors = 1 % control plotting of priors
% , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters % , prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters
% , TeX % print TeX tables and graphics % , TeX % print TeX tables and graphics
%
% Data and model options % Data and model options
%, first_obs = 501 % number of first observation % , first_obs = 501 % number of first observation
% , logdata % if loglinear is set, this option is necessary if the user provides data already in logs, otherwise the log transformation will be applied twice (this may result in complex data) % , logdata % if data is already in logs
% , loglinear % computes a log-linear approximation of the model instead of a linear approximation , nobs = 250 % number of observations
%, nobs = 500 % number of observations , prefilter=1 % demean each data series by its empirical mean and use centered moments
% , xls_sheet = willi % name of sheet with data in Excel
% , xls_sheet = data % name/number of sheet with data in Excel
% , xls_range = B2:D200 % range of data in Excel sheet % , xls_range = B2:D200 % range of data in Excel sheet
%
% Optimization options that can be set by the user in the mod file, otherwise default values are provided % Optimization options that can be set by the user in the mod file, otherwise default values are provided
% , analytic_derivation % uses analytic derivatives to compute standard errors for GMM % , huge_number=1e7 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons , mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance
, mode_compute = @{optimizer} % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer %, additional_optimizer_steps = [7] % vector of additional mode-finders run after mode_compute
%, optim = ('TolFun', 1e-3 % optim: a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute, some exemplary common options:
% ,'TolX', 1e-5 % , optim = ('TolFun' , 1e-6 % termination tolerance on the function value, a positive scalar
% ) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute % ,'TolX' , 1e-6 % termination tolerance on x, a positive scalar
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between % ,'MaxIter' , 3000 % maximum number of iterations allowed, a positive integer
% % ,'MaxFunEvals' , 1D6 % maximum number of function evaluations allowed, a positive integer
% % Numerical algorithms options % ,'UseParallel' , 1 % when true (and supported by optimizer) solver estimates gradients in parallel (using Matlab/Octave's parallel toolbox)
% ,'Jacobian' , 'off' % when 'off' gradient-based solvers approximate Jacobian using finite differences; for GMM we can also pass the analytical Jacobian to gradient-based solvers by setting this 'on'
% )
% , silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
% Numerical algorithms options
% , aim_solver % Use AIM algorithm to compute perturbation approximation % , aim_solver % Use AIM algorithm to compute perturbation approximation
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION % , dr=default % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION
% , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm % , dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm % , dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm
% , dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the cycle reduction algorithm
% , k_order_solver % use k_order_solver in higher order perturbation approximations
% , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER % , lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER
% , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver % , lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver
% , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver % , lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver
% , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm % , lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm
% , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT % , sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT
% , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver % , sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver
% , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems [IS THIS CORRET @wmutschl] % , qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems
% , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition % , qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
% , schur_vec_tol=1e-11 % tolerance level used to find nonstationary variables in Schur decomposition of the transition matrix
% , mode_check % plot the target function for values around the computed minimum for each estimated parameter in turn
% , mode_check_neighbourhood_size = 5 % width of the window (expressed in percentage deviation) around the computed minimum to be displayed on the diagnostic plots
% , mode_check_symmetric_plots=1 % ensure that the check plots are symmetric around the minimum
% , mode_check_number_of_points = 20 % number of points around the minimum where the target function is evaluated (for each parameter)
); );
@#endfor @#endfor

40
tests/estimation/method_of_moments/RBC/RBC_MoM_steady_helper.m Normal file
View File

@ -0,0 +1,40 @@
% =========================================================================
% Copyright (C) 2020-2021 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/>.
% =========================================================================
function [N, info]= RBC_MoM_steady_helper(THETA,ETAl,ETAc,BETTA,B,C_O_N,W)
info=0;
if ~isreal(C_O_N)
info=1;
N=NaN;
return;
end
if ETAc == 1 && ETAl == 1
N = (1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA/(1+(1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA);
else
% No closed-form solution use a fixed-point algorithm
N0 = 1/3;
try
[N, ~, exitflag] = fsolve(@(N) THETA*(1-N)^(-ETAl)*N^ETAc - (1-BETTA*B)*(C_O_N*(1-B))^(-ETAc)*W, N0,optimset('Display','off','TolX',1e-12,'TolFun',1e-12));
if exitflag<1
info=1;
end
catch
N=NaN;
info=1;
end
end

227
tests/estimation/method_of_moments/RBC_MoM_Andreasen.mod
View File

@ -1,227 +0,0 @@
% Tests SMM and GMM routines
%
% Copyright (C) 2020 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/>.
% =========================================================================
% Define testscenario
@#define orderApp = 2
@#define estimParams = 1
% Note that we will set the numerical optimization tolerance levels very large to speed up the testsuite
@#define optimizer = 13
@#include "RBC_MoM_common.inc"
shocks;
var u_a; stderr 0.0072;
end;
varobs c iv n;
@#if estimParams == 0
estimated_params;
DELTA, 0.025;
BETTA, 0.984;
B, 0.5;
ETAc, 2;
ALFA, 0.667;
RHOA, 0.979;
stderr u_a, 0.0072;
end;
@#endif
@#if estimParams == 1
estimated_params;
DELTA, , 0, 1;
BETTA, , 0, 1;
B, , 0, 1;
ETAc, , 0, 10;
ALFA, , 0, 1;
RHOA, , 0, 1;
stderr u_a, , 0, 1;
end;
@#endif
@#if estimParams == 2
estimated_params;
DELTA, 0.025, 0, 1, normal_pdf, 0.02, 0.5;
BETTA, 0.98, 0, 1, beta_pdf, 0.90, 0.25;
B, 0.45, 0, 1, normal_pdf, 0.40, 0.5;
%ETAl, 1, 0, 10, normal_pdf, 0.25, 0.0.1;
ETAc, 1.8, 0, 10, normal_pdf, 1.80, 0.5;
ALFA, 0.65, 0, 1, normal_pdf, 0.60, 0.5;
RHOA, 0.95, 0, 1, normal_pdf, 0.90, 0.5;
stderr u_a, 0.01, 0, 1, normal_pdf, 0.01, 0.5;
%THETA, 3.48, 0, 10, normal_pdf, 0.25, 0.0.1;
end;
@#endif
% Simulate data
%stoch_simul(order=@{orderApp},pruning,nodisplay,nomoments,periods=500);
%save('RBC_MoM_data_@{orderApp}.mat', options_.varobs{:} );
%pause(1);
estimated_params_init(use_calibration);
end;
%--------------------------------------------------------------------------
% Method of Moments Estimation
%--------------------------------------------------------------------------
matched_moments;
c;
n;
iv;
c*c;
c*iv;
iv*n;
iv*iv;
n*c;
n*n;
c*c(-1);
n*n(-1);
iv*iv(-1);
c*c(-3);
n*n(-3);
iv*iv(-3);
c*c(-5);
n*n(-5);
iv*iv(-5);
end;
% get indices in declaration order
ic = strmatch('c', M_.endo_names,'exact');
iiv = strmatch('iv', M_.endo_names,'exact');
in = strmatch('n', M_.endo_names,'exact');
% first entry: number of variable in declaration order
% second entry: lag
% third entry: power
matched_moments_ = {
[ic ] [0 ], [1 ];
[in ] [0 ], [1 ];
[iiv ] [0 ], [1 ];
[ic ic ] [0 0], [1 1];
[ic iiv] [0 0], [1 1];
%[ic in ] [0 0], [1 1];
%[iiv ic ] [0 0], [1 1];
[in iiv] [0 0], [1 1];
[iiv iiv] [0 0], [1 1];
[ic in] [0 0], [1 1];
%[in iiv] [0 0], [1 1];
[in in ] [0 0], [1 1];
[ic ic ] [0 -1], [1 1];
[in in ] [0 -1], [1 1];
[iiv iiv] [0 -1], [1 1];
[ic ic ] [0 -3], [1 1];
[in in ] [0 -3], [1 1];
[iiv iiv] [0 -3], [1 1];
[ic ic ] [0 -5], [1 1];
[in in ] [0 -5], [1 1];
[iiv iiv] [0 -5], [1 1];
};
if ~isequal(M_.matched_moments,matched_moments_)
error('Translation to matched_moments-block failed')
end
method_of_moments(
% Necessery options
mom_method = GMM % method of moments method; possible values: GMM|SMM
, datafile = 'RBC_Andreasen_Data_2.mat' % name of filename with data
% Options for both GMM and SMM
%, bartlett_kernel_lag = 20 % bandwith in optimal weighting matrix
, order = @{orderApp} % order of Taylor approximation in perturbation
%, penalized_estimator % use penalized optimization
%, pruning % use pruned state space system at higher-order
%, verbose % display and store intermediate estimation results
, weighting_matrix = ['DIAGONAL','OPTIMAL'] % weighting matrix in moments distance objective function; possible values: OPTIMAL|IDENTITY_MATRIX|DIAGONAL|filename
%, weighting_matrix_scaling_factor=1
, additional_optimizer_steps = [13] % vector of additional mode-finders run after mode_compute
%, prefilter=0 % demean each data series by its empirical mean and use centered moments
%
% Options for SMM
%, bounded_shock_support % trim shocks in simulation to +- 2 stdev
%, drop = 500 % number of periods dropped at beginning of simulation
%, seed = 24051986 % seed used in simulations
%, simulation_multiple = 5 % multiple of the data length used for simulation
%, burnin = 200
%
% General options
%, dirname = 'MM' % directory in which to store estimation output
%, graph_format = EPS % specify the file format(s) for graphs saved to disk
%, nodisplay % do not display the graphs, but still save them to disk
%, nograph % do not create graphs (which implies that they are not saved to the disk nor displayed)
%, noprint % do not print stuff to console
%, plot_priors = 1 % control plotting of priors
%, prior_trunc = 1e-10 % probability of extreme values of the prior density that is ignored when computing bounds for the parameters
, TeX % print TeX tables and graphics
%
% Data and model options
%, first_obs = 501 % number of first observation
%, logdata % if loglinear is set, this option is necessary if the user provides data already in logs, otherwise the log transformation will be applied twice (this may result in complex data)
%, loglinear % computes a log-linear approximation of the model instead of a linear approximation
%, nobs = 50 % number of observations
% , xls_sheet = willi % name of sheet with data in Excel
% , xls_range = B2:D200 % range of data in Excel sheet
%
% Optimization options that can be set by the user in the mod file, otherwise default values are provided
%, analytic_derivation % uses analytic derivatives to compute standard errors for GMM
%, huge_number=1D10 % value for replacing the infinite bounds on parameters by finite numbers. Used by some optimizers for numerical reasons
, mode_compute = 13 % specifies the optimizer for minimization of moments distance, note that by default there is a new optimizer
, optim = ('TolFun', 1D-6
,'TolX', 1D-6
) % a list of NAME and VALUE pairs to set options for the optimization routines. Available options depend on mode_compute
%, silent_optimizer % run minimization of moments distance silently without displaying results or saving files in between
, se_tolx = 1e-6 % convergence criterion on funciton input for numerical differentiation
%
% % Numerical algorithms options
%, aim_solver % Use AIM algorithm to compute perturbation approximation
%, dr=DEFAULT % method used to compute the decision rule; possible values are DEFAULT, CYCLE_REDUCTION, LOGARITHMIC_REDUCTION
%, dr_cycle_reduction_tol = 1e-7 % convergence criterion used in the cycle reduction algorithm
%, dr_logarithmic_reduction_maxiter = 100 % maximum number of iterations used in the logarithmic reduction algorithm
%, dr_logarithmic_reduction_tol = 1e-12 % convergence criterion used in the cycle reduction algorithm
%, k_order_solver % use k_order_solver in higher order perturbation approximations
%, lyapunov = DEFAULT % algorithm used to solve lyapunov equations; possible values are DEFAULT, FIXED_POINT, DOUBLING, SQUARE_ROOT_SOLVER
%, lyapunov_complex_threshold = 1e-15 % complex block threshold for the upper triangular matrix in symmetric Lyapunov equation solver
%, lyapunov_fixed_point_tol = 1e-10 % convergence criterion used in the fixed point Lyapunov solver
%, lyapunov_doubling_tol = 1e-16 % convergence criterion used in the doubling algorithm
%, sylvester = default % algorithm to solve Sylvester equation; possible values are DEFAULT, FIXED_POINT
%, sylvester_fixed_point_tol = 1e-12 % convergence criterion used in the fixed point Sylvester solver
%, qz_criterium = 0.999999 % value used to split stable from unstable eigenvalues in reordering the Generalized Schur decomposition used for solving first order problems [IS THIS CORRET @wmutschl]
%, qz_zero_threshold = 1e-6 % value used to test if a generalized eigenvalue is 0/0 in the generalized Schur decomposition
, mode_check
%, mode_check_neighbourhood_size=0.5
%, mode_check_symmetric_plots=0
%, mode_check_number_of_points=25
);

22
tests/estimation/method_of_moments/RBC_MoM_steady_helper.m
View File

@ -1,22 +0,0 @@
function [N, info]= RBC_MoM_steady_helper(THETA,ETAl,ETAc,BETTA,B,C_O_N,W)
info=0;
if ~isreal(C_O_N)
info=1;
N=NaN;
return;
end
if ETAc == 1 && ETAl == 1
N = (1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA/(1+(1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA);
else
% No closed-form solution use a fixed-point algorithm
N0 = 1/3;
try
[N, ~, exitflag] = fsolve(@(N) THETA*(1-N)^(-ETAl)*N^ETAc - (1-BETTA*B)*(C_O_N*(1-B))^(-ETAc)*W, N0,optimset('Display','off','TolX',1e-12,'TolFun',1e-12));
if exitflag<1
info=1;
end
catch
N=NaN;
info=1;
end
end

74
tests/estimation/method_of_moments/RBC_MoM_steadystate.m
View File

@ -1,74 +0,0 @@
% By Willi Mutschler, September 26, 2016. Email: willi@mutschler.eu
function [ys,params,check] = RBCmodel_steadystate(ys,exo,M_,options_)
%% Step 0: initialize indicator and set options for numerical solver
check = 0;
options = optimset('Display','off','TolX',1e-12,'TolFun',1e-12);
params = M_.params;
%% Step 1: read out parameters to access them with their name
for ii = 1:M_.param_nbr
eval([ M_.param_names{ii} ' = M_.params(' int2str(ii) ');']);
end
%% Step 2: Check parameter restrictions
if ETAc*ETAl<1 % parameter violates restriction (here it is artifical)
check=1; %set failure indicator
return; %return without updating steady states
end
%% Step 3: Enter model equations here
A = 1;
RK = 1/BETTA - (1-DELTA);
K_O_N = (RK/(A*(1-ALFA)))^(-1/ALFA);
if K_O_N <= 0
check = 1; % set failure indicator
return; % return without updating steady states
end
W = A*ALFA*(K_O_N)^(1-ALFA);
IV_O_N = DELTA*K_O_N;
Y_O_N = A*K_O_N^(1-ALFA);
C_O_N = Y_O_N - IV_O_N;
if C_O_N <= 0
check = 1; % set failure indicator
return; % return without updating steady states
end
% The labor level
if ETAc == 1 && ETAl == 1
N = (1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA/(1+(1-BETTA*B)*(C_O_N*(1-B))^-1*W/THETA);
else
% No closed-form solution use a fixed-point algorithm
N0 = 1/3;
[N,~,exitflag] = fsolve(@(N) THETA*(1-N)^(-ETAl)*N^ETAc - (1-BETTA*B)*(C_O_N*(1-B))^(-ETAc)*W, N0,options);
if exitflag <= 0
check = 1; % set failure indicator
return % return without updating steady states
end
end
C=C_O_N*N;
Y=Y_O_N*N;
IV=IV_O_N*N;
K=K_O_N*N;
LA = (C-B*C)^(-ETAc)-BETTA*B*(C-B*C)^(-ETAc);
k=log(K);
c=log(C);
a=log(A);
iv=log(IV);
y=log(Y);
la=log(LA);
n=log(N);
rk=log(RK);
w=log(W);
%% Step 4: Update parameters and variables
params=NaN(M_.param_nbr,1);
for iter = 1:M_.param_nbr %update parameters set in the file
eval([ 'params(' num2str(iter) ') = ' M_.param_names{iter} ';' ])
end
for ii = 1:M_.orig_endo_nbr %auxiliary variables are set automatically
eval(['ys(' int2str(ii) ') = ' M_.endo_names{ii} ';']);
end
end