Chi squared inverse left tailed (#1964)
* Implementation of the CHISQ.INV() method for ChiSquared distribution left-tail
This commit is contained in:
parent
e2ff14fe89
commit
e68978f1c7
|
|
@ -9,7 +9,7 @@ and this project adheres to [Semantic Versioning](https://semver.org).
|
||||||
|
|
||||||
### Added
|
### 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)
|
- Support for ActiveSheet and SelectedCells in the ODS Reader and Writer. [PR #1908](https://github.com/PHPOffice/PhpSpreadsheet/pull/1908)
|
||||||
|
|
||||||
### Changed
|
### Changed
|
||||||
|
|
|
||||||
|
|
@ -508,7 +508,7 @@ class Calculation
|
||||||
],
|
],
|
||||||
'CHISQ.INV' => [
|
'CHISQ.INV' => [
|
||||||
'category' => Category::CATEGORY_STATISTICAL,
|
'category' => Category::CATEGORY_STATISTICAL,
|
||||||
'functionCall' => [Functions::class, 'DUMMY'],
|
'functionCall' => [Statistical\Distributions\ChiSquared::class, 'inverseLeftTail'],
|
||||||
'argumentCount' => '2',
|
'argumentCount' => '2',
|
||||||
],
|
],
|
||||||
'CHISQ.INV.RT' => [
|
'CHISQ.INV.RT' => [
|
||||||
|
|
|
||||||
|
|
@ -11,6 +11,8 @@ class ChiSquared
|
||||||
|
|
||||||
private const MAX_ITERATIONS = 256;
|
private const MAX_ITERATIONS = 256;
|
||||||
|
|
||||||
|
private const EPS = 2.22e-16;
|
||||||
|
|
||||||
/**
|
/**
|
||||||
* CHIDIST.
|
* CHIDIST.
|
||||||
*
|
*
|
||||||
|
|
@ -127,6 +129,35 @@ class ChiSquared
|
||||||
return $newtonRaphson->execute($probability);
|
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)
|
public static function test($actual, $expected)
|
||||||
{
|
{
|
||||||
$rows = count($actual);
|
$rows = count($actual);
|
||||||
|
|
@ -167,4 +198,104 @@ class ChiSquared
|
||||||
|
|
||||||
return ($columns - 1) * ($rows - 1);
|
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);
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
|
||||||
|
|
@ -0,0 +1,37 @@
|
||||||
|
<?php
|
||||||
|
|
||||||
|
namespace PhpOffice\PhpSpreadsheetTests\Calculation\Functions\Statistical;
|
||||||
|
|
||||||
|
use PhpOffice\PhpSpreadsheet\Calculation\Functions;
|
||||||
|
use PhpOffice\PhpSpreadsheet\Calculation\Statistical;
|
||||||
|
use PHPUnit\Framework\TestCase;
|
||||||
|
|
||||||
|
class ChiInvLeftTailTest extends TestCase
|
||||||
|
{
|
||||||
|
protected function setUp(): void
|
||||||
|
{
|
||||||
|
Functions::setCompatibilityMode(Functions::COMPATIBILITY_EXCEL);
|
||||||
|
}
|
||||||
|
|
||||||
|
/**
|
||||||
|
* @dataProvider providerCHIINV
|
||||||
|
*
|
||||||
|
* @param mixed $expectedResult
|
||||||
|
* @param mixed $probability
|
||||||
|
* @param mixed $degrees
|
||||||
|
*/
|
||||||
|
public function testCHIINV($expectedResult, $probability, $degrees): void
|
||||||
|
{
|
||||||
|
$result = Statistical\Distributions\ChiSquared::inverseLeftTail($probability, $degrees);
|
||||||
|
if (!is_string($expectedResult)) {
|
||||||
|
$reverse = Statistical\Distributions\ChiSquared::distributionLeftTail($result, $degrees, true);
|
||||||
|
self::assertEqualsWithDelta($probability, $reverse, 1E-12);
|
||||||
|
}
|
||||||
|
self::assertEqualsWithDelta($expectedResult, $result, 1E-12);
|
||||||
|
}
|
||||||
|
|
||||||
|
public function providerCHIINV()
|
||||||
|
{
|
||||||
|
return require 'tests/data/Calculation/Statistical/CHIINVLeftTail.php';
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
@ -0,0 +1,64 @@
|
||||||
|
<?php
|
||||||
|
|
||||||
|
return [
|
||||||
|
[
|
||||||
|
4.671330448981,
|
||||||
|
0.3, 7,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
12.548861396889,
|
||||||
|
0.75, 10,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
3.283020286760,
|
||||||
|
0.93, 1,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
1.832581463748,
|
||||||
|
0.6, 2,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
0.454936423120,
|
||||||
|
0.5, 1,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
4.351460191096,
|
||||||
|
0.5, 5,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
1.323303696932,
|
||||||
|
0.75, 1,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
0.210721031316,
|
||||||
|
0.1, 2,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
3.218875824869,
|
||||||
|
0.8, 2,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
1.212532903046,
|
||||||
|
0.25, 3,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
'#VALUE!',
|
||||||
|
'NaN', 3,
|
||||||
|
],
|
||||||
|
[
|
||||||
|
'#VALUE!',
|
||||||
|
0.25, 'NaN',
|
||||||
|
],
|
||||||
|
'Probability < 0' => [
|
||||||
|
'#NUM!',
|
||||||
|
-0.1, 3,
|
||||||
|
],
|
||||||
|
'Probability > 1' => [
|
||||||
|
'#NUM!',
|
||||||
|
1.1, 3,
|
||||||
|
],
|
||||||
|
'Freedom > 1' => [
|
||||||
|
'#NUM!',
|
||||||
|
0.1, 0.5,
|
||||||
|
],
|
||||||
|
];
|
||||||
Loading…
Reference in New Issue