You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by er...@apache.org on 2013/07/30 17:04:23 UTC
svn commit: r1508481 [3/5] - in /commons/proper/math/trunk/src:
main/java/org/apache/commons/math3/fitting/leastsquares/
main/java/org/apache/commons/math3/optim/
test/java/org/apache/commons/math3/fitting/leastsquares/
test/resources/org/apache/common...
Added: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java?rev=1508481&view=auto
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java (added)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java Tue Jul 30 15:04:22 2013
@@ -0,0 +1,109 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math3.fitting.leastsquares;
+
+import java.io.IOException;
+import org.apache.commons.math3.exception.ConvergenceException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.exception.MathUnsupportedOperationException;
+import org.apache.commons.math3.optim.SimpleVectorValueChecker;
+import org.apache.commons.math3.linear.DiagonalMatrix;
+import org.junit.Test;
+
+/**
+ * <p>Some of the unit tests are re-implementations of the MINPACK <a
+ * href="http://www.netlib.org/minpack/ex/file17">file17</a> and <a
+ * href="http://www.netlib.org/minpack/ex/file22">file22</a> test files.
+ * The redistribution policy for MINPACK is available <a
+ * href="http://www.netlib.org/minpack/disclaimer">here</a>/
+ *
+ * @version $Id$
+ */
+public class GaussNewtonOptimizerTest
+ extends AbstractLeastSquaresOptimizerAbstractTest<GaussNewtonOptimizer> {
+ @Override
+ public GaussNewtonOptimizer createOptimizer() {
+ return GaussNewtonOptimizer.create()
+ .withConvergenceChecker(new SimpleVectorValueChecker(1e-6, 1e-6));
+ }
+
+ @Override
+ public int getMaxIterations() {
+ return 1000;
+ }
+
+ @Override
+ @Test(expected=ConvergenceException.class)
+ public void testMoreEstimatedParametersSimple() {
+ /*
+ * Exception is expected with this optimizer
+ */
+ super.testMoreEstimatedParametersSimple();
+ }
+
+ @Override
+ @Test(expected=ConvergenceException.class)
+ public void testMoreEstimatedParametersUnsorted() {
+ /*
+ * Exception is expected with this optimizer
+ */
+ super.testMoreEstimatedParametersUnsorted();
+ }
+
+ @Test(expected=TooManyEvaluationsException.class)
+ public void testMaxEvaluations() throws Exception {
+ CircleVectorial circle = new CircleVectorial();
+ circle.addPoint( 30.0, 68.0);
+ circle.addPoint( 50.0, -6.0);
+ circle.addPoint(110.0, -20.0);
+ circle.addPoint( 35.0, 15.0);
+ circle.addPoint( 45.0, 97.0);
+
+ GaussNewtonOptimizer optimizer = createOptimizer()
+ .withConvergenceChecker(new SimpleVectorValueChecker(1e-30, 1e-30))
+ .withMaxIterations(Integer.MAX_VALUE)
+ .withMaxEvaluations(100)
+ .withModelAndJacobian(circle.getModelFunction(),
+ circle.getModelFunctionJacobian())
+ .withTarget(new double[] { 0, 0, 0, 0, 0 })
+ .withWeight(new DiagonalMatrix(new double[] { 1, 1, 1, 1, 1 }))
+ .withStartPoint(new double[] { 98.680, 47.345 });
+
+ optimizer.optimize();
+ }
+
+ @Override
+ @Test(expected=ConvergenceException.class)
+ public void testCircleFittingBadInit() {
+ /*
+ * This test does not converge with this optimizer.
+ */
+ super.testCircleFittingBadInit();
+ }
+
+ @Override
+ @Test(expected=ConvergenceException.class)
+ public void testHahn1()
+ throws IOException {
+ /*
+ * TODO This test leads to a singular problem with the Gauss-Newton
+ * optimizer. This should be inquired.
+ */
+ super.testHahn1();
+ }
+}
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerTest.java
------------------------------------------------------------------------------
svn:keywords = Id Revision
Added: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java?rev=1508481&view=auto
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java (added)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java Tue Jul 30 15:04:22 2013
@@ -0,0 +1,356 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math3.fitting.leastsquares;
+
+import java.io.Serializable;
+import java.util.ArrayList;
+import java.util.List;
+import org.apache.commons.math3.optim.PointVectorValuePair;
+import org.apache.commons.math3.analysis.MultivariateVectorFunction;
+import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
+import org.apache.commons.math3.exception.ConvergenceException;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.TooManyEvaluationsException;
+import org.apache.commons.math3.exception.MathUnsupportedOperationException;
+import org.apache.commons.math3.geometry.euclidean.twod.Vector2D;
+import org.apache.commons.math3.linear.SingularMatrixException;
+import org.apache.commons.math3.linear.DiagonalMatrix;
+import org.apache.commons.math3.util.FastMath;
+import org.apache.commons.math3.util.Precision;
+import org.junit.Assert;
+import org.junit.Test;
+import org.junit.Ignore;
+
+/**
+ * <p>Some of the unit tests are re-implementations of the MINPACK <a
+ * href="http://www.netlib.org/minpack/ex/file17">file17</a> and <a
+ * href="http://www.netlib.org/minpack/ex/file22">file22</a> test files.
+ * The redistribution policy for MINPACK is available <a
+ * href="http://www.netlib.org/minpack/disclaimer">here</a>.
+ *
+ * @version $Id$
+ */
+public class LevenbergMarquardtOptimizerTest
+ extends AbstractLeastSquaresOptimizerAbstractTest<LevenbergMarquardtOptimizer> {
+ @Override
+ public LevenbergMarquardtOptimizer createOptimizer() {
+ return LevenbergMarquardtOptimizer.create();
+ }
+
+ @Override
+ public int getMaxIterations() {
+ return 25;
+ }
+
+ @Override
+ @Test(expected=SingularMatrixException.class)
+ public void testNonInvertible() {
+ /*
+ * Overrides the method from parent class, since the default singularity
+ * threshold (1e-14) does not trigger the expected exception.
+ */
+ LinearProblem problem = new LinearProblem(new double[][] {
+ { 1, 2, -3 },
+ { 2, 1, 3 },
+ { -3, 0, -9 }
+ }, new double[] { 1, 1, 1 });
+
+ final LevenbergMarquardtOptimizer optimizer = createOptimizer()
+ .withMaxEvaluations(100)
+ .withMaxIterations(20)
+ .withModelAndJacobian(problem.getModelFunction(),
+ problem.getModelFunctionJacobian())
+ .withTarget(problem.getTarget())
+ .withWeight(new DiagonalMatrix(new double[] { 1, 1, 1 }))
+ .withStartPoint(new double[] { 0, 0, 0 });
+
+ final double[] optimum = optimizer.optimize().getPoint();
+ Assert.assertTrue(FastMath.sqrt(optimizer.getTarget().length) * optimizer.computeRMS(optimum) > 0.6);
+
+ optimizer.computeCovariances(optimum, 1.5e-14);
+ }
+
+ @Test
+ public void testControlParameters() {
+ CircleVectorial circle = new CircleVectorial();
+ circle.addPoint( 30.0, 68.0);
+ circle.addPoint( 50.0, -6.0);
+ circle.addPoint(110.0, -20.0);
+ circle.addPoint( 35.0, 15.0);
+ circle.addPoint( 45.0, 97.0);
+ checkEstimate(circle.getModelFunction(),
+ circle.getModelFunctionJacobian(),
+ 0.1, 10, 1.0e-14, 1.0e-16, 1.0e-10, false);
+ checkEstimate(circle.getModelFunction(),
+ circle.getModelFunctionJacobian(),
+ 0.1, 10, 1.0e-15, 1.0e-17, 1.0e-10, true);
+ checkEstimate(circle.getModelFunction(),
+ circle.getModelFunctionJacobian(),
+ 0.1, 5, 1.0e-15, 1.0e-16, 1.0e-10, true);
+ circle.addPoint(300, -300);
+ checkEstimate(circle.getModelFunction(),
+ circle.getModelFunctionJacobian(),
+ 0.1, 20, 1.0e-18, 1.0e-16, 1.0e-10, true);
+ }
+
+ private void checkEstimate(MultivariateVectorFunction problem,
+ MultivariateMatrixFunction problemJacobian,
+ double initialStepBoundFactor, int maxCostEval,
+ double costRelativeTolerance, double parRelativeTolerance,
+ double orthoTolerance, boolean shouldFail) {
+ try {
+ final LevenbergMarquardtOptimizer optimizer = LevenbergMarquardtOptimizer.create()
+ .withTuningParameters(initialStepBoundFactor,
+ costRelativeTolerance,
+ parRelativeTolerance,
+ orthoTolerance,
+ Precision.SAFE_MIN)
+ .withMaxEvaluations(maxCostEval)
+ .withMaxIterations(100)
+ .withModelAndJacobian(problem, problemJacobian)
+ .withTarget(new double[] { 0, 0, 0, 0, 0 })
+ .withWeight(new DiagonalMatrix(new double[] { 1, 1, 1, 1, 1 }))
+ .withStartPoint(new double[] { 98.680, 47.345 });
+
+ optimizer.optimize();
+
+ Assert.assertTrue(!shouldFail);
+ } catch (DimensionMismatchException ee) {
+ Assert.assertTrue(shouldFail);
+ } catch (TooManyEvaluationsException ee) {
+ Assert.assertTrue(shouldFail);
+ }
+ }
+
+ /**
+ * Non-linear test case: fitting of decay curve (from Chapter 8 of
+ * Bevington's textbook, "Data reduction and analysis for the physical sciences").
+ * XXX The expected ("reference") values may not be accurate and the tolerance too
+ * relaxed for this test to be currently really useful (the issue is under
+ * investigation).
+ */
+ @Test
+ public void testBevington() {
+ final double[][] dataPoints = {
+ // column 1 = times
+ { 15, 30, 45, 60, 75, 90, 105, 120, 135, 150,
+ 165, 180, 195, 210, 225, 240, 255, 270, 285, 300,
+ 315, 330, 345, 360, 375, 390, 405, 420, 435, 450,
+ 465, 480, 495, 510, 525, 540, 555, 570, 585, 600,
+ 615, 630, 645, 660, 675, 690, 705, 720, 735, 750,
+ 765, 780, 795, 810, 825, 840, 855, 870, 885, },
+ // column 2 = measured counts
+ { 775, 479, 380, 302, 185, 157, 137, 119, 110, 89,
+ 74, 61, 66, 68, 48, 54, 51, 46, 55, 29,
+ 28, 37, 49, 26, 35, 29, 31, 24, 25, 35,
+ 24, 30, 26, 28, 21, 18, 20, 27, 17, 17,
+ 14, 17, 24, 11, 22, 17, 12, 10, 13, 16,
+ 9, 9, 14, 21, 17, 13, 12, 18, 10, },
+ };
+
+ final BevingtonProblem problem = new BevingtonProblem();
+
+ final int len = dataPoints[0].length;
+ final double[] weights = new double[len];
+ for (int i = 0; i < len; i++) {
+ problem.addPoint(dataPoints[0][i],
+ dataPoints[1][i]);
+
+ weights[i] = 1 / dataPoints[1][i];
+ }
+
+ final LevenbergMarquardtOptimizer optimizer = LevenbergMarquardtOptimizer.create()
+ .withMaxEvaluations(100)
+ .withMaxIterations(20)
+ .withModelAndJacobian(problem.getModelFunction(),
+ problem.getModelFunctionJacobian())
+ .withTarget(dataPoints[1])
+ .withWeight(new DiagonalMatrix(weights))
+ .withStartPoint(new double[] { 10, 900, 80, 27, 225 });
+
+ final PointVectorValuePair optimum = optimizer.optimize();
+ final double[] solution = optimum.getPoint();
+ final double[] expectedSolution = { 10.4, 958.3, 131.4, 33.9, 205.0 };
+
+ final double[][] covarMatrix = optimizer.computeCovariances(solution, 1e-14);
+ final double[][] expectedCovarMatrix = {
+ { 3.38, -3.69, 27.98, -2.34, -49.24 },
+ { -3.69, 2492.26, 81.89, -69.21, -8.9 },
+ { 27.98, 81.89, 468.99, -44.22, -615.44 },
+ { -2.34, -69.21, -44.22, 6.39, 53.80 },
+ { -49.24, -8.9, -615.44, 53.8, 929.45 }
+ };
+
+ final int numParams = expectedSolution.length;
+
+ // Check that the computed solution is within the reference error range.
+ for (int i = 0; i < numParams; i++) {
+ final double error = FastMath.sqrt(expectedCovarMatrix[i][i]);
+ Assert.assertEquals("Parameter " + i, expectedSolution[i], solution[i], error);
+ }
+
+ // Check that each entry of the computed covariance matrix is within 10%
+ // of the reference matrix entry.
+ for (int i = 0; i < numParams; i++) {
+ for (int j = 0; j < numParams; j++) {
+ Assert.assertEquals("Covariance matrix [" + i + "][" + j + "]",
+ expectedCovarMatrix[i][j],
+ covarMatrix[i][j],
+ FastMath.abs(0.1 * expectedCovarMatrix[i][j]));
+ }
+ }
+ }
+
+ @Test
+ public void testCircleFitting2() {
+ final double xCenter = 123.456;
+ final double yCenter = 654.321;
+ final double xSigma = 10;
+ final double ySigma = 15;
+ final double radius = 111.111;
+ // The test is extremely sensitive to the seed.
+ final long seed = 59421061L;
+ final RandomCirclePointGenerator factory
+ = new RandomCirclePointGenerator(xCenter, yCenter, radius,
+ xSigma, ySigma,
+ seed);
+ final CircleProblem circle = new CircleProblem(xSigma, ySigma);
+
+ final int numPoints = 10;
+ for (Vector2D p : factory.generate(numPoints)) {
+ circle.addPoint(p.getX(), p.getY());
+ }
+
+ // First guess for the center's coordinates and radius.
+ final double[] init = { 90, 659, 115 };
+
+ final LevenbergMarquardtOptimizer optimizer = LevenbergMarquardtOptimizer.create()
+ .withMaxEvaluations(100)
+ .withMaxIterations(50)
+ .withModelAndJacobian(circle.getModelFunction(),
+ circle.getModelFunctionJacobian())
+ .withTarget(circle.target())
+ .withWeight(new DiagonalMatrix(circle.weight()))
+ .withStartPoint(init);
+
+ final PointVectorValuePair optimum = optimizer.optimize();
+ final double[] paramFound = optimum.getPoint();
+
+ // Retrieve errors estimation.
+ final double[] asymptoticStandardErrorFound = optimizer.computeSigma(paramFound, 1e-14);
+
+ // Check that the parameters are found within the assumed error bars.
+ Assert.assertEquals(xCenter, paramFound[0], asymptoticStandardErrorFound[0]);
+ Assert.assertEquals(yCenter, paramFound[1], asymptoticStandardErrorFound[1]);
+ Assert.assertEquals(radius, paramFound[2], asymptoticStandardErrorFound[2]);
+ }
+
+ private static class QuadraticProblem {
+ private List<Double> x;
+ private List<Double> y;
+
+ public QuadraticProblem() {
+ x = new ArrayList<Double>();
+ y = new ArrayList<Double>();
+ }
+
+ public void addPoint(double x, double y) {
+ this.x.add(x);
+ this.y.add(y);
+ }
+
+ public MultivariateVectorFunction getModelFunction() {
+ return new MultivariateVectorFunction() {
+ public double[] value(double[] variables) {
+ double[] values = new double[x.size()];
+ for (int i = 0; i < values.length; ++i) {
+ values[i] = (variables[0] * x.get(i) + variables[1]) * x.get(i) + variables[2];
+ }
+ return values;
+ }
+ };
+ }
+
+ public MultivariateMatrixFunction getModelFunctionJacobian() {
+ return new MultivariateMatrixFunction() {
+ public double[][] value(double[] params) {
+ double[][] jacobian = new double[x.size()][3];
+ for (int i = 0; i < jacobian.length; ++i) {
+ jacobian[i][0] = x.get(i) * x.get(i);
+ jacobian[i][1] = x.get(i);
+ jacobian[i][2] = 1.0;
+ }
+ return jacobian;
+ }
+ };
+ }
+ }
+
+ private static class BevingtonProblem {
+ private List<Double> time;
+ private List<Double> count;
+
+ public BevingtonProblem() {
+ time = new ArrayList<Double>();
+ count = new ArrayList<Double>();
+ }
+
+ public void addPoint(double t, double c) {
+ time.add(t);
+ count.add(c);
+ }
+
+ public MultivariateVectorFunction getModelFunction() {
+ return new MultivariateVectorFunction() {
+ public double[] value(double[] params) {
+ double[] values = new double[time.size()];
+ for (int i = 0; i < values.length; ++i) {
+ final double t = time.get(i);
+ values[i] = params[0] +
+ params[1] * Math.exp(-t / params[3]) +
+ params[2] * Math.exp(-t / params[4]);
+ }
+ return values;
+ }
+ };
+ }
+
+ public MultivariateMatrixFunction getModelFunctionJacobian() {
+ return new MultivariateMatrixFunction() {
+ public double[][] value(double[] params) {
+ double[][] jacobian = new double[time.size()][5];
+
+ for (int i = 0; i < jacobian.length; ++i) {
+ final double t = time.get(i);
+ jacobian[i][0] = 1;
+
+ final double p3 = params[3];
+ final double p4 = params[4];
+ final double tOp3 = t / p3;
+ final double tOp4 = t / p4;
+ jacobian[i][1] = Math.exp(-tOp3);
+ jacobian[i][2] = Math.exp(-tOp4);
+ jacobian[i][3] = params[1] * Math.exp(-tOp3) * tOp3 / p3;
+ jacobian[i][4] = params[2] * Math.exp(-tOp4) * tOp4 / p4;
+ }
+ return jacobian;
+ }
+ };
+ }
+ }
+}
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/fitting/leastsquares/LevenbergMarquardtOptimizerTest.java
------------------------------------------------------------------------------
svn:keywords = Id Revision