001 /*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements. See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License. You may obtain a copy of the License at
008 *
009 * http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017
018 package org.apache.commons.math3.optim.nonlinear.scalar;
019
020 import org.apache.commons.math3.analysis.MultivariateFunction;
021 import org.apache.commons.math3.analysis.MultivariateVectorFunction;
022 import org.apache.commons.math3.exception.DimensionMismatchException;
023 import org.apache.commons.math3.linear.RealMatrix;
024
025 /**
026 * This class converts
027 * {@link MultivariateVectorFunction vectorial objective functions} to
028 * {@link MultivariateFunction scalar objective functions}
029 * when the goal is to minimize them.
030 * <br/>
031 * This class is mostly used when the vectorial objective function represents
032 * a theoretical result computed from a point set applied to a model and
033 * the models point must be adjusted to fit the theoretical result to some
034 * reference observations. The observations may be obtained for example from
035 * physical measurements whether the model is built from theoretical
036 * considerations.
037 * <br/>
038 * This class computes a possibly weighted squared sum of the residuals, which is
039 * a scalar value. The residuals are the difference between the theoretical model
040 * (i.e. the output of the vectorial objective function) and the observations. The
041 * class implements the {@link MultivariateFunction} interface and can therefore be
042 * minimized by any optimizer supporting scalar objectives functions.This is one way
043 * to perform a least square estimation. There are other ways to do this without using
044 * this converter, as some optimization algorithms directly support vectorial objective
045 * functions.
046 * <br/>
047 * This class support combination of residuals with or without weights and correlations.
048 *
049 * @see MultivariateFunction
050 * @see MultivariateVectorFunction
051 * @version $Id: LeastSquaresConverter.java 1416643 2012-12-03 19:37:14Z tn $
052 * @since 2.0
053 */
054
055 public class LeastSquaresConverter implements MultivariateFunction {
056 /** Underlying vectorial function. */
057 private final MultivariateVectorFunction function;
058 /** Observations to be compared to objective function to compute residuals. */
059 private final double[] observations;
060 /** Optional weights for the residuals. */
061 private final double[] weights;
062 /** Optional scaling matrix (weight and correlations) for the residuals. */
063 private final RealMatrix scale;
064
065 /**
066 * Builds a simple converter for uncorrelated residuals with identical
067 * weights.
068 *
069 * @param function vectorial residuals function to wrap
070 * @param observations observations to be compared to objective function to compute residuals
071 */
072 public LeastSquaresConverter(final MultivariateVectorFunction function,
073 final double[] observations) {
074 this.function = function;
075 this.observations = observations.clone();
076 this.weights = null;
077 this.scale = null;
078 }
079
080 /**
081 * Builds a simple converter for uncorrelated residuals with the
082 * specified weights.
083 * <p>
084 * The scalar objective function value is computed as:
085 * <pre>
086 * objective = ∑weight<sub>i</sub>(observation<sub>i</sub>-objective<sub>i</sub>)<sup>2</sup>
087 * </pre>
088 * </p>
089 * <p>
090 * Weights can be used for example to combine residuals with different standard
091 * deviations. As an example, consider a residuals array in which even elements
092 * are angular measurements in degrees with a 0.01° standard deviation and
093 * odd elements are distance measurements in meters with a 15m standard deviation.
094 * In this case, the weights array should be initialized with value
095 * 1.0/(0.01<sup>2</sup>) in the even elements and 1.0/(15.0<sup>2</sup>) in the
096 * odd elements (i.e. reciprocals of variances).
097 * </p>
098 * <p>
099 * The array computed by the objective function, the observations array and the
100 * weights array must have consistent sizes or a {@link DimensionMismatchException}
101 * will be triggered while computing the scalar objective.
102 * </p>
103 *
104 * @param function vectorial residuals function to wrap
105 * @param observations observations to be compared to objective function to compute residuals
106 * @param weights weights to apply to the residuals
107 * @throws DimensionMismatchException if the observations vector and the weights
108 * vector dimensions do not match (objective function dimension is checked only when
109 * the {@link #value(double[])} method is called)
110 */
111 public LeastSquaresConverter(final MultivariateVectorFunction function,
112 final double[] observations,
113 final double[] weights) {
114 if (observations.length != weights.length) {
115 throw new DimensionMismatchException(observations.length, weights.length);
116 }
117 this.function = function;
118 this.observations = observations.clone();
119 this.weights = weights.clone();
120 this.scale = null;
121 }
122
123 /**
124 * Builds a simple converter for correlated residuals with the
125 * specified weights.
126 * <p>
127 * The scalar objective function value is computed as:
128 * <pre>
129 * objective = y<sup>T</sup>y with y = scale×(observation-objective)
130 * </pre>
131 * </p>
132 * <p>
133 * The array computed by the objective function, the observations array and the
134 * the scaling matrix must have consistent sizes or a {@link DimensionMismatchException}
135 * will be triggered while computing the scalar objective.
136 * </p>
137 *
138 * @param function vectorial residuals function to wrap
139 * @param observations observations to be compared to objective function to compute residuals
140 * @param scale scaling matrix
141 * @throws DimensionMismatchException if the observations vector and the scale
142 * matrix dimensions do not match (objective function dimension is checked only when
143 * the {@link #value(double[])} method is called)
144 */
145 public LeastSquaresConverter(final MultivariateVectorFunction function,
146 final double[] observations,
147 final RealMatrix scale) {
148 if (observations.length != scale.getColumnDimension()) {
149 throw new DimensionMismatchException(observations.length, scale.getColumnDimension());
150 }
151 this.function = function;
152 this.observations = observations.clone();
153 this.weights = null;
154 this.scale = scale.copy();
155 }
156
157 /** {@inheritDoc} */
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 }
167
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 }
184
185 return sumSquares;
186 }
187 }