From 282f7175a4d8944e4744cdc8307e58b1efd19e72 Mon Sep 17 00:00:00 2001 From: Sebastien Brisard Date: Wed, 16 May 2012 05:36:40 +0000 Subject: [PATCH] Unit tests for GammaDistribution, based on reference data generated with Maxima. Solves MATH-753. git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1339014 13f79535-47bb-0310-9956-ffa450edef68 --- .../distribution/GammaDistributionTest.java | 159 ++++++++++++++++++ 1 file changed, 159 insertions(+) diff --git a/src/test/java/org/apache/commons/math3/distribution/GammaDistributionTest.java b/src/test/java/org/apache/commons/math3/distribution/GammaDistributionTest.java index 6b5c8fdfa..ddf4302e1 100644 --- a/src/test/java/org/apache/commons/math3/distribution/GammaDistributionTest.java +++ b/src/test/java/org/apache/commons/math3/distribution/GammaDistributionTest.java @@ -17,7 +17,15 @@ package org.apache.commons.math3.distribution; +import java.io.BufferedReader; +import java.io.IOException; +import java.io.InputStream; +import java.io.InputStreamReader; + import org.apache.commons.math3.exception.NotStrictlyPositiveException; +import org.apache.commons.math3.special.Gamma; +import org.apache.commons.math3.stat.descriptive.SummaryStatistics; +import org.apache.commons.math3.util.FastMath; import org.junit.Assert; import org.junit.Test; @@ -167,4 +175,155 @@ public class GammaDistributionTest extends RealDistributionAbstractTest { Assert.assertEquals(dist.getNumericalMean(), 1.1d * 4.2d, tol); Assert.assertEquals(dist.getNumericalVariance(), 1.1d * 4.2d * 4.2d, tol); } + + public static double density(final double x, final double shape, + final double scale) { + /* + * This is a copy of + * double GammaDistribution.density(double) + * prior to r1338548. + */ + if (x < 0) { + return 0; + } + return FastMath.pow(x / scale, shape - 1) / scale * + FastMath.exp(-x / scale) / FastMath.exp(Gamma.logGamma(shape)); + } + + /* + * MATH-753: large values of x or shape parameter cause density(double) to + * overflow. Reference data is generated with the Maxima script + * gamma-distribution.mac, which can be found in + * src/test/resources/org/apache/commons/math3/distribution. + */ + + private void doTestMath753(final double shape, + final double meanNoOF, final double sdNoOF, + final double meanOF, final double sdOF, + final String resourceName) throws IOException { + final GammaDistribution distribution = new GammaDistribution(shape, 1.0); + final SummaryStatistics statOld = new SummaryStatistics(); + final SummaryStatistics statNewNoOF = new SummaryStatistics(); + final SummaryStatistics statNewOF = new SummaryStatistics(); + + final InputStream resourceAsStream; + resourceAsStream = this.getClass().getResourceAsStream(resourceName); + Assert.assertNotNull("Could not find resource " + resourceName, + resourceAsStream); + final BufferedReader in; + in = new BufferedReader(new InputStreamReader(resourceAsStream)); + + try { + for (String line = in.readLine(); line != null; line = in + .readLine()) { + final String[] tokens = line.split(", "); + Assert.assertTrue("expected two floating-point values", + tokens.length == 2); + final double x = Double.parseDouble(tokens[0]); + final String msg = "x = " + x + ", shape = " + shape + + ", scale = 1.0"; + final double expected = Double.parseDouble(tokens[1]); + final double ulp = FastMath.ulp(expected); + final double actualOld = density(x, shape, 1.0); + final double actualNew = distribution.density(x); + final double errOld, errNew; + errOld = FastMath.abs((actualOld - expected) / ulp); + errNew = FastMath.abs((actualNew - expected) / ulp); + + if (Double.isNaN(actualOld) || Double.isInfinite(actualOld)) { + Assert.assertFalse(msg, Double.isNaN(actualNew)); + Assert.assertFalse(msg, Double.isInfinite(actualNew)); + statNewOF.addValue(errNew); + } else { + statOld.addValue(errOld); + statNewNoOF.addValue(errNew); + } + } + if (statOld.getN() != 0) { + /* + * If no overflow occurs, check that new implementation is + * better than old one. + */ + final StringBuilder sb = new StringBuilder("shape = "); + sb.append(shape); + sb.append(", scale = 1.0\n"); + sb.append("Old implementation\n"); + sb.append("------------------\n"); + sb.append(statOld.toString()); + sb.append("New implementation\n"); + sb.append("------------------\n"); + sb.append(statNewNoOF.toString()); + final String msg = sb.toString(); + + final double oldMin = statOld.getMin(); + final double newMin = statNewNoOF.getMin(); + Assert.assertTrue(msg, newMin <= oldMin); + + final double oldMax = statOld.getMax(); + final double newMax = statNewNoOF.getMax(); + Assert.assertTrue(msg, newMax <= oldMax); + + final double oldMean = statOld.getMean(); + final double newMean = statNewNoOF.getMean(); + Assert.assertTrue(msg, newMean <= oldMean); + + final double oldSd = statOld.getStandardDeviation(); + final double newSd = statNewNoOF.getStandardDeviation(); + Assert.assertTrue(msg, newSd <= oldSd); + + Assert.assertTrue(msg, newMean <= meanNoOF); + Assert.assertTrue(msg, newSd <= sdNoOF); + } + if (statNewOF.getN() != 0) { + final double newMean = statNewOF.getMean(); + final double newSd = statNewOF.getStandardDeviation(); + + final StringBuilder sb = new StringBuilder("shape = "); + sb.append(shape); + sb.append(", scale = 1.0"); + sb.append(", max. mean error (ulps) = "); + sb.append(meanOF); + sb.append(", actual mean error (ulps) = "); + sb.append(newMean); + sb.append(", max. sd of error (ulps) = "); + sb.append(sdOF); + sb.append(", actual sd of error (ulps) = "); + sb.append(newSd); + final String msg = sb.toString(); + + Assert.assertTrue(msg, newMean <= meanOF); + Assert.assertTrue(msg, newSd <= sdOF); + } + } catch (IOException e) { + Assert.fail(e.getMessage()); + } finally { + in.close(); + } + } + + + @Test + public void testMath753Shape1() throws IOException { + doTestMath753(1.0, 1.5, 0.5, 0.0, 0.0, "gamma-distribution-shape-1.csv"); + } + + @Test + public void testMath753Shape10() throws IOException { + doTestMath753(10.0, 1.0, 1.0, 0.0, 0.0, "gamma-distribution-shape-10.csv"); + } + + @Test + public void testMath753Shape100() throws IOException { + doTestMath753(100.0, 1.5, 1.0, 0.0, 0.0, "gamma-distribution-shape-100.csv"); + } + + @Test + public void testMath753Shape142() throws IOException { + doTestMath753(142.0, 0.5, 1.5, 40.0, 40.0, "gamma-distribution-shape-142.csv"); + } + + @Test + public void testMath753Shape1000() throws IOException { + doTestMath753(1000.0, 1.0, 1.0, 160.0, 220.0, "gamma-distribution-shape-1000.csv"); + } }