[commons-statistics] 02/04: Correct inverse CDF for Geometric distribution.

[ah...@apache.org]

This is an automated email from the ASF dual-hosted git repository.

aherbert pushed a commit to branch master
in repository https://gitbox.apache.org/repos/asf/commons-statistics.git

commit cc74870d70349f3f6e4325a6f8caec18fcd4d5e7
Author: aherbert <ah...@apache.org>
AuthorDate: Tue Sep 21 13:09:12 2021 +0100

    Correct inverse CDF for Geometric distribution.
    
    This requires a check of the result using the cdf and then correction to
    round down.
---
 .../distribution/GeometricDistribution.java          |  9 ++++++++-
 .../distribution/GeometricDistributionTest.java      | 20 ++++++++++++++++++++
 2 files changed, 28 insertions(+), 1 deletion(-)

diff --git a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/GeometricDistribution.java b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/GeometricDistribution.java
index b7be11f..0724715 100644
--- a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/GeometricDistribution.java
+++ b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/GeometricDistribution.java
@@ -108,7 +108,14 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
         if (p == 0) {
             return 0;
         }
-        return Math.max(0, (int) Math.ceil(Math.log1p(-p) / log1mProbabilityOfSuccess - 1));
+        final int x = (int) Math.ceil(Math.log1p(-p) / log1mProbabilityOfSuccess - 1);
+        // Note: x may be too high due to floating-point error and rounding up with ceil.
+        // Return the next value down if that is also above the input cumulative probability.
+        // This ensures x == icdf(cdf(x))
+        if (x <= 0) {
+            return 0;
+        }
+        return cumulativeProbability(x - 1) >= p ? x - 1 : x;
     }
 
     /**
diff --git a/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/GeometricDistributionTest.java b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/GeometricDistributionTest.java
index 41e8e06..9081408 100644
--- a/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/GeometricDistributionTest.java
+++ b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/GeometricDistributionTest.java
@@ -16,6 +16,7 @@
  */
 package org.apache.commons.statistics.distribution;
 
+import java.util.Arrays;
 import org.apache.commons.rng.simple.RandomSource;
 import org.junit.jupiter.api.Assertions;
 import org.junit.jupiter.api.BeforeEach;
@@ -185,6 +186,25 @@ class GeometricDistributionTest extends DiscreteDistributionAbstractTest {
         }
     }
 
+    /**
+     * Test the inverse CDF returns the correct x from the CDF result.
+     * This case was identified using various probabilities to discover a mismatch
+     * of x != icdf(cdf(x)). This occurs due to rounding errors on the inversion.
+     *
+     * @param p Probability of success
+     */
+    @ParameterizedTest
+    @ValueSource(doubles = {0.2, 0.8})
+    void testInverseCDF(double p) {
+        final GeometricDistribution dist = new GeometricDistribution(p);
+        final int[] x = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10};
+        final double[] values = Arrays.stream(x).mapToDouble(dist::cumulativeProbability).toArray();
+        setDistribution(dist);
+        setInverseCumulativeTestPoints(values);
+        setInverseCumulativeTestValues(x);
+        verifyInverseCumulativeProbabilities();
+    }
+
     @ParameterizedTest
     @ValueSource(doubles = {0.1, 0.456, 0.999})
     void testParameterAccessors(double p) {