From 1e27af9cf577478a2b6b582e9882022f91e36610 Mon Sep 17 00:00:00 2001 From: mrFendel Date: Wed, 19 Apr 2023 17:13:47 +0300 Subject: [PATCH] - Zelen-Severo CDF aproximation - p-value for varianceRatioTest --- .../kmath/series/VarianceRatioTest.kt | 11 +++++++-- .../kscience/kmath/series/zSNormalCdf.kt | 24 +++++++++++++++++++ .../kmath/series/TestVarianceRatioTest.kt | 5 ++++ 3 files changed, 38 insertions(+), 2 deletions(-) create mode 100644 kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/zSNormalCdf.kt diff --git a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/VarianceRatioTest.kt b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/VarianceRatioTest.kt index d81a6c317..280a9db64 100644 --- a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/VarianceRatioTest.kt +++ b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/VarianceRatioTest.kt @@ -9,10 +9,11 @@ import space.kscience.kmath.operations.DoubleField.pow import space.kscience.kmath.operations.algebra import space.kscience.kmath.operations.bufferAlgebra import space.kscience.kmath.operations.fold +import kotlin.math.absoluteValue // TODO: add p-value with formula: 2*(1 - cdf(|zScore|)) -public data class VarianceRatioTestResult(val varianceRatio: Double=1.0, val zScore: Double=0.0) +public data class VarianceRatioTestResult(val varianceRatio: Double=1.0, val zScore: Double=0.0, val pValue: Double=0.5) /** * Container class for Variance Ratio Test result: * ratio itself, corresponding Z-score, also it's p-value @@ -23,6 +24,7 @@ public fun varianceRatioTest(series: Series, shift: Int, homoscedastic: /** * Calculates the Z-statistic and the p-value for the Lo and MacKinlay's Variance Ratio test (1987) * under Homoscedastic or Heteroscedstic assumptions + * with two-sided p-value test * https://ssrn.com/abstract=346975 * **/ @@ -63,6 +65,11 @@ public fun varianceRatioTest(series: Series, shift: Int, homoscedastic: } val zScore = (varianceRatio - 1) / phi.pow(0.5) - return VarianceRatioTestResult(varianceRatio, zScore) + val pValue = 2*(1 - zSNormalCDF(zScore.absoluteValue)) + return VarianceRatioTestResult(varianceRatio, zScore, pValue) } } + + + + diff --git a/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/zSNormalCdf.kt b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/zSNormalCdf.kt new file mode 100644 index 000000000..1eb733ea9 --- /dev/null +++ b/kmath-stat/src/commonMain/kotlin/space/kscience/kmath/series/zSNormalCdf.kt @@ -0,0 +1,24 @@ +/* + * Copyright 2018-2023 KMath contributors. + * Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file. + */ + +package space.kscience.kmath.series + +import space.kscience.kmath.operations.DoubleField.pow +import kotlin.math.PI +import kotlin.math.absoluteValue +import kotlin.math.exp + +public fun zSNormalCDF(x: Double): Double { + + /** + * Zelen & Severo approximation for the standard normal CDF. + * The error is bounded by 7.5 * 10e-8. + * */ + + val t = 1 / (1 + 0.2316419 * x.absoluteValue) + val summ = 0.319381530*t - 0.356563782*t.pow(2) + 1.781477937*t.pow(3) - 1.821255978*t.pow(4) + 1.330274429*t.pow(5) + val temp = summ * exp(-x.absoluteValue.pow(2) / 2) / (2 * PI).pow(0.5) + return if (x >= 0) 1 - temp else temp +} \ No newline at end of file diff --git a/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/TestVarianceRatioTest.kt b/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/TestVarianceRatioTest.kt index 8e559276c..afc0d541d 100644 --- a/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/TestVarianceRatioTest.kt +++ b/kmath-stat/src/commonTest/kotlin/space/kscience/kmath/series/TestVarianceRatioTest.kt @@ -21,9 +21,11 @@ class TestVarianceRatioTest { assertEquals(1.818181, resultHomo.varianceRatio, 1e-6) // homoscedastic zScore assertEquals(2.587318, resultHomo.zScore, 1e-6) + assertEquals(.0096, resultHomo.pValue, 1e-4) val resultHetero = varianceRatioTest(monotonicData, 2, homoscedastic = false) // heteroscedastic zScore assertEquals(0.819424, resultHetero.zScore, 1e-6) + assertEquals(.4125, resultHetero.pValue, 1e-4) } } @@ -35,9 +37,11 @@ class TestVarianceRatioTest { assertEquals(0.0, resultHomo.varianceRatio, 1e-6) // homoscedastic zScore assertEquals(-3.162277, resultHomo.zScore, 1e-6) + assertEquals(.0015, resultHomo.pValue, 1e-4) val resultHetero = varianceRatioTest(volatileData, 2, homoscedastic = false) // heteroscedastic zScore assertEquals(-1.0540925, resultHetero.zScore, 1e-6) + assertEquals(.2918, resultHetero.pValue, 1e-4) } } @@ -62,6 +66,7 @@ class TestVarianceRatioTest { val result = varianceRatioTest(zeroVolData, 4) assertEquals(1.0, result.varianceRatio, 1e-6) assertEquals(0.0, result.zScore, 1e-6) + assertEquals(0.5, result.pValue, 1e-4) } } } \ No newline at end of file