MATH-1128

Introducing a "LazyUnweightedEvaluation": the computation of the model
and Jacobian are deferred until one or the other is actually accessed.
Class "LocalMultivariateJacobianFunction" replaces the anonymous class
that was created when calling the "model" method.
The "evaluate" method of "LocalLeastSquaresProblem" instantiates either
"UnweightedEvaluation" or "LazyUnweightedEvaluation", as requested by
the user of the "create" factory method.


git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1603219 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Gilles Sadowski 2014-06-17 16:48:12 +00:00
parent d4f978ddd5
commit 0714c7cbe6
3 changed files with 273 additions and 79 deletions

View File

@ -73,7 +73,15 @@ Users are encouraged to upgrade to this version as this release not
2. A few methods in the FastMath class are in fact slower that their 2. A few methods in the FastMath class are in fact slower that their
counterpart in either Math or StrictMath (cf. MATH-740 and MATH-901). counterpart in either Math or StrictMath (cf. MATH-740 and MATH-901).
"> ">
<action dev="luc" type="fix" issue="MATH-1127"> <action dev="erans" type="fix" issue="MATH-1129">
"Percentile": wrong sorting in the presence of NaN.
</action>
<action dev="erans" type="update" issue="MATH-1128">
Added lazy evaluation to "LeastSquaresFactory" (in "o.a.c.m.fitting.leastsquares")
to avoid evaluating the model when the optimization algorithm does not actually
require it.
</action>
<action dev="luc" type="fix" issue="MATH-1127">
Fixed overflow in Precision.equals with ulps (both double and float versions). Fixed overflow in Precision.equals with ulps (both double and float versions).
</action> </action>
<action dev="tn" type="fix" issue="MATH-1125" due-to="Ajo Fod"> <action dev="tn" type="fix" issue="MATH-1125" due-to="Ajo Fod">

View File

