This is an implementation of Powell's conjugate gradient descent
method that might go in the minimization section of optim.
Best,
Nir
## Copyright (C) 2011 Nir Krakauer
##
## This program 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.
##
## This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
## -*- texinfo -*-
## @deftypefn {Function File} [@var{p}, @var{obj_value}, @var{convergence}, @var{iters}, @var{nevs}] = powell (@var{f}, @var{args}, @var{control})
##powell: implements a direction-set (Powell's) method for minimizing a function of several variables without calculation of the gradient [1, 2]
##
## @subheading Arguments
##
## @itemize @bullet
## @item
## @var{f}: name of function to minimize (string or handle)
##
## @item
## @var{args}: a cell array that holds all arguments of the function
## The argument with respect to which minimization is done
## must be numeric
##
## @item
## @var{control}: an optional cell array of 1-10 elements. If a cell array shorter than 10 elements is provided, the trailing elements are provided with default values. (currently, only elements 1, 4, 6, 9, 10 are used)
## @itemize @minus
## @item
## elem 1: maximum iterations (positive integer, or -1 or Inf for unlimited (default))
## @item
## elem 2: verbosity
## @itemize
## @item
## 0 = no screen output (default)
## @item
## 1 = only final results
## @item
## 2 = summary every iteration
## @item
## 3 = detailed information
## @end itemize
## @item
## elem 3: convergence criterion
## @itemize
## @item
## 1 = strict (function, gradient and param change) (default)
## @item
## 0 = weak - only function convergence required
## @end itemize
## @item
## elem 4: order of argument in args with respect to which minimization is done (default is 1, first argument to the objective function)
## @item
## elem 5: (optional) Memory limit
## @item
## elem 6: function change tolerance, default 1e-8
## @item
## elem 7: parameter change tolerance, default 1e-6
## @item
## elem 8: gradient tolerance, default 1e-5
## @item
## elem 9: maximum function evaluations (positive integer, or -1 or Inf for unlimited (default))
## @item
## elem 10: an n*n matrix containing the initial set of (presumably orthogonal) directions to minimize along, where n is the number of elements in the argument to be minimized for; or an n*1 vector of magnitudes for the initial directions (defaults to the set of unit direction vectors)
## @end itemize
## @end itemize
##
## @subheading Examples
##
## @example
## @group
## y = @@(x, s) x(1) .^ 2 + x(2) .^ 2 + s;
## [x_optim, y_min, conv, iters, nevs] = powell(y, @{[1 0.5], 1@});
## %should return something like x_optim = [5E-14 1E-14], y_min = 1, conv = 1, iters = 2, nevs = 24
## @end group
##
## @end example
##
## @subheading Returns:
##
## @itemize @bullet
## @item
## @var{p}: the minimizing value of the function argument
## @item
## @var{obj_value}: the value of @var{f}() at @var{p}
## @item
## @var{convergence}: 1 if normal convergence, 0 if not
## @item
## @var{iters}: number of iterations performed
## @item
## @var{nevs}: number of function evaluations
## @end itemize
##
## @subheading References
##
## @enumerate
## @item
## Powell MJD (1964), An efficient method for finding the minimum of a function of several variables without calculating derivatives, @cite{Computer Journal}, 7 :155-162
##
## @item
## Press, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP (1992). @cite{Numerical Recipes in Fortran: The Art of Scientific Computing} (2nd Ed.). New York: Cambridge University Press (Section 10.5)
## @end enumerate
## @end deftypefn
## Author: Nir Krakauer <nkraka...@ccny.cuny.edu>
## Description: Multidimensional minimization (direction-set method)
function [p, obj_value, convergence, iters, nevs] = powell(f, args, control);
# check number and types of arguments
if ((nargin < 2) || (nargin > 3))
usage('powell: you must supply 2 or 3 arguments');
endif
if (~iscell(args))
usage('powell: the second argument should be a cell array of function input arguments');
end
%default values
narg = 1; %position of minimized variable
ftol = 1E-8;
maxiter = Inf;
maxev = Inf;
xi_set = 0;
%read in user settings, if provided
if (nargin == 3) && iscell(control)
if (numel(control) >= 1) && ~isempty(control{1}) && isnumeric(control{1}) && (control{1} >= 0)
maxiter = control{1};
end
if (numel(control) >= 4) && ~isempty(control{4}) && isnumeric(control{4}) && (control{4} >= 1)
narg = control{4};
end
if (numel(control) >= 6) && ~isempty(control{6}) && isnumeric(control{6}) && (control{6} >= 0)
ftol = control{6};
end
if (numel(control) >= 9) && ~isempty(control{9}) && isnumeric(control{9}) && (control{9} >= 0)
maxev = control{9};
end
if (numel(control) >= 10) && ~isempty(control{10}) && isnumeric(control{10})
if isvector(control{10}) %control{10} is n*1 or 1*n
xi = diag(control{10});
else %control{10} is n*n
xi = control{10};
end
xi_set = 1;
end
end
nevs = 0;
iters = 0;
convergence = 0;
try
obj_value = feval(f, args{:});
catch
error("function does not exist or cannot be evaluated");
end
nevs++;
n = numel(args{narg}); %number of dimensions to minimize over
xit = zeros(n, 1);
if ~xi_set
xi = eye(n);
end
p = args{narg}; %initial value of the argument being minimized
%do an iteration
while (iters <= maxiter) && (nevs <= maxev) && (~convergence)
iters++;
pt = p; %best point as iteration begins
fp = obj_value; %value of the objective function as iteration begins
ibig = 0; %will hold direction along which the objective function decreased the most in this iteration
dlt = 0; %will hold decrease in objective function value in this iteration
for i = 1:n
xit = reshape(xi(:, i), size(p));
fptt = obj_value;
%%args = splice(args, narg, 1, {p});
args{narg} = p;
[a, obj_value, nev] = line_min(f, xit, args, narg);
%[a, obj_value, nev] = linmin(f, args, xit, narg, ftol, narg);
nevs = nevs + nev;
p = p + a*xit;
change = fptt - obj_value;
if (change > dlt)
dlt = change;
ibig = i;
end
end
if ( 2*abs(fp-obj_value) <= ftol*(abs(fp) + abs(obj_value)) )
convergence = 1;
return
end
if iters == maxiter
disp('iteration maximum exceeded')
return
end
%attempt parabolic extrapolation
ptt = 2*p - pt;
xit = p - pt;
%%args = splice(args, narg, 1, {ptt});
args{narg} = ptt;
fptt = feval(f, args{:});
nevs++;
if fptt < fp %check whether the extrapolation actually makes the objective function smaller
t = 2 * (fp - 2*obj_value + fptt) * (fp-obj_value-dlt)^2 - dlt * (fp-fptt)^2;
if t < 0
p = ptt;
%%args = splice(args, narg, 1, {p});
args{narg} = p;
[a,obj_value,nev] = line_min(f, xit, args, narg);
nevs = nevs + nev;
p = p + a*xit;
%add the net direction from this iteration to the direction set
xi(:, ibig) = xi(:, n);
xi(:, n) = xit(:);
end
end
end
------------------------------------------------------------------------------
All the data continuously generated in your IT infrastructure
contains a definitive record of customers, application performance,
security threats, fraudulent activity, and more. Splunk takes this
data and makes sense of it. IT sense. And common sense.
http://p.sf.net/sfu/splunk-novd2d
_______________________________________________
Octave-dev mailing list
Octave-dev@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/octave-dev
