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
9
* http://www.apache.org/licenses/LICENSE-2.0
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.
18
package org.apache.commons.math.estimation;
20
import java.util.ArrayList;
21
import java.util.HashSet;
22
import java.util.Iterator;
24
import org.apache.commons.math.estimation.EstimatedParameter;
25
import org.apache.commons.math.estimation.EstimationException;
26
import org.apache.commons.math.estimation.EstimationProblem;
27
import org.apache.commons.math.estimation.GaussNewtonEstimator;
28
import org.apache.commons.math.estimation.WeightedMeasurement;
30
import junit.framework.*;
33
* <p>Some of the unit tests are re-implementations of the MINPACK <a
34
* href="http://www.netlib.org/minpack/ex/file17">file17</a> and <a
35
* href="http://www.netlib.org/minpack/ex/file22">file22</a> test files.
36
* The redistribution policy for MINPACK is available <a
37
* href="http://www.netlib.org/minpack/disclaimer">here</a>, for
38
* convenience, it is reproduced below.</p>
40
* <table border="0" width="80%" cellpadding="10" align="center" bgcolor="#E0E0E0">
42
* Minpack Copyright Notice (1999) University of Chicago.
46
* Redistribution and use in source and binary forms, with or without
47
* modification, are permitted provided that the following conditions
50
* <li>Redistributions of source code must retain the above copyright
51
* notice, this list of conditions and the following disclaimer.</li>
52
* <li>Redistributions in binary form must reproduce the above
53
* copyright notice, this list of conditions and the following
54
* disclaimer in the documentation and/or other materials provided
55
* with the distribution.</li>
56
* <li>The end-user documentation included with the redistribution, if any,
57
* must include the following acknowledgment:
58
* <code>This product includes software developed by the University of
59
* Chicago, as Operator of Argonne National Laboratory.</code>
60
* Alternately, this acknowledgment may appear in the software itself,
61
* if and wherever such third-party acknowledgments normally appear.</li>
62
* <li><strong>WARRANTY DISCLAIMER. THE SOFTWARE IS SUPPLIED "AS IS"
63
* WITHOUT WARRANTY OF ANY KIND. THE COPYRIGHT HOLDER, THE
64
* UNITED STATES, THE UNITED STATES DEPARTMENT OF ENERGY, AND
65
* THEIR EMPLOYEES: (1) DISCLAIM ANY WARRANTIES, EXPRESS OR
66
* IMPLIED, INCLUDING BUT NOT LIMITED TO ANY IMPLIED WARRANTIES
67
* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, TITLE
68
* OR NON-INFRINGEMENT, (2) DO NOT ASSUME ANY LEGAL LIABILITY
69
* OR RESPONSIBILITY FOR THE ACCURACY, COMPLETENESS, OR
70
* USEFULNESS OF THE SOFTWARE, (3) DO NOT REPRESENT THAT USE OF
71
* THE SOFTWARE WOULD NOT INFRINGE PRIVATELY OWNED RIGHTS, (4)
72
* DO NOT WARRANT THAT THE SOFTWARE WILL FUNCTION
73
* UNINTERRUPTED, THAT IT IS ERROR-FREE OR THAT ANY ERRORS WILL
74
* BE CORRECTED.</strong></li>
75
* <li><strong>LIMITATION OF LIABILITY. IN NO EVENT WILL THE COPYRIGHT
76
* HOLDER, THE UNITED STATES, THE UNITED STATES DEPARTMENT OF
77
* ENERGY, OR THEIR EMPLOYEES: BE LIABLE FOR ANY INDIRECT,
78
* INCIDENTAL, CONSEQUENTIAL, SPECIAL OR PUNITIVE DAMAGES OF
79
* ANY KIND OR NATURE, INCLUDING BUT NOT LIMITED TO LOSS OF
80
* PROFITS OR LOSS OF DATA, FOR ANY REASON WHATSOEVER, WHETHER
81
* SUCH LIABILITY IS ASSERTED ON THE BASIS OF CONTRACT, TORT
82
* (INCLUDING NEGLIGENCE OR STRICT LIABILITY), OR OTHERWISE,
83
* EVEN IF ANY OF SAID PARTIES HAS BEEN WARNED OF THE
84
* POSSIBILITY OF SUCH LOSS OR DAMAGES.</strong></li>
88
* @author Argonne National Laboratory. MINPACK project. March 1980 (original fortran minpack tests)
89
* @author Burton S. Garbow (original fortran minpack tests)
90
* @author Kenneth E. Hillstrom (original fortran minpack tests)
91
* @author Jorge J. More (original fortran minpack tests)
92
* @author Luc Maisonobe (non-minpack tests and minpack tests Java translation)
94
public class GaussNewtonEstimatorTest
97
public GaussNewtonEstimatorTest(String name) {
101
public void testTrivial() throws EstimationException {
102
LinearProblem problem =
103
new LinearProblem(new LinearMeasurement[] {
104
new LinearMeasurement(new double[] {2},
105
new EstimatedParameter[] {
106
new EstimatedParameter("p0", 0)
109
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
110
estimator.estimate(problem);
111
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
113
problem.getUnboundParameters()[0].getEstimate(),
117
public void testQRColumnsPermutation() throws EstimationException {
119
EstimatedParameter[] x = {
120
new EstimatedParameter("p0", 0), new EstimatedParameter("p1", 0)
122
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
123
new LinearMeasurement(new double[] { 1.0, -1.0 },
124
new EstimatedParameter[] { x[0], x[1] },
126
new LinearMeasurement(new double[] { 2.0 },
127
new EstimatedParameter[] { x[1] },
129
new LinearMeasurement(new double[] { 1.0, -2.0 },
130
new EstimatedParameter[] { x[0], x[1] },
134
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
135
estimator.estimate(problem);
136
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
137
assertEquals(7.0, x[0].getEstimate(), 1.0e-10);
138
assertEquals(3.0, x[1].getEstimate(), 1.0e-10);
142
public void testNoDependency() throws EstimationException {
143
EstimatedParameter[] p = new EstimatedParameter[] {
144
new EstimatedParameter("p0", 0),
145
new EstimatedParameter("p1", 0),
146
new EstimatedParameter("p2", 0),
147
new EstimatedParameter("p3", 0),
148
new EstimatedParameter("p4", 0),
149
new EstimatedParameter("p5", 0)
151
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
152
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[0] }, 0.0),
153
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[1] }, 1.1),
154
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[2] }, 2.2),
155
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[3] }, 3.3),
156
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[4] }, 4.4),
157
new LinearMeasurement(new double[] {2}, new EstimatedParameter[] { p[5] }, 5.5)
159
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
160
estimator.estimate(problem);
161
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
162
for (int i = 0; i < p.length; ++i) {
163
assertEquals(0.55 * i, p[i].getEstimate(), 1.0e-10);
167
public void testOneSet() throws EstimationException {
169
EstimatedParameter[] p = {
170
new EstimatedParameter("p0", 0),
171
new EstimatedParameter("p1", 0),
172
new EstimatedParameter("p2", 0)
174
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
175
new LinearMeasurement(new double[] { 1.0 },
176
new EstimatedParameter[] { p[0] },
178
new LinearMeasurement(new double[] { -1.0, 1.0 },
179
new EstimatedParameter[] { p[0], p[1] },
181
new LinearMeasurement(new double[] { -1.0, 1.0 },
182
new EstimatedParameter[] { p[1], p[2] },
186
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
187
estimator.estimate(problem);
188
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
189
assertEquals(1.0, p[0].getEstimate(), 1.0e-10);
190
assertEquals(2.0, p[1].getEstimate(), 1.0e-10);
191
assertEquals(3.0, p[2].getEstimate(), 1.0e-10);
195
public void testTwoSets() throws EstimationException {
196
EstimatedParameter[] p = {
197
new EstimatedParameter("p0", 0),
198
new EstimatedParameter("p1", 1),
199
new EstimatedParameter("p2", 2),
200
new EstimatedParameter("p3", 3),
201
new EstimatedParameter("p4", 4),
202
new EstimatedParameter("p5", 5)
205
double epsilon = 1.0e-7;
206
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
208
// 4 elements sub-problem
209
new LinearMeasurement(new double[] { 2.0, 1.0, 4.0 },
210
new EstimatedParameter[] { p[0], p[1], p[3] },
212
new LinearMeasurement(new double[] { -4.0, -2.0, 3.0, -7.0 },
213
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
215
new LinearMeasurement(new double[] { 4.0, 1.0, -2.0, 8.0 },
216
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
218
new LinearMeasurement(new double[] { -3.0, -12.0, -1.0 },
219
new EstimatedParameter[] { p[1], p[2], p[3] },
222
// 2 elements sub-problem
223
new LinearMeasurement(new double[] { epsilon, 1.0 },
224
new EstimatedParameter[] { p[4], p[5] },
225
1.0 + epsilon * epsilon),
226
new LinearMeasurement(new double[] { 1.0, 1.0 },
227
new EstimatedParameter[] { p[4], p[5] },
232
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
233
estimator.estimate(problem);
234
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
235
assertEquals( 3.0, p[0].getEstimate(), 1.0e-10);
236
assertEquals( 4.0, p[1].getEstimate(), 1.0e-10);
237
assertEquals(-1.0, p[2].getEstimate(), 1.0e-10);
238
assertEquals(-2.0, p[3].getEstimate(), 1.0e-10);
239
assertEquals( 1.0 + epsilon, p[4].getEstimate(), 1.0e-10);
240
assertEquals( 1.0 - epsilon, p[5].getEstimate(), 1.0e-10);
244
public void testNonInversible() throws EstimationException {
246
EstimatedParameter[] p = {
247
new EstimatedParameter("p0", 0),
248
new EstimatedParameter("p1", 0),
249
new EstimatedParameter("p2", 0)
251
LinearMeasurement[] m = new LinearMeasurement[] {
252
new LinearMeasurement(new double[] { 1.0, 2.0, -3.0 },
253
new EstimatedParameter[] { p[0], p[1], p[2] },
255
new LinearMeasurement(new double[] { 2.0, 1.0, 3.0 },
256
new EstimatedParameter[] { p[0], p[1], p[2] },
258
new LinearMeasurement(new double[] { -3.0, -9.0 },
259
new EstimatedParameter[] { p[0], p[2] },
262
LinearProblem problem = new LinearProblem(m);
264
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
266
estimator.estimate(problem);
267
fail("an exception should have been caught");
268
} catch (EstimationException ee) {
270
} catch (Exception e) {
271
fail("wrong exception type caught");
275
public void testIllConditioned() throws EstimationException {
276
EstimatedParameter[] p = {
277
new EstimatedParameter("p0", 0),
278
new EstimatedParameter("p1", 1),
279
new EstimatedParameter("p2", 2),
280
new EstimatedParameter("p3", 3)
283
LinearProblem problem1 = new LinearProblem(new LinearMeasurement[] {
284
new LinearMeasurement(new double[] { 10.0, 7.0, 8.0, 7.0 },
285
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
287
new LinearMeasurement(new double[] { 7.0, 5.0, 6.0, 5.0 },
288
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
290
new LinearMeasurement(new double[] { 8.0, 6.0, 10.0, 9.0 },
291
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
293
new LinearMeasurement(new double[] { 7.0, 5.0, 9.0, 10.0 },
294
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
297
GaussNewtonEstimator estimator1 = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
298
estimator1.estimate(problem1);
299
assertEquals(0, estimator1.getRMS(problem1), 1.0e-10);
300
assertEquals(1.0, p[0].getEstimate(), 1.0e-10);
301
assertEquals(1.0, p[1].getEstimate(), 1.0e-10);
302
assertEquals(1.0, p[2].getEstimate(), 1.0e-10);
303
assertEquals(1.0, p[3].getEstimate(), 1.0e-10);
305
LinearProblem problem2 = new LinearProblem(new LinearMeasurement[] {
306
new LinearMeasurement(new double[] { 10.0, 7.0, 8.1, 7.2 },
307
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
309
new LinearMeasurement(new double[] { 7.08, 5.04, 6.0, 5.0 },
310
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
312
new LinearMeasurement(new double[] { 8.0, 5.98, 9.89, 9.0 },
313
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
315
new LinearMeasurement(new double[] { 6.99, 4.99, 9.0, 9.98 },
316
new EstimatedParameter[] { p[0], p[1], p[2], p[3] },
319
GaussNewtonEstimator estimator2 = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
320
estimator2.estimate(problem2);
321
assertEquals(0, estimator2.getRMS(problem2), 1.0e-10);
322
assertEquals(-81.0, p[0].getEstimate(), 1.0e-8);
323
assertEquals(137.0, p[1].getEstimate(), 1.0e-8);
324
assertEquals(-34.0, p[2].getEstimate(), 1.0e-8);
325
assertEquals( 22.0, p[3].getEstimate(), 1.0e-8);
329
public void testMoreEstimatedParametersSimple() throws EstimationException {
331
EstimatedParameter[] p = {
332
new EstimatedParameter("p0", 7),
333
new EstimatedParameter("p1", 6),
334
new EstimatedParameter("p2", 5),
335
new EstimatedParameter("p3", 4)
337
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
338
new LinearMeasurement(new double[] { 3.0, 2.0 },
339
new EstimatedParameter[] { p[0], p[1] },
341
new LinearMeasurement(new double[] { 1.0, -1.0, 1.0 },
342
new EstimatedParameter[] { p[1], p[2], p[3] },
344
new LinearMeasurement(new double[] { 2.0, 1.0 },
345
new EstimatedParameter[] { p[0], p[2] },
349
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
351
estimator.estimate(problem);
352
fail("an exception should have been caught");
353
} catch (EstimationException ee) {
355
} catch (Exception e) {
356
fail("wrong exception type caught");
361
public void testMoreEstimatedParametersUnsorted() throws EstimationException {
362
EstimatedParameter[] p = {
363
new EstimatedParameter("p0", 2),
364
new EstimatedParameter("p1", 2),
365
new EstimatedParameter("p2", 2),
366
new EstimatedParameter("p3", 2),
367
new EstimatedParameter("p4", 2),
368
new EstimatedParameter("p5", 2)
370
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
371
new LinearMeasurement(new double[] { 1.0, 1.0 },
372
new EstimatedParameter[] { p[0], p[1] },
374
new LinearMeasurement(new double[] { 1.0, 1.0, 1.0 },
375
new EstimatedParameter[] { p[2], p[3], p[4] },
377
new LinearMeasurement(new double[] { 1.0, -1.0 },
378
new EstimatedParameter[] { p[4], p[5] },
380
new LinearMeasurement(new double[] { 1.0, -1.0, 1.0 },
381
new EstimatedParameter[] { p[3], p[2], p[5] },
383
new LinearMeasurement(new double[] { 1.0, -1.0 },
384
new EstimatedParameter[] { p[4], p[3] },
388
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
390
estimator.estimate(problem);
391
fail("an exception should have been caught");
392
} catch (EstimationException ee) {
394
} catch (Exception e) {
395
fail("wrong exception type caught");
400
public void testRedundantEquations() throws EstimationException {
401
EstimatedParameter[] p = {
402
new EstimatedParameter("p0", 1),
403
new EstimatedParameter("p1", 1)
405
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
406
new LinearMeasurement(new double[] { 1.0, 1.0 },
407
new EstimatedParameter[] { p[0], p[1] },
409
new LinearMeasurement(new double[] { 1.0, -1.0 },
410
new EstimatedParameter[] { p[0], p[1] },
412
new LinearMeasurement(new double[] { 1.0, 3.0 },
413
new EstimatedParameter[] { p[0], p[1] },
417
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
418
estimator.estimate(problem);
419
assertEquals(0, estimator.getRMS(problem), 1.0e-10);
420
EstimatedParameter[] all = problem.getAllParameters();
421
for (int i = 0; i < all.length; ++i) {
422
assertEquals(all[i].getName().equals("p0") ? 2.0 : 1.0,
423
all[i].getEstimate(), 1.0e-10);
428
public void testInconsistentEquations() throws EstimationException {
429
EstimatedParameter[] p = {
430
new EstimatedParameter("p0", 1),
431
new EstimatedParameter("p1", 1)
433
LinearProblem problem = new LinearProblem(new LinearMeasurement[] {
434
new LinearMeasurement(new double[] { 1.0, 1.0 },
435
new EstimatedParameter[] { p[0], p[1] },
437
new LinearMeasurement(new double[] { 1.0, -1.0 },
438
new EstimatedParameter[] { p[0], p[1] },
440
new LinearMeasurement(new double[] { 1.0, 3.0 },
441
new EstimatedParameter[] { p[0], p[1] },
445
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
446
estimator.estimate(problem);
447
assertTrue(estimator.getRMS(problem) > 0.1);
451
public void testMaxIterations() {
452
Circle circle = new Circle(98.680, 47.345);
453
circle.addPoint( 30.0, 68.0);
454
circle.addPoint( 50.0, -6.0);
455
circle.addPoint(110.0, -20.0);
456
circle.addPoint( 35.0, 15.0);
457
circle.addPoint( 45.0, 97.0);
459
GaussNewtonEstimator estimator = new GaussNewtonEstimator(4, 1.0e-14, 1.0e-14);
460
estimator.estimate(circle);
461
fail("an exception should have been caught");
462
} catch (EstimationException ee) {
464
} catch (Exception e) {
465
fail("wrong exception type caught");
469
public void testCircleFitting() throws EstimationException {
470
Circle circle = new Circle(98.680, 47.345);
471
circle.addPoint( 30.0, 68.0);
472
circle.addPoint( 50.0, -6.0);
473
circle.addPoint(110.0, -20.0);
474
circle.addPoint( 35.0, 15.0);
475
circle.addPoint( 45.0, 97.0);
476
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-10, 1.0e-10);
477
estimator.estimate(circle);
478
double rms = estimator.getRMS(circle);
479
assertEquals(1.768262623567235, Math.sqrt(circle.getM()) * rms, 1.0e-10);
480
assertEquals(69.96016176931406, circle.getRadius(), 1.0e-10);
481
assertEquals(96.07590211815305, circle.getX(), 1.0e-10);
482
assertEquals(48.13516790438953, circle.getY(), 1.0e-10);
485
public void testCircleFittingBadInit() throws EstimationException {
486
Circle circle = new Circle(-12, -12);
487
double[][] points = new double[][] {
488
{-0.312967, 0.072366}, {-0.339248, 0.132965}, {-0.379780, 0.202724},
489
{-0.390426, 0.260487}, {-0.361212, 0.328325}, {-0.346039, 0.392619},
490
{-0.280579, 0.444306}, {-0.216035, 0.470009}, {-0.149127, 0.493832},
491
{-0.075133, 0.483271}, {-0.007759, 0.452680}, { 0.060071, 0.410235},
492
{ 0.103037, 0.341076}, { 0.118438, 0.273884}, { 0.131293, 0.192201},
493
{ 0.115869, 0.129797}, { 0.072223, 0.058396}, { 0.022884, 0.000718},
494
{-0.053355, -0.020405}, {-0.123584, -0.032451}, {-0.216248, -0.032862},
495
{-0.278592, -0.005008}, {-0.337655, 0.056658}, {-0.385899, 0.112526},
496
{-0.405517, 0.186957}, {-0.415374, 0.262071}, {-0.387482, 0.343398},
497
{-0.347322, 0.397943}, {-0.287623, 0.458425}, {-0.223502, 0.475513},
498
{-0.135352, 0.478186}, {-0.061221, 0.483371}, { 0.003711, 0.422737},
499
{ 0.065054, 0.375830}, { 0.108108, 0.297099}, { 0.123882, 0.222850},
500
{ 0.117729, 0.134382}, { 0.085195, 0.056820}, { 0.029800, -0.019138},
501
{-0.027520, -0.072374}, {-0.102268, -0.091555}, {-0.200299, -0.106578},
502
{-0.292731, -0.091473}, {-0.356288, -0.051108}, {-0.420561, 0.014926},
503
{-0.471036, 0.074716}, {-0.488638, 0.182508}, {-0.485990, 0.254068},
504
{-0.463943, 0.338438}, {-0.406453, 0.404704}, {-0.334287, 0.466119},
505
{-0.254244, 0.503188}, {-0.161548, 0.495769}, {-0.075733, 0.495560},
506
{ 0.001375, 0.434937}, { 0.082787, 0.385806}, { 0.115490, 0.323807},
507
{ 0.141089, 0.223450}, { 0.138693, 0.131703}, { 0.126415, 0.049174},
508
{ 0.066518, -0.010217}, {-0.005184, -0.070647}, {-0.080985, -0.103635},
509
{-0.177377, -0.116887}, {-0.260628, -0.100258}, {-0.335756, -0.056251},
510
{-0.405195, -0.000895}, {-0.444937, 0.085456}, {-0.484357, 0.175597},
511
{-0.472453, 0.248681}, {-0.438580, 0.347463}, {-0.402304, 0.422428},
512
{-0.326777, 0.479438}, {-0.247797, 0.505581}, {-0.152676, 0.519380},
513
{-0.071754, 0.516264}, { 0.015942, 0.472802}, { 0.076608, 0.419077},
514
{ 0.127673, 0.330264}, { 0.159951, 0.262150}, { 0.153530, 0.172681},
515
{ 0.140653, 0.089229}, { 0.078666, 0.024981}, { 0.023807, -0.037022},
516
{-0.048837, -0.077056}, {-0.127729, -0.075338}, {-0.221271, -0.067526}
518
for (int i = 0; i < points.length; ++i) {
519
circle.addPoint(points[i][0], points[i][1]);
521
GaussNewtonEstimator estimator = new GaussNewtonEstimator(100, 1.0e-6, 1.0e-6);
523
estimator.estimate(circle);
524
fail("an exception should have been caught");
525
} catch (EstimationException ee) {
527
} catch (Exception e) {
528
fail("wrong exception type caught");
532
private static class LinearProblem extends SimpleEstimationProblem {
534
public LinearProblem(LinearMeasurement[] measurements) {
535
HashSet set = new HashSet();
536
for (int i = 0; i < measurements.length; ++i) {
537
addMeasurement(measurements[i]);
538
EstimatedParameter[] parameters = measurements[i].getParameters();
539
for (int j = 0; j < parameters.length; ++j) {
540
set.add(parameters[j]);
543
for (Iterator iterator = set.iterator(); iterator.hasNext();) {
544
addParameter((EstimatedParameter) iterator.next());
550
private static class LinearMeasurement extends WeightedMeasurement {
552
public LinearMeasurement(double[] factors, EstimatedParameter[] parameters,
554
super(1.0, setPoint, true);
555
this.factors = factors;
556
this.parameters = parameters;
560
public double getTheoreticalValue() {
562
for (int i = 0; i < factors.length; ++i) {
563
v += factors[i] * parameters[i].getEstimate();
568
public double getPartial(EstimatedParameter parameter) {
569
for (int i = 0; i < parameters.length; ++i) {
570
if (parameters[i] == parameter) {
577
public EstimatedParameter[] getParameters() {
581
private double[] factors;
582
private EstimatedParameter[] parameters;
583
private static final long serialVersionUID = -3922448707008868580L;
587
private static class Circle implements EstimationProblem {
589
public Circle(double cx, double cy) {
590
this.cx = new EstimatedParameter("cx", cx);
591
this.cy = new EstimatedParameter(new EstimatedParameter("cy", cy));
592
points = new ArrayList();
595
public void addPoint(double px, double py) {
596
points.add(new PointModel(px, py));
600
return points.size();
603
public WeightedMeasurement[] getMeasurements() {
604
return (WeightedMeasurement[]) points.toArray(new PointModel[points.size()]);
607
public EstimatedParameter[] getAllParameters() {
608
return new EstimatedParameter[] { cx, cy };
611
public EstimatedParameter[] getUnboundParameters() {
612
return new EstimatedParameter[] { cx, cy };
615
public double getPartialRadiusX() {
617
for (Iterator iterator = points.iterator(); iterator.hasNext();) {
618
dRdX += ((PointModel) iterator.next()).getPartialDiX();
620
return dRdX / points.size();
623
public double getPartialRadiusY() {
625
for (Iterator iterator = points.iterator(); iterator.hasNext();) {
626
dRdY += ((PointModel) iterator.next()).getPartialDiY();
628
return dRdY / points.size();
631
public double getRadius() {
633
for (Iterator iterator = points.iterator(); iterator.hasNext();) {
634
r += ((PointModel) iterator.next()).getCenterDistance();
636
return r / points.size();
639
public double getX() {
640
return cx.getEstimate();
643
public double getY() {
644
return cy.getEstimate();
647
private class PointModel extends WeightedMeasurement {
649
public PointModel(double px, double py) {
655
public double getPartial(EstimatedParameter parameter) {
656
if (parameter == cx) {
657
return getPartialDiX() - getPartialRadiusX();
658
} else if (parameter == cy) {
659
return getPartialDiY() - getPartialRadiusY();
664
public double getCenterDistance() {
665
double dx = px - cx.getEstimate();
666
double dy = py - cy.getEstimate();
667
return Math.sqrt(dx * dx + dy * dy);
670
public double getPartialDiX() {
671
return (cx.getEstimate() - px) / getCenterDistance();
674
public double getPartialDiY() {
675
return (cy.getEstimate() - py) / getCenterDistance();
678
public double getTheoreticalValue() {
679
return getCenterDistance() - getRadius();
684
private static final long serialVersionUID = 1L;
688
private EstimatedParameter cx;
689
private EstimatedParameter cy;
690
private ArrayList points;
694
public static Test suite() {
695
return new TestSuite(GaussNewtonEstimatorTest.class);
added
removed
src/test/org/apache/commons/math/estimation/GaussNewtonEstimatorTest.java
