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 2022/11/17 15:40:11 UTC
[commons-statistics] branch master updated: STATISTICS-57: Add a trapezoidal distribution
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
The following commit(s) were added to refs/heads/master by this push:
new 0de1692 STATISTICS-57: Add a trapezoidal distribution
0de1692 is described below
commit 0de1692d7b328b1393ac2948c08cda19f1160ab1
Author: aherbert <ah...@apache.org>
AuthorDate: Fri Nov 11 14:28:14 2022 +0000
STATISTICS-57: Add a trapezoidal distribution
---
.../distribution/TrapezoidalDistribution.java | 469 +++++++++++++++++++++
.../distribution/TrapezoidalDistributionTest.java | 161 +++++++
.../distribution/test.trapezoidal.1.properties | 43 ++
.../distribution/test.trapezoidal.10.properties | 29 ++
.../distribution/test.trapezoidal.11.properties | 55 +++
.../distribution/test.trapezoidal.12.properties | 52 +++
.../distribution/test.trapezoidal.13.properties | 41 ++
.../distribution/test.trapezoidal.2.properties | 50 +++
.../distribution/test.trapezoidal.3.properties | 53 +++
.../distribution/test.trapezoidal.4.properties | 47 +++
.../distribution/test.trapezoidal.5.properties | 46 ++
.../distribution/test.trapezoidal.6.properties | 53 +++
.../distribution/test.trapezoidal.7.properties | 44 ++
.../distribution/test.trapezoidal.8.properties | 44 ++
.../distribution/test.trapezoidal.9.properties | 29 ++
.../distribution/DistributionsApplication.java | 1 +
.../examples/distribution/TrapezoidalCommand.java | 176 ++++++++
17 files changed, 1393 insertions(+)
diff --git a/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/TrapezoidalDistribution.java b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/TrapezoidalDistribution.java
new file mode 100644
index 0000000..98dec47
--- /dev/null
+++ b/commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/TrapezoidalDistribution.java
@@ -0,0 +1,469 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.statistics.distribution;
+
+import org.apache.commons.rng.UniformRandomProvider;
+
+/**
+ * Implementation of the trapezoidal distribution.
+ *
+ * <p>The probability density function of \( X \) is:
+ *
+ * <p>\[ f(x; a, b, c, d) = \begin{cases}
+ * \frac{2}{d+c-a-b}\frac{x-a}{b-a} & \text{for } a\le x \lt b \\
+ * \frac{2}{d+c-a-b} & \text{for } b\le x \lt c \\
+ * \frac{2}{d+c-a-b}\frac{d-x}{d-c} & \text{for } c\le x \le d
+ * \end{cases} \]
+ *
+ * <p>for \( -\infty \lt a \le b \le c \le d \lt \infty \) and
+ * \( x \in [a, d] \).
+ *
+ * @see <a href="https://en.wikipedia.org/wiki/Trapezoidal_distribution">Trapezoidal distribution (Wikipedia)</a>
+ */
+public abstract class TrapezoidalDistribution extends AbstractContinuousDistribution {
+ /** Lower limit of this distribution (inclusive). */
+ protected final double a;
+ /** Start of the trapezoid constant density. */
+ protected final double b;
+ /** End of the trapezoid constant density. */
+ protected final double c;
+ /** Upper limit of this distribution (inclusive). */
+ protected final double d;
+
+ /**
+ * Specialisation of the trapezoidal distribution used when the distribution simplifies
+ * to an alternative distribution.
+ */
+ private static class DelegatedTrapezoidalDistribution extends TrapezoidalDistribution {
+ /** Distribution delegate. */
+ private final ContinuousDistribution delegate;
+
+ /**
+ * @param a Lower limit of this distribution (inclusive).
+ * @param b Start of the trapezoid constant density.
+ * @param c End of the trapezoid constant density.
+ * @param d Upper limit of this distribution (inclusive).
+ * @param delegate Distribution delegate.
+ */
+ DelegatedTrapezoidalDistribution(double a, double b, double c, double d,
+ ContinuousDistribution delegate) {
+ super(a, b, c, d);
+ this.delegate = delegate;
+ }
+
+ @Override
+ public double density(double x) {
+ return delegate.density(x);
+ }
+
+ @Override
+ public double probability(double x0, double x1) {
+ return delegate.probability(x0, x1);
+ }
+
+ @Override
+ public double logDensity(double x) {
+ return delegate.logDensity(x);
+ }
+
+ @Override
+ public double cumulativeProbability(double x) {
+ return delegate.cumulativeProbability(x);
+ }
+
+ @Override
+ public double inverseCumulativeProbability(double p) {
+ return delegate.inverseCumulativeProbability(p);
+ }
+
+ @Override
+ public double survivalProbability(double x) {
+ return delegate.survivalProbability(x);
+ }
+
+ @Override
+ public double inverseSurvivalProbability(double p) {
+ return delegate.inverseSurvivalProbability(p);
+ }
+
+ @Override
+ public double getMean() {
+ return delegate.getMean();
+ }
+
+ @Override
+ public double getVariance() {
+ return delegate.getVariance();
+ }
+
+ @Override
+ public Sampler createSampler(UniformRandomProvider rng) {
+ return delegate.createSampler(rng);
+ }
+ }
+
+ /**
+ * Specialisation of the trapezoidal distribution used when {@code b == c}.
+ *
+ * <p>This delegates all methods to the triangular distribution.
+ */
+ private static class TriangularTrapezoidalDistribution extends DelegatedTrapezoidalDistribution {
+ /**
+ * @param a Lower limit of this distribution (inclusive).
+ * @param b Start/end of the trapezoid constant density (mode).
+ * @param d Upper limit of this distribution (inclusive).
+ */
+ TriangularTrapezoidalDistribution(double a, double b, double d) {
+ super(a, b, b, d, TriangularDistribution.of(a, b, d));
+ }
+ }
+
+ /**
+ * Specialisation of the trapezoidal distribution used when {@code a == b} and {@code c == d}.
+ *
+ * <p>This delegates all methods to the uniform distribution.
+ */
+ private static class UniformTrapezoidalDistribution extends DelegatedTrapezoidalDistribution {
+ /**
+ * @param a Lower limit of this distribution (inclusive).
+ * @param d Upper limit of this distribution (inclusive).
+ */
+ UniformTrapezoidalDistribution(double a, double d) {
+ super(a, a, d, d, UniformContinuousDistribution.of(a, d));
+ }
+ }
+
+ /**
+ * Regular implementation of the trapezoidal distribution.
+ */
+ private static class RegularTrapezoidalDistribution extends TrapezoidalDistribution {
+ /** Cached value (d + c - a - b). */
+ private final double divisor;
+ /** Cached value (b - a). */
+ private final double bma;
+ /** Cached value (d - c). */
+ private final double dmc;
+ /** Cumulative probability at b. */
+ private final double cdfB;
+ /** Cumulative probability at c. */
+ private final double cdfC;
+ /** Survival probability at b. */
+ private final double sfB;
+ /** Survival probability at c. */
+ private final double sfC;
+
+ /**
+ * @param a Lower limit of this distribution (inclusive).
+ * @param b Start of the trapezoid constant density.
+ * @param c End of the trapezoid constant density.
+ * @param d Upper limit of this distribution (inclusive).
+ */
+ RegularTrapezoidalDistribution(double a, double b, double c, double d) {
+ super(a, b, c, d);
+
+ // Sum positive terms
+ divisor = (d - a) + (c - b);
+ bma = b - a;
+ dmc = d - c;
+
+ cdfB = bma / divisor;
+ sfB = 1 - cdfB;
+ sfC = dmc / divisor;
+ cdfC = 1 - sfC;
+ }
+
+ @Override
+ public double density(double x) {
+ // Note: x < a allows correct density where a == b
+ if (x < a) {
+ return 0;
+ }
+ if (x < b) {
+ final double divident = (x - a) / bma;
+ return 2 * (divident / divisor);
+ }
+ if (x < c) {
+ return 2 / divisor;
+ }
+ if (x < d) {
+ final double divident = (d - x) / dmc;
+ return 2 * (divident / divisor);
+ }
+ return 0;
+ }
+
+ @Override
+ public double cumulativeProbability(double x) {
+ if (x <= a) {
+ return 0;
+ }
+ if (x < b) {
+ final double divident = (x - a) * (x - a) / bma;
+ return divident / divisor;
+ }
+ if (x < c) {
+ final double divident = 2 * x - b - a;
+ return divident / divisor;
+ }
+ if (x < d) {
+ final double divident = (d - x) * (d - x) / dmc;
+ return 1 - divident / divisor;
+ }
+ return 1;
+ }
+
+ @Override
+ public double survivalProbability(double x) {
+ // By symmetry:
+ if (x <= a) {
+ return 1;
+ }
+ if (x < b) {
+ final double divident = (x - a) * (x - a) / bma;
+ return 1 - divident / divisor;
+ }
+ if (x < c) {
+ final double divident = 2 * x - b - a;
+ return 1 - divident / divisor;
+ }
+ if (x < d) {
+ final double divident = (d - x) * (d - x) / dmc;
+ return divident / divisor;
+ }
+ return 0;
+ }
+
+ @Override
+ public double inverseCumulativeProbability(double p) {
+ ArgumentUtils.checkProbability(p);
+ if (p == 0) {
+ return a;
+ }
+ if (p == 1) {
+ return d;
+ }
+ if (p < cdfB) {
+ return a + Math.sqrt(p * divisor * bma);
+ }
+ if (p < cdfC) {
+ return 0.5 * ((p * divisor) + a + b);
+ }
+ return d - Math.sqrt((1 - p) * divisor * dmc);
+ }
+
+ @Override
+ public double inverseSurvivalProbability(double p) {
+ // By symmetry:
+ ArgumentUtils.checkProbability(p);
+ if (p == 1) {
+ return a;
+ }
+ if (p == 0) {
+ return d;
+ }
+ if (p > sfB) {
+ return a + Math.sqrt((1 - p) * divisor * bma);
+ }
+ if (p > sfC) {
+ return 0.5 * (((1 - p) * divisor) + a + b);
+ }
+ return d - Math.sqrt(p * divisor * dmc);
+ }
+
+ @Override
+ public double getMean() {
+ // Compute using a standardized distribution
+ // b' = (b-a) / (d-a)
+ // c' = (c-a) / (d-a)
+ final double scale = d - a;
+ final double bp = bma / scale;
+ final double cp = (c - a) / scale;
+ return nonCentralMoment(1, bp, cp) * scale + a;
+ }
+
+ @Override
+ public double getVariance() {
+ // Compute using a standardized distribution
+ // b' = (b-a) / (d-a)
+ // c' = (c-a) / (d-a)
+ final double scale = d - a;
+ final double bp = bma / scale;
+ final double cp = (c - a) / scale;
+ final double mu = nonCentralMoment(1, bp, cp);
+ return (nonCentralMoment(2, bp, cp) - mu * mu) * scale * scale;
+ }
+
+ /**
+ * Compute the {@code k}-th non-central moment of the standardized trapezoidal
+ * distribution.
+ *
+ * <p>Shifting the distribution by scale {@code (d - a)} and location {@code a}
+ * creates a standardized trapezoidal distribution. This has a simplified
+ * non-central moment as {@code a = 0, d = 1, 0 <= b < c <= 1}.
+ * <pre>
+ * 2 1 ( 1 - c^(k+2) )
+ * E[X^k] = ----------- -------------- ( ----------- - b^(k+1) )
+ * (1 + c - b) (k + 1)(k + 2) ( 1 - c )
+ * </pre>
+ *
+ * <p>Simplification eliminates issues computing the moments when {@code a == b}
+ * or {@code c == d} in the original (non-standardized) distribution.
+ *
+ * @param k Moment to compute
+ * @param b Start of the trapezoid constant density (shape parameter in [0, 1]).
+ * @param c End of the trapezoid constant density (shape parameter in [0, 1]).
+ * @return the moment
+ */
+ private static double nonCentralMoment(int k, double b, double c) {
+ // As c -> 1 then (1 - c^(k+2)) loses precision
+ // 1 - x^y == -(x^y - 1) [high precision powm1]
+ // == -(exp(y * log(x)) - 1)
+ // Note: avoid log(1) using the limit:
+ // (1 - c^(k+2)) / (1-c) -> (k+2) as c -> 1
+ final double term1 = c == 1 ? k + 2 : Math.expm1((k + 2) * Math.log(c)) / (c - 1);
+ final double term2 = Math.pow(b, k + 1);
+ return 2 * ((term1 - term2) / (c - b + 1) / ((k + 1) * (k + 2)));
+ }
+ }
+
+ /**
+ * @param a Lower limit of this distribution (inclusive).
+ * @param b Start of the trapezoid constant density.
+ * @param c End of the trapezoid constant density.
+ * @param d Upper limit of this distribution (inclusive).
+ */
+ private TrapezoidalDistribution(double a, double b, double c, double d) {
+ this.a = a;
+ this.b = b;
+ this.c = c;
+ this.d = d;
+ }
+
+ /**
+ * Creates a trapezoidal distribution.
+ *
+ * <p>The distribution density is represented as an up sloping line from
+ * {@code a} to {@code b}, constant from {@code b} to {@code c}, and then a down
+ * sloping line from {@code c} to {@code d}.
+ *
+ * @param a Lower limit of this distribution (inclusive).
+ * @param b Start of the trapezoid constant density (first shape parameter).
+ * @param c End of the trapezoid constant density (second shape parameter).
+ * @param d Upper limit of this distribution (inclusive).
+ * @return the distribution
+ * @throws IllegalArgumentException if {@code a >= d}, if {@code b < a}, if
+ * {@code c < b} or if {@code c > d}.
+ */
+ public static TrapezoidalDistribution of(double a, double b, double c, double d) {
+ if (a >= d) {
+ throw new DistributionException(DistributionException.INVALID_RANGE_LOW_GTE_HIGH,
+ a, d);
+ }
+ if (b < a) {
+ throw new DistributionException(DistributionException.TOO_SMALL,
+ b, a);
+ }
+ if (c < b) {
+ throw new DistributionException(DistributionException.TOO_SMALL,
+ c, b);
+ }
+ if (c > d) {
+ throw new DistributionException(DistributionException.TOO_LARGE,
+ c, d);
+ }
+ // For consistency, delegate to the appropriate simplified distribution.
+ // Note: Floating-point equality comparison is intentional.
+ if (b == c) {
+ return new TriangularTrapezoidalDistribution(a, b, d);
+ }
+ if (d - a == c - b) {
+ return new UniformTrapezoidalDistribution(a, d);
+ }
+ return new RegularTrapezoidalDistribution(a, b, c, d);
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * <p>For lower limit {@code a}, start of the density constant region {@code b},
+ * end of the density constant region {@code c} and upper limit {@code d}, the
+ * mean is:
+ *
+ * <p>\[ \frac{1}{3(d+c-b-a)}\left(\frac{d^3-c^3}{d-c}-\frac{b^3-a^3}{b-a}\right) \]
+ */
+ @Override
+ public abstract double getMean();
+
+ /**
+ * {@inheritDoc}
+ *
+ * <p>For lower limit {@code a}, start of the density constant region {@code b},
+ * end of the density constant region {@code c} and upper limit {@code d}, the
+ * variance is:
+ *
+ * <p>\[ \frac{1}{6(d+c-b-a)}\left(\frac{d^4-c^4}{d-c}-\frac{b^4-a^4}{b-a}\right) - \mu^2 \]
+ *
+ * <p>where \( \mu \) is the mean.
+ */
+ @Override
+ public abstract double getVariance();
+
+ /**
+ * Gets the start of the constant region of the density function.
+ *
+ * <p>This is the first shape parameter {@code b} of the distribution.
+ *
+ * @return the first shape parameter {@code b}
+ */
+ public double getB() {
+ return b;
+ }
+
+ /**
+ * Gets the end of the constant region of the density function.
+ *
+ * <p>This is the second shape parameter {@code c} of the distribution.
+ *
+ * @return the second shape parameter {@code c}
+ */
+ public double getC() {
+ return c;
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * <p>The lower bound of the support is equal to the lower limit parameter
+ * {@code a} of the distribution.
+ */
+ @Override
+ public double getSupportLowerBound() {
+ return a;
+ }
+
+ /**
+ * {@inheritDoc}
+ *
+ * <p>The upper bound of the support is equal to the upper limit parameter
+ * {@code d} of the distribution.
+ */
+ @Override
+ public double getSupportUpperBound() {
+ return d;
+ }
+}
diff --git a/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/TrapezoidalDistributionTest.java b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/TrapezoidalDistributionTest.java
new file mode 100644
index 0000000..5578cf4
--- /dev/null
+++ b/commons-statistics-distribution/src/test/java/org/apache/commons/statistics/distribution/TrapezoidalDistributionTest.java
@@ -0,0 +1,161 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.statistics.distribution;
+
+import java.math.BigDecimal;
+import java.math.MathContext;
+import java.util.stream.Stream;
+import org.junit.jupiter.api.Assertions;
+import org.junit.jupiter.api.Test;
+import org.junit.jupiter.params.ParameterizedTest;
+import org.junit.jupiter.params.provider.Arguments;
+import org.junit.jupiter.params.provider.MethodSource;
+
+/**
+ * Test cases for {@link TrapezoidalDistribution}.
+ * Extends {@link BaseContinuousDistributionTest}. See javadoc of that class for details.
+ */
+class TrapezoidalDistributionTest extends BaseContinuousDistributionTest {
+ @Override
+ ContinuousDistribution makeDistribution(Object... parameters) {
+ final double a = (Double) parameters[0];
+ final double b = (Double) parameters[1];
+ final double c = (Double) parameters[2];
+ final double d = (Double) parameters[3];
+ return TrapezoidalDistribution.of(a, b, c, d);
+ }
+
+ @Override
+ Object[][] makeInvalidParameters() {
+ return new Object[][] {
+ {0.0, 0.0, 0.0, 0.0},
+ // 1.0, 2.0, 3.0, 4.0 is OK - move points to incorrect locations
+ {5.0, 2.0, 3.0, 4.0}, // a > d
+ {1.0, 5.0, 3.0, 4.0}, // b > d
+ {1.0, 2.0, 5.0, 4.0}, // c > d
+ {3.5, 2.0, 3.0, 4.0}, // a > c
+ {1.0, 3.5, 3.0, 4.0}, // b > c
+ {2.5, 2.0, 3.0, 4.0}, // a > b
+ {1.0, 2.0, 3.0, 0.0}, // d < a
+ {1.0, 2.0, 0.0, 4.0}, // c < a
+ {1.0, 0.0, 3.0, 4.0}, // b < a
+ {1.0, 2.0, 3.0, 1.5}, // d < b
+ {1.0, 2.0, 1.5, 4.0}, // c < b
+ {1.0, 2.0, 3.0, 2.5}, // d < c
+ };
+ }
+
+ @Override
+ String[] getParameterNames() {
+ return new String[] {"SupportLowerBound", "B", "C", "SupportUpperBound"};
+ }
+
+ @Override
+ protected double getRelativeTolerance() {
+ // Tolerance is 4.440892098500626E-15.
+ return 20 * RELATIVE_EPS;
+ }
+
+ //-------------------- Additional test cases -------------------------------
+
+ @ParameterizedTest
+ @MethodSource
+ void testAdditionalMoments(double a, double b, double c, double d, double mean, double var) {
+ final TrapezoidalDistribution dist = TrapezoidalDistribution.of(a, b, c, d);
+ final DoubleTolerance tol = DoubleTolerances.ulps(8);
+ TestUtils.assertEquals(mean, dist.getMean(), tol);
+ TestUtils.assertEquals(var, dist.getVariance(), tol);
+ }
+
+ static Stream<Arguments> testAdditionalMoments() {
+ return Stream.of(
+ // Computed using scipy.stats.trapezoid
+ // Up slope, then flat
+ Arguments.of(0, 0.1, 1, 1, 0.5245614035087719, 0.07562480763311791),
+ Arguments.of(0, 1e-3, 1, 1, 0.5002499583124894, 0.08325006249999839),
+ Arguments.of(0, 1e-6, 1, 1, 0.5000002499999582, 0.08333325000006259),
+ Arguments.of(0, 1e-9, 1, 1, 0.50000000025, 0.08333333324999997),
+ Arguments.of(0, 1e-12, 1, 1, 0.50000000000025, 0.08333333333324999),
+ Arguments.of(0, 1e-15, 1, 1, 0.5000000000000003, 0.0833333333333332),
+ Arguments.of(0, 0, 1, 1, 0.5, 0.08333333333333331),
+ // Flat, then down slope
+ Arguments.of(0, 0, 0.9, 1, 0.47543859649122816, 0.07562480763311777),
+ Arguments.of(0, 0, 0.999, 1, 0.49975004168751025, 0.08325006249999842),
+ Arguments.of(0, 0, 0.999999, 1, 0.4999997500000417, 0.08333325000006248),
+ Arguments.of(0, 0, 0.999999999, 1, 0.49999999975000003, 0.08333333325000003),
+ Arguments.of(0, 0, 0.999999999999, 1, 0.49999999999975003, 0.08333333333324999),
+ Arguments.of(0, 0, 0.999999999999999, 1, 0.4999999999999998, 0.08333333333333326),
+ Arguments.of(0, 0, 1, 1, 0.5, 0.08333333333333331)
+ );
+ }
+
+ /**
+ * Create a trapezoid with a very long upper tail to explicitly test the survival
+ * probability is high precision.
+ */
+ @Test
+ void testHighPrecisionSurvivalProbabilities() {
+ final double a = 0;
+ final double b = 0;
+ final double c = 1;
+ final double d = 1 << 14;
+ final double x1 = Math.nextDown(d);
+ final double x2 = Math.nextDown(x1);
+ final double p1 = survivalProbability(a, b, c, d, x1);
+ final double p2 = survivalProbability(a, b, c, d, x2);
+ final TrapezoidalDistribution dist = TrapezoidalDistribution.of(a, b, c, d);
+ final double[] points = {x1, x2};
+ final double[] probabilities = {p1, p2};
+
+ // This fails if the sf(x) = 1 - cdf(x)
+ testSurvivalProbabilityHighPrecision(
+ dist,
+ points,
+ probabilities,
+ DoubleTolerances.relative(1e-15));
+
+ // This fails if the isf(p) = icdf(1 - p)
+ testInverseSurvivalProbability(
+ dist,
+ probabilities,
+ points,
+ DoubleTolerances.ulps(0));
+ }
+
+ /**
+ * Compute the trapezoid distribution survival probability for the value {@code x}
+ * in the region {@code [c, d]}.
+ *
+ * @param a Lower limit of the distribution (inclusive).
+ * @param b Start of the trapezoid constant density.
+ * @param c End of the trapezoid constant density.
+ * @param d Upper limit of the distribution (inclusive).
+ * @param x Value in [c, d].
+ * @return the probability
+ */
+ private static double survivalProbability(double a, double b, double c, double d, double x) {
+ Assertions.assertTrue(c <= x && x <= d, "Domain error");
+ final BigDecimal aa = new BigDecimal(a);
+ final BigDecimal bb = new BigDecimal(b);
+ final BigDecimal cc = new BigDecimal(c);
+ final BigDecimal dd = new BigDecimal(d);
+ final BigDecimal divisor = dd.add(cc).subtract(aa).subtract(bb).multiply(dd.subtract(cc));
+ return dd.subtract(new BigDecimal(x)).pow(2)
+ .divide(divisor, MathContext.DECIMAL128).doubleValue();
+ }
+}
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.1.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.1.properties
new file mode 100644
index 0000000..4d8c7ff
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.1.properties
@@ -0,0 +1,43 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using R 3.4.4: require(trapezoid)
+# dtrapezoid(x, min=0, mode1=0.375, mode2=0.75, max=1)
+parameters = 0.0, 0.375, 0.75, 1.0
+# Computed manually using Maxima
+mean = 5.265151515151514916e-1
+variance = 4.7829143709825576357e-2
+lower = 0
+upper = 1
+
+cdf.points = 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
+
+cdf.values = \
+ 0.000000000000000000000 0.019393939393939393923 0.077575757575757575690 \
+ 0.174545454545454586937 0.309090909090909138346 0.454545454545454530315 \
+ 0.600000000000000088818 0.745454545454545591809 0.883636363636363664220 \
+ 0.970909090909090943811 1.000000000000000000000
+
+sf.values = \
+ 1.000000000000000000000 0.980606060606060592200 0.922424242424242368799 \
+ 0.825454545454545440819 0.690909090909090917165 0.545454545454545414174 \
+ 0.399999999999999911182 0.254545454545454408191 0.116363636363636335780 \
+ 0.029090909090909056189 0.000000000000000000000
+
+pdf.values = \
+ 0.00000000000000000000 0.38787878787878787845 0.77575757575757575690 \
+ 1.16363636363636380189 1.45454545454545458583 1.45454545454545458583 \
+ 1.45454545454545458583 1.45454545454545458583 1.16363636363636335780 \
+ 0.58181818181818167890 0.00000000000000000000
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.10.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.10.properties
new file mode 100644
index 0000000..ac4189b
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.10.properties
@@ -0,0 +1,29 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per triangular distribution (1, 1, 2)
+# From [test.triangular.4.properties]
+parameters = 1.0, 1.0, 1.0, 2.0
+# Computed manually
+mean = 1.3333333333333333
+variance = 0.05555555555555555
+lower = 1
+upper = 2
+
+cdf.points = 0, 0.5, 1, 1.25, 1.5, 1.75, 2, 3
+# cdf(x) = 1 - (2 - x)^2
+cdf.values = 0, 0, 0, 0.4375, 0.75, 0.9375, 1, 1
+# pdf(x) = 2 - 2 (x - 1)
+pdf.values = 0, 0, 2, 1.5, 1, 0.5, 0, 0
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.11.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.11.properties
new file mode 100644
index 0000000..bcfa2d2
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.11.properties
@@ -0,0 +1,55 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per uniform distribution (-0.5, 1.25)
+# From [test.uniformcontinuous.1.properties]
+parameters = -0.5 -0.5 1.25 1.25
+# Lower tolerance to pass survival function test.
+# scipy uniform does not appear to compute sf to high precision.
+tolerance.relative = 1e-14
+# probability(x0, x1) test expects this to be cdf(x1) - cdf(x0).
+# The computation is 0.5 ulp exact using (x1 - x0) / (upper - lower).
+# Configure the absolute error to allow this.
+tolerance.absolute = 5e-17
+
+# Computed using scipy.stats uniform(-0.5, 1.75)
+mean = 0.375
+variance = 0.2552083333333333
+lower = -0.5
+upper = 1.25
+cdf.points = \
+ -0.5001, -0.5, -0.4999, -0.25, -0.0001, 0.0,\
+ 0.0001, 0.25, 1.0, 1.2499, 1.25, 1.2501
+cdf.values = \
+ 0.0000000000000000e+00, 0.0000000000000000e+00,\
+ 5.7142857142850847e-05, 1.4285714285714285e-01,\
+ 2.8565714285714289e-01, 2.8571428571428570e-01,\
+ 2.8577142857142857e-01, 4.2857142857142855e-01,\
+ 8.5714285714285710e-01, 9.9994285714285713e-01,\
+ 1.0000000000000000e+00, 1.0000000000000000e+00
+pdf.values = \
+ 0. , 0.5714285714285714, 0.5714285714285714,\
+ 0.5714285714285714, 0.5714285714285714, 0.5714285714285714,\
+ 0.5714285714285714, 0.5714285714285714, 0.5714285714285714,\
+ 0.5714285714285714, 0.5714285714285714, 0.
+# 1-cdf is inaccurate.
+# Create sf(1.2499) manually: (1.25-1.2499) / (1.25- -0.5)
+sf.values = \
+ 1.000000000000000e+00, 1.000000000000000e+00,\
+ 9.999428571428571e-01, 8.571428571428572e-01,\
+ 7.143428571428572e-01, 7.142857142857143e-01,\
+ 7.142285714285714e-01, 5.714285714285714e-01,\
+ 1.428571428571429e-01, 5.7142857142850847e-05,\
+ 0.000000000000000e+00, 0.000000000000000e+00
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.12.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.12.properties
new file mode 100644
index 0000000..815c70d
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.12.properties
@@ -0,0 +1,52 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per uniform distribution (-4.6, 12.2)
+# From [test.uniformcontinuous.2.properties]
+parameters = -4.6 -4.6 12.2 12.2
+# Lower tolerance to pass survival function test.
+# scipy uniform does not appear to compute sf to high precision.
+tolerance.relative = 1e-14
+# Computed using scipy.stats uniform(-4.6, 16.8)
+mean = 3.8
+variance = 23.52
+lower = -4.6
+upper = 12.2
+cdf.points = \
+ -4.6, -4, -3, -1.5, 0, 1.1, 1.9, 2.7, 5.5, 9.1, 11.3, 12.2
+cdf.values = \
+ 0. , 0.03571428571428569, 0.09523809523809522, \
+ 0.1845238095238095 , 0.2738095238095238 , 0.33928571428571425, \
+ 0.38690476190476186, 0.4345238095238095 , 0.6011904761904762 , \
+ 0.8154761904761904 , 0.9464285714285714 , 1.0
+pdf.values = \
+ 0.05952380952380952, 0.05952380952380952, 0.05952380952380952, \
+ 0.05952380952380952, 0.05952380952380952, 0.05952380952380952, \
+ 0.05952380952380952, 0.05952380952380952, 0.05952380952380952, \
+ 0.05952380952380952, 0.05952380952380952, 0.05952380952380952
+# 1-cdf is inaccurate
+sf.values = \
+ 1.0000000000000000e+00, 9.6428571428571430e-01,\
+ 9.0476190476190477e-01, 8.1547619047619047e-01,\
+ 7.2619047619047628e-01, 6.6071428571428581e-01,\
+ 6.1309523809523814e-01, 5.6547619047619047e-01,\
+ 3.9880952380952384e-01, 1.8452380952380965e-01,\
+ 5.3571428571428603e-02, 0.0
+cdf.hp.points = -4.599999999999999
+cdf.hp.values = 5.2867763077388404e-17
+# Created manually as scipy.stats does not have this to high precision.
+# (12.2 - Math.nextDown(12.2)) / (12.2 - -4.6)
+sf.hp.points = 12.199999999999998
+sf.hp.values = 1.0573552615477683E-16
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.13.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.13.properties
new file mode 100644
index 0000000..93e64a6
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.13.properties
@@ -0,0 +1,41 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per uniform distribution (10.0, 100000000.0)
+# From [test.uniformcontinuous.3.properties]
+# Large range to require high precision survival function
+parameters = 10.0 10.0 100000000.0 100000000.0
+# Computed manually using JShell:
+# double u = 100000000.0
+# double l = 10.0
+# pdf = 1.0 / (u-l)
+# cdf = (x-l) / (u-l)
+# sf = (u-x) / (u-l)
+mean = 50000005.0
+variance = 833333166666675.0
+lower = 10
+upper = 100000000
+cdf.points = \
+ 11, 1000, 10000, 10000000
+cdf.values = \
+ 1.00000010000001e-08, 9.9000009900001e-06, 9.9900009990001e-05, 0.099999909999991
+pdf.values = \
+ 1.00000010000001e-08, 1.00000010000001e-08 \
+ 1.00000010000001e-08, 1.00000010000001e-08
+sf.values = \
+ 0.999999989999999, 0.99999009999901, 0.99990009999001, 0.900000090000009
+# (100000000.0 - Math.nextDown(100000000.0)) / (100000000.0 - 10)
+sf.hp.points = 9.999999999999999E7
+sf.hp.values = 1.4901162683963924e-16
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.2.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.2.properties
new file mode 100644
index 0000000..16e7604
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.2.properties
@@ -0,0 +1,50 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using R 3.4.4: require(trapezoid)
+# dtrapezoid(x, min=2, mode1=5, mode2=6, max=6)
+parameters = 2.0 5.0 6.0 6.0
+# Computed manually using scipy.stats trapezoid
+mean = 4.6
+variance = 0.8733333333333331
+lower = 2
+upper = 6
+
+cdf.points = 2.00 2.25 2.50 2.75 3.00 3.25 3.50 3.75 4.00 4.25 4.50 4.75 5.00 5.25 5.50 \
+ 5.75 6.00
+
+cdf.values = \
+ 0.0000000000000000000000 0.0041666666666666666088 0.0166666666666666664354 \
+ 0.0374999999999999986122 0.0666666666666666657415 0.1041666666666666712926 \
+ 0.1499999999999999944489 0.2041666666666666907215 0.2666666666666666629659 \
+ 0.3374999999999999666933 0.4166666666666666851704 0.5041666666666665408414 \
+ 0.5999999999999999777955 0.6999999999999999555911 0.8000000000000000444089 \
+ 0.9000000000000000222045 1.0000000000000000000000
+
+sf.values = \
+ 1.000000000000000000000 0.995833333333333348136 0.983333333333333281523 \
+ 0.962500000000000022204 0.933333333333333348136 0.895833333333333370341 \
+ 0.849999999999999977796 0.795833333333333281523 0.733333333333333392545 \
+ 0.662500000000000088818 0.583333333333333259318 0.495833333333333459159 \
+ 0.400000000000000022204 0.300000000000000044409 0.199999999999999955591 \
+ 0.099999999999999977796 0.000000000000000000000
+
+pdf.values = \
+ 0.000000000000000000000 0.033333333333333332871 0.066666666666666665741 \
+ 0.100000000000000005551 0.133333333333333331483 0.166666666666666685170 \
+ 0.200000000000000011102 0.233333333333333364790 0.266666666666666662966 \
+ 0.300000000000000044409 0.333333333333333370341 0.366666666666666696273 \
+ 0.400000000000000022204 0.400000000000000022204 0.400000000000000022204 \
+ 0.400000000000000022204 0.00000000000000000000
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.3.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.3.properties
new file mode 100644
index 0000000..6281740
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.3.properties
@@ -0,0 +1,53 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using R 3.4.4: require(trapezoid)
+# dtrapezoid(x, min=2, mode1=2, mode2=4, max=6)
+parameters = 2.0 2.0 4.0 6.0
+# Limited by SF
+tolerance.relative = 8e-15
+
+# Computed manually using scipy.stats trapezoid
+mean = 3.5555555555555554
+variance = 0.9135802469135812
+lower = 2
+upper = 6
+
+cdf.points = 2.00 2.25 2.50 2.75 3.00 3.25 3.50 3.75 4.00 4.25 4.50 4.75 5.00 5.25 5.50 \
+ 5.75 6.00
+
+cdf.values = \
+ 0.000000000000000000000 0.083333333333333328707 0.166666666666666657415 \
+ 0.250000000000000000000 0.333333333333333314830 0.416666666666666685170 \
+ 0.500000000000000000000 0.583333333333333370341 0.666666666666666740682 \
+ 0.744791666666666740682 0.812500000000000000000 0.869791666666666740682 \
+ 0.916666666666666629659 0.953125000000000000000 0.979166666666666629659 \
+ 0.994791666666666629659 1.000000000000000000000
+
+sf.values = \
+ 1.0000000000000000000000 0.9166666666666666296592 0.8333333333333333703408 \
+ 0.7500000000000000000000 0.6666666666666667406815 0.5833333333333332593185 \
+ 0.5000000000000000000000 0.4166666666666666296592 0.3333333333333332593185 \
+ 0.2552083333333332593185 0.1875000000000000000000 0.1302083333333332593185 \
+ 0.0833333333333333703408 0.0468750000000000000000 0.0208333333333333703408 \
+ 0.0052083333333333703408 0.0000000000000000000000
+
+pdf.values = \
+ 0.333333333333333314830 0.333333333333333314830 0.333333333333333314830 \
+ 0.333333333333333314830 0.333333333333333314830 0.333333333333333314830 \
+ 0.333333333333333314830 0.333333333333333314830 0.333333333333333314830 \
+ 0.291666666666666629659 0.250000000000000000000 0.208333333333333314830 \
+ 0.166666666666666657415 0.125000000000000000000 0.083333333333333328707 \
+ 0.041666666666666664354 0.000000000000000000000
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.4.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.4.properties
new file mode 100644
index 0000000..fe5869d
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.4.properties
@@ -0,0 +1,47 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using scipy.stats 1.7.3:
+# (a, b, c, d) = (-5, -1, 0, 2)
+# scipy.stats.trapezoid((b-a)/(d-a), (c-a)/(d-a), a, d-a)
+parameters = -5.0, -1.0, 0.0, 2.0
+mean = -1.1250000000000004
+variance = 2.1510416666666643
+lower = -5
+upper = 2
+
+cdf.points = -5. , -4.5, -4. , -3.5, -3. , -2.5, -2. , -1.5, -1. , -0.5, 0. , \
+ 0.5, 1. , 1.5, 2.
+
+cdf.values = \
+ 0. , 0.0078125 , 0.03125 , \
+ 0.0703125 , 0.125 , 0.19531250000000003, \
+ 0.28125 , 0.3828125 , 0.5 , \
+ 0.6250000000000001 , 0.7500000000000001 , 0.859375 , \
+ 0.9375 , 0.984375 , 1.
+
+sf.values = \
+ 1. , 0.9921875 , 0.96875 , \
+ 0.9296875 , 0.875 , 0.8046875 , \
+ 0.71875 , 0.6171875 , 0.5 , \
+ 0.3749999999999999, 0.2499999999999999, 0.140625 , \
+ 0.0625 , 0.015625 , 0.
+
+pdf.values = \
+ 0. , 0.03125 , 0.0625 , \
+ 0.09375 , 0.125 , 0.15625 , \
+ 0.1875 , 0.21875 , 0.25 , \
+ 0.25 , 0.25 , 0.1875 , \
+ 0.12500000000000006, 0.06249999999999998, 0.
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.5.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.5.properties
new file mode 100644
index 0000000..f8b39d4
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.5.properties
@@ -0,0 +1,46 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using scipy.stats 1.7.3:
+# (a, b, c, d) = (3, 4, 5, 6)
+# scipy.stats.trapezoid((b-a)/(d-a), (c-a)/(d-a), a, d-a)
+parameters = 3.0, 4.0, 5.0, 6.0
+mean = 4.5
+variance = 0.4166666666666675
+lower = 3
+upper = 6
+
+cdf.points = 3. , 3.25, 3.5 , 3.75, 4. , 4.25, 4.5 , 4.75, 5. , 5.25, 5.5 , 5.75, 6.
+
+cdf.values = \
+ 0. , 0.015625 , 0.0625 , \
+ 0.140625 , 0.25 , 0.375 , \
+ 0.5000000000000001, 0.6250000000000001, 0.75 , \
+ 0.859375 , 0.9375 , 0.984375 , \
+ 1.
+
+sf.values = \
+ 1. , 0.984375 , 0.9375 , \
+ 0.859375 , 0.75 , 0.625 , \
+ 0.4999999999999999, 0.3749999999999999, 0.25 , \
+ 0.140625 , 0.0625 , 0.015625 , \
+ 0.
+
+pdf.values = \
+ 0. , 0.125 , 0.25 , \
+ 0.375 , 0.5 , 0.5 , \
+ 0.5 , 0.5 , 0.5 , \
+ 0.37499999999999994, 0.24999999999999992, 0.12500000000000003, \
+ 0.
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.6.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.6.properties
new file mode 100644
index 0000000..f982291
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.6.properties
@@ -0,0 +1,53 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# computed using scipy.stats 1.7.3:
+# (a, b, c, d) = (-5, -3.5, 2, 4)
+# scipy.stats.trapezoid((b-a)/(d-a), (c-a)/(d-a), a, d-a)
+parameters = -5.0, -3.5, 2.0, 4.0
+mean = -0.6149425287356323
+variance = 4.6405238472717665
+lower = -5
+upper = 4
+
+cdf.points = -5. , -4.5, -4. , -3.5, -3. , -2.5, -2. , -1.5, -1. , -0.5, 0. , \
+ 0.5, 1. , 1.5, 2. , 2.5, 3. , 3.5, 4.
+
+cdf.values = \
+ 0. , 0.011494252873563218, 0.04597701149425287 , \
+ 0.10344827586206895 , 0.1724137931034483 , 0.24137931034482762 , \
+ 0.3103448275862069 , 0.37931034482758624 , 0.44827586206896547 , \
+ 0.5172413793103449 , 0.5862068965517242 , 0.6551724137931035 , \
+ 0.7241379310344828 , 0.7931034482758621 , 0.8620689655172414 , \
+ 0.9224137931034483 , 0.9655172413793103 , 0.9913793103448276 , \
+ 1.
+
+sf.values = \
+ 1. , 0.9885057471264368 , 0.9540229885057472 , \
+ 0.896551724137931 , 0.8275862068965517 , 0.7586206896551724 , \
+ 0.6896551724137931 , 0.6206896551724137 , 0.5517241379310345 , \
+ 0.48275862068965514 , 0.4137931034482758 , 0.34482758620689646 , \
+ 0.27586206896551724 , 0.2068965517241379 , 0.13793103448275856 , \
+ 0.07758620689655171 , 0.034482758620689724, 0.008620689655172376, \
+ 0.
+
+pdf.values = \
+ 0. , 0.04597701149425287, 0.09195402298850575, \
+ 0.13793103448275862, 0.13793103448275862, 0.13793103448275862, \
+ 0.13793103448275862, 0.13793103448275862, 0.13793103448275862, \
+ 0.13793103448275862, 0.13793103448275862, 0.13793103448275862, \
+ 0.13793103448275862, 0.13793103448275862, 0.13793103448275862, \
+ 0.10344827586206895, 0.06896551724137934, 0.03448275862068967, \
+ 0.
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.7.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.7.properties
new file mode 100644
index 0000000..36c25a6
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.7.properties
@@ -0,0 +1,44 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per triangular distribution (-3, 2, 12)
+# From [test.triangular.1.properties]
+parameters = -3.0, 2.0, 2.0, 12.0
+# Limited by SF
+tolerance.relative = 6e-14
+
+mean = 3.6666666666666665
+variance = 9.722222222222221
+lower = -3
+upper = 12
+cdf.points = \
+ -3.0, \
+ -2.0, -1.0, 0.0, 1.0,\
+ 2.0,\
+ 3.0, 4.0, 10.0, 11.0,\
+ 12.0
+cdf.values = \
+ 0.0,\
+ 0.013333333333333, 0.053333333333333, 0.12, 0.213333333333333,\
+ 0.333333333333333,\
+ 0.46, 0.573333333333333,\
+ 0.973333333333333, 0.993333333333333,\
+ 1.0,
+pdf.values = \
+ 0,\
+ 0.02666666666666667, 0.05333333333333334, 0.08, 0.10666666666666667,\
+ 0.13333333333333333,\
+ 0.12, 0.10666666666666667, 0.02666666666666667, 0.013333333333333334,\
+ 0,
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.8.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.8.properties
new file mode 100644
index 0000000..50a08ef
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.8.properties
@@ -0,0 +1,44 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per triangular distribution (1, 2, 5)
+# From [test.triangular.2.properties]
+parameters = 1.0 2.0 2.0 5.0
+# Limited by SF
+tolerance.relative = 1e-14
+
+mean = 2.6666666666666667
+variance = 0.722222222222222
+lower = 1
+upper = 5
+cdf.points = \
+ 1.0, 1.25, 1.5, 1.75, \
+ 2.0, 2.25, 2.5, 2.75, \
+ 3.0, 3.25, 3.5, 3.75, \
+ 4.0, 4.25, 4.5, 4.75, 5.0
+cdf.values = \
+ 0. , 0.015625 , 0.0625 , \
+ 0.140625 , 0.25 , 0.3697916666666667, \
+ 0.4791666666666667, 0.578125 , 0.6666666666666666, \
+ 0.7447916666666666, 0.8125 , 0.8697916666666666, \
+ 0.9166666666666666, 0.953125 , 0.9791666666666666, \
+ 0.9947916666666666, 1.
+pdf.values = \
+ 0. , 0.125 , 0.25 , \
+ 0.375 , 0.5 , 0.4583333333333333 , \
+ 0.4166666666666667 , 0.375 , 0.3333333333333333 , \
+ 0.2916666666666667 , 0.25 , 0.20833333333333334 , \
+ 0.16666666666666666 , 0.125 , 0.08333333333333333 , \
+ 0.041666666666666664, 0.
diff --git a/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.9.properties b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.9.properties
new file mode 100644
index 0000000..bd97975
--- /dev/null
+++ b/commons-statistics-distribution/src/test/resources/org/apache/commons/statistics/distribution/test.trapezoidal.9.properties
@@ -0,0 +1,29 @@
+# Licensed to the Apache Software Foundation (ASF) under one or more
+# contributor license agreements. See the NOTICE file distributed with
+# this work for additional information regarding copyright ownership.
+# The ASF licenses this file to You under the Apache License, Version 2.0
+# (the "License"); you may not use this file except in compliance with
+# the License. You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+# As per triangular distribution (1, 2, 2)
+# From [test.triangular.3.properties]
+parameters = 1.0, 2.0, 2.0, 2.0
+# Computed manually
+mean = 1.6666666666666667
+variance = 0.05555555555555555
+lower = 1
+upper = 2
+
+cdf.points = 0, 0.5, 1, 1.25, 1.5, 1.75, 2, 3
+# cdf(x) = (x - 1)^2
+cdf.values = 0, 0, 0, 0.0625, 0.25, 0.5625, 1, 1
+# pdf(x) = 2 (x - 1)
+pdf.values = 0, 0, 0, 0.5, 1, 1.5, 2, 0
diff --git a/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/DistributionsApplication.java b/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/DistributionsApplication.java
index c1f003e..758d483 100644
--- a/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/DistributionsApplication.java
+++ b/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/DistributionsApplication.java
@@ -88,6 +88,7 @@ public final class DistributionsApplication {
PascalCommand.class,
PoissonCommand.class,
TCommand.class,
+ TrapezoidalCommand.class,
TriangularCommand.class,
TruncatedNormalCommand.class,
UniformContinuousCommand.class,
diff --git a/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/TrapezoidalCommand.java b/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/TrapezoidalCommand.java
new file mode 100644
index 0000000..25422d2
--- /dev/null
+++ b/commons-statistics-examples/examples-distribution/src/main/java/org/apache/commons/statistics/examples/distribution/TrapezoidalCommand.java
@@ -0,0 +1,176 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.statistics.examples.distribution;
+
+import java.util.ArrayList;
+import java.util.List;
+import org.apache.commons.statistics.distribution.ContinuousDistribution;
+import org.apache.commons.statistics.distribution.TrapezoidalDistribution;
+import picocli.CommandLine.ArgGroup;
+import picocli.CommandLine.Command;
+import picocli.CommandLine.Option;
+
+/**
+ * Command for the {@link TrapezoidalDistribution}.
+ */
+@Command(name = "trapezoidal",
+ aliases = {"trap"},
+ description = "Trapezoidal distribution.",
+ subcommands = {
+ TrapezoidalCommand.Check.class,
+ TrapezoidalCommand.PDF.class,
+ TrapezoidalCommand.LPDF.class,
+ TrapezoidalCommand.CDF.class,
+ TrapezoidalCommand.SF.class,
+ TrapezoidalCommand.ICDF.class,
+ TrapezoidalCommand.ISF.class,
+ })
+class TrapezoidalCommand extends AbstractDistributionCommand {
+
+ /** Base command for the distribution that defines the parameters. */
+ private abstract static class BaseCommand extends AbstractContinuousDistributionCommand {
+ /** Distribution parameters. */
+ @ArgGroup(validate = false, heading = "Distribution parameters:%n", order = 1)
+ private Params params = new Params();
+
+ /** Parameters class. */
+ static class Params {
+ /** The distribution lower limit. */
+ @Option(names = {"-a", "--lower"},
+ paramLabel = "a",
+ arity = "1..*",
+ split = ",",
+ description = {"lower bound (default: ${DEFAULT-VALUE})."})
+ private double[] lower = {-4, 0, 0};
+
+ /** The start of the density constant region. */
+ @Option(names = {"-b", "--shape1"},
+ paramLabel = "b",
+ arity = "1..*",
+ split = ",",
+ description = {"constant start (default: ${DEFAULT-VALUE})."})
+ private double[] start = {-1.5, 0, 1};
+
+ /** The end of the density constant region. */
+ @Option(names = {"-c", "--shape2"},
+ paramLabel = "c",
+ arity = "1..*",
+ split = ",",
+ description = {"constant end (default: ${DEFAULT-VALUE})."})
+ private double[] end = {2.5, 2, 3};
+
+ /** The distribution upper limit. */
+ @Option(names = {"-d", "--upper"},
+ paramLabel = "d",
+ arity = "1..*",
+ split = ",",
+ description = {"upper bound (default: ${DEFAULT-VALUE})."})
+ private double[] upper = {4, 3, 3};
+ }
+
+ /** Extend the options to set the default values for this distribution. */
+ static final class Options extends ContinuousDistributionOptions {
+ /** Set defaults. */
+ private Options() {
+ min = -5;
+ max = 5;
+ }
+ }
+
+ @Override
+ protected List<Distribution<ContinuousDistribution>> getDistributions() {
+ double[] a = params.lower;
+ double[] b = params.start;
+ double[] c = params.end;
+ double[] d = params.upper;
+ final int n = DistributionUtils.validateLengths(a.length, b.length, c.length, d.length);
+
+ a = DistributionUtils.expandToLength(a, n);
+ b = DistributionUtils.expandToLength(b, n);
+ c = DistributionUtils.expandToLength(c, n);
+ d = DistributionUtils.expandToLength(d, n);
+
+ // Create distributions
+ final ArrayList<Distribution<ContinuousDistribution>> list = new ArrayList<>();
+ for (int i = 0; i < n; i++) {
+ final ContinuousDistribution dist = TrapezoidalDistribution.of(a[i], b[i], c[i], d[i]);
+ list.add(new Distribution<>(dist, "a=" + a[i] + ",b=" + b[i] + ",c=" + c[i] + ",d=" + d[i]));
+ }
+ return list;
+ }
+ }
+
+ /** Base command for the distribution that defines the parameters. */
+ private abstract static class ProbabilityCommand extends BaseCommand {
+ /** The distribution options. */
+ @ArgGroup(validate = false, heading = "Evaluation options:%n", order = 2)
+ private Options distributionOptions = new Options();
+
+ @Override
+ protected DistributionOptions getDistributionOptions() {
+ return distributionOptions;
+ }
+ }
+
+ /** Base command for the distribution that defines the parameters for inverse probability functions. */
+ private abstract static class InverseProbabilityCommand extends BaseCommand {
+ /** The distribution options. */
+ @ArgGroup(validate = false, heading = "Evaluation options:%n", order = 2)
+ private InverseContinuousDistributionOptions distributionOptions = new InverseContinuousDistributionOptions();
+
+ @Override
+ protected DistributionOptions getDistributionOptions() {
+ return distributionOptions;
+ }
+ }
+
+ /** Verification checks command. */
+ @Command(name = "check",
+ hidden = true,
+ description = "Trapezoidal distribution verification checks.")
+ static class Check extends ProbabilityCommand {}
+
+ /** PDF command. */
+ @Command(name = "pdf",
+ description = "Trapezoidal distribution PDF.")
+ static class PDF extends ProbabilityCommand {}
+
+ /** LPDF command. */
+ @Command(name = "lpdf",
+ description = "Trapezoidal distribution natural logarithm of the PDF.")
+ static class LPDF extends ProbabilityCommand {}
+
+ /** CDF command. */
+ @Command(name = "cdf",
+ description = "Trapezoidal distribution CDF.")
+ static class CDF extends ProbabilityCommand {}
+
+ /** SF command. */
+ @Command(name = "sf",
+ description = "Trapezoidal distribution survival probability.")
+ static class SF extends ProbabilityCommand {}
+
+ /** ICDF command. */
+ @Command(name = "icdf",
+ description = "Trapezoidal distribution inverse CDF.")
+ static class ICDF extends InverseProbabilityCommand {}
+
+ /** ISF command. */
+ @Command(name = "isf",
+ description = "Trapezoidal distribution inverse SF.")
+ static class ISF extends InverseProbabilityCommand {}
+}