completely rewrote estimation package documentation
with downloadable example and explanation diagrams git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@618726 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Luc Maisonobe
parent
1d7fd3318b
commit
225e5c703c
|
|
@ -0,0 +1,187 @@
|
|||
/*
|
||||
* 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.
|
||||
*/
|
||||
|
||||
import org.apache.commons.math.estimation.EstimationException;
|
||||
import org.apache.commons.math.estimation.EstimatedParameter;
|
||||
import org.apache.commons.math.estimation.EstimationProblem;
|
||||
import org.apache.commons.math.estimation.LevenbergMarquardtEstimator;
|
||||
import org.apache.commons.math.estimation.SimpleEstimationProblem;
|
||||
import org.apache.commons.math.estimation.WeightedMeasurement;
|
||||
|
||||
public class TrajectoryDeterminationProblem extends SimpleEstimationProblem {
|
||||
|
||||
public static void main(String[] args) {
|
||||
try {
|
||||
TrajectoryDeterminationProblem problem =
|
||||
new TrajectoryDeterminationProblem(0.0, 100.0, 800.0, 1.0, 0.0);
|
||||
|
||||
double[][] distances = {
|
||||
{ 0.0, 806.5849 }, { 20.0, 796.8148 }, { 40.0, 791.0833 }, { 60.0, 789.6712 },
|
||||
{ 80.0, 793.1334 }, { 100.0, 797.7248 }, { 120.0, 803.2785 }, { 140.0, 813.4939 },
|
||||
{ 160.0, 826.9295 }, { 180.0, 844.0640 }, { 200.0, 863.3829 }, { 220.0, 883.3143 },
|
||||
{ 240.0, 908.6867 }, { 260.0, 934.8561 }, { 280.0, 964.0730 }, { 300.0, 992.1033 },
|
||||
{ 320.0, 1023.998 }, { 340.0, 1057.439 }, { 360.0, 1091.912 }, { 380.0, 1125.968 },
|
||||
{ 400.0, 1162.789 }, { 420.0, 1201.517 }, { 440.0, 1239.176 }, { 460.0, 1279.347 } };
|
||||
for (int i = 0; i < distances.length; ++i) {
|
||||
problem.addDistanceMeasurement(1.0, distances[i][0], distances[i][1]);
|
||||
};
|
||||
|
||||
double[][] angles = {
|
||||
{ 10.0, 1.415423 }, { 30.0, 1.352643 }, { 50.0, 1.289290 }, { 70.0, 1.225249 },
|
||||
{ 90.0, 1.161203 }, {110.0, 1.098538 }, {130.0, 1.036263 }, {150.0, 0.976052 },
|
||||
{170.0, 0.917921 }, {190.0, 0.861830 }, {210.0, 0.808237 }, {230.0, 0.757043 },
|
||||
{250.0, 0.708650 }, {270.0, 0.662949 }, {290.0, 0.619903 }, {310.0, 0.579160 },
|
||||
{330.0, 0.541033 }, {350.0, 0.505590 }, {370.0, 0.471746 }, {390.0, 0.440155 },
|
||||
{410.0, 0.410522 }, {430.0, 0.382701 }, {450.0, 0.356957 }, {470.0, 0.332400 } };
|
||||
for (int i = 0; i < angles.length; ++i) {
|
||||
problem.addAngularMeasurement(3.0e7, angles[i][0], angles[i][1]);
|
||||
};
|
||||
|
||||
LevenbergMarquardtEstimator estimator = new LevenbergMarquardtEstimator();
|
||||
estimator.estimate(problem);
|
||||
System.out.println("initial position: " + problem.getX0() + " " + problem.getY0());
|
||||
System.out.println("velocity: " + problem.getVx0() + " " + problem.getVy0());
|
||||
|
||||
} catch (EstimationException ee) {
|
||||
System.err.println(ee.getMessage());
|
||||
}
|
||||
}
|
||||
|
||||
public TrajectoryDeterminationProblem(double t0,
|
||||
double x0Guess, double y0Guess,
|
||||
double vx0Guess, double vy0Guess) {
|
||||
this.t0 = t0;
|
||||
x0 = new EstimatedParameter( "x0", x0Guess);
|
||||
y0 = new EstimatedParameter( "y0", y0Guess);
|
||||
vx0 = new EstimatedParameter("vx0", vx0Guess);
|
||||
vy0 = new EstimatedParameter("vy0", vy0Guess);
|
||||
|
||||
// inform the base class about the parameters
|
||||
addParameter(x0);
|
||||
addParameter(y0);
|
||||
addParameter(vx0);
|
||||
addParameter(vy0);
|
||||
|
||||
}
|
||||
|
||||
public double getX0() {
|
||||
return x0.getEstimate();
|
||||
}
|
||||
|
||||
public double getY0() {
|
||||
return y0.getEstimate();
|
||||
}
|
||||
|
||||
public double getVx0() {
|
||||
return vx0.getEstimate();
|
||||
}
|
||||
|
||||
public double getVy0() {
|
||||
return vy0.getEstimate();
|
||||
}
|
||||
|
||||
public void addAngularMeasurement(double wi, double ti, double ai) {
|
||||
// let the base class handle the measurement
|
||||
addMeasurement(new AngularMeasurement(wi, ti, ai));
|
||||
}
|
||||
|
||||
public void addDistanceMeasurement(double wi, double ti, double di) {
|
||||
// let the base class handle the measurement
|
||||
addMeasurement(new DistanceMeasurement(wi, ti, di));
|
||||
}
|
||||
|
||||
public double x(double t) {
|
||||
return x0.getEstimate() + (t - t0) * vx0.getEstimate();
|
||||
}
|
||||
|
||||
public double y(double t) {
|
||||
return y0.getEstimate() + (t - t0) * vy0.getEstimate();
|
||||
}
|
||||
|
||||
private class AngularMeasurement extends WeightedMeasurement {
|
||||
|
||||
public AngularMeasurement(double weight, double t, double angle) {
|
||||
super(weight, angle);
|
||||
this.t = t;
|
||||
}
|
||||
|
||||
public double getTheoreticalValue() {
|
||||
return Math.atan2(y(t), x(t));
|
||||
}
|
||||
|
||||
public double getPartial(EstimatedParameter parameter) {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
double r = Math.sqrt(xt * xt + yt * yt);
|
||||
double u = yt / (r + xt);
|
||||
double c = 2 * u / (1 + u * u);
|
||||
if (parameter == x0) {
|
||||
return -c;
|
||||
} else if (parameter == vx0) {
|
||||
return -c * t;
|
||||
} else if (parameter == y0) {
|
||||
return c * xt / yt;
|
||||
} else {
|
||||
return c * t * xt / yt;
|
||||
}
|
||||
}
|
||||
|
||||
private final double t;
|
||||
private static final long serialVersionUID = -5990040582592763282L;
|
||||
|
||||
}
|
||||
|
||||
private class DistanceMeasurement extends WeightedMeasurement {
|
||||
|
||||
public DistanceMeasurement(double weight, double t, double angle) {
|
||||
super(weight, angle);
|
||||
this.t = t;
|
||||
}
|
||||
|
||||
public double getTheoreticalValue() {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
return Math.sqrt(xt * xt + yt * yt);
|
||||
}
|
||||
|
||||
public double getPartial(EstimatedParameter parameter) {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
double r = Math.sqrt(xt * xt + yt * yt);
|
||||
if (parameter == x0) {
|
||||
return xt / r;
|
||||
} else if (parameter == vx0) {
|
||||
return xt * t / r;
|
||||
} else if (parameter == y0) {
|
||||
return yt / r;
|
||||
} else {
|
||||
return yt * t / r;
|
||||
}
|
||||
}
|
||||
|
||||
private final double t;
|
||||
private static final long serialVersionUID = 3257286197740459503L;
|
||||
|
||||
}
|
||||
|
||||
private double t0;
|
||||
private EstimatedParameter x0;
|
||||
private EstimatedParameter y0;
|
||||
private EstimatedParameter vx0;
|
||||
private EstimatedParameter vy0;
|
||||
|
||||
}
|
||||
Binary file not shown.
|
After Width: | Height: | Size: 20 KiB |
Binary file not shown.
|
After Width: | Height: | Size: 16 KiB |
|
|
@ -1,21 +1,21 @@
|
|||
<?xml version="1.0"?>
|
||||
|
||||
<!--
|
||||
Licensed to the Apache Software Foundation (ASF) under one or more
|
||||
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
|
||||
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.
|
||||
-->
|
||||
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.
|
||||
-->
|
||||
|
||||
<?xml-stylesheet type="text/xsl" href="./xdoc.xsl"?>
|
||||
<!-- $Revision: 480435 $ $Date: 2006-11-29 08:06:35 +0100 (mer., 29 nov. 2006) $ -->
|
||||
|
|
@ -29,86 +29,336 @@
|
|||
<section name="12 Parametric Estimation">
|
||||
<subsection name="12.1 Overview" href="overview">
|
||||
<p>
|
||||
The estimation package provides classes to fit some non-linear model
|
||||
to available observations depending on it. These problems are commonly
|
||||
called estimation problems.
|
||||
The estimation package provides classes to fit some non-linear
|
||||
model to available observations depending on it. These
|
||||
problems are commonly called estimation problems.
|
||||
</p>
|
||||
<p>
|
||||
The estimation problems considered here are parametric problems where
|
||||
a user-provided model depends on initially unknown scalar parameters and
|
||||
several measurements made on values that depend on the model are available.
|
||||
As examples, one can consider the center and radius of a circle given
|
||||
points approximately lying on a ring, or a satellite orbit given range,
|
||||
range-rate and angular measurements from various ground stations.
|
||||
The estimation problems considered here are parametric
|
||||
problems where a user-provided model depends on initially
|
||||
unknown scalar parameters and several measurements made on
|
||||
values that depend on the model are available. As examples,
|
||||
one can consider the center and radius of a circle given
|
||||
points approximately lying on a ring, or a satellite orbit
|
||||
given range, range-rate and angular measurements from various
|
||||
ground stations.
|
||||
</p>
|
||||
<p>
|
||||
One important class of estimation problems is weighted least squares problems.
|
||||
They basically consist in finding the values for some parameters p<sub>k</sub>
|
||||
such that a cost function J =
|
||||
<!-- TODO: get entity for summation imported -->
|
||||
sum(w<sub>i</sub> r<sub>i</sub><sup>2</sup>)
|
||||
is minimized. The various r<sub>i</sub> terms represent the deviation
|
||||
r<sub>i</sub> = mes<sub>i</sub> - mod<sub>i</sub> between the measurements and
|
||||
the parameterized models. The w<sub>i</sub> factors are the measurements weights,
|
||||
they are often chosen either all equal to 1.0 or proportional to the inverse of
|
||||
the variance of the measurement type. The solver adjusts the values of the
|
||||
estimated parameters p<sub>k</sub> which are not bound (i.e. the free parameters).
|
||||
It does not touch the parameters which have been put in a bound state by the user.
|
||||
One important class of estimation problems is weighted least
|
||||
squares problems. They basically consist in finding the values
|
||||
for some parameters p<sub>k</sub> such that a cost function
|
||||
J = sum(w<sub>i</sub>r<sub>i</sub><sup>2</sup>) is minimized.
|
||||
The various r<sub>i</sub> terms represent the deviation
|
||||
r<sub>i</sub> = mes<sub>i</sub> - mod<sub>i</sub>
|
||||
between the measurements and the parameterized models. The
|
||||
w<sub>i</sub> factors are the measurements weights, they are often
|
||||
chosen either all equal to 1.0 or proportional to the inverse of the
|
||||
variance of the measurement type. The solver adjusts the values of
|
||||
the estimated parameters p<sub>k</sub> which are not bound (i.e. the
|
||||
free parameters). It does not touch the parameters which have been
|
||||
put in a bound state by the user.
|
||||
</p>
|
||||
<p>
|
||||
The aim of this package is similar to the aim of the optimization package, but the
|
||||
algorithms are entirely differents as:
|
||||
The aim of this package is similar to the aim of the
|
||||
optimization package, but the algorithms are entirely
|
||||
different as:
|
||||
<ul>
|
||||
<li>
|
||||
they need the partial derivatives of the measurements
|
||||
with respect to the free parameters
|
||||
they need the partial derivatives of the measurements with
|
||||
respect to the free parameters
|
||||
</li>
|
||||
<li>
|
||||
they are residuals based instead of generic cost functions based
|
||||
they are residuals based instead of generic cost functions
|
||||
based
|
||||
</li>
|
||||
</ul>
|
||||
</p>
|
||||
</subsection>
|
||||
|
||||
</subsection>
|
||||
<subsection name="12.2 Models" href="models">
|
||||
The <a href="../apidocs/org/apache/commons/math/estimation/EstimationProblem.html">
|
||||
org.apache.commons.math.estimation.EstimationProblem</a> interface is a very
|
||||
simple container packing together parameters and measurements.
|
||||
</subsection>
|
||||
<subsection name="12.3 Parameters" href="parameters">
|
||||
<subsection name="12.2 Problem modeling" href="problem">
|
||||
<p>
|
||||
The <a href="../apidocs/org/apache/commons/math/estimation/EstimatedParameter.html">
|
||||
org.apache.commons.math.estimation.EstimatedParameter</a> class to represent each
|
||||
estimated parameter. The parameters are set up with a guessed value which will be
|
||||
adjusted by the solver until a best fit is achieved. It is possible to change which
|
||||
parameters are modified and which are preserved thanks to a bound property. Such
|
||||
settings are often needed by expert users to handle contingency cases with very
|
||||
low observability.
|
||||
The problem modeling is the most important part for the
|
||||
user. Understanding it is the key to proper use of the
|
||||
package. One interface and two classes are provided for this
|
||||
purpose: <a href="../apidocs/org/apache/commons/math/estimation/EstimationProblem.html">
|
||||
EstimationProblem</a>, <a href="../apidocs/org/apache/commons/math/estimation/EstimatedParameter.html">
|
||||
EstimatedParameter</a> and <a href="../apidocs/org/apache/commons/math/estimation/WeightedMeasurement.html">
|
||||
WeightedMeasurement</a>.
|
||||
</p>
|
||||
</subsection>
|
||||
<subsection name="12.4 Measurements" href="measurements">
|
||||
<p>
|
||||
The user extends the <a href="../apidocs/org/apache/commons/math/estimation/WeightedMeasurement.html">
|
||||
org.apache.commons.math.estimation.WeightedMeasurement</a> abstract class to define its
|
||||
own measurements. Each measurement types should have its own implementing class, for
|
||||
example in the satellite example above , the user should define three classes, one
|
||||
for range measurements, one for range-rates measurements and one for angular measurements.
|
||||
Each measurement would correspond to an instance of the appropriate class, set up with
|
||||
the date, a reference to the ground station, the weight and the measured value.
|
||||
</p>
|
||||
</subsection>
|
||||
<subsection name="12.5 Solvers" href="solvers">
|
||||
<p>
|
||||
The package provides two common <a href="../apidocs/org/apache/commons/math/estimation/Estimator.html">
|
||||
org.apache.commons.math.estimation.Estimator</a> implementations to solve weighted
|
||||
least squares problems. The first one is based on the
|
||||
<a href="../apidocs/org/apache/commons/math/estimation/GaussNewtonEstimator.html">Gauss-Newton</a> method.
|
||||
The second one is based on the
|
||||
<a href="../apidocs/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.html">Levenberg-Marquardt</a>
|
||||
method. The first one is best suited when a good approximation of the parameters is known while the second one
|
||||
is more robust and can handle starting points far from the solution.
|
||||
Consider the following example problem: we want to determine the
|
||||
linear trajectory of a sailing ship by performing angular and
|
||||
distance measurements from an observing spot on the shore. The
|
||||
problem model is represented by two equations:
|
||||
</p>
|
||||
</subsection>
|
||||
<p>
|
||||
x(t) = x<sub>0</sub>+(t-t<sub>0</sub>)vx<sub>0</sub><br/>
|
||||
y(t) = y<sub>0</sub>+(t-t<sub>0</sub>)vy<sub>0</sub>
|
||||
</p>
|
||||
<p>
|
||||
These two equations depend on four parameters (x<sub>0</sub>, y<sub>0</sub>,
|
||||
vx<sub>0</sub> and vy<sub>0</sub>). We want to determine these four parameters.
|
||||
</p>
|
||||
<p>
|
||||
Assuming the observing spot is located at the origin of the coordinates
|
||||
system and that the angular measurements correspond to the angle between
|
||||
the x axis and the line of sight, the theoretical values of the angular
|
||||
measurements at t<sub>i</sub> and of the distance measurements at
|
||||
t<sub>j</sub> are modeled as follows:
|
||||
</p>
|
||||
<p>
|
||||
angle<sub>i,theo</sub> = atan2(y(t<sub>i</sub>), x(t<sub>i</sub>))<br/>
|
||||
distance<sub>j,theo</sub> = sqrt(x(t<sub>j</sub>)<sup>2</sup>+y(t<sub>j</sub>)<sup>2</sup>)
|
||||
</p>
|
||||
<p>
|
||||
The real observations generate a set of measurements values angle<sub>i,meas</sub>
|
||||
and distance<sub>j,meas</sub>.
|
||||
</p>
|
||||
<p>
|
||||
The following class diagram shows one way to solve this problem using the
|
||||
estimation package. The grey elements are already provided by the package
|
||||
whereas the purple elements are developed by the user.
|
||||
</p>
|
||||
<img src="./estimation-class-diagram.png"/>
|
||||
<p>
|
||||
The <code>TrajectoryDeterminationProblem</code> class holds the linear model
|
||||
equations x(t) and y(t). It delegate storage of the four parameters x<sub>0</sub>,
|
||||
y<sub>0</sub>, vx<sub>0</sub> and vy<sub>0</sub> and of the various measurements
|
||||
angle<sub>i,meas</sub> and distance<sub>j,meas</sub> to its base class
|
||||
<code>SimpleEstimationProblem</code>. Since the theoretical values of the measurements
|
||||
angle<sub>i,theo</sub> and distance<sub>j,theo</sub> depend on the linear model,
|
||||
the two classes <code>AngularMeasurement</code> and <code>DistanceMeasurement</code>
|
||||
are implemented as internal classes, thus having access to the equations of the
|
||||
linear model and to the parameters.
|
||||
</p>
|
||||
<p>
|
||||
Here are the various parts of the <code>TrajectoryDeterminationProblem.java</code>
|
||||
source file. This example, with an additional <code>main</code> method is
|
||||
available <a href="./TrajectoryDeterminationProblem.java">here</a>.
|
||||
</p>
|
||||
<dd>First, the general setup of the class: declarations, fields, constructor, setters and getters:
|
||||
<source>
|
||||
public class TrajectoryDeterminationProblem extends SimpleEstimationProblem {
|
||||
public TrajectoryDeterminationProblem(double t0,
|
||||
double x0Guess, double y0Guess,
|
||||
double vx0Guess, double vy0Guess) {
|
||||
this.t0 = t0;
|
||||
x0 = new EstimatedParameter( "x0", x0Guess);
|
||||
y0 = new EstimatedParameter( "y0", y0Guess);
|
||||
vx0 = new EstimatedParameter("vx0", vx0Guess);
|
||||
vy0 = new EstimatedParameter("vy0", vy0Guess);
|
||||
|
||||
// inform the base class about the parameters
|
||||
addParameter(x0);
|
||||
addParameter(y0);
|
||||
addParameter(vx0);
|
||||
addParameter(vy0);
|
||||
|
||||
}
|
||||
|
||||
public double getX0() {
|
||||
return x0.getEstimate();
|
||||
}
|
||||
|
||||
public double getY0() {
|
||||
return y0.getEstimate();
|
||||
}
|
||||
|
||||
public double getVx0() {
|
||||
return vx0.getEstimate();
|
||||
}
|
||||
|
||||
public double getVy0() {
|
||||
return vy0.getEstimate();
|
||||
}
|
||||
|
||||
public void addAngularMeasurement(double wi, double ti, double ai) {
|
||||
// let the base class handle the measurement
|
||||
addMeasurement(new AngularMeasurement(wi, ti, ai));
|
||||
}
|
||||
|
||||
public void addDistanceMeasurement(double wi, double ti, double di) {
|
||||
// let the base class handle the measurement
|
||||
addMeasurement(new DistanceMeasurement(wi, ti, di));
|
||||
}
|
||||
|
||||
public double x(double t) {
|
||||
return x0.getEstimate() + (t - t0) * vx0.getEstimate();
|
||||
}
|
||||
|
||||
public double y(double t) {
|
||||
return y0.getEstimate() + (t - t0) * vy0.getEstimate();
|
||||
}
|
||||
|
||||
// measurements internal classes go here
|
||||
|
||||
private double t0;
|
||||
private EstimatedParameter x0;
|
||||
private EstimatedParameter y0;
|
||||
private EstimatedParameter vx0;
|
||||
private EstimatedParameter vy0;
|
||||
|
||||
}
|
||||
</source>
|
||||
</dd>
|
||||
<dd>The two specialized measurements class are simple internal classes that
|
||||
implement the equation for their respective measurement type, using the
|
||||
enclosing class to get the parameters references and the linear models x(t)
|
||||
and y(t). The <code>serialVersionUID</code> static fields are present because
|
||||
the <code>WeightedMeasurement</code> class implements the
|
||||
<code>Serializable</code> interface.
|
||||
<source>
|
||||
private class AngularMeasurement extends WeightedMeasurement {
|
||||
|
||||
public AngularMeasurement(double weight, double t, double angle) {
|
||||
super(weight, angle);
|
||||
this.t = t;
|
||||
}
|
||||
|
||||
public double getTheoreticalValue() {
|
||||
return Math.atan2(y(t), x(t));
|
||||
}
|
||||
|
||||
public double getPartial(EstimatedParameter parameter) {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
double r = Math.sqrt(xt * xt + yt * yt);
|
||||
double u = yt / (r + xt);
|
||||
double c = 2 * u / (1 + u * u);
|
||||
if (parameter == x0) {
|
||||
return -c;
|
||||
} else if (parameter == vx0) {
|
||||
return -c * t;
|
||||
} else if (parameter == y0) {
|
||||
return c * xt / yt;
|
||||
} else {
|
||||
return c * t * xt / yt;
|
||||
}
|
||||
}
|
||||
|
||||
private final double t;
|
||||
private static final long serialVersionUID = -5990040582592763282L;
|
||||
|
||||
}
|
||||
</source>
|
||||
<source>
|
||||
private class DistanceMeasurement extends WeightedMeasurement {
|
||||
|
||||
public DistanceMeasurement(double weight, double t, double angle) {
|
||||
super(weight, angle);
|
||||
this.t = t;
|
||||
}
|
||||
|
||||
public double getTheoreticalValue() {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
return Math.sqrt(xt * xt + yt * yt);
|
||||
}
|
||||
|
||||
public double getPartial(EstimatedParameter parameter) {
|
||||
double xt = x(t);
|
||||
double yt = y(t);
|
||||
double r = Math.sqrt(xt * xt + yt * yt);
|
||||
if (parameter == x0) {
|
||||
return xt / r;
|
||||
} else if (parameter == vx0) {
|
||||
return xt * t / r;
|
||||
} else if (parameter == y0) {
|
||||
return yt / r;
|
||||
} else {
|
||||
return yt * t / r;
|
||||
}
|
||||
}
|
||||
|
||||
private final double t;
|
||||
private static final long serialVersionUID = 3257286197740459503L;
|
||||
|
||||
}
|
||||
</source>
|
||||
</dd>
|
||||
</subsection>
|
||||
<subsection name="12.3 Problem solving" href="solving">
|
||||
<p>
|
||||
Solving the problem is simply a matter of choosing an implementation
|
||||
of the <a href="../apidocs/org/apache/commons/math/estimation/Estimator.html">
|
||||
Estimator</a> interface and to pass the problem instance to its <code>estimate</code>
|
||||
method. Two implementations are already provided by the library: <a
|
||||
href="../apidocs/org/apache/commons/math/estimation/GaussNewtonEstimator.html">
|
||||
GaussNewtonEstimator</a> and <a
|
||||
href="../apidocs/org/apache/commons/math/estimation/LevenbergMarquardtEstimator.html">
|
||||
LevenbergMarquardtEstimator</a>. The first one implements a simple Gauss-Newton
|
||||
algorithm, which is sufficient when the starting point (initial guess) is close
|
||||
enough to the solution. The second one implements a more complex Levenberg-Marquardt
|
||||
algorithm which is more robust when the initial guess is far from the solution.
|
||||
</p>
|
||||
<p>
|
||||
The following sequence diagram explains roughly what occurs under the hood
|
||||
in the <code>estimate</code> method.
|
||||
</p>
|
||||
<img src="./estimation-sequence-diagram.png"/>
|
||||
<p>
|
||||
Basically, the estimator first retrieves the parameters and the measurements.
|
||||
The estimation loop is based on the gradient of the sum of the squares of the
|
||||
residuals, hence, the estimators get the various partial derivatives of all
|
||||
measurements with respect to all parameters. A new state hopefully globally
|
||||
reducing the residuals is estimated, and the parameters value are updated.
|
||||
This estimation loops stops when either the convergence conditions are met
|
||||
or the maximal number of iterations is exceeded.
|
||||
</p>
|
||||
</subsection>
|
||||
<subsection name="12.4 Fine tuning" href="tuning">
|
||||
<p>
|
||||
One important tuning parameter for weighted least-squares solving is the
|
||||
weight attributed to each measurement. This weights has two purposes:
|
||||
</p>
|
||||
<ul>
|
||||
<li>fixing unit problems when combining different types of measurements</li>
|
||||
<li>adjusting the influence of good or bad measurements on the solution</li>
|
||||
</ul>
|
||||
<p>
|
||||
The weight is a multiplicative factor for the <em>square</em> of the residuals.
|
||||
A common choice is to use the inverse of the variance of the measurements error
|
||||
as the weighting factor for all measurements for one type. On our sailing ship
|
||||
example, we may have a range measurements accuracy of about 1 meter and an angular
|
||||
measurements accuracy of about 0.01 degree, or 1.7 10<sup>-4</sup> radians. So we
|
||||
would use w=1.0 for distance measurements weight and w=3 10<sup>7</sup> for
|
||||
angular measurements weight. If we knew that the measurements quality is bad
|
||||
at tracking start because of measurement system warm-up delay for example, then
|
||||
we would reduce the weight for the first measurements and use for example
|
||||
w=0.1 and w=3 10<sup>6</sup> respectively, depending on the type.
|
||||
</p>
|
||||
<p>
|
||||
After a problem has been set up, it is possible to fine tune the
|
||||
way it will be solved. For example, it may appear the measurements are not
|
||||
sufficient to get some parameters with sufficient confidence due to observability
|
||||
problems. It is possible to fix some parameters in order to prevent the solver
|
||||
from changing them. This is realized by passing <code>true</code> to the
|
||||
<code>setBound</code> method of the parameter.
|
||||
</p>
|
||||
<p>
|
||||
It is also possible to ignore some measurements by passing <code>true</code> to the
|
||||
<code>setIgnored</code> method of the measurement. A typical use is to
|
||||
<ol>
|
||||
<li>
|
||||
perform a first determination with all parameters, to check each measurement
|
||||
residual after convergence (i.e. to compute the difference between the
|
||||
measurement and its theoretical value as computed from the estimated parameters),
|
||||
</li>
|
||||
<li>
|
||||
compute standard deviation for the measurements samples (one sample for each
|
||||
measurements type)
|
||||
</li>
|
||||
<li>
|
||||
ignore measurements whose residual are above some threshold (for example three
|
||||
time the standard deviation on the residuals) assuming they correspond to
|
||||
bad measurements,
|
||||
</li>
|
||||
<li>
|
||||
perform another determination on the reduced measurements set.
|
||||
</li>
|
||||
</ol>
|
||||
</p>
|
||||
</subsection>
|
||||
</section>
|
||||
</body>
|
||||
</document>
|
||||
|
|
|
|||
Loading…
Reference in New Issue