diff --git a/matlab/+mom/Jtest.m b/matlab/+mom/Jtest.m index 254d55cec..5c5001a3d 100644 --- a/matlab/+mom/Jtest.m +++ b/matlab/+mom/Jtest.m @@ -1,5 +1,5 @@ -function oo_ = Jtest(xparam, objective_function, Woptflag, oo_, options_mom_, bayestopt_, Bounds, estim_params_, M_, nobs) -% function oo_ = Jtest(xparam, objective_function, Woptflag, oo_, options_mom_, bayestopt_, Bounds, estim_params_, M_, nobs) +function oo_ = Jtest(xparam, objective_function, Woptflag, oo_, options_mom_, bayestopt_, BoundsInfo, estim_params_, M_, nobs) +% function oo_ = Jtest(xparam, objective_function, Woptflag, oo_, options_mom_, bayestopt_, BoundsInfo, estim_params_, M_, nobs) % ------------------------------------------------------------------------- % Computes the J-test statistic and p-value for a GMM/SMM estimation % ========================================================================= @@ -10,7 +10,7 @@ function oo_ = Jtest(xparam, objective_function, Woptflag, oo_, options_mom_, ba % oo_: [struct] results % options_mom_: [struct] options % bayestopt_: [struct] information on priors -% Bounds: [struct] bounds on parameters +% BoundsInfo: [struct] bounds on parameters % estim_params_: [struct] information on estimated parameters % M_: [struct] information on the model % nobs: [scalar] number of observations @@ -49,7 +49,7 @@ if options_mom_.mom.mom_nbr > length(xparam) W_opt = mom.optimal_weighting_matrix(oo_.mom.m_data, oo_.mom.model_moments, options_mom_.mom.bartlett_kernel_lag); oo_J = oo_; oo_J.mom.Sw = chol(W_opt); - fval = feval(objective_function, xparam, Bounds, oo_J, estim_params_, M_, options_mom_); + fval = feval(objective_function, xparam, BoundsInfo, oo_J, estim_params_, M_, options_mom_); else fval = oo_.mom.Q; end diff --git a/matlab/+mom/mode_compute_gmm_smm.m b/matlab/+mom/mode_compute_gmm_smm.m index 2c1fc5ef6..5514b0bd6 100644 --- a/matlab/+mom/mode_compute_gmm_smm.m +++ b/matlab/+mom/mode_compute_gmm_smm.m @@ -1,5 +1,5 @@ -function [xparam1, oo_, Woptflag] = mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, Bounds) -% function [xparam1, oo_, Woptflag] = mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, Bounds) +function [xparam1, oo_, Woptflag] = mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, BoundsInfo) +% function [xparam1, oo_, Woptflag] = mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, BoundsInfo) % ------------------------------------------------------------------------- % Iterated method of moments for GMM and SMM, computes the minimum of the % objective function (distance between data moments and model moments) @@ -14,7 +14,7 @@ function [xparam1, oo_, Woptflag] = mode_compute_gmm_smm(xparam0, objective_func % options_mom_: [structure] options % estim_params_: [structure] information on estimated parameters % bayestopt_: [structure] information on priors -% Bounds: [structure] bounds for optimization +% BoundsInfo: [structure] bounds for optimization % ------------------------------------------------------------------------- % OUTPUT % xparam1: [vector] mode of objective function @@ -100,7 +100,7 @@ for stage_iter = 1:size(options_mom_.mom.weighting_matrix,1) options_mom_.current_optimizer = options_mom_.optimizer_vec{optim_iter}; if options_mom_.optimizer_vec{optim_iter} == 0 xparam1 = xparam0; % no minimization, evaluate objective at current values - fval = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_); + fval = feval(objective_function, xparam1, BoundsInfo, oo_, estim_params_, M_, options_mom_); else if options_mom_.optimizer_vec{optim_iter} == 13 options_mom_.mom.vector_output = true; @@ -112,8 +112,8 @@ for stage_iter = 1:size(options_mom_.mom.weighting_matrix,1) else options_mom_.mom.compute_derivs = false; end - [xparam1, fval] = dynare_minimize_objective(objective_function, xparam0, options_mom_.optimizer_vec{optim_iter}, options_mom_, [Bounds.lb Bounds.ub], bayestopt_.name, bayestopt_, [],... - Bounds, oo_, estim_params_, M_, options_mom_); + [xparam1, fval] = dynare_minimize_objective(objective_function, xparam0, options_mom_.optimizer_vec{optim_iter}, options_mom_, [BoundsInfo.lb BoundsInfo.ub], bayestopt_.name, bayestopt_, [],... + BoundsInfo, oo_, estim_params_, M_, options_mom_); if options_mom_.mom.vector_output fval = fval'*fval; end @@ -125,5 +125,5 @@ for stage_iter = 1:size(options_mom_.mom.weighting_matrix,1) xparam0 = xparam1; end options_mom_.vector_output = false; - [~, ~, ~,~,~, oo_] = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_); % get oo_.mom.model_moments for iterated GMM/SMM to compute optimal weighting matrix + [~, ~, ~,~,~, oo_] = feval(objective_function, xparam1, BoundsInfo, oo_, estim_params_, M_, options_mom_); % get oo_.mom.model_moments for iterated GMM/SMM to compute optimal weighting matrix end diff --git a/matlab/+mom/objective_function.m b/matlab/+mom/objective_function.m index 414467353..c6b06948e 100644 --- a/matlab/+mom/objective_function.m +++ b/matlab/+mom/objective_function.m @@ -1,11 +1,11 @@ -function [fval, info, exit_flag, df, junkHessian, oo_, M_] = objective_function(xparam, Bounds, oo_, estim_params_, M_, options_mom_) -% [fval, info, exit_flag, df, junk1, oo_, M_] = objective_function(xparam, Bounds, oo_, estim_params_, M_, options_mom_) +function [fval, info, exit_flag, df, junkHessian, oo_, M_] = objective_function(xparam, BoundsInfo, oo_, estim_params_, M_, options_mom_) +% [fval, info, exit_flag, df, junk1, oo_, M_] = objective_function(xparam, BoundsInfo, oo_, estim_params_, M_, options_mom_) % ------------------------------------------------------------------------- % This function evaluates the objective function for method of moments estimation % ========================================================================= % INPUTS % o xparam: [vector] current value of estimated parameters as returned by set_prior() -% o Bounds: [structure] containing parameter bounds +% o BoundsInfo: [structure] containing parameter bounds % o oo_: [structure] for results % o estim_params_: [structure] describing the estimated_parameters % o M_ [structure] describing the model @@ -86,7 +86,7 @@ end % Ensure that xparam1 is a column vector; particleswarm.m requires this. xparam = xparam(:); M_ = set_all_parameters(xparam, estim_params_, M_); -[fval,info,exit_flag] = check_bounds_and_definiteness_estimation(xparam, M_, estim_params_, Bounds); +[fval,info,exit_flag] = check_bounds_and_definiteness_estimation(xparam, M_, estim_params_, BoundsInfo); if info(1) if options_mom_.mom.vector_output == 1 % lsqnonlin requires vector output fval = ones(size(oo_.mom.data_moments,1),1)*options_mom_.huge_number; diff --git a/matlab/+mom/run.m b/matlab/+mom/run.m index b5579aeb4..57464756a 100644 --- a/matlab/+mom/run.m +++ b/matlab/+mom/run.m @@ -333,14 +333,14 @@ if ~isempty(bayestopt_) && doBayesianEstimation end end % set prior bounds - Bounds = prior_bounds(bayestopt_, options_mom_.prior_trunc); - Bounds.lb = max(Bounds.lb,lb); - Bounds.ub = min(Bounds.ub,ub); + BoundsInfo = prior_bounds(bayestopt_, options_mom_.prior_trunc); + BoundsInfo.lb = max(BoundsInfo.lb,lb); + BoundsInfo.ub = min(BoundsInfo.ub,ub); else % no priors are declared so Dynare will estimate the parameters with % Frequentist methods using inequality constraints for the parameters - Bounds.lb = lb; - Bounds.ub = ub; + BoundsInfo.lb = lb; + BoundsInfo.ub = ub; if options_mom_.mom.penalized_estimator fprintf('Penalized estimation turned off as you did not declare priors\n') options_mom_.mom.penalized_estimator = 0; @@ -348,18 +348,18 @@ else end % set correct bounds for standard deviations and correlations -Bounds = mom.set_correct_bounds_for_stderr_corr(estim_params_,Bounds); +BoundsInfo = mom.set_correct_bounds_for_stderr_corr(estim_params_,BoundsInfo); % test if initial values of the estimated parameters are all between the prior lower and upper bounds if options_mom_.use_calibration_initialization try - check_prior_bounds(xparam0,Bounds,M_,estim_params_,options_mom_,bayestopt_); + check_prior_bounds(xparam0,BoundsInfo,M_,estim_params_,options_mom_,bayestopt_); catch last_error fprintf('Cannot use parameter values from calibration as they violate the prior bounds.') rethrow(last_error); end else - check_prior_bounds(xparam0,Bounds,M_,estim_params_,options_mom_,bayestopt_); + check_prior_bounds(xparam0,BoundsInfo,M_,estim_params_,options_mom_,bayestopt_); end % check for positive definiteness @@ -493,7 +493,7 @@ try oo_.mom.Sw = eye(options_mom_.mom.mom_nbr); % initialize with identity weighting matrix end tic_id = tic; - [fval, info, ~, ~, ~, oo_, M_] = feval(objective_function, xparam0, Bounds, oo_, estim_params_, M_, options_mom_); + [fval, info, ~, ~, ~, oo_, M_] = feval(objective_function, xparam0, BoundsInfo, oo_, estim_params_, M_, options_mom_); elapsed_time = toc(tic_id); if isnan(fval) error('method_of_moments: The initial value of the objective function with identity weighting matrix is NaN!') @@ -535,16 +535,16 @@ mom.print_info_on_estimation_settings(options_mom_, number_of_estimated_paramete % ------------------------------------------------------------------------- if strcmp(options_mom_.mom.mom_method,'GMM') || strcmp(options_mom_.mom.mom_method,'SMM') % compute mode - [xparam1, oo_, Woptflag] = mom.mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, Bounds); + [xparam1, oo_, Woptflag] = mom.mode_compute_gmm_smm(xparam0, objective_function, oo_, M_, options_mom_, estim_params_, bayestopt_, BoundsInfo); % compute standard errors at mode options_mom_.mom.vector_output = false; % make sure flag is reset M_ = set_all_parameters(xparam1,estim_params_,M_); % update M_ and oo_ (in particular to get oo_.mom.model_moments) if strcmp(options_mom_.mom.mom_method,'GMM') && options_mom_.mom.analytic_standard_errors options_mom_.mom.compute_derivs = true; % for GMM we compute derivatives analytically in the objective function with this flag end - [~, ~, ~,~,~, oo_] = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_); % compute model moments and oo_.mom.model_moments_params_derivs + [~, ~, ~,~,~, oo_] = feval(objective_function, xparam1, BoundsInfo, oo_, estim_params_, M_, options_mom_); % compute model moments and oo_.mom.model_moments_params_derivs options_mom_.mom.compute_derivs = false; % reset to not compute derivatives in objective function during optimization - [stdh,hessian_xparam1] = mom.standard_errors(xparam1, objective_function, Bounds, oo_, estim_params_, M_, options_mom_, Woptflag); + [stdh,hessian_xparam1] = mom.standard_errors(xparam1, objective_function, BoundsInfo, oo_, estim_params_, M_, options_mom_, Woptflag); end @@ -552,8 +552,8 @@ end % checks for mode and hessian at the mode % ------------------------------------------------------------------------- if options_mom_.mode_check.status - mode_check(objective_function, xparam1, hessian_xparam1, options_mom_, M_, estim_params_, bayestopt_, Bounds, true,... - Bounds, oo_, estim_params_, M_, options_mom_); + mode_check(objective_function, xparam1, hessian_xparam1, options_mom_, M_, estim_params_, bayestopt_, BoundsInfo, true,... + BoundsInfo, oo_, estim_params_, M_, options_mom_); end @@ -564,7 +564,7 @@ if strcmp(options_mom_.mom.mom_method,'SMM') || strcmp(options_mom_.mom.mom_meth % Store results in output structure oo_.mom = display_estimation_results_table(xparam1,stdh,M_,options_mom_,estim_params_,bayestopt_,oo_.mom,prior_dist_names,options_mom_.mom.mom_method,lower(options_mom_.mom.mom_method)); % J test - oo_ = mom.Jtest(xparam1, objective_function, Woptflag, oo_, options_mom_, bayestopt_, Bounds, estim_params_, M_, dataset_.nobs); + oo_ = mom.Jtest(xparam1, objective_function, Woptflag, oo_, options_mom_, bayestopt_, BoundsInfo, estim_params_, M_, dataset_.nobs); % display comparison of model moments and data moments mom.display_comparison_moments(M_, options_mom_, oo_.mom.data_moments, oo_.mom.model_moments); end diff --git a/matlab/+mom/set_correct_bounds_for_stderr_corr.m b/matlab/+mom/set_correct_bounds_for_stderr_corr.m index d65e9473b..1bfac50a8 100644 --- a/matlab/+mom/set_correct_bounds_for_stderr_corr.m +++ b/matlab/+mom/set_correct_bounds_for_stderr_corr.m @@ -1,14 +1,14 @@ -function Bounds = set_correct_bounds_for_stderr_corr(estim_params_,Bounds) -% function Bounds = set_correct_bounds_for_stderr_corr(estim_params_,Bounds) +function BoundsInfo = set_correct_bounds_for_stderr_corr(estim_params_,BoundsInfo) +% function BoundsInfo = set_correct_bounds_for_stderr_corr(estim_params_,BoundsInfo) % ------------------------------------------------------------------------- % Set correct bounds for standard deviation and corrrelation parameters % ========================================================================= % INPUTS % o estim_params_ [struct] information on estimated parameters -% o Bounds [struct] information on bounds +% o BoundsInfo [struct] information on bounds % ------------------------------------------------------------------------- % OUTPUT -% o Bounds [struct] updated bounds +% o BoundsInfo [struct] updated bounds % ------------------------------------------------------------------------- % This function is called by % o mom.run @@ -34,10 +34,10 @@ function Bounds = set_correct_bounds_for_stderr_corr(estim_params_,Bounds) number_of_estimated_parameters = estim_params_.nvx+estim_params_.nvn+estim_params_.ncx+estim_params_.ncn+estim_params_.np; % set correct bounds for standard deviations and corrrelations param_of_interest = (1:number_of_estimated_parameters)'<=estim_params_.nvx+estim_params_.nvn; -LB_below_0 = (Bounds.lb<0 & param_of_interest); -Bounds.lb(LB_below_0) = 0; +LB_below_0 = (BoundsInfo.lb<0 & param_of_interest); +BoundsInfo.lb(LB_below_0) = 0; param_of_interest = (1:number_of_estimated_parameters)'> estim_params_.nvx+estim_params_.nvn & (1:number_of_estimated_parameters)'1 & param_of_interest); -Bounds.lb(LB_below_minus_1) = -1; -Bounds.ub(UB_above_1) = 1; \ No newline at end of file +LB_below_minus_1 = (BoundsInfo.lb<-1 & param_of_interest); +UB_above_1 = (BoundsInfo.ub>1 & param_of_interest); +BoundsInfo.lb(LB_below_minus_1) = -1; +BoundsInfo.ub(UB_above_1) = 1; \ No newline at end of file diff --git a/matlab/+mom/standard_errors.m b/matlab/+mom/standard_errors.m index 6bfec1e72..0fac91d94 100644 --- a/matlab/+mom/standard_errors.m +++ b/matlab/+mom/standard_errors.m @@ -1,5 +1,5 @@ -function [SE_values, Asympt_Var] = standard_errors(xparam, objective_function, Bounds, oo_, estim_params_, M_, options_mom_, Wopt_flag) -% [SE_values, Asympt_Var] = standard_errors(xparam, objective_function, Bounds, oo_, estim_params_, M_, options_mom_, Wopt_flag) +function [SE_values, Asympt_Var] = standard_errors(xparam, objective_function, BoundsInfo, oo_, estim_params_, M_, options_mom_, Wopt_flag) +% [SE_values, Asympt_Var] = standard_errors(xparam, objective_function, BoundsInfo, oo_, estim_params_, M_, options_mom_, Wopt_flag) % ------------------------------------------------------------------------- % This function computes standard errors to the method of moments estimates % Adapted from replication codes of @@ -8,7 +8,7 @@ function [SE_values, Asympt_Var] = standard_errors(xparam, objective_function, B % INPUTS % o xparam: value of estimated parameters as returned by set_prior() % o objective_function string of objective function -% o Bounds: structure containing parameter bounds +% o BoundsInfo: structure containing parameter bounds % o oo_: structure for results % o estim_params_: structure describing the estimated_parameters % o M_ structure describing the model @@ -75,12 +75,12 @@ else %Positive step xparam_eps_p = xparam; xparam_eps_p(i,1) = xparam_eps_p(i) + eps_value; - [~, info_p, ~, ~,~, oo__p] = feval(objective_function, xparam_eps_p, Bounds, oo_, estim_params_, M_, options_mom_); + [~, info_p, ~, ~,~, oo__p] = feval(objective_function, xparam_eps_p, BoundsInfo, oo_, estim_params_, M_, options_mom_); % Negative step xparam_eps_m = xparam; xparam_eps_m(i,1) = xparam_eps_m(i) - eps_value; - [~, info_m, ~, ~,~, oo__m] = feval(objective_function, xparam_eps_m, Bounds, oo_, estim_params_, M_, options_mom_); + [~, info_m, ~, ~,~, oo__m] = feval(objective_function, xparam_eps_m, BoundsInfo, oo_, estim_params_, M_, options_mom_); % The Jacobian: if nnz(info_p)==0 && nnz(info_m)==0