public class LevenbergMarquardtOptimizer extends java.lang.Object implements LeastSquaresOptimizer
This implementation should work even for over-determined systems (i.e. systems having more point than equations). Over-determined systems are solved by ignoring the point which have the smallest impact according to their jacobian column norm. Only the rank of the matrix and some loop bounds are changed to implement this.
The resolution engine is a simple translation of the MINPACK lmder routine with minor changes. The changes include the over-determined resolution, the use of inherited convergence checker and the Q.R. decomposition which has been rewritten following the algorithm described in the P. Lascaux and R. Theodor book Analyse numérique matricielle appliquée à l'art de l'ingénieur, Masson 1986.
The authors of the original fortran version are:
Minpack Copyright Notice (1999) University of Chicago. All rights reserved |
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions
are met:
|
Modifier and Type | Class and Description |
---|---|
private static class |
LevenbergMarquardtOptimizer.InternalData
Holds internal data.
|
LeastSquaresOptimizer.Optimum
Modifier and Type | Field and Description |
---|---|
private double |
costRelativeTolerance
Desired relative error in the sum of squares.
|
private double |
initialStepBoundFactor
Positive input variable used in determining the initial step bound.
|
private double |
orthoTolerance
Desired max cosine on the orthogonality between the function vector
and the columns of the jacobian.
|
private double |
parRelativeTolerance
Desired relative error in the approximate solution parameters.
|
private double |
qrRankingThreshold
Threshold for QR ranking.
|
private static double |
TWO_EPS
Twice the "epsilon machine".
|
Constructor and Description |
---|
LevenbergMarquardtOptimizer()
Default constructor.
|
LevenbergMarquardtOptimizer(double initialStepBoundFactor,
double costRelativeTolerance,
double parRelativeTolerance,
double orthoTolerance,
double qrRankingThreshold)
Construct an instance with all parameters specified.
|
Modifier and Type | Method and Description |
---|---|
private void |
determineLMDirection(double[] qy,
double[] diag,
double[] lmDiag,
LevenbergMarquardtOptimizer.InternalData internalData,
int solvedCols,
double[] work,
double[] lmDir)
Solve a*x = b and d*x = 0 in the least squares sense.
|
private double |
determineLMParameter(double[] qy,
double delta,
double[] diag,
LevenbergMarquardtOptimizer.InternalData internalData,
int solvedCols,
double[] work1,
double[] work2,
double[] work3,
double[] lmDir,
double lmPar)
Determines the Levenberg-Marquardt parameter.
|
double |
getCostRelativeTolerance()
Gets the value of a tuning parameter.
|
double |
getInitialStepBoundFactor()
Gets the value of a tuning parameter.
|
double |
getOrthoTolerance()
Gets the value of a tuning parameter.
|
double |
getParameterRelativeTolerance()
Gets the value of a tuning parameter.
|
double |
getRankingThreshold()
Gets the value of a tuning parameter.
|
LeastSquaresOptimizer.Optimum |
optimize(LeastSquaresProblem problem)
Solve the non-linear least squares problem.
|
private LevenbergMarquardtOptimizer.InternalData |
qrDecomposition(RealMatrix jacobian,
int solvedCols)
Decompose a matrix A as A.P = Q.R using Householder transforms.
|
private void |
qTy(double[] y,
LevenbergMarquardtOptimizer.InternalData internalData)
Compute the product Qt.y for some Q.R.
|
LevenbergMarquardtOptimizer |
withCostRelativeTolerance(double newCostRelativeTolerance) |
LevenbergMarquardtOptimizer |
withInitialStepBoundFactor(double newInitialStepBoundFactor) |
LevenbergMarquardtOptimizer |
withOrthoTolerance(double newOrthoTolerance)
Modifies the given parameter.
|
LevenbergMarquardtOptimizer |
withParameterRelativeTolerance(double newParRelativeTolerance) |
LevenbergMarquardtOptimizer |
withRankingThreshold(double newQRRankingThreshold) |
private static final double TWO_EPS
private final double initialStepBoundFactor
private final double costRelativeTolerance
private final double parRelativeTolerance
private final double orthoTolerance
private final double qrRankingThreshold
public LevenbergMarquardtOptimizer()
The default values for the algorithm settings are:
Precision.SAFE_MIN
public LevenbergMarquardtOptimizer(double initialStepBoundFactor, double costRelativeTolerance, double parRelativeTolerance, double orthoTolerance, double qrRankingThreshold)
initialStepBoundFactor
- initial step bound factorcostRelativeTolerance
- cost relative toleranceparRelativeTolerance
- parameters relative toleranceorthoTolerance
- orthogonality toleranceqrRankingThreshold
- threshold in the QR decomposition. Columns with a 2
norm less than this threshold are considered to be
all 0s.public LevenbergMarquardtOptimizer withInitialStepBoundFactor(double newInitialStepBoundFactor)
newInitialStepBoundFactor
- Positive input variable used in
determining the initial step bound. This bound is set to the
product of initialStepBoundFactor and the euclidean norm of
diag * x
if non-zero, or else to newInitialStepBoundFactor
itself. In most cases factor should lie in the interval
(0.1, 100.0)
. 100
is a generally recommended value.
of the matrix is reduced.public LevenbergMarquardtOptimizer withCostRelativeTolerance(double newCostRelativeTolerance)
newCostRelativeTolerance
- Desired relative error in the sum of squares.public LevenbergMarquardtOptimizer withParameterRelativeTolerance(double newParRelativeTolerance)
newParRelativeTolerance
- Desired relative error in the approximate solution
parameters.public LevenbergMarquardtOptimizer withOrthoTolerance(double newOrthoTolerance)
newOrthoTolerance
- Desired max cosine on the orthogonality between
the function vector and the columns of the Jacobian.public LevenbergMarquardtOptimizer withRankingThreshold(double newQRRankingThreshold)
newQRRankingThreshold
- Desired threshold for QR ranking.
If the squared norm of a column vector is smaller or equal to this
threshold during QR decomposition, it is considered to be a zero vector
and hence the rank of the matrix is reduced.public double getInitialStepBoundFactor()
withInitialStepBoundFactor(double)
public double getCostRelativeTolerance()
withCostRelativeTolerance(double)
public double getParameterRelativeTolerance()
withParameterRelativeTolerance(double)
public double getOrthoTolerance()
withOrthoTolerance(double)
public double getRankingThreshold()
withRankingThreshold(double)
public LeastSquaresOptimizer.Optimum optimize(LeastSquaresProblem problem)
optimize
in interface LeastSquaresOptimizer
problem
- the problem definition, including model function and
convergence criteria.private double determineLMParameter(double[] qy, double delta, double[] diag, LevenbergMarquardtOptimizer.InternalData internalData, int solvedCols, double[] work1, double[] work2, double[] work3, double[] lmDir, double lmPar)
This implementation is a translation in Java of the MINPACK lmpar routine.
This method sets the lmPar and lmDir attributes.
The authors of the original fortran function are:
Luc Maisonobe did the Java translation.
qy
- Array containing qTy.delta
- Upper bound on the euclidean norm of diagR * lmDir.diag
- Diagonal matrix.internalData
- Data (modified in-place in this method).solvedCols
- Number of solved point.work1
- work arraywork2
- work arraywork3
- work arraylmDir
- the "returned" LM direction will be stored in this array.lmPar
- the value of the LM parameter from the previous iteration.private void determineLMDirection(double[] qy, double[] diag, double[] lmDiag, LevenbergMarquardtOptimizer.InternalData internalData, int solvedCols, double[] work, double[] lmDir)
This implementation is a translation in Java of the MINPACK qrsolv routine.
This method sets the lmDir and lmDiag attributes.
The authors of the original fortran function are:
Luc Maisonobe did the Java translation.
qy
- array containing qTydiag
- diagonal matrixlmDiag
- diagonal elements associated with lmDirinternalData
- Data (modified in-place in this method).solvedCols
- Number of sloved point.work
- work arraylmDir
- the "returned" LM direction is stored in this arrayprivate LevenbergMarquardtOptimizer.InternalData qrDecomposition(RealMatrix jacobian, int solvedCols) throws ConvergenceException
As suggested in the P. Lascaux and R. Theodor book Analyse numérique matricielle appliquée à l'art de l'ingénieur (Masson, 1986), instead of representing the Householder transforms with uk unit vectors such that:
Hk = I - 2uk.uktwe use k non-unit vectors such that:
Hk = I - betakvk.vktwhere vk = ak - alphak ek. The betak coefficients are provided upon exit as recomputing them from the vk vectors would be costly.
This decomposition handles rank deficient cases since the tranformations are performed in non-increasing columns norms order thanks to columns pivoting. The diagonal elements of the R matrix are therefore also in non-increasing absolute values order.
jacobian
- Weighted Jacobian matrix at the current point.solvedCols
- Number of solved point.ConvergenceException
- if the decomposition cannot be performed.private void qTy(double[] y, LevenbergMarquardtOptimizer.InternalData internalData)
y
- vector to multiply (will be overwritten with the result)internalData
- Data.Copyright (c) 2003-2016 Apache Software Foundation