ILIAS  release_5-4 Revision v5.4.26-12-gabc799a52e6
Beta.php
Go to the documentation of this file.
1<?php
2
4
7
8class Beta
9{
10 private const MAX_ITERATIONS = 256;
11
12 private const LOG_GAMMA_X_MAX_VALUE = 2.55e305;
13
14 private const XMININ = 2.23e-308;
15
29 public static function distribution($value, $alpha, $beta, $rMin = 0.0, $rMax = 1.0)
30 {
31 $value = Functions::flattenSingleValue($value);
32 $alpha = Functions::flattenSingleValue($alpha);
33 $beta = Functions::flattenSingleValue($beta);
34 $rMin = ($rMin === null) ? 0.0 : Functions::flattenSingleValue($rMin);
35 $rMax = ($rMax === null) ? 1.0 : Functions::flattenSingleValue($rMax);
36
37 try {
43 } catch (Exception $e) {
44 return $e->getMessage();
45 }
46
47 if ($rMin > $rMax) {
48 $tmp = $rMin;
49 $rMin = $rMax;
50 $rMax = $tmp;
51 }
52 if (($value < $rMin) || ($value > $rMax) || ($alpha <= 0) || ($beta <= 0) || ($rMin == $rMax)) {
53 return Functions::NAN();
54 }
55
56 $value -= $rMin;
57 $value /= ($rMax - $rMin);
58
59 return self::incompleteBeta($value, $alpha, $beta);
60 }
61
75 public static function inverse($probability, $alpha, $beta, $rMin = 0.0, $rMax = 1.0)
76 {
77 $probability = Functions::flattenSingleValue($probability);
78 $alpha = Functions::flattenSingleValue($alpha);
79 $beta = Functions::flattenSingleValue($beta);
80 $rMin = ($rMin === null) ? 0.0 : Functions::flattenSingleValue($rMin);
81 $rMax = ($rMax === null) ? 1.0 : Functions::flattenSingleValue($rMax);
82
83 try {
84 $probability = DistributionValidations::validateProbability($probability);
89 } catch (Exception $e) {
90 return $e->getMessage();
91 }
92
93 if ($rMin > $rMax) {
94 $tmp = $rMin;
95 $rMin = $rMax;
96 $rMax = $tmp;
97 }
98 if (($alpha <= 0) || ($beta <= 0) || ($rMin == $rMax) || ($probability <= 0.0)) {
99 return Functions::NAN();
100 }
101
102 return self::calculateInverse($probability, $alpha, $beta, $rMin, $rMax);
103 }
104
108 private static function calculateInverse(float $probability, float $alpha, float $beta, float $rMin, float $rMax)
109 {
110 $a = 0;
111 $b = 2;
112
113 $i = 0;
114 while ((($b - $a) > Functions::PRECISION) && (++$i <= self::MAX_ITERATIONS)) {
115 $guess = ($a + $b) / 2;
116 $result = self::distribution($guess, $alpha, $beta);
117 if (($result === $probability) || ($result === 0.0)) {
118 $b = $a;
119 } elseif ($result > $probability) {
120 $b = $guess;
121 } else {
122 $a = $guess;
123 }
124 }
125
126 if ($i === self::MAX_ITERATIONS) {
127 return Functions::NA();
128 }
129
130 return round($rMin + $guess * ($rMax - $rMin), 12);
131 }
132
147 public static function incompleteBeta(float $x, float $p, float $q): float
148 {
149 if ($x <= 0.0) {
150 return 0.0;
151 } elseif ($x >= 1.0) {
152 return 1.0;
153 } elseif (($p <= 0.0) || ($q <= 0.0) || (($p + $q) > self::LOG_GAMMA_X_MAX_VALUE)) {
154 return 0.0;
155 }
156
157 $beta_gam = exp((0 - self::logBeta($p, $q)) + $p * log($x) + $q * log(1.0 - $x));
158 if ($x < ($p + 1.0) / ($p + $q + 2.0)) {
159 return $beta_gam * self::betaFraction($x, $p, $q) / $p;
160 }
161
162 return 1.0 - ($beta_gam * self::betaFraction(1 - $x, $q, $p) / $q);
163 }
164
165 // Function cache for logBeta function
166 private static $logBetaCacheP = 0.0;
167
168 private static $logBetaCacheQ = 0.0;
169
170 private static $logBetaCacheResult = 0.0;
171
182 private static function logBeta(float $p, float $q): float
183 {
184 if ($p != self::$logBetaCacheP || $q != self::$logBetaCacheQ) {
185 self::$logBetaCacheP = $p;
186 self::$logBetaCacheQ = $q;
187 if (($p <= 0.0) || ($q <= 0.0) || (($p + $q) > self::LOG_GAMMA_X_MAX_VALUE)) {
188 self::$logBetaCacheResult = 0.0;
189 } else {
190 self::$logBetaCacheResult = Gamma::logGamma($p) + Gamma::logGamma($q) - Gamma::logGamma($p + $q);
191 }
192 }
193
195 }
196
203 private static function betaFraction(float $x, float $p, float $q): float
204 {
205 $c = 1.0;
206 $sum_pq = $p + $q;
207 $p_plus = $p + 1.0;
208 $p_minus = $p - 1.0;
209 $h = 1.0 - $sum_pq * $x / $p_plus;
210 if (abs($h) < self::XMININ) {
212 }
213 $h = 1.0 / $h;
214 $frac = $h;
215 $m = 1;
216 $delta = 0.0;
217 while ($m <= self::MAX_ITERATIONS && abs($delta - 1.0) > Functions::PRECISION) {
218 $m2 = 2 * $m;
219 // even index for d
220 $d = $m * ($q - $m) * $x / (($p_minus + $m2) * ($p + $m2));
221 $h = 1.0 + $d * $h;
222 if (abs($h) < self::XMININ) {
224 }
225 $h = 1.0 / $h;
226 $c = 1.0 + $d / $c;
227 if (abs($c) < self::XMININ) {
229 }
230 $frac *= $h * $c;
231 // odd index for d
232 $d = -($p + $m) * ($sum_pq + $m) * $x / (($p + $m2) * ($p_plus + $m2));
233 $h = 1.0 + $d * $h;
234 if (abs($h) < self::XMININ) {
236 }
237 $h = 1.0 / $h;
238 $c = 1.0 + $d / $c;
239 if (abs($c) < self::XMININ) {
241 }
242 $delta = $h * $c;
243 $frac *= $delta;
244 ++$m;
245 }
246
247 return $frac;
248 }
249
250 private static function betaValue(float $a, float $b): float
251 {
252 return (Gamma::gammaValue($a) * Gamma::gammaValue($b)) /
253 Gamma::gammaValue($a + $b);
254 }
255
256 private static function regularizedIncompleteBeta(float $value, float $a, float $b): float
257 {
258 return self::incompleteBeta($value, $a, $b) / self::betaValue($a, $b);
259 }
260}
$result
An exception for terminatinating execution or to throw for unit testing.
static flattenSingleValue($value='')
Convert an array to a single scalar value by extracting the first element.
Definition: Functions.php:649
static regularizedIncompleteBeta(float $value, float $a, float $b)
Definition: Beta.php:256
static distribution($value, $alpha, $beta, $rMin=0.0, $rMax=1.0)
BETADIST.
Definition: Beta.php:29
static calculateInverse(float $probability, float $alpha, float $beta, float $rMin, float $rMax)
Definition: Beta.php:108
static inverse($probability, $alpha, $beta, $rMin=0.0, $rMax=1.0)
BETAINV.
Definition: Beta.php:75
static logBeta(float $p, float $q)
The natural logarithm of the beta function.
Definition: Beta.php:182
static incompleteBeta(float $x, float $p, float $q)
Incomplete beta function.
Definition: Beta.php:147
static betaFraction(float $x, float $p, float $q)
Evaluates of continued fraction part of incomplete beta function.
Definition: Beta.php:203
$x
Definition: complexTest.php:9
for( $i=6;$i< 13;$i++) for($i=1; $i< 13; $i++) $d
Definition: date.php:296
$i
Definition: disco.tpl.php:19
$h