001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017 018package org.apache.commons.statistics.distribution; 019 020/** 021 * Implementation of the triangular distribution. 022 * 023 * <p>The probability density function of \( X \) is: 024 * 025 * <p>\[ f(x; a, b, c) = \begin{cases} 026 * \frac{2(x-a)}{(b-a)(c-a)} & \text{for } a \le x \lt c \\ 027 * \frac{2}{b-a} & \text{for } x = c \\ 028 * \frac{2(b-x)}{(b-a)(b-c)} & \text{for } c \lt x \le b \\ 029 * \end{cases} \] 030 * 031 * <p>for \( -\infty \lt a \le c \le b \lt \infty \) and 032 * \( x \in [a, b] \). 033 * 034 * @see <a href="https://en.wikipedia.org/wiki/Triangular_distribution">Triangular distribution (Wikipedia)</a> 035 * @see <a href="https://mathworld.wolfram.com/TriangularDistribution.html">Triangular distribution (MathWorld)</a> 036 */ 037public final class TriangularDistribution extends AbstractContinuousDistribution { 038 /** Lower limit of this distribution (inclusive). */ 039 private final double a; 040 /** Upper limit of this distribution (inclusive). */ 041 private final double b; 042 /** Mode of this distribution. */ 043 private final double c; 044 /** Cached value ((b - a) * (c - a). */ 045 private final double divisor1; 046 /** Cached value ((b - a) * (b - c)). */ 047 private final double divisor2; 048 /** Cumulative probability at the mode. */ 049 private final double cdfMode; 050 /** Survival probability at the mode. */ 051 private final double sfMode; 052 053 /** 054 * @param a Lower limit of this distribution (inclusive). 055 * @param c Mode of this distribution. 056 * @param b Upper limit of this distribution (inclusive). 057 */ 058 private TriangularDistribution(double a, 059 double c, 060 double b) { 061 this.a = a; 062 this.c = c; 063 this.b = b; 064 divisor1 = (b - a) * (c - a); 065 divisor2 = (b - a) * (b - c); 066 cdfMode = (c - a) / (b - a); 067 sfMode = (b - c) / (b - a); 068 } 069 070 /** 071 * Creates a triangular distribution. 072 * 073 * @param a Lower limit of this distribution (inclusive). 074 * @param c Mode of this distribution. 075 * @param b Upper limit of this distribution (inclusive). 076 * @return the distribution 077 * @throws IllegalArgumentException if {@code a >= b}, if {@code c > b} or if 078 * {@code c < a}. 079 */ 080 public static TriangularDistribution of(double a, 081 double c, 082 double b) { 083 if (a >= b) { 084 throw new DistributionException(DistributionException.INVALID_RANGE_LOW_GTE_HIGH, 085 a, b); 086 } 087 if (c < a) { 088 throw new DistributionException(DistributionException.TOO_SMALL, 089 c, a); 090 } 091 if (c > b) { 092 throw new DistributionException(DistributionException.TOO_LARGE, 093 c, b); 094 } 095 return new TriangularDistribution(a, c, b); 096 } 097 098 /** 099 * Gets the mode parameter of this distribution. 100 * 101 * @return the mode. 102 */ 103 public double getMode() { 104 return c; 105 } 106 107 /** {@inheritDoc} */ 108 @Override 109 public double density(double x) { 110 if (x < a) { 111 return 0; 112 } 113 if (x < c) { 114 final double divident = 2 * (x - a); 115 return divident / divisor1; 116 } 117 if (x == c) { 118 return 2 / (b - a); 119 } 120 if (x <= b) { 121 final double divident = 2 * (b - x); 122 return divident / divisor2; 123 } 124 return 0; 125 } 126 127 /** {@inheritDoc} */ 128 @Override 129 public double cumulativeProbability(double x) { 130 if (x <= a) { 131 return 0; 132 } 133 if (x < c) { 134 final double divident = (x - a) * (x - a); 135 return divident / divisor1; 136 } 137 if (x == c) { 138 return cdfMode; 139 } 140 if (x < b) { 141 final double divident = (b - x) * (b - x); 142 return 1 - (divident / divisor2); 143 } 144 return 1; 145 } 146 147 148 /** {@inheritDoc} */ 149 @Override 150 public double survivalProbability(double x) { 151 // By symmetry: 152 if (x <= a) { 153 return 1; 154 } 155 if (x < c) { 156 final double divident = (x - a) * (x - a); 157 return 1 - (divident / divisor1); 158 } 159 if (x == c) { 160 return sfMode; 161 } 162 if (x < b) { 163 final double divident = (b - x) * (b - x); 164 return divident / divisor2; 165 } 166 return 0; 167 } 168 169 /** {@inheritDoc} */ 170 @Override 171 public double inverseCumulativeProbability(double p) { 172 ArgumentUtils.checkProbability(p); 173 if (p == 0) { 174 return a; 175 } 176 if (p == 1) { 177 return b; 178 } 179 if (p < cdfMode) { 180 return a + Math.sqrt(p * divisor1); 181 } 182 return b - Math.sqrt((1 - p) * divisor2); 183 } 184 185 /** {@inheritDoc} */ 186 @Override 187 public double inverseSurvivalProbability(double p) { 188 // By symmetry: 189 ArgumentUtils.checkProbability(p); 190 if (p == 1) { 191 return a; 192 } 193 if (p == 0) { 194 return b; 195 } 196 if (p >= sfMode) { 197 return a + Math.sqrt((1 - p) * divisor1); 198 } 199 return b - Math.sqrt(p * divisor2); 200 } 201 202 /** 203 * {@inheritDoc} 204 * 205 * <p>For lower limit \( a \), upper limit \( b \), and mode \( c \), 206 * the mean is \( (a + b + c) / 3 \). 207 */ 208 @Override 209 public double getMean() { 210 return (a + b + c) / 3; 211 } 212 213 /** 214 * {@inheritDoc} 215 * 216 * <p>For lower limit \( a \), upper limit \( b \), and mode \( c \), 217 * the variance is \( (a^2 + b^2 + c^2 - ab - ac - bc) / 18 \). 218 */ 219 @Override 220 public double getVariance() { 221 return (a * a + b * b + c * c - a * b - a * c - b * c) / 18; 222 } 223 224 /** 225 * {@inheritDoc} 226 * 227 * <p>The lower bound of the support is equal to the lower limit parameter 228 * {@code a} of the distribution. 229 */ 230 @Override 231 public double getSupportLowerBound() { 232 return a; 233 } 234 235 /** 236 * {@inheritDoc} 237 * 238 * <p>The upper bound of the support is equal to the upper limit parameter 239 * {@code b} of the distribution. 240 */ 241 @Override 242 public double getSupportUpperBound() { 243 return b; 244 } 245}