src/main/java/org/apache/commons/math/optimization/LeastSquaresConverter.java - platform/external/apache-commons-math - Git at Google

 /*
  * 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.math.optimization;

 import org.apache.commons.math.FunctionEvaluationException;
 import org.apache.commons.math.MathRuntimeException;
 import org.apache.commons.math.analysis.MultivariateRealFunction;
 import org.apache.commons.math.analysis.MultivariateVectorialFunction;
 import org.apache.commons.math.exception.util.LocalizedFormats;
 import org.apache.commons.math.linear.RealMatrix;

 /** This class converts {@link MultivariateVectorialFunction vectorial
  * objective functions} to {@link MultivariateRealFunction scalar objective functions}
  * when the goal is to minimize them.
  * <p>
  * This class is mostly used when the vectorial objective function represents
  * a theoretical result computed from a point set applied to a model and
  * the models point must be adjusted to fit the theoretical result to some
  * reference observations. The observations may be obtained for example from
  * physical measurements whether the model is built from theoretical
  * considerations.
  * </p>
  * <p>
  * This class computes a possibly weighted squared sum of the residuals, which is
  * a scalar value. The residuals are the difference between the theoretical model
  * (i.e. the output of the vectorial objective function) and the observations. The
  * class implements the {@link MultivariateRealFunction} interface and can therefore be
  * minimized by any optimizer supporting scalar objectives functions.This is one way
  * to perform a least square estimation. There are other ways to do this without using
  * this converter, as some optimization algorithms directly support vectorial objective
  * functions.
  * </p>
  * <p>
  * This class support combination of residuals with or without weights and correlations.
  * </p>
   *
  * @see MultivariateRealFunction
  * @see MultivariateVectorialFunction
  * @version $Revision: 1070725 $ $Date: 2011-02-15 02:31:12 +0100 (mar. 15 févr. 2011) $
  * @since 2.0
  */

 public class LeastSquaresConverter implements MultivariateRealFunction {

     /** Underlying vectorial function. */
     private final MultivariateVectorialFunction function;

     /** Observations to be compared to objective function to compute residuals. */
     private final double[] observations;

     /** Optional weights for the residuals. */
     private final double[] weights;

     /** Optional scaling matrix (weight and correlations) for the residuals. */
     private final RealMatrix scale;

     /** Build a simple converter for uncorrelated residuals with the same weight.
      * @param function vectorial residuals function to wrap
      * @param observations observations to be compared to objective function to compute residuals
      */
     public LeastSquaresConverter(final MultivariateVectorialFunction function,
                                  final double[] observations) {
         this.function     = function;
         this.observations = observations.clone();
         this.weights      = null;
         this.scale        = null;
     }

     /** Build a simple converter for uncorrelated residuals with the specific weights.
      * <p>
      * The scalar objective function value is computed as:
      * <pre>
      * objective = &sum;weight<sub>i</sub>(observation<sub>i</sub>-objective<sub>i</sub>)<sup>2</sup>
      * </pre>
      * </p>
      * <p>
      * Weights can be used for example to combine residuals with different standard
      * deviations. As an example, consider a residuals array in which even elements
      * are angular measurements in degrees with a 0.01&deg; standard deviation and
      * odd elements are distance measurements in meters with a 15m standard deviation.
      * In this case, the weights array should be initialized with value
      * 1.0/(0.01<sup>2</sup>) in the even elements and 1.0/(15.0<sup>2</sup>) in the
      * odd elements (i.e. reciprocals of variances).
      * </p>
      * <p>
      * The array computed by the objective function, the observations array and the
      * weights array must have consistent sizes or a {@link FunctionEvaluationException} will be
      * triggered while computing the scalar objective.
      * </p>
      * @param function vectorial residuals function to wrap
      * @param observations observations to be compared to objective function to compute residuals
      * @param weights weights to apply to the residuals
      * @exception IllegalArgumentException if the observations vector and the weights
      * vector dimensions don't match (objective function dimension is checked only when
      * the {@link #value(double[])} method is called)
      */
     public LeastSquaresConverter(final MultivariateVectorialFunction function,
                                  final double[] observations, final double[] weights)
         throws IllegalArgumentException {
         if (observations.length != weights.length) {
             throw MathRuntimeException.createIllegalArgumentException(
                     LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
                     observations.length, weights.length);
         }
         this.function     = function;
         this.observations = observations.clone();
         this.weights      = weights.clone();
         this.scale        = null;
     }

     /** Build a simple converter for correlated residuals with the specific weights.
      * <p>
      * The scalar objective function value is computed as:
      * <pre>
      * objective = y<sup>T</sup>y with y = scale&times;(observation-objective)
      * </pre>
      * </p>
      * <p>
      * The array computed by the objective function, the observations array and the
      * the scaling matrix must have consistent sizes or a {@link FunctionEvaluationException}
      * will be triggered while computing the scalar objective.
      * </p>
      * @param function vectorial residuals function to wrap
      * @param observations observations to be compared to objective function to compute residuals
      * @param scale scaling matrix
      * @exception IllegalArgumentException if the observations vector and the scale
      * matrix dimensions don't match (objective function dimension is checked only when
      * the {@link #value(double[])} method is called)
      */
     public LeastSquaresConverter(final MultivariateVectorialFunction function,
                                  final double[] observations, final RealMatrix scale)
         throws IllegalArgumentException {
         if (observations.length != scale.getColumnDimension()) {
             throw MathRuntimeException.createIllegalArgumentException(
                     LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
                     observations.length, scale.getColumnDimension());
         }
         this.function     = function;
         this.observations = observations.clone();
         this.weights      = null;
         this.scale        = scale.copy();
     }

     /** {@inheritDoc} */
     public double value(final double[] point) throws FunctionEvaluationException {

         // compute residuals
         final double[] residuals = function.value(point);
         if (residuals.length != observations.length) {
             throw new FunctionEvaluationException(point,LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
                                         residuals.length, observations.length);
         }
         for (int i = 0; i < residuals.length; ++i) {
             residuals[i] -= observations[i];
         }

         // compute sum of squares
         double sumSquares = 0;
         if (weights != null) {
             for (int i = 0; i < residuals.length; ++i) {
                 final double ri = residuals[i];
                 sumSquares +=  weights[i] * ri * ri;
             }
         } else if (scale != null) {
             for (final double yi : scale.operate(residuals)) {
                 sumSquares += yi * yi;
             }
         } else {
             for (final double ri : residuals) {
                 sumSquares += ri * ri;
             }
         }

         return sumSquares;

     }

 }
	/*
	* 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.math.optimization;

	import org.apache.commons.math.FunctionEvaluationException;
	import org.apache.commons.math.MathRuntimeException;
	import org.apache.commons.math.analysis.MultivariateRealFunction;
	import org.apache.commons.math.analysis.MultivariateVectorialFunction;
	import org.apache.commons.math.exception.util.LocalizedFormats;
	import org.apache.commons.math.linear.RealMatrix;

	/** This class converts {@link MultivariateVectorialFunction vectorial
	* objective functions} to {@link MultivariateRealFunction scalar objective functions}
	* when the goal is to minimize them.
	* <p>
	* This class is mostly used when the vectorial objective function represents
	* a theoretical result computed from a point set applied to a model and
	* the models point must be adjusted to fit the theoretical result to some
	* reference observations. The observations may be obtained for example from
	* physical measurements whether the model is built from theoretical
	* considerations.
	* </p>
	* <p>
	* This class computes a possibly weighted squared sum of the residuals, which is
	* a scalar value. The residuals are the difference between the theoretical model
	* (i.e. the output of the vectorial objective function) and the observations. The
	* class implements the {@link MultivariateRealFunction} interface and can therefore be
	* minimized by any optimizer supporting scalar objectives functions.This is one way
	* to perform a least square estimation. There are other ways to do this without using
	* this converter, as some optimization algorithms directly support vectorial objective
	* functions.
	* </p>
	* <p>
	* This class support combination of residuals with or without weights and correlations.
	* </p>
	*
	* @see MultivariateRealFunction
	* @see MultivariateVectorialFunction
	* @version $Revision: 1070725 $ $Date: 2011-02-15 02:31:12 +0100 (mar. 15 févr. 2011) $
	* @since 2.0
	*/

	public class LeastSquaresConverter implements MultivariateRealFunction {

	/** Underlying vectorial function. */
	private final MultivariateVectorialFunction function;

	/** Observations to be compared to objective function to compute residuals. */
	private final double[] observations;

	/** Optional weights for the residuals. */
	private final double[] weights;

	/** Optional scaling matrix (weight and correlations) for the residuals. */
	private final RealMatrix scale;

	/** Build a simple converter for uncorrelated residuals with the same weight.
	* @param function vectorial residuals function to wrap
	* @param observations observations to be compared to objective function to compute residuals
	*/
	public LeastSquaresConverter(final MultivariateVectorialFunction function,
	final double[] observations) {
	this.function = function;
	this.observations = observations.clone();
	this.weights = null;
	this.scale = null;
	}

	/** Build a simple converter for uncorrelated residuals with the specific weights.
	* <p>
	* The scalar objective function value is computed as:
	* <pre>
	* objective = ∑weight<sub>i</sub>(observation<sub>i</sub>-objective<sub>i</sub>)<sup>2</sup>
	* </pre>
	* </p>
	* <p>
	* Weights can be used for example to combine residuals with different standard
	* deviations. As an example, consider a residuals array in which even elements
	* are angular measurements in degrees with a 0.01° standard deviation and
	* odd elements are distance measurements in meters with a 15m standard deviation.
	* In this case, the weights array should be initialized with value
	* 1.0/(0.01<sup>2</sup>) in the even elements and 1.0/(15.0<sup>2</sup>) in the
	* odd elements (i.e. reciprocals of variances).
	* </p>
	* <p>
	* The array computed by the objective function, the observations array and the
	* weights array must have consistent sizes or a {@link FunctionEvaluationException} will be
	* triggered while computing the scalar objective.
	* </p>
	* @param function vectorial residuals function to wrap
	* @param observations observations to be compared to objective function to compute residuals
	* @param weights weights to apply to the residuals
	* @exception IllegalArgumentException if the observations vector and the weights
	* vector dimensions don't match (objective function dimension is checked only when
	* the {@link #value(double[])} method is called)
	*/
	public LeastSquaresConverter(final MultivariateVectorialFunction function,
	final double[] observations, final double[] weights)
	throws IllegalArgumentException {
	if (observations.length != weights.length) {
	throw MathRuntimeException.createIllegalArgumentException(
	LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
	observations.length, weights.length);
	}
	this.function = function;
	this.observations = observations.clone();
	this.weights = weights.clone();
	this.scale = null;
	}

	/** Build a simple converter for correlated residuals with the specific weights.
	* <p>
	* The scalar objective function value is computed as:
	* <pre>
	* objective = y<sup>T</sup>y with y = scale×(observation-objective)
	* </pre>
	* </p>
	* <p>
	* The array computed by the objective function, the observations array and the
	* the scaling matrix must have consistent sizes or a {@link FunctionEvaluationException}
	* will be triggered while computing the scalar objective.
	* </p>
	* @param function vectorial residuals function to wrap
	* @param observations observations to be compared to objective function to compute residuals
	* @param scale scaling matrix
	* @exception IllegalArgumentException if the observations vector and the scale
	* matrix dimensions don't match (objective function dimension is checked only when
	* the {@link #value(double[])} method is called)
	*/
	public LeastSquaresConverter(final MultivariateVectorialFunction function,
	final double[] observations, final RealMatrix scale)
	throws IllegalArgumentException {
	if (observations.length != scale.getColumnDimension()) {
	throw MathRuntimeException.createIllegalArgumentException(
	LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
	observations.length, scale.getColumnDimension());
	}
	this.function = function;
	this.observations = observations.clone();
	this.weights = null;
	this.scale = scale.copy();
	}

	/** {@inheritDoc} */
	public double value(final double[] point) throws FunctionEvaluationException {

	// compute residuals
	final double[] residuals = function.value(point);
	if (residuals.length != observations.length) {
	throw new FunctionEvaluationException(point,LocalizedFormats.DIMENSIONS_MISMATCH_SIMPLE,
	residuals.length, observations.length);
	}
	for (int i = 0; i < residuals.length; ++i) {
	residuals[i] -= observations[i];
	}

	// compute sum of squares
	double sumSquares = 0;
	if (weights != null) {
	for (int i = 0; i < residuals.length; ++i) {
	final double ri = residuals[i];
	sumSquares += weights[i] * ri * ri;
	}
	} else if (scale != null) {
	for (final double yi : scale.operate(residuals)) {
	sumSquares += yi * yi;
	}
	} else {
	for (final double ri : residuals) {
	sumSquares += ri * ri;
	}
	}

	return sumSquares;

	}

	}