diff --git a/CHANGELOG.md b/CHANGELOG.md index d7b123ef..eb3ab111 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,7 +9,7 @@ and this project adheres to [Semantic Versioning](https://semver.org). ### Added -- Implemented the CHITEST() and CHISQ.DIST() Statistical function. +- Implemented the CHITEST(), CHISQ.DIST() and CHISQ.INV() and equivalent Statistical functions, for both left- and right-tailed distributions. - Support for ActiveSheet and SelectedCells in the ODS Reader and Writer. [PR #1908](https://github.com/PHPOffice/PhpSpreadsheet/pull/1908) ### Changed diff --git a/src/PhpSpreadsheet/Calculation/Calculation.php b/src/PhpSpreadsheet/Calculation/Calculation.php index c4a88bc3..4dfcf9ab 100644 --- a/src/PhpSpreadsheet/Calculation/Calculation.php +++ b/src/PhpSpreadsheet/Calculation/Calculation.php @@ -508,7 +508,7 @@ class Calculation ], 'CHISQ.INV' => [ 'category' => Category::CATEGORY_STATISTICAL, - 'functionCall' => [Functions::class, 'DUMMY'], + 'functionCall' => [Statistical\Distributions\ChiSquared::class, 'inverseLeftTail'], 'argumentCount' => '2', ], 'CHISQ.INV.RT' => [ diff --git a/src/PhpSpreadsheet/Calculation/Statistical/Distributions/ChiSquared.php b/src/PhpSpreadsheet/Calculation/Statistical/Distributions/ChiSquared.php index 409e5883..df3451cd 100644 --- a/src/PhpSpreadsheet/Calculation/Statistical/Distributions/ChiSquared.php +++ b/src/PhpSpreadsheet/Calculation/Statistical/Distributions/ChiSquared.php @@ -11,6 +11,8 @@ class ChiSquared private const MAX_ITERATIONS = 256; + private const EPS = 2.22e-16; + /** * CHIDIST. * @@ -127,6 +129,35 @@ class ChiSquared return $newtonRaphson->execute($probability); } + /** + * CHIINV. + * + * Returns the inverse of the left-tailed probability of the chi-squared distribution. + * + * @param mixed (float) $probability Probability for the function + * @param mixed (int) $degrees degrees of freedom + * + * @return float|string + */ + public static function inverseLeftTail($probability, $degrees) + { + $probability = Functions::flattenSingleValue($probability); + $degrees = Functions::flattenSingleValue($degrees); + + try { + $probability = self::validateFloat($probability); + $degrees = self::validateInt($degrees); + } catch (Exception $e) { + return $e->getMessage(); + } + + if ($probability < 0.0 || $probability > 1.0 || $degrees < 1) { + return Functions::NAN(); + } + + return self::inverseLeftTailCalculation($probability, $degrees); + } + public static function test($actual, $expected) { $rows = count($actual); @@ -167,4 +198,104 @@ class ChiSquared return ($columns - 1) * ($rows - 1); } + + private static function inverseLeftTailCalculation($probability, $degrees) + { + // bracket the root + $min = 0; + $sd = sqrt(2.0 * $degrees); + $max = 2 * $sd; + $s = -1; + + while ($s * self::pchisq($max, $degrees) > $probability * $s) { + $min = $max; + $max += 2 * $sd; + } + + // Find root using bisection + $chi2 = 0.5 * ($min + $max); + + while (($max - $min) > self::EPS * $chi2) { + if ($s * self::pchisq($chi2, $degrees) > $probability * $s) { + $min = $chi2; + } else { + $max = $chi2; + } + $chi2 = 0.5 * ($min + $max); + } + + return $chi2; + } + + private static function pchisq($chi2, $degrees) + { + return self::gammp($degrees, 0.5 * $chi2); + } + + private static function gammp($n, $x) + { + if ($x < 0.5 * $n + 1) { + return self::gser($n, $x); + } + + return 1 - self::gcf($n, $x); + } + + // Return the incomplete gamma function P(n/2,x) evaluated by + // series representation. Algorithm from numerical recipe. + // Assume that n is a positive integer and x>0, won't check arguments. + // Relative error controlled by the eps parameter + private static function gser($n, $x) + { + $gln = Gamma::ln($n / 2); + $a = 0.5 * $n; + $ap = $a; + $sum = 1.0 / $a; + $del = $sum; + for ($i = 1; $i < 101; ++$i) { + ++$ap; + $del = $del * $x / $ap; + $sum += $del; + if ($del < $sum * self::EPS) { + break; + } + } + + return $sum * exp(-$x + $a * log($x) - $gln); + } + + // Return the incomplete gamma function Q(n/2,x) evaluated by + // its continued fraction representation. Algorithm from numerical recipe. + // Assume that n is a postive integer and x>0, won't check arguments. + // Relative error controlled by the eps parameter + private static function gcf($n, $x) + { + $gln = Gamma::ln($n / 2); + $a = 0.5 * $n; + $b = $x + 1 - $a; + $fpmin = 1.e-300; + $c = 1 / $fpmin; + $d = 1 / $b; + $h = $d; + for ($i = 1; $i < 101; ++$i) { + $an = -$i * ($i - $a); + $b += 2; + $d = $an * $d + $b; + if (abs($d) < $fpmin) { + $d = $fpmin; + } + $c = $b + $an / $c; + if (abs($c) < $fpmin) { + $c = $fpmin; + } + $d = 1 / $d; + $del = $d * $c; + $h = $h * $del; + if (abs($del - 1) < self::EPS) { + break; + } + } + + return $h * exp(-$x + $a * log($x) - $gln); + } } diff --git a/tests/PhpSpreadsheetTests/Calculation/Functions/Statistical/ChiInvLeftTailTest.php b/tests/PhpSpreadsheetTests/Calculation/Functions/Statistical/ChiInvLeftTailTest.php new file mode 100644 index 00000000..962e20ad --- /dev/null +++ b/tests/PhpSpreadsheetTests/Calculation/Functions/Statistical/ChiInvLeftTailTest.php @@ -0,0 +1,37 @@ + [ + '#NUM!', + -0.1, 3, + ], + 'Probability > 1' => [ + '#NUM!', + 1.1, 3, + ], + 'Freedom > 1' => [ + '#NUM!', + 0.1, 0.5, + ], +];