@ -41,12 +41,41 @@ import org.apache.commons.math3.util.Pair;
public class LeastSquaresFactory { public class LeastSquaresFactory {
/** Prevent instantiation. */ /** Prevent instantiation. */
private LeastSquaresFactory() { private LeastSquaresFactory() {}
/**
* Create a {@link org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem}
* from the given elements. There will be no weights applied (unit weights).
*
* @param model the model function. Produces the computed values.
* @param observed the observed (target) values
* @param start the initial guess.
* @param checker convergence checker
* @param maxEvaluations the maximum number of times to evaluate the model
* @param maxIterations the maximum number to times to iterate in the algorithm
* @param lazyEvaluation Whether the call to {@link Evaluation#evaluate(RealVector)}
* will defer the evaluation until access to the value is requested.
* @return the specified General Least Squares problem.
*/
public static LeastSquaresProblem create(final MultivariateJacobianFunction model,
final RealVector observed,
final RealVector start,
final ConvergenceChecker<Evaluation> checker,
final int maxEvaluations,
final int maxIterations,
final boolean lazyEvaluation) {
return new LocalLeastSquaresProblem(model,
observed,
start,
checker,
maxEvaluations,
maxIterations,
lazyEvaluation);
} }
/** /**
* Create a {@link org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem} * Create a {@link org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem}
* from the given elements. There will be no weights applied (Identity weights). * from the given elements. There will be no weights applied (unit weights).
* *
* @param model the model function. Produces the computed values. * @param model the model function. Produces the computed values.
* @param observed the observed (target) values * @param observed the observed (target) values
@ -62,14 +91,13 @@ public class LeastSquaresFactory {
final ConvergenceChecker<Evaluation> checker, final ConvergenceChecker<Evaluation> checker,
final int maxEvaluations, final int maxEvaluations,
final int maxIterations) { final int maxIterations) {
return new LocalLeastSquaresProblem( return new LocalLeastSquaresProblem(model,
model, observed,
observed, start,
start, checker,
checker, maxEvaluations,
maxEvaluations, maxIterations,
maxIterations false);
);
} }
/** /**
@ -92,22 +120,19 @@ public class LeastSquaresFactory {
final ConvergenceChecker<Evaluation> checker, final ConvergenceChecker<Evaluation> checker,
final int maxEvaluations, final int maxEvaluations,
final int maxIterations) { final int maxIterations) {
return weightMatrix( return weightMatrix(create(model,
create( observed,
model, start,
observed, checker,
start, maxEvaluations,
checker, maxIterations),
maxEvaluations, weight);
maxIterations
),
weight);
} }
/** /**
* Create a {@link org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem} * Create a {@link org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem}
* from the given elements. * from the given elements.
* <p/> * <p>
* This factory method is provided for continuity with previous interfaces. Newer * This factory method is provided for continuity with previous interfaces. Newer
* applications should use {@link #create(MultivariateJacobianFunction, RealVector, * applications should use {@link #create(MultivariateJacobianFunction, RealVector,
* RealVector, ConvergenceChecker, int, int)}, or {@link #create(MultivariateJacobianFunction, * RealVector, ConvergenceChecker, int, int)}, or {@link #create(MultivariateJacobianFunction,
@ -131,15 +156,13 @@ public class LeastSquaresFactory {
final ConvergenceChecker<Evaluation> checker, final ConvergenceChecker<Evaluation> checker,
final int maxEvaluations, final int maxEvaluations,
final int maxIterations) { final int maxIterations) {
return create( return create(model(model, jacobian),
model(model, jacobian), new ArrayRealVector(observed, false),
new ArrayRealVector(observed, false), new ArrayRealVector(start, false),
new ArrayRealVector(start, false), weight,
weight, checker,
checker, maxEvaluations,
maxEvaluations, maxIterations);
maxIterations
);
} }
/** /**
@ -171,7 +194,7 @@ public class LeastSquaresFactory {
*/ */
public static LeastSquaresProblem weightDiagonal(final LeastSquaresProblem problem, public static LeastSquaresProblem weightDiagonal(final LeastSquaresProblem problem,
final RealVector weights) { final RealVector weights) {
//TODO more efficient implementation // TODO more efficient implementation
return weightMatrix(problem, new DiagonalMatrix(weights.toArray())); return weightMatrix(problem, new DiagonalMatrix(weights.toArray()));
} }
@ -193,8 +216,7 @@ public class LeastSquaresFactory {
return super.evaluate(point); return super.evaluate(point);
} }
/* delegate the rest */ // Delegate the rest.
}; };
} }
@ -205,9 +227,7 @@ public class LeastSquaresFactory {
* @param checker the convergence checker to adapt. * @param checker the convergence checker to adapt.
* @return a convergence checker that delegates to {@code checker}. * @return a convergence checker that delegates to {@code checker}.
*/ */
public static ConvergenceChecker<Evaluation> evaluationChecker( public static ConvergenceChecker<Evaluation> evaluationChecker(final ConvergenceChecker<PointVectorValuePair> checker) {
final ConvergenceChecker<PointVectorValuePair> checker
) {
return new ConvergenceChecker<Evaluation>() { return new ConvergenceChecker<Evaluation>() {
public boolean converged(final int iteration, public boolean converged(final int iteration,
final Evaluation previous, final Evaluation previous,
@ -255,22 +275,69 @@ public class LeastSquaresFactory {
* @param jacobian the Jacobian function * @param jacobian the Jacobian function
* @return a function that computes both at the same time * @return a function that computes both at the same time
*/ */
public static MultivariateJacobianFunction model( public static MultivariateJacobianFunction model(final MultivariateVectorFunction value,
final MultivariateVectorFunction value, final MultivariateMatrixFunction jacobian) {
final MultivariateMatrixFunction jacobian return new LocalMultivariateJacobianFunction(value, jacobian);
) {
return new MultivariateJacobianFunction() {
public Pair<RealVector, RealMatrix> value(final RealVector point) {
//TODO get array from RealVector without copying?
final double[] pointArray = point.toArray();
//evaluate and return data without copying
return new Pair<RealVector, RealMatrix>(
new ArrayRealVector(value.value(pointArray), false),
new Array2DRowRealMatrix(jacobian.value(pointArray), false));
}
};
} }
/**
* Combine a {@link MultivariateVectorFunction} with a {@link
* MultivariateMatrixFunction} to produce a {@link MultivariateJacobianFunction}.
*
* @param value the vector value function
* @param jacobian the Jacobian function
* @return a function that computes both at the same time
*/
private static class LocalMultivariateJacobianFunction
implements MultivariateJacobianFunction {
/** Model. */
private final MultivariateVectorFunction value;
/** Model's Jacobian. */
private final MultivariateMatrixFunction jacobian;
/**
* @param value Model function.
* @param jacobian Model's Jacobian function.
*/
LocalMultivariateJacobianFunction(final MultivariateVectorFunction value,
final MultivariateMatrixFunction jacobian) {
this.value = value;
this.jacobian = jacobian;
}
/** {@inheritDoc} */
public Pair<RealVector, RealMatrix> value(final RealVector point) {
//TODO get array from RealVector without copying?
final double[] p = point.toArray();
// Evaluate.
return new Pair<RealVector, RealMatrix>(computeValue(p),
computeJacobian(p));
}
/**
* Compute the value.
*
* @param params Point.
* @return the Jacobian at the given point.
*/
public RealVector computeValue(final double[] params) {
return new ArrayRealVector(value.value(params), false);
}
/**
* Compute the Jacobian.
*
* @param params Point.
* @return the Jacobian at the given point.
*/
public RealMatrix computeJacobian(final double[] params) {
return new Array2DRowRealMatrix(jacobian.value(params), false);
}
}
/** /**
* A private, "field" immutable (not "real" immutable) implementation of {@link * A private, "field" immutable (not "real" immutable) implementation of {@link
* LeastSquaresProblem}. * LeastSquaresProblem}.
@ -281,11 +348,13 @@ public class LeastSquaresFactory {
implements LeastSquaresProblem { implements LeastSquaresProblem {
/** Target values for the model function at optimum. */ /** Target values for the model function at optimum. */
private RealVector target; private final RealVector target;
/** Model function. */ /** Model function. */
private MultivariateJacobianFunction model; private final MultivariateJacobianFunction model;
/** Initial guess. */ /** Initial guess. */
private RealVector start; private final RealVector start;
/** Whether to use lazy evaluation. */
private final boolean lazyEvaluation;
/** /**
* Create a {@link LeastSquaresProblem} from the given data. * Create a {@link LeastSquaresProblem} from the given data.
@ -296,17 +365,21 @@ public class LeastSquaresFactory {
* @param checker the convergence checker * @param checker the convergence checker
* @param maxEvaluations the allowed evaluations * @param maxEvaluations the allowed evaluations
* @param maxIterations the allowed iterations * @param maxIterations the allowed iterations
* @param lazyEvaluation Whether the call to {@link Evaluation#evaluate(RealVector)}
* will defer the evaluation until access to the value is requested.
*/ */
LocalLeastSquaresProblem(final MultivariateJacobianFunction model, LocalLeastSquaresProblem(final MultivariateJacobianFunction model,
final RealVector target, final RealVector target,
final RealVector start, final RealVector start,
final ConvergenceChecker<Evaluation> checker, final ConvergenceChecker<Evaluation> checker,
final int maxEvaluations, final int maxEvaluations,
final int maxIterations) { final int maxIterations,
boolean lazyEvaluation) {
super(maxEvaluations, maxIterations, checker); super(maxEvaluations, maxIterations, checker);
this.target = target; this.target = target;
this.model = model; this.model = model;
this.start = start; this.start = start;
this.lazyEvaluation = lazyEvaluation;
} }
/** {@inheritDoc} */ /** {@inheritDoc} */
@ -326,28 +399,32 @@ public class LeastSquaresFactory {
/** {@inheritDoc} */ /** {@inheritDoc} */
public Evaluation evaluate(final RealVector point) { public Evaluation evaluate(final RealVector point) {
//evaluate value and jacobian in one function call // Copy so optimizer can change point without changing our instance.
final Pair<RealVector, RealMatrix> value = this.model.value(point); final RealVector p = point.copy();
return new UnweightedEvaluation(
value.getFirst(), if (lazyEvaluation) {
value.getSecond(), return new LazyUnweightedEvaluation(model,
this.target, target,
// copy so optimizer can change point without changing our instance p);
point.copy()); } else {
// Evaluate value and jacobian in one function call.
final Pair<RealVector, RealMatrix> value = model.value(p);
return new UnweightedEvaluation(value.getFirst(),
value.getSecond(),
target,
p);
}
} }
/** /**
* Container with the model evaluation at a particular point. * Container with the model evaluation at a particular point.
* <p/>
* TODO revisit lazy evaluation
*/ */
private static class UnweightedEvaluation extends AbstractEvaluation { private static class UnweightedEvaluation extends AbstractEvaluation {
/** Point of evaluation. */
/** the point of evaluation */
private final RealVector point; private final RealVector point;
/** deriviative at point */ /** Derivative at point. */
private final RealMatrix jacobian; private final RealMatrix jacobian;
/** the computed residuals. */ /** Computed residuals. */
private final RealVector residuals; private final RealVector residuals;
/** /**
@ -370,22 +447,63 @@ public class LeastSquaresFactory {
/** {@inheritDoc} */ /** {@inheritDoc} */
public RealMatrix getJacobian() { public RealMatrix getJacobian() {
return this.jacobian; return jacobian;
} }
/** {@inheritDoc} */ /** {@inheritDoc} */
public RealVector getPoint() { public RealVector getPoint() {
return this.point; return point;
} }
/** {@inheritDoc} */ /** {@inheritDoc} */
public RealVector getResiduals() { public RealVector getResiduals() {
return this.residuals; return residuals;
} }
} }
} /**
* Container with the model <em>lazy</em> evaluation at a particular point.
*/
private static class LazyUnweightedEvaluation extends AbstractEvaluation {
/** Point of evaluation. */
private final RealVector point;
/** Model and Jacobian functions. */
private final LocalMultivariateJacobianFunction model;
/** Target values for the model function at optimum. */
private final RealVector target;
/**
* Create an {@link Evaluation} with no weights.
*
* @param model the model function
* @param target the observed values
* @param point the abscissa
*/
private LazyUnweightedEvaluation(final MultivariateJacobianFunction model,
final RealVector target,
final RealVector point) {
super(target.getDimension());
// Safe to cast as long as we control usage of this class.
this.model = (LocalMultivariateJacobianFunction) model;
this.point = point;
this.target = target;
}
/** {@inheritDoc} */
public RealMatrix getJacobian() {
return model.computeJacobian(point.toArray());
}
/** {@inheritDoc} */
public RealVector getPoint() {
return point;
}
/** {@inheritDoc} */
public RealVector getResiduals() {
return target.subtract(model.computeValue(point.toArray()));
}
}
}
} }

View File

@ -14,6 +14,8 @@
package org.apache.commons.math3.fitting.leastsquares; package org.apache.commons.math3.fitting.leastsquares;
import org.apache.commons.math3.TestUtils; import org.apache.commons.math3.TestUtils;
import org.apache.commons.math3.analysis.MultivariateMatrixFunction;
import org.apache.commons.math3.analysis.MultivariateVectorFunction;
import org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem.Evaluation; import org.apache.commons.math3.fitting.leastsquares.LeastSquaresProblem.Evaluation;
import org.apache.commons.math3.linear.ArrayRealVector; import org.apache.commons.math3.linear.ArrayRealVector;
import org.apache.commons.math3.linear.DiagonalMatrix; import org.apache.commons.math3.linear.DiagonalMatrix;
@ -217,4 +219,70 @@ public class EvaluationTest {
Assert.assertEquals(evaluation.getPoint().getEntry(0), 0, 0); Assert.assertEquals(evaluation.getPoint().getEntry(0), 0, 0);
} }
@Test
public void testLazyEvaluation() {
final RealVector dummy = new ArrayRealVector(new double[] { 0 });
final LeastSquaresProblem p
= LeastSquaresFactory.create(LeastSquaresFactory.model(dummyModel(), dummyJacobian()),
dummy, dummy, null, 0, 0, true);
// Should not throw because actual evaluation is deferred.
final Evaluation eval = p.evaluate(dummy);
try {
eval.getResiduals();
Assert.fail("Exception expected");
} catch (RuntimeException e) {
// Expecting exception.
Assert.assertEquals("dummyModel", e.getMessage());
}
try {
eval.getJacobian();
Assert.fail("Exception expected");
} catch (RuntimeException e) {
// Expecting exception.
Assert.assertEquals("dummyJacobian", e.getMessage());
}
}
@Test
public void testDirectEvaluation() {
final RealVector dummy = new ArrayRealVector(new double[] { 0 });
final LeastSquaresProblem p
= LeastSquaresFactory.create(LeastSquaresFactory.model(dummyModel(), dummyJacobian()),
dummy, dummy, null, 0, 0, false);
try {
// Should throw.
final Evaluation eval = p.evaluate(dummy);
Assert.fail("Exception expected");
} catch (RuntimeException e) {
// Expecting exception.
// Whether it is model of Jacobian that caused it is not significant.
final String msg = e.getMessage();
Assert.assertTrue(msg.equals("dummyModel") ||
msg.equals("dummyJacobian"));
}
}
/** Used for testing direct vs lazy evaluation. */
private MultivariateVectorFunction dummyModel() {
return new MultivariateVectorFunction() {
public double[] value(double[] p) {
throw new RuntimeException("dummyModel");
}
};
}
/** Used for testing direct vs lazy evaluation. */
private MultivariateMatrixFunction dummyJacobian() {
return new MultivariateMatrixFunction() {
public double[][] value(double[] p) {
throw new RuntimeException("dummyJacobian");
}
};
}
} }