LoessInterpolator implements UnivariateInterpolator, Serializable { public final double[] smooth(final double[] xval, final double[] yval, final double[] weights) throws NonMonotonicSequenceException, DimensionMismatchException, NoDataException, NotFiniteNumberException, NumberIsTooSmallException { if (xval.length != yval.length) { throw new DimensionMismatchException(xval.length, yval.length); } final int n = xval.length; if (n == 0) { throw new NoDataException(); } checkAllFiniteReal(xval); checkAllFiniteReal(yval); checkAllFiniteReal(weights); MathArrays.checkOrder(xval); if (n == 1) { return new double[]{yval[0]}; } if (n == 2) { return new double[]{yval[0], yval[1]}; } int bandwidthInPoints = (int) (bandwidth * n); if (bandwidthInPoints < 2) { throw new NumberIsTooSmallException(LocalizedFormats.BANDWIDTH, bandwidthInPoints, 2, true); } final double[] res = new double[n]; final double[] residuals = new double[n]; final double[] sortedResiduals = new double[n]; final double[] robustnessWeights = new double[n]; Arrays.fill(robustnessWeights, 1); for (int iter = 0; iter <= robustnessIters; ++iter) { final int[] bandwidthInterval = {0, bandwidthInPoints - 1}; for (int i = 0; i < n; ++i) { final double x = xval[i]; if (i > 0) { updateBandwidthInterval(xval, weights, i, bandwidthInterval); } final int ileft = bandwidthInterval[0]; final int iright = bandwidthInterval[1]; final int edge; if (xval[i] - xval[ileft] > xval[iright] - xval[i]) { edge = ileft; } else { edge = iright; } double sumWeights = 0; double sumX = 0; double sumXSquared = 0; double sumY = 0; double sumXY = 0; double denom = FastMath.abs(1.0 / (xval[edge] - x)); for (int k = ileft; k <= iright; ++k) { final double xk = xval[k]; final double yk = yval[k]; final double dist = (k < i) ? x - xk : xk - x; final double w = tricube(dist * denom) * robustnessWeights[k] * weights[k]; final double xkw = xk * w; sumWeights += w; sumX += xkw; sumXSquared += xk * xkw; sumY += yk * w; sumXY += yk * xkw; } final double meanX = sumX / sumWeights; final double meanY = sumY / sumWeights; final double meanXY = sumXY / sumWeights; final double meanXSquared = sumXSquared / sumWeights; final double beta; if (FastMath.sqrt(FastMath.abs(meanXSquared - meanX * meanX)) < accuracy) { beta = 0; } else { beta = (meanXY - meanX * meanY) / (meanXSquared - meanX * meanX); } final double alpha = meanY - beta * meanX; res[i] = beta * x + alpha; residuals[i] = FastMath.abs(yval[i] - res[i]); } if (iter == robustnessIters) { break; } System.arraycopy(residuals, 0, sortedResiduals, 0, n); Arrays.sort(sortedResiduals); final double medianResidual = sortedResiduals[n / 2]; if (FastMath.abs(medianResidual) < accuracy) { break; } for (int i = 0; i < n; ++i) { final double arg = residuals[i] / (6 * medianResidual); if (arg >= 1) { robustnessWeights[i] = 0; } else { final double w = 1 - arg * arg; robustnessWeights[i] = w * w; } } } return res; } LoessInterpolator();  LoessInterpolator(double bandwidth, int robustnessIters);  LoessInterpolator(double bandwidth, int robustnessIters, double accuracy); final PolynomialSplineFunction interpolate(final double[] xval,
                                                      final double[] yval); final double[] smooth(final double[] xval, final double[] yval,
                                 final double[] weights); final double[] smooth(final double[] xval, final double[] yval); static final double DEFAULT_BANDWIDTH; static final int DEFAULT_ROBUSTNESS_ITERS; static final double DEFAULT_ACCURACY; }
