LeastSquaresConverter.java

  1. /*
  2.  * Licensed to the Apache Software Foundation (ASF) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * The ASF licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *      http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */

  18. package org.apache.commons.math4.legacy.optim.nonlinear.scalar;

  20. import org.apache.commons.math4.legacy.analysis.MultivariateFunction;
  21. import org.apache.commons.math4.legacy.analysis.MultivariateVectorFunction;
  22. import org.apache.commons.math4.legacy.exception.DimensionMismatchException;
  23. import org.apache.commons.math4.legacy.linear.RealMatrix;

  25. /**
  26.  * This class converts
  27.  * {@link MultivariateVectorFunction vectorial objective functions} to
  28.  * {@link MultivariateFunction scalar objective functions}
  29.  * when the goal is to minimize them.
  30.  * <br>
  31.  * This class is mostly used when the vectorial objective function represents
  32.  * a theoretical result computed from a point set applied to a model and
  33.  * the models point must be adjusted to fit the theoretical result to some
  34.  * reference observations. The observations may be obtained for example from
  35.  * physical measurements whether the model is built from theoretical
  36.  * considerations.
  37.  * <br>
  38.  * This class computes a possibly weighted squared sum of the residuals, which is
  39.  * a scalar value. The residuals are the difference between the theoretical model
  40.  * (i.e. the output of the vectorial objective function) and the observations. The
  41.  * class implements the {@link MultivariateFunction} interface and can therefore be
  42.  * minimized by any optimizer supporting scalar objectives functions.This is one way
  43.  * to perform a least square estimation. There are other ways to do this without using
  44.  * this converter, as some optimization algorithms directly support vectorial objective
  45.  * functions.
  46.  * <br>
  47.  * This class support combination of residuals with or without weights and correlations.
  48.   *
  49.  * @see MultivariateFunction
  50.  * @see MultivariateVectorFunction
  51.  * @since 2.0
  52.  */

  54. public class LeastSquaresConverter implements MultivariateFunction {
  55.     /** Underlying vectorial function. */
  56.     private final MultivariateVectorFunction function;
  57.     /** Observations to be compared to objective function to compute residuals. */
  58.     private final double[] observations;
  59.     /** Optional weights for the residuals. */
  60.     private final double[] weights;
  61.     /** Optional scaling matrix (weight and correlations) for the residuals. */
  62.     private final RealMatrix scale;

  64.     /**
  65.      * Builds a simple converter for uncorrelated residuals with identical
  66.      * weights.
  67.      *
  68.      * @param function vectorial residuals function to wrap
  69.      * @param observations observations to be compared to objective function to compute residuals
  70.      */
  71.     public LeastSquaresConverter(final MultivariateVectorFunction function,
  72.                                  final double[] observations) {
  73.         this.function     = function;
  74.         this.observations = observations.clone();
  75.         this.weights      = null;
  76.         this.scale        = null;
  77.     }

  79.     /**
  80.      * Builds a simple converter for uncorrelated residuals with the
  81.      * specified weights.
  82.      * <p>
  83.      * The scalar objective function value is computed as:
  84.      * <div style="white-space: pre"><code>
  85.      * objective = &sum;weight<sub>i</sub>(observation<sub>i</sub>-objective<sub>i</sub>)<sup>2</sup>
  86.      * </code></div>
  87.      *
  88.      * <p>
  89.      * Weights can be used for example to combine residuals with different standard
  90.      * deviations. As an example, consider a residuals array in which even elements
  91.      * are angular measurements in degrees with a 0.01&deg; standard deviation and
  92.      * odd elements are distance measurements in meters with a 15m standard deviation.
  93.      * In this case, the weights array should be initialized with value
  94.      * 1.0/(0.01<sup>2</sup>) in the even elements and 1.0/(15.0<sup>2</sup>) in the
  95.      * odd elements (i.e. reciprocals of variances).
  96.      * </p>
  97.      * <p>
  98.      * The array computed by the objective function, the observations array and the
  99.      * weights array must have consistent sizes or a {@link DimensionMismatchException}
  100.      * will be triggered while computing the scalar objective.
  101.      * </p>
  102.      *
  103.      * @param function vectorial residuals function to wrap
  104.      * @param observations observations to be compared to objective function to compute residuals
  105.      * @param weights weights to apply to the residuals
  106.      * @throws DimensionMismatchException if the observations vector and the weights
  107.      * vector dimensions do not match (objective function dimension is checked only when
  108.      * the {@link #value(double[])} method is called)
  109.      */
  110.     public LeastSquaresConverter(final MultivariateVectorFunction function,
  111.                                  final double[] observations,
  112.                                  final double[] weights) {
  113.         if (observations.length != weights.length) {
  114.             throw new DimensionMismatchException(observations.length, weights.length);
  115.         }
  116.         this.function     = function;
  117.         this.observations = observations.clone();
  118.         this.weights      = weights.clone();
  119.         this.scale        = null;
  120.     }

  122.     /**
  123.      * Builds a simple converter for correlated residuals with the
  124.      * specified weights.
  125.      * <p>
  126.      * The scalar objective function value is computed as:
  127.      * <div style="white-space: pre"><code>
  128.      * objective = y<sup>T</sup>y with y = scale&times;(observation-objective)
  129.      * </code></div>
  130.      *
  131.      * <p>
  132.      * The array computed by the objective function, the observations array and the
  133.      * the scaling matrix must have consistent sizes or a {@link DimensionMismatchException}
  134.      * will be triggered while computing the scalar objective.
  135.      * </p>
  136.      *
  137.      * @param function vectorial residuals function to wrap
  138.      * @param observations observations to be compared to objective function to compute residuals
  139.      * @param scale scaling matrix
  140.      * @throws DimensionMismatchException if the observations vector and the scale
  141.      * matrix dimensions do not match (objective function dimension is checked only when
  142.      * the {@link #value(double[])} method is called)
  143.      */
  144.     public LeastSquaresConverter(final MultivariateVectorFunction function,
  145.                                  final double[] observations,
  146.                                  final RealMatrix scale) {
  147.         if (observations.length != scale.getColumnDimension()) {
  148.             throw new DimensionMismatchException(observations.length, scale.getColumnDimension());
  149.         }
  150.         this.function     = function;
  151.         this.observations = observations.clone();
  152.         this.weights      = null;
  153.         this.scale        = scale.copy();
  154.     }

  156.     /** {@inheritDoc} */
  157.     @Override
  158.     public double value(final double[] point) {
  159.         // compute residuals
  160.         final double[] residuals = function.value(point);
  161.         if (residuals.length != observations.length) {
  162.             throw new DimensionMismatchException(residuals.length, observations.length);
  163.         }
  164.         for (int i = 0; i < residuals.length; ++i) {
  165.             residuals[i] -= observations[i];
  166.         }

  168.         // compute sum of squares
  169.         double sumSquares = 0;
  170.         if (weights != null) {
  171.             for (int i = 0; i < residuals.length; ++i) {
  172.                 final double ri = residuals[i];
  173.                 sumSquares +=  weights[i] * ri * ri;
  174.             }
  175.         } else if (scale != null) {
  176.             for (final double yi : scale.operate(residuals)) {
  177.                 sumSquares += yi * yi;
  178.             }
  179.         } else {
  180.             for (final double ri : residuals) {
  181.                 sumSquares += ri * ri;
  182.             }
  183.         }

  185.         return sumSquares;
  186.     }
  187. }
Created with JaCoCo 0.8.10.202304240956Code Coverage Report for Apache Commons Math 4.0-SNAPSHOT