LevMar
- class sherpa.optmethods.LevMar(name='levmar')[source] [edit on github]
Bases:
OptMethod
Levenberg-Marquardt optimization method.
The Levenberg-Marquardt method is an interface to the MINPACK subroutine lmdif to find the local minimum of nonlinear least squares functions of several variables by a modification of the Levenberg-Marquardt algorithm 1.
- ftol
The function tolerance to terminate the search for the minimum; the default is FLT_EPSILON ~ 1.19209289551e-07, where FLT_EPSILON is the smallest number x such that
1.0 != 1.0 + x
. The conditions are satisfied when both the actual and predicted relative reductions in the sum of squares are, at most, ftol.- Type
number
- xtol
The relative error desired in the approximate solution; default is FLT_EPSILON ~ 1.19209289551e-07, where FLT_EPSILON is the smallest number x such that
1.0 != 1.0 + x
. The conditions are satisfied when the relative error between two consecutive iterates is, at most,xtol
.- Type
number
- gtol
The orthogonality desired between the function vector and the columns of the jacobian; default is FLT_EPSILON ~ 1.19209289551e-07, where FLT_EPSILON is the smallest number x such that
1.0 != 1.0 + x
. The conditions are satisfied when the cosine of the angle between fvec and any column of the jacobian is, at most,gtol
in absolute value.- Type
number
- maxfev
The maximum number of function evaluations; the default value of
None
means to use1024 * n
, wheren
is the number of free parameters.- Type
int or
None
- epsfcn
This is used in determining a suitable step length for the forward-difference approximation; default is FLT_EPSILON ~ 1.19209289551e-07, where FLT_EPSILON is the smallest number x such that
1.0 != 1.0 + x
. This approximation assumes that the relative errors in the functions are of the order ofepsfcn
. Ifepsfcn
is less than the machine precision, it is assumed that the relative errors in the functions are of the order of the machine precision.- Type
number
- factor
Used in determining the initial step bound; default is 100. The initial step bound is set to the product of
factor
and the euclidean norm of diag*x if nonzero, or else to factor itself. In most cases,factor
should be from the interval (.1,100.).- Type
- verbose
The amount of information to print during the fit. The default is
0
, which means no output.- Type
References
- 1
J.J. More, “The Levenberg Marquardt algorithm: implementation and theory,” in Lecture Notes in Mathematics 630: Numerical Analysis, G.A. Watson (Ed.), Springer-Verlag: Berlin, 1978, pp.105-116.
Attributes Summary
The default settings for the optimiser.
Methods Summary
fit
(statfunc, pars, parmins, parmaxes[, ...])Run the optimiser.
Attributes Documentation
- default_config
The default settings for the optimiser.
Methods Documentation
- fit(statfunc, pars, parmins, parmaxes, statargs=(), statkwargs={}) [edit on github]
Run the optimiser.
- Parameters
statfunc (function) – Given a list of parameter values as the first argument and, as the remaining positional arguments,
statargs
andstatkwargs
as keyword arguments, return the statistic value.pars (sequence) – The start position of the model parameter values.
parmins (sequence) – The minimum allowed values for each model parameter. This must match the length of
pars
.parmaxes (sequence) – The maximum allowed values for each model parameter. This must match the length of
pars
.statargs (optional) – Additional positional arguments to send to
statfunc
.statkwargs (optional) – Additional keyword arguments to send to
statfunc
.
- Returns
newpars – The tuple contains: boolean indicating whether the optimization succeeded or not, the best fit parameters as a NumPy array, the statistic value at the best-fit location, a string message indicating the status, and a dictionary containing information about the optimisation (this depends on the optimiser).
- Return type