@Test public void testOnOnePoint() { double[] xval = {0.5}; double[] yval = {0.7}; double[] res = new LoessInterpolator().smooth(xval, yval); Assert.assertEquals(1, res.length); Assert.assertEquals(0.7, res[0], 0.0); }
@Test public void testOnTwoPoints() { double[] xval = {0.5, 0.6}; double[] yval = {0.7, 0.8}; double[] res = new LoessInterpolator().smooth(xval, yval); Assert.assertEquals(2, res.length); Assert.assertEquals(0.7, res[0], 0.0); Assert.assertEquals(0.8, res[1], 0.0); }
@Test public void testOnStraightLine() { double[] xval = {1,2,3,4,5}; double[] yval = {2,4,6,8,10}; LoessInterpolator li = new LoessInterpolator(0.6, 2, 1e-12); double[] res = li.smooth(xval, yval); Assert.assertEquals(5, res.length); for(int i = 0; i < 5; ++i) { Assert.assertEquals(yval[i], res[i], 1e-8); } }
@Test public void testOnDistortedSine() { int numPoints = 100; double[] xval = new double[numPoints]; double[] yval = new double[numPoints]; double xnoise = 0.1; double ynoise = 0.2; generateSineData(xval, yval, xnoise, ynoise); LoessInterpolator li = new LoessInterpolator(0.3, 4, 1e-12); double[] res = li.smooth(xval, yval); double noisyResidualSum = 0; double fitResidualSum = 0; for(int i = 0; i < numPoints; ++i) { double expected = FastMath.sin(xval[i]); double noisy = yval[i]; double fit = res[i]; noisyResidualSum += FastMath.pow(noisy - expected, 2); fitResidualSum += FastMath.pow(fit - expected, 2); } Assert.assertTrue(fitResidualSum < noisyResidualSum); }
@Test public void testIncreasingBandwidthIncreasesSmoothness() { int numPoints = 100; double[] xval = new double[numPoints]; double[] yval = new double[numPoints]; double xnoise = 0.1; double ynoise = 0.1; generateSineData(xval, yval, xnoise, ynoise); double[] bandwidths = {0.1, 0.5, 1.0}; double[] variances = new double[bandwidths.length]; for (int i = 0; i < bandwidths.length; i++) { double bw = bandwidths[i]; LoessInterpolator li = new LoessInterpolator(bw, 4, 1e-12); double[] res = li.smooth(xval, yval); for (int j = 1; j < res.length; ++j) { variances[i] += FastMath.pow(res[j] - res[j-1], 2); } } for(int i = 1; i < variances.length; ++i) { Assert.assertTrue(variances[i] < variances[i-1]); } }
@Test public void testIncreasingRobustnessItersIncreasesSmoothnessWithOutliers() { int numPoints = 100; double[] xval = new double[numPoints]; double[] yval = new double[numPoints]; double xnoise = 0.1; double ynoise = 0.1; generateSineData(xval, yval, xnoise, ynoise); yval[numPoints/3] *= 100; yval[2 * numPoints/3] *= -100; double[] variances = new double[4]; for (int i = 0; i < 4; i++) { LoessInterpolator li = new LoessInterpolator(0.3, i, 1e-12); double[] res = li.smooth(xval, yval); for (int j = 1; j < res.length; ++j) { variances[i] += FastMath.abs(res[j] - res[j-1]); } } for(int i = 1; i < variances.length; ++i) { Assert.assertTrue(variances[i] < variances[i-1]); } }
@Test(expected=DimensionMismatchException.class) public void testUnequalSizeArguments() { new LoessInterpolator().smooth(new double[] {1,2,3}, new double[] {1,2,3,4}); }
@Test(expected=NoDataException.class) public void testEmptyData() { new LoessInterpolator().smooth(new double[] {}, new double[] {}); }
@Test(expected=NonMonotonicSequenceException.class) public void testNonStrictlyIncreasing1() { new LoessInterpolator().smooth(new double[] {4,3,1,2}, new double[] {3,4,5,6}); }
@Test(expected=NonMonotonicSequenceException.class) public void testNonStrictlyIncreasing2() { new LoessInterpolator().smooth(new double[] {1,2,2,3}, new double[] {3,4,5,6}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal1() { new LoessInterpolator().smooth(new double[] {1,2,Double.NaN}, new double[] {3,4,5}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal2() { new LoessInterpolator().smooth(new double[] {1,2,Double.POSITIVE_INFINITY}, new double[] {3,4,5}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal3() { new LoessInterpolator().smooth(new double[] {1,2,Double.NEGATIVE_INFINITY}, new double[] {3,4,5}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal4() { new LoessInterpolator().smooth(new double[] {3,4,5}, new double[] {1,2,Double.NaN}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal5() { new LoessInterpolator().smooth(new double[] {3,4,5}, new double[] {1,2,Double.POSITIVE_INFINITY}); }
@Test(expected=NotFiniteNumberException.class) public void testNotAllFiniteReal6() { new LoessInterpolator().smooth(new double[] {3,4,5}, new double[] {1,2,Double.NEGATIVE_INFINITY}); }
@Test(expected=NumberIsTooSmallException.class) public void testInsufficientBandwidth() { LoessInterpolator li = new LoessInterpolator(0.1, 3, 1e-12); li.smooth(new double[] {1,2,3,4,5,6,7,8,9,10,11,12}, new double[] {1,2,3,4,5,6,7,8,9,10,11,12}); }
@Test public void testMath296withoutWeights() { double[] xval = { 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0}; double[] yval = { 0.47, 0.48, 0.55, 0.56, -0.08, -0.04, -0.07, -0.07, -0.56, -0.46, -0.56, -0.52, -3.03, -3.08, -3.09, -3.04, 3.54, 3.46, 3.36, 3.35}; double[] yref = { 0.461, 0.499, 0.541, 0.308, 0.175, -0.042, -0.072, -0.196, -0.311, -0.446, -0.557, -1.497, -2.133, -3.08, -3.09, -0.621, 0.982, 3.449, 3.389, 3.336 }; LoessInterpolator li = new LoessInterpolator(0.3, 4, 1e-12); double[] res = li.smooth(xval, yval); Assert.assertEquals(xval.length, res.length); for(int i = 0; i < res.length; ++i) { Assert.assertEquals(yref[i], res[i], 0.02); } }