This class solve a least-square problem by + * Gauss-Newton least-squares solver. + *
This class solve a least-square problem by * solving the normal equations of the linearized problem at each iteration. Either LU - * decomposition or QR decomposition can be used to solve the normal equations. LU - * decomposition is faster but QR decomposition is more robust for difficult problems. + * decomposition or Cholesky decomposition can be used to solve the normal equations, + * or QR decomposition or SVD decomposition can be used to solve the linear system. LU + * decomposition is faster but QR decomposition is more robust for difficult problems, + * and SVD can compute a solution for rank-deficient problems. *
* * @version $Id$ @@ -119,6 +123,22 @@ public class GaussNewtonOptimizer implements LeastSquaresOptimizer { throw new ConvergenceException(LocalizedFormats.UNABLE_TO_SOLVE_SINGULAR_PROBLEM, e); } } + }, + /** + * Solve the linear least squares problem using the {@link + * SingularValueDecomposition}. + * + *This method is slower, but can provide a solution for rank deficient and + * nearly singular systems. + */ + SVD { + @Override + protected RealVector solve(final RealMatrix jacobian, + final RealVector residuals) { + return new SingularValueDecomposition(jacobian) + .getSolver() + .solve(residuals); + } }; /** diff --git a/src/test/java/org/apache/commons/math3/fitting/leastsquares/AbstractLeastSquaresOptimizerAbstractTest.java b/src/test/java/org/apache/commons/math3/fitting/leastsquares/AbstractLeastSquaresOptimizerAbstractTest.java index 256925902..86c3ab37a 100644 --- a/src/test/java/org/apache/commons/math3/fitting/leastsquares/AbstractLeastSquaresOptimizerAbstractTest.java +++ b/src/test/java/org/apache/commons/math3/fitting/leastsquares/AbstractLeastSquaresOptimizerAbstractTest.java @@ -303,8 +303,12 @@ public abstract class AbstractLeastSquaresOptimizerAbstractTest { problem.getBuilder().start(new double[]{2, 2, 2, 2, 2, 2}).build()); Assert.assertEquals(0, optimum.getRMS(), TOl); - //TODO the first two elements of point were not previously checked - assertEquals(TOl, optimum.getPoint(), 2, 1, 3, 4, 5, 6); + RealVector point = optimum.getPoint(); + //the first two elements are under constrained + //check first two elements obey the constraint: sum to 3 + Assert.assertEquals(3, point.getEntry(0) + point.getEntry(1), TOl); + //#constrains = #states fro the last 4 elements + assertEquals(TOl, point.getSubVector(2, 4), 3, 4, 5, 6); } @Test diff --git a/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerWithSVDTest.java b/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerWithSVDTest.java new file mode 100644 index 000000000..a205d5398 --- /dev/null +++ b/src/test/java/org/apache/commons/math3/fitting/leastsquares/GaussNewtonOptimizerWithSVDTest.java @@ -0,0 +1,149 @@ +/* + * 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 org.apache.commons.math3.exception.ConvergenceException; +import org.apache.commons.math3.exception.TooManyEvaluationsException; +import org.apache.commons.math3.fitting.leastsquares.GaussNewtonOptimizer.Decomposition; +import org.apache.commons.math3.fitting.leastsquares.LeastSquaresOptimizer.Optimum; +import org.apache.commons.math3.geometry.euclidean.threed.Plane; +import org.apache.commons.math3.geometry.euclidean.threed.Vector3D; +import org.apache.commons.math3.optim.SimpleVectorValueChecker; +import org.apache.commons.math3.util.FastMath; +import org.junit.Assert; +import org.junit.Test; + +import java.io.IOException; + +/** + *
Some of the unit tests are re-implementations of the MINPACK file17 and file22 test files. + * The redistribution policy for MINPACK is available here/ + * + * @version $Id$ + */ +public class GaussNewtonOptimizerWithSVDTest + extends AbstractLeastSquaresOptimizerAbstractTest { + + @Override + public int getMaxIterations() { + return 1000; + } + + @Override + public LeastSquaresOptimizer getOptimizer() { + return new GaussNewtonOptimizer(Decomposition.SVD); + } + + @Test + public void testMaxEvaluations() throws Exception { + try{ + 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); + + LeastSquaresProblem lsp = builder(circle) + .checkerPair(new SimpleVectorValueChecker(1e-30, 1e-30)) + .maxIterations(Integer.MAX_VALUE) + .start(new double[]{98.680, 47.345}) + .build(); + + optimizer.optimize(lsp); + + fail(optimizer); + }catch (TooManyEvaluationsException e){ + //expected + } + } + + @Override + @Test + public void testCircleFittingBadInit() { + /* + * This test converged to the wrong solution with this optimizer. + * It seems that the state becomes so large that the convergence + * checker's relative tolerance test passes. + */ + try { + super.testCircleFittingBadInit(); + fail(optimizer); + } catch (AssertionError e) { + //expected + } + } + + @Override + @Test + public void testHahn1() + throws IOException { + /* + * TODO This test leads to a singular problem with the Gauss-Newton + * optimizer. This should be inquired. + */ + try{ + super.testHahn1(); + fail(optimizer); + } catch (ConvergenceException e){ + //expected for LU + } catch (TooManyEvaluationsException e){ + //expected for QR + } + } + + @Test + @Override + public void testGetIterations() { + /* this diverges with SVD */ + try { + super.testGetIterations(); + fail(optimizer); + } catch (TooManyEvaluationsException e) { + //expected + } + } + + @Test + @Override + public void testNonInvertible() throws Exception { + /* SVD can compute a solution to singular problems. + * In this case the target vector, b, is not in the + * span of the jacobian matrix, A. The closes point + * to b on the plane spanned by A is computed. + */ + LinearProblem problem = new LinearProblem(new double[][]{ + {1, 2, -3}, + {2, 1, 3}, + {-3, 0, -9} + }, new double[]{1, 1, 1}); + + Optimum optimum = optimizer.optimize(problem.getBuilder().build()); + + Plane span = new Plane(Vector3D.ZERO, new Vector3D(1, 2, -3), new Vector3D(2, 1, 0), TOl); + double expected = FastMath.abs(span.getOffset(new Vector3D(1, 1, 1))); + double actual = optimum.getResiduals().getNorm(); + + //verify + Assert.assertEquals(expected, actual, TOl); + } + +}