You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by ah...@apache.org on 2021/09/21 16:32:46 UTC
[commons-statistics] branch master updated (5822059 -> c1edfd4)
This is an automated email from the ASF dual-hosted git repository.
aherbert pushed a change to branch master
in repository https://gitbox.apache.org/repos/asf/commons-statistics.git.
from 5822059 Fix geometric distribution for p=1
new de9832c Consistent testing of p=0 or p=1 in Geometric and Pascal distribution
new cc74870 Correct inverse CDF for Geometric distribution.
new 5f44ae3 STATISTICS-34: Switch the Geometric PMF based on the probability
new c1edfd4 Add error message for 0 < p <= 1
The 4 revisions listed above as "new" are entirely new to this
repository and will be described in separate emails. The revisions
listed as "add" were already present in the repository and have only
been added to this reference.
Summary of changes:
.../distribution/DistributionException.java | 2 +
.../distribution/GeometricDistribution.java | 37 +++++++++++--
.../distribution/PascalDistribution.java | 25 ++++-----
.../distribution/GeometricDistributionTest.java | 62 +++++++++++++++++++---
.../distribution/PascalDistributionTest.java | 61 +++++++++++++++------
5 files changed, 144 insertions(+), 43 deletions(-)
[commons-statistics] 01/04: Consistent testing of p=0 or p=1 in
Geometric and Pascal distribution
Posted by 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 de9832c69f094fb1d9050693122dead2e73ad81b
Author: aherbert <ah...@apache.org>
AuthorDate: Tue Sep 21 11:43:33 2021 +0100
Consistent testing of p=0 or p=1 in Geometric and Pascal distribution
Fix Pascal distribution upper bound when p=1.
---
.../distribution/GeometricDistribution.java | 6 ++-
.../distribution/PascalDistribution.java | 25 ++++-----
.../distribution/GeometricDistributionTest.java | 13 ++---
.../distribution/PascalDistributionTest.java | 61 ++++++++++++++++------
4 files changed, 65 insertions(+), 40 deletions(-)
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 3f16e05..b7be11f 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
@@ -59,7 +59,8 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
@Override
public double probability(int x) {
if (x <= 0) {
- return x < 0 ? 0.0 : probabilityOfSuccess;
+ // Special case of x=0 exploiting cancellation.
+ return x == 0 ? probabilityOfSuccess : 0.0;
}
return Math.exp(log1mProbabilityOfSuccess * x) * probabilityOfSuccess;
}
@@ -68,7 +69,8 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
@Override
public double logProbability(int x) {
if (x <= 0) {
- return x < 0 ? Double.NEGATIVE_INFINITY : logProbabilityOfSuccess;
+ // Special case of x=0 exploiting cancellation.
+ return x == 0 ? logProbabilityOfSuccess : Double.NEGATIVE_INFINITY;
}
return x * log1mProbabilityOfSuccess + logProbabilityOfSuccess;
}
diff --git a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/PascalDistribution.java b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/PascalDistribution.java
index b75e55e..164cec1 100644
--- a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/PascalDistribution.java
+++ b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/PascalDistribution.java
@@ -108,11 +108,9 @@ public class PascalDistribution extends AbstractDiscreteDistribution {
/** {@inheritDoc} */
@Override
public double probability(int x) {
- if (x < 0) {
- return 0.0;
- } else if (x == 0) {
- // Special case exploiting cancellation.
- return probabilityOfSuccessPowNumOfSuccesses;
+ if (x <= 0) {
+ // Special case of x=0 exploiting cancellation.
+ return x == 0 ? probabilityOfSuccessPowNumOfSuccesses : 0.0;
}
final int n = x + numberOfSuccesses - 1;
if (n < 0) {
@@ -127,11 +125,9 @@ public class PascalDistribution extends AbstractDiscreteDistribution {
/** {@inheritDoc} */
@Override
public double logProbability(int x) {
- if (x < 0) {
- return Double.NEGATIVE_INFINITY;
- } else if (x == 0) {
- // Special case exploiting cancellation.
- return logProbabilityOfSuccessByNumOfSuccesses;
+ if (x <= 0) {
+ // Special case of x=0 exploiting cancellation.
+ return x == 0 ? logProbabilityOfSuccessByNumOfSuccesses : Double.NEGATIVE_INFINITY;
}
final int n = x + numberOfSuccesses - 1;
if (n < 0) {
@@ -206,15 +202,14 @@ public class PascalDistribution extends AbstractDiscreteDistribution {
/**
* {@inheritDoc}
*
- * <p>The upper bound of the support is always positive infinity no matter the
- * parameters. Positive infinity is symbolized by {@code Integer.MAX_VALUE}.
+ * <p>The upper bound of the support is positive infinity except for the
+ * probability parameter {@code p = 1.0}.
*
- * @return upper bound of the support (always {@code Integer.MAX_VALUE}
- * for positive infinity)
+ * @return upper bound of the support ({@code Integer.MAX_VALUE} or 0)
*/
@Override
public int getSupportUpperBound() {
- return Integer.MAX_VALUE;
+ return probabilityOfSuccess < 1 ? Integer.MAX_VALUE : 0;
}
/**
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 25b528a..41e8e06 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
@@ -154,19 +154,20 @@ class GeometricDistributionTest extends DiscreteDistributionAbstractTest {
//-------------------- Additional test cases -------------------------------
+ /** Test degenerate case p = 1. */
@Test
- void testProbabilityOfSuccessOne() {
+ void testDegenerate1() {
final GeometricDistribution dist = new GeometricDistribution(1.0);
Assertions.assertEquals(1.0, dist.getProbabilityOfSuccess());
Assertions.assertEquals(0.0, dist.getMean());
Assertions.assertEquals(0.0, dist.getVariance());
setDistribution(dist);
- setProbabilityTestPoints(new int[] {0, 1, 2});
- setProbabilityTestValues(new double[] {1.0, 0.0, 0.0});
- setCumulativeTestPoints(new int[] {0, 1, 2});
- setCumulativeTestValues(new double[] {1.0, 1.0, 1.0});
- setInverseCumulativeTestPoints(new double[] {0, 0.5, 1.0});
+ setProbabilityTestPoints(new int[] {-1, 0, 1, 2, 5, 10});
+ setProbabilityTestValues(new double[] {0.0, 1.0, 0.0, 0.0, 0.0, 0.0});
+ setCumulativeTestPoints(new int[] {-1, 0, 1, 2, 5, 10 });
+ setCumulativeTestValues(new double[] {0.0, 1.0, 1.0, 1.0, 1.0, 1.0});
+ setInverseCumulativeTestPoints(new double[] {0.1, 0.5, 1.0});
setInverseCumulativeTestValues(new int[] {0, 0, 0});
verifyProbabilities();
verifyLogProbabilities();
diff --git a/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/PascalDistributionTest.java b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/PascalDistributionTest.java
index 6df69ee..22ddd2c 100644
--- a/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/PascalDistributionTest.java
+++ b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/PascalDistributionTest.java
@@ -16,6 +16,7 @@
*/
package org.apache.commons.statistics.distribution;
+import org.apache.commons.rng.simple.RandomSource;
import org.junit.jupiter.api.Assertions;
import org.junit.jupiter.api.BeforeEach;
import org.junit.jupiter.api.Test;
@@ -93,40 +94,66 @@ class PascalDistributionTest extends DiscreteDistributionAbstractTest {
//-------------------- Additional test cases -------------------------------
- /** Test degenerate case p = 0 */
+ /** Test degenerate case p = 0. */
@Test
void testDegenerate0() {
- setDistribution(new PascalDistribution(5, 0.0d));
- setCumulativeTestPoints(new int[] {-1, 0, 1, 5, 10 });
- setCumulativeTestValues(new double[] {0d, 0d, 0d, 0d, 0d});
- setProbabilityTestPoints(new int[] {-1, 0, 1, 10, 11});
- setProbabilityTestValues(new double[] {0d, 0d, 0d, 0d, 0d});
- setInverseCumulativeTestPoints(new double[] {0.1d, 0.5d});
- setInverseCumulativeTestValues(new int[] {Integer.MAX_VALUE, Integer.MAX_VALUE});
+ final PascalDistribution dist = new PascalDistribution(5, 0.0d);
+ Assertions.assertEquals(0.0, dist.getProbabilityOfSuccess());
+ Assertions.assertEquals(Double.POSITIVE_INFINITY, dist.getMean());
+ Assertions.assertEquals(Double.POSITIVE_INFINITY, dist.getVariance());
+
+ setDistribution(dist);
+ setProbabilityTestPoints(new int[] {-1, 0, 1, 2, 5, 10});
+ setProbabilityTestValues(new double[] {0.0, 0.0, 0.0, 0.0, 0.0, 0.0});
+ setCumulativeTestPoints(new int[] {-1, 0, 1, 2, 5, 10});
+ setCumulativeTestValues(new double[] {0.0, 0.0, 0.0, 0.0, 0.0, 0.0});
+ setInverseCumulativeTestPoints(new double[] {0.1, 0.5, 1.0});
+ setInverseCumulativeTestValues(new int[] {Integer.MAX_VALUE, Integer.MAX_VALUE, Integer.MAX_VALUE});
verifyProbabilities();
verifyLogProbabilities();
verifyCumulativeProbabilities();
verifySurvivalProbability();
verifySurvivalAndCumulativeProbabilityComplement();
verifyInverseCumulativeProbabilities();
+
+ Assertions.assertEquals(0, dist.getSupportLowerBound());
+ Assertions.assertEquals(Integer.MAX_VALUE, dist.getSupportUpperBound());
+
+ final DiscreteDistribution.Sampler s = dist.createSampler(RandomSource.SPLIT_MIX_64.create());
+ for (int i = 0; i < 5; i++) {
+ Assertions.assertEquals(Integer.MAX_VALUE, s.sample());
+ }
}
- /** Test degenerate case p = 1 */
+ /** Test degenerate case p = 1. */
@Test
void testDegenerate1() {
- setDistribution(new PascalDistribution(5, 1.0d));
- setCumulativeTestPoints(new int[] {-1, 0, 1, 2, 5, 10 });
- setCumulativeTestValues(new double[] {0d, 1d, 1d, 1d, 1d, 1d});
+ final PascalDistribution dist = new PascalDistribution(5, 1.0d);
+ Assertions.assertEquals(1.0, dist.getProbabilityOfSuccess());
+ Assertions.assertEquals(0.0, dist.getMean());
+ Assertions.assertEquals(0.0, dist.getVariance());
+
+ setDistribution(dist);
setProbabilityTestPoints(new int[] {-1, 0, 1, 2, 5, 10});
- setProbabilityTestValues(new double[] {0d, 1d, 0d, 0d, 0d, 0d});
- setInverseCumulativeTestPoints(new double[] {0.1d, 0.5d});
- setInverseCumulativeTestValues(new int[] {0, 0});
+ setProbabilityTestValues(new double[] {0.0, 1.0, 0.0, 0.0, 0.0, 0.0});
+ setCumulativeTestPoints(new int[] {-1, 0, 1, 2, 5, 10});
+ setCumulativeTestValues(new double[] {0.0, 1.0, 1.0, 1.0, 1.0, 1.0});
+ setInverseCumulativeTestPoints(new double[] {0.1, 0.5, 1.0});
+ setInverseCumulativeTestValues(new int[] {0, 0, 0});
verifyProbabilities();
verifyLogProbabilities();
verifyCumulativeProbabilities();
verifySurvivalProbability();
verifySurvivalAndCumulativeProbabilityComplement();
verifyInverseCumulativeProbabilities();
+
+ Assertions.assertEquals(0, dist.getSupportLowerBound());
+ Assertions.assertEquals(0, dist.getSupportUpperBound());
+
+ final DiscreteDistribution.Sampler s = dist.createSampler(RandomSource.SPLIT_MIX_64.create());
+ for (int i = 0; i < 5; i++) {
+ Assertions.assertEquals(0, s.sample());
+ }
}
@ParameterizedTest
@@ -157,11 +184,11 @@ class PascalDistributionTest extends DiscreteDistributionAbstractTest {
PascalDistribution dist;
dist = new PascalDistribution(10, 0.5);
- Assertions.assertEquals((10d * 0.5d) / 0.5d, dist.getMean(), tol);
+ Assertions.assertEquals((10d * 0.5d) / 0.5, dist.getMean(), tol);
Assertions.assertEquals((10d * 0.5d) / (0.5d * 0.5d), dist.getVariance(), tol);
dist = new PascalDistribution(25, 0.7);
- Assertions.assertEquals((25d * 0.3d) / 0.7d, dist.getMean(), tol);
+ Assertions.assertEquals((25d * 0.3d) / 0.7, dist.getMean(), tol);
Assertions.assertEquals((25d * 0.3d) / (0.7d * 0.7d), dist.getVariance(), tol);
}
}
[commons-statistics] 02/04: Correct inverse CDF for Geometric
distribution.
Posted by 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) {
[commons-statistics] 04/04: Add error message for 0 < p <= 1
Posted by 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 c1edfd47243f0c1dcf79455893b649845d8e4966
Author: aherbert <ah...@apache.org>
AuthorDate: Tue Sep 21 17:32:43 2021 +0100
Add error message for 0 < p <= 1
---
.../apache/commons/statistics/distribution/DistributionException.java | 2 ++
.../apache/commons/statistics/distribution/GeometricDistribution.java | 2 +-
2 files changed, 3 insertions(+), 1 deletion(-)
diff --git a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/DistributionException.java b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/DistributionException.java
index 42416d0..5d28eaf 100644
--- a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/DistributionException.java
+++ b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/DistributionException.java
@@ -34,6 +34,8 @@ class DistributionException extends IllegalArgumentException {
static final String INVALID_RANGE_LOW_GT_HIGH = "Lower bound %s > upper bound %s";
/** Error message for "invalid probability" condition when "x not in [0, 1]". */
static final String INVALID_PROBABILITY = "Not a probability: %s is out of range [0, 1]";
+ /** Error message for "invalid non-zero probability" condition when "x not in (0, 1]". */
+ static final String INVALID_NON_ZERO_PROBABILITY = "Not a non-zero probability: %s is out of range (0, 1]";
/** Error message for "negative" condition when "x < 0". */
static final String NEGATIVE = "Number %s is negative";
/** Error message for "not strictly positive" condition when "x <= 0". */
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 a0c6aaa..8164ddf 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
@@ -44,7 +44,7 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
*/
public GeometricDistribution(double p) {
if (p <= 0 || p > 1) {
- throw new DistributionException(DistributionException.INVALID_PROBABILITY, p);
+ throw new DistributionException(DistributionException.INVALID_NON_ZERO_PROBABILITY, p);
}
probabilityOfSuccess = p;
[commons-statistics] 03/04: STATISTICS-34: Switch the Geometric PMF
based on the probability
Posted by 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 5f44ae372410b1994933875bd1ea79b73464b0df
Author: aherbert <ah...@apache.org>
AuthorDate: Tue Sep 21 17:28:04 2021 +0100
STATISTICS-34: Switch the Geometric PMF based on the probability
Using the direct implementation of the PMF is more accurate when 1-p is
exact.
---
.../distribution/GeometricDistribution.java | 20 ++++++++++++++-
.../distribution/GeometricDistributionTest.java | 29 ++++++++++++++++++++++
2 files changed, 48 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 0724715..a0c6aaa 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
@@ -16,6 +16,7 @@
*/
package org.apache.commons.statistics.distribution;
+import java.util.function.IntToDoubleFunction;
import org.apache.commons.rng.UniformRandomProvider;
import org.apache.commons.rng.sampling.distribution.GeometricSampler;
@@ -23,12 +24,17 @@ import org.apache.commons.rng.sampling.distribution.GeometricSampler;
* Implementation of the <a href="http://en.wikipedia.org/wiki/Geometric_distribution">geometric distribution</a>.
*/
public class GeometricDistribution extends AbstractDiscreteDistribution {
+ /** 1/2. */
+ private static final double HALF = 0.5;
+
/** The probability of success. */
private final double probabilityOfSuccess;
/** {@code log(p)} where p is the probability of success. */
private final double logProbabilityOfSuccess;
/** {@code log(1 - p)} where p is the probability of success. */
private final double log1mProbabilityOfSuccess;
+ /** Implementation of PMF(x). Assumes that {@code x > 0}. */
+ private final IntToDoubleFunction pmf;
/**
* Creates a geometric distribution.
@@ -44,6 +50,18 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
probabilityOfSuccess = p;
logProbabilityOfSuccess = Math.log(p);
log1mProbabilityOfSuccess = Math.log1p(-p);
+
+ // Choose the PMF implementation.
+ // When p >= 0.5 then 1 - p is exact and using the power function
+ // is consistently more accurate than the use of the exponential function.
+ // When p -> 0 then the exponential function avoids large error propagation
+ // of the power function used with an inexact 1 - p.
+ if (p >= HALF) {
+ final double q = 1 - p;
+ pmf = x -> Math.pow(q, x) * probabilityOfSuccess;
+ } else {
+ pmf = x -> Math.exp(log1mProbabilityOfSuccess * x) * probabilityOfSuccess;
+ }
}
/**
@@ -62,7 +80,7 @@ public class GeometricDistribution extends AbstractDiscreteDistribution {
// Special case of x=0 exploiting cancellation.
return x == 0 ? probabilityOfSuccess : 0.0;
}
- return Math.exp(log1mProbabilityOfSuccess * x) * probabilityOfSuccess;
+ return pmf.applyAsDouble(x);
}
/** {@inheritDoc} */
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 9081408..1f41e00 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
@@ -187,6 +187,35 @@ class GeometricDistributionTest extends DiscreteDistributionAbstractTest {
}
/**
+ * Test the PMF is computed using the power function when p is above 0.5.
+ * <p>Note: The geometric distribution PMF is defined as:
+ * <pre>
+ * pmf(x) = (1-p)^x * p
+ * </pre>
+ * <p>As {@code p -> 0} use of the power function should be avoided as it will
+ * propagate the inexact computation of {@code 1 - p}. The implementation can
+ * switch to using a rearrangement with the exponential function which avoid
+ * computing {@code 1 - p}.
+ * <p>See STATISTICS-34.
+ *
+ * @param p Probability of success
+ */
+ @ParameterizedTest
+ @ValueSource(doubles = {0.5, 0.6658665, 0.75, 0.8125347, 0.9, 0.95, 0.99})
+ void testPMF(double p) {
+ final GeometricDistribution dist = new GeometricDistribution(p);
+ setDistribution(dist);
+
+ // The PMF should be an exact match to the direct implementation with Math.pow.
+ setTolerance(0);
+ final int[] x = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 20, 30, 40};
+ final double[] values = Arrays.stream(x).mapToDouble(k -> p * Math.pow(1 - p, k)).toArray();
+ setProbabilityTestPoints(x);
+ setProbabilityTestValues(values);
+ verifyProbabilities();
+ }
+
+ /**
* 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.