From 24218a5278ca3d99c209efc990d16048e29c8536 Mon Sep 17 00:00:00 2001 From: Luc Maisonobe Date: Sat, 29 May 2010 18:15:50 +0000 Subject: [PATCH] Fixed Levenberg-Marquardt optimizer that did not use the vectorial convergence checker. Now this optimizer can use either the general vectorial convergence checker or its own specialized convergence settings. Minor changes had to be introduced in the test data, they have been validated JIRA: MATH-362 git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@949433 13f79535-47bb-0310-9956-ffa450edef68 --- .../general/LevenbergMarquardtOptimizer.java | 41 ++++++++++++++----- src/site/xdoc/changes.xml | 5 +++ .../LevenbergMarquardtOptimizerTest.java | 2 +- .../optimization/general/MinpackTest.java | 25 +++++------ 4 files changed, 48 insertions(+), 25 deletions(-) diff --git a/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java b/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java index ea0f20da1..28aee4348 100644 --- a/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java +++ b/src/main/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizer.java @@ -34,8 +34,9 @@ import org.apache.commons.math.optimization.VectorialPointValuePair; * *

The resolution engine is a simple translation of the MINPACK lmder routine with minor - * changes. The changes include the over-determined resolution and the Q.R. - * decomposition which has been rewritten following the algorithm described in the + * 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: @@ -143,6 +144,7 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { * Build an optimizer for least squares problems. *

The default values for the algorithm settings are: *

*

+ *

These default values may be overridden after construction. If the {@link + * #setConvergenceChecker vectorial convergence checker} is set to a non-null value, it + * will be used instead of the {@link #setCostRelativeTolerance cost relative tolerance} + * and {@link #setParRelativeTolerance parameters relative tolerance} settings. */ public LevenbergMarquardtOptimizer() { @@ -157,6 +163,7 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { setMaxIterations(1000); // default values for the tuning parameters + setConvergenceChecker(null); setInitialStepBoundFactor(100.0); setCostRelativeTolerance(1.0e-10); setParRelativeTolerance(1.0e-10); @@ -179,7 +186,8 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { /** * Set the desired relative error in the sum of squares. - * + *

This setting is used only if the {@link #setConvergenceChecker vectorial + * convergence checker} is set to null.

* @param costRelativeTolerance desired relative error in the sum of squares */ public void setCostRelativeTolerance(double costRelativeTolerance) { @@ -188,7 +196,8 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { /** * Set the desired relative error in the approximate solution parameters. - * + *

This setting is used only if the {@link #setConvergenceChecker vectorial + * convergence checker} is set to null.

* @param parRelativeTolerance desired relative error * in the approximate solution parameters */ @@ -198,7 +207,8 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { /** * Set the desired max cosine on the orthogonality. - * + *

This setting is always used, regardless of the {@link #setConvergenceChecker + * vectorial convergence checker} being null or non-null.

* @param orthoTolerance desired max cosine on the orthogonality * between the function vector and the columns of the jacobian */ @@ -235,11 +245,13 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { // outer loop lmPar = 0; boolean firstIteration = true; + VectorialPointValuePair current = new VectorialPointValuePair(point, objective); while (true) { incrementIterationsCounter(); // compute the Q.R. decomposition of the jacobian matrix + VectorialPointValuePair previous = current; updateJacobian(); qrDecomposition(); @@ -291,7 +303,7 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { } if (maxCosine <= orthoTolerance) { // convergence has been reached - return new VectorialPointValuePair(point, objective); + return current; } // rescale if necessary @@ -333,6 +345,7 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { // evaluate the function at x + p and calculate its norm updateResidualsAndCost(); + current = new VectorialPointValuePair(point, objective); // compute the scaled actual reduction double actRed = -1.0; @@ -401,11 +414,19 @@ public class LevenbergMarquardtOptimizer extends AbstractLeastSquaresOptimizer { } // tests for convergence. - if (((Math.abs(actRed) <= costRelativeTolerance) && - (preRed <= costRelativeTolerance) && - (ratio <= 2.0)) || + if (checker != null) { + // we use the vectorial convergence checker + if (checker.converged(getIterations(), previous, current)) { + return current; + } + } else { + // we use the Levenberg-Marquardt specific convergence parameters + if (((Math.abs(actRed) <= costRelativeTolerance) && + (preRed <= costRelativeTolerance) && + (ratio <= 2.0)) || (delta <= parRelativeTolerance * xNorm)) { - return new VectorialPointValuePair(point, objective); + return current; + } } // tests for termination and stringent tolerances diff --git a/src/site/xdoc/changes.xml b/src/site/xdoc/changes.xml index 69b7f604b..e6611f35d 100644 --- a/src/site/xdoc/changes.xml +++ b/src/site/xdoc/changes.xml @@ -52,6 +52,11 @@ The type attribute can be add,update,fix,remove. If the output is not quite correct, check for invisible trailing spaces! --> + + Fixed Levenberg-Marquardt optimizer that did not use the vectorial convergence checker. + Now this optimizer can use either the general vectorial convergence checker or its own + specialized convergence settings. + Fixed loss of significance error in PersonsCorrelation p-value computation causing p-values smaller than the machine epsilon (~1E-16) to be reported as 0. diff --git a/src/test/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizerTest.java b/src/test/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizerTest.java index aabeab6bc..06f1721fb 100644 --- a/src/test/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizerTest.java +++ b/src/test/java/org/apache/commons/math/optimization/general/LevenbergMarquardtOptimizerTest.java @@ -485,7 +485,7 @@ public class LevenbergMarquardtOptimizerTest circle.addPoint(points[i][0], points[i][1]); } LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer(); - optimizer.setConvergenceChecker(new SimpleVectorialValueChecker(1.0e-10, 1.0e-10)); + optimizer.setConvergenceChecker(new SimpleVectorialValueChecker(1.0e-8, 1.0e-8)); VectorialPointValuePair optimum = optimizer.optimize(circle, target, weights, new double[] { -12, -12 }); Point2D.Double center = new Point2D.Double(optimum.getPointRef()[0], optimum.getPointRef()[1]); diff --git a/src/test/java/org/apache/commons/math/optimization/general/MinpackTest.java b/src/test/java/org/apache/commons/math/optimization/general/MinpackTest.java index 7ff333113..59467b224 100644 --- a/src/test/java/org/apache/commons/math/optimization/general/MinpackTest.java +++ b/src/test/java/org/apache/commons/math/optimization/general/MinpackTest.java @@ -152,14 +152,14 @@ public class MinpackTest extends TestCase { minpackTest(new FreudensteinRothFunction(new double[] { 5.0, -20.0 }, 12432.833948863, 6.9988751744895, new double[] { - 11.4130046614746, - -0.896796038685958 + 11.4121122022341, + -0.8968550851268697 }), false); minpackTest(new FreudensteinRothFunction(new double[] { 50.0, -200.0 }, 11426454.595762, 6.99887517242903, new double[] { - 11.4127817857886, - -0.89680510749204 + 11.412069435091231, + -0.8968582807605691 }), false); } @@ -325,7 +325,7 @@ public class MinpackTest extends TestCase { minpackTest(new JennrichSampsonFunction(10, new double[] { 0.3, 0.4 }, 64.5856498144943, 11.1517793413499, new double[] { - 0.257819926636811, 0.257829976764542 + 0.2578330049, 0.257829976764542 }), false); } @@ -499,8 +499,8 @@ public class MinpackTest extends TestCase { function.getTarget(), function.getWeight(), function.getStartPoint()); assertFalse(exceptionExpected); - assertTrue(function.checkTheoreticalMinCost(optimizer.getRMS())); - assertTrue(function.checkTheoreticalMinParams(optimum)); + function.checkTheoreticalMinCost(optimizer.getRMS()); + function.checkTheoreticalMinParams(optimum); } catch (OptimizationException lsse) { assertTrue(exceptionExpected); } catch (FunctionEvaluationException fe) { @@ -561,23 +561,20 @@ public class MinpackTest extends TestCase { return startParams.length; } - public boolean checkTheoreticalMinCost(double rms) { + public void checkTheoreticalMinCost(double rms) { double threshold = costAccuracy * (1.0 + theoreticalMinCost); - return Math.abs(Math.sqrt(m) * rms - theoreticalMinCost) <= threshold; + assertEquals(theoreticalMinCost, Math.sqrt(m) * rms, threshold); } - public boolean checkTheoreticalMinParams(VectorialPointValuePair optimum) { + public void checkTheoreticalMinParams(VectorialPointValuePair optimum) { double[] params = optimum.getPointRef(); if (theoreticalMinParams != null) { for (int i = 0; i < theoreticalMinParams.length; ++i) { double mi = theoreticalMinParams[i]; double vi = params[i]; - if (Math.abs(mi - vi) > (paramsAccuracy * (1.0 + Math.abs(mi)))) { - return false; - } + assertEquals(mi, vi, paramsAccuracy * (1.0 + Math.abs(mi))); } } - return true; } public MultivariateMatrixFunction jacobian() {