Commit 6a627122 authored by Sébastien Villemot's avatar Sébastien Villemot

Merge branch 'minimize_objective' into 'master'

Minimize objective

See merge request Dynare/dynare!1737
parents 37362723 4a9f9263
......@@ -513,7 +513,7 @@ newrat.tolerance.f=1e-5;
newrat.tolerance.f_analytic=1e-7;
newrat.maxiter=1000;
newrat.verbosity=1;
newrat.Save_files=1;
newrat.Save_files=0;
options_.newrat=newrat;
......
......@@ -265,7 +265,7 @@ if ~isequal(options_.mode_compute,0) && ~options_.mh_posterior_mode_estimation
if compute_hessian
crit = options_.newrat.tolerance.f;
newratflag = newratflag>0;
hh = reshape(mr_hessian(xparam1,objective_function,fval,newratflag,crit,new_rat_hess_info,dataset_, dataset_info, options_,M_,estim_params_,bayestopt_,bounds,oo_), nx, nx);
hh = reshape(mr_hessian(xparam1,objective_function,fval,newratflag,crit,new_rat_hess_info,[bounds.lb bounds.ub],bayestopt_.p2,dataset_, dataset_info, options_,M_,estim_params_,bayestopt_,bounds,oo_), nx, nx);
end
options_.kalman_algo = kalman_algo0;
end
......
......@@ -222,6 +222,9 @@ switch minimizer_algorithm
csminwel1(objective_function, start_par_value, H0, analytic_grad, crit, nit, numgrad, epsilon, Verbose, Save_files, varargin{:});
hessian_mat=inv(inverse_hessian_mat);
case 5
if isempty(prior_information) %mr_hessian requires it, but can be NaN
prior_information.p2=NaN(n_params,1);
end
if options_.analytic_derivation==-1 %set outside as code for use of analytic derivation
analytic_grad=1;
crit = options_.newrat.tolerance.f_analytic;
......@@ -265,9 +268,11 @@ switch minimizer_algorithm
hess_info.gstep=options_.gstep;
hess_info.htol = 1.e-4;
hess_info.h1=options_.gradient_epsilon*ones(n_params,1);
[opt_par_values,hessian_mat,gg,fval,invhess,new_rat_hess_info] = newrat(objective_function,start_par_value,bounds,analytic_grad,crit,nit,0,Verbose, Save_files,hess_info,varargin{:});
%hessian_mat is the plain outer product gradient Hessian
[opt_par_values,hessian_mat,gg,fval,invhess,new_rat_hess_info] = newrat(objective_function,start_par_value,bounds,analytic_grad,crit,nit,0,Verbose, Save_files,hess_info,prior_information.p2,varargin{:}); %hessian_mat is the plain outer product gradient Hessian
case 6
if isempty(prior_information) %Inf will be reset
prior_information.p2=Inf(n_params,1);
end
[opt_par_values, hessian_mat, Scale, fval] = gmhmaxlik(objective_function, start_par_value, ...
Initial_Hessian, options_.mh_jscale, bounds, prior_information.p2, options_.gmhmaxlik, options_.optim_opt, varargin{:});
case 7
......
function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,func,penalty,hflag,htol0,hess_info,varargin)
% function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,func,penalty,hflag,htol0,hess_info,varargin)
function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,func,penalty,hflag,htol0,hess_info,bounds,prior_std,varargin)
% function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,func,penalty,hflag,htol0,hess_info,bounds,prior_std,varargin)
% numerical gradient and Hessian, with 'automatic' check of numerical
% error
%
......@@ -24,7 +24,10 @@ function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,f
% derivatives
% - hess_info structure storing the step sizes for
% computation of Hessian
% - varargin other inputs:
% - bounds prior bounds of parameters
% - prior_std prior standard devation of parameters (can be NaN)
% - varargin other inputs
% e.g. in dsge_likelihood
% varargin{1} --> DynareDataset
% varargin{2} --> DatasetInfo
% varargin{3} --> DynareOptions
......@@ -64,9 +67,9 @@ function [hessian_mat, gg, htol1, ihh, hh_mat0, hh1, hess_info] = mr_hessian(x,f
n=size(x,1);
[f0,exit_flag, ff0]=penalty_objective_function(x,func,penalty,varargin{:});
h2=varargin{7}.ub-varargin{7}.lb;
hmax=varargin{7}.ub-x;
hmax=min(hmax,x-varargin{7}.lb);
h2=bounds(:,2)-bounds(:,1);
hmax=bounds(:,2)-x;
hmax=min(hmax,x-bounds(:,1));
if isempty(ff0)
outer_product_gradient=0;
else
......@@ -240,7 +243,7 @@ if outer_product_gradient
sd=sqrt(diag(ihh)); %standard errors
sdh=sqrt(1./diag(hh)); %diagonal standard errors
for j=1:length(sd)
sd0(j,1)=min(varargin{6}.p2(j), sd(j)); %prior std
sd0(j,1)=min(prior_std(j), sd(j)); %prior std
sd0(j,1)=10^(0.5*(log10(sd0(j,1))+log10(sdh(j,1))));
end
ihh=ihh./(sd*sd').*(sd0*sd0'); %inverse outer product with modified std's
......
function [xparam1, hh, gg, fval, igg, hess_info] = newrat(func0, x, bounds, analytic_derivation, ftol0, nit, flagg, Verbose, Save_files, hess_info, varargin)
function [xparam1, hh, gg, fval, igg, hess_info] = newrat(func0, x, bounds, analytic_derivation, ftol0, nit, flagg, Verbose, Save_files, hess_info, prior_std, varargin)
% [xparam1, hh, gg, fval, igg, hess_info] = newrat(func0, x, bounds, analytic_derivation, ftol0, nit, flagg, Verbose, Save_files, hess_info, varargin)
%
% Optimiser with outer product gradient and with sequences of univariate steps
......@@ -8,6 +8,7 @@ function [xparam1, hh, gg, fval, igg, hess_info] = newrat(func0, x, bounds, anal
% - func0 name of the function that also outputs the single contributions at times t=1,...,T
% of the log-likelihood to compute outer product gradient
% - x starting guess
% - bounds prior bounds of parameters
% - analytic_derivation 1 if analytic derivatives, 0 otherwise
% - ftol0 termination criterion for function change
% - nit maximum number of iterations
......@@ -21,7 +22,10 @@ function [xparam1, hh, gg, fval, igg, hess_info] = newrat(func0, x, bounds, anal
% - Save_files 1 if intermediate output is to be saved
% - hess_info structure storing the step sizes for
% computation of Hessian
% - varargin other inputs:
% - prior_std prior standard devation of parameters (can be NaN);
% passed to mr_hessian
% - varargin other inputs
% e.g. in dsge_likelihood and others:
% varargin{1} --> DynareDataset
% varargin{2} --> DatasetInfo
% varargin{3} --> DynareOptions
......@@ -69,7 +73,6 @@ ftol=ftol0;
gtol=1.e-3;
htol=htol_base;
htol0=htol_base;
gibbstol=length(varargin{6}.pshape)/50; %25;
% force fcn, grad to function handle
if ischar(func0)
......@@ -85,7 +88,7 @@ fval=fval0;
outer_product_gradient=1;
if isempty(hh)
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(x,func0,penalty,flagit,htol,hess_info,varargin{:});
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(x,func0,penalty,flagit,htol,hess_info,bounds,prior_std,varargin{:});
if isempty(dum)
outer_product_gradient=0;
igg = 1e-4*eye(nx);
......@@ -203,7 +206,7 @@ while norm(gg)>gtol && check==0 && jit<nit
if flagit==2
hh=hh0;
elseif flagg>0
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(xparam1,func0,penalty,flagg,ftol0,hess_info,varargin{:});
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(xparam1,func0,penalty,flagg,ftol0,hess_info,bounds,prior_std,varargin{:});
if flagg==2
hh = reshape(dum,nx,nx);
ee=eig(hh);
......@@ -243,7 +246,7 @@ while norm(gg)>gtol && check==0 && jit<nit
save('m1.mat','x','fval0','nig')
end
end
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(xparam1,func0,penalty,flagit,htol,hess_info,varargin{:});
[dum, gg, htol0, igg, hhg, h1, hess_info]=mr_hessian(xparam1,func0,penalty,flagit,htol,hess_info,bounds,prior_std,varargin{:});
if isempty(dum)
outer_product_gradient=0;
end
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment