Files
2022-07-02 22:37:44 +02:00

121 lines
3.6 KiB
Python

# -*- coding: utf-8 -*-
from numpy import array, infty, log, nan, pi, sqrt
from pandas_ta._typing import Array, IntFloat, Number, Union
from pandas_ta.maps import Imports
from pandas_ta.utils import hpoly
__all__ = [
'inv_norm',
]
def _gaussian_poly_coefficients() -> Array:
"""Three pairs of Polynomial Approximation Coefficients
for the Gaussian Normal CDF"""
p0 = array([
-5.99633501014107895267E1, 9.80010754185999661536E1,
-5.66762857469070293439E1, 1.39312609387279679503E1,
-1.23916583867381258016E0
])
q0 = array([
1.00000000000000000000E0, 1.95448858338141759834E0,
4.67627912898881538453E0, 8.63602421390890590575E1,
-2.25462687854119370527E2, 2.00260212380060660359E2,
-8.20372256168333339912E1, 1.59056225126211695515E1,
-1.18331621121330003142E0
])
p1 = array([
4.05544892305962419923E0, 3.15251094599893866154E1,
5.71628192246421288162E1, 4.40805073893200834700E1,
1.46849561928858024014E1, 2.18663306850790267539E0,
-1.40256079171354495875E-1, -3.50424626827848203418E-2,
-8.57456785154685413611E-4
])
q1 = array([
1.00000000000000000000E0, 1.57799883256466749731E1,
4.53907635128879210584E1, 4.13172038254672030440E1,
1.50425385692907503408E1, 2.50464946208309415979E0,
-1.42182922854787788574E-1, -3.80806407691578277194E-2,
-9.33259480895457427372E-4
])
p2 = array([
3.23774891776946035970E0, 6.91522889068984211695E0,
3.93881025292474443415E0, 1.33303460815807542389E0,
2.01485389549179081538E-1, 1.23716634817820021358E-2,
3.01581553508235416007E-4, 2.65806974686737550832E-6,
6.23974539184983293730E-9
])
q2 = array([
1.00000000000000000000E0, 6.02427039364742014255E0,
3.67983563856160859403E0, 1.37702099489081330271E0,
2.16236993594496635890E-1, 1.34204006088543189037E-2,
3.28014464682127739104E-4, 2.89247864745380683936E-6,
6.79019408009981274425E-9
])
return [p0, q0, p1, q1, p2, q2]
def inv_norm(value: IntFloat) -> Union[None, Number]:
"""Inverse Normal (inv_norm)
Calculates the 'x' in which the area under the Gaussian PDF is
equal to value.
If the user has package "statsmodels" installed, the method will call and
return norm().ppf(value)
Source: https://github.com/scipy/scipy/blob/701ffcc8a6f04509d115aac5e5681c538b5265a2/scipy/special/cephes/ndtri.c
"""
if Imports["statsmodels"]:
from scipy.stats import norm
return norm().ppf(value)
negate = True
v = value
# if v == 0.0: return -npInfty
if v == 0.0:
return -infty
if v == 1.0:
return infty
if v < 0.0 or value > 1.0:
return nan
p0, q0, p1, q1, p2, q2 = _gaussian_poly_coefficients()
sqrt2pi = sqrt(2 * pi)
threshold = 0.13533528323661269189
if v > 1.0 - threshold:
v, negate = 1.0 - v, False
# 0 <= |x0 - 0.5| <= 3/8
if v > threshold:
v -= 0.5
v2 = v * v
y = v + v * (v2 * hpoly(p0, v2) / hpoly(q0, v2))
y *= sqrt2pi
return y
y = sqrt(-2.0 * log(v))
y0 = y - log(y) / y
z = 1.0 / y
if y < 8.0:
# Approximation for interval z = sqrt(-2 log y ) between 2 and 8
# i.e., x between exp(-2) = .135 and exp(-32) = 1.27e-14.
y1 = z * hpoly(p1, z) / hpoly(q1, z)
else:
# Approximation for interval z = sqrt(-2 log y ) between 8 and 64
# i.e., x between exp(-32) = 1.27e-14 and exp(-2048) = 3.67e-890.
y1 = z * hpoly(p2, z) / hpoly(q2, z)
y = y0 - y1
if negate:
y = -y
return y