Normal Probability Function

[Albert Chan]

(03-30-2022 08:02 AM)trojdor Wrote:  Additionally, modern CDF function libraries tend to replace the Sqrt(2*PI) in P(x) denominator with a pre-calculated constant, 2.50662827463, for speed / efficiency.

It may be better to remove factor altogether, using P(x) = erfc(x/√(2))/2

Below formula practically matched erfc (A&S 7.1.26) based P(x) formula.
In other words, same max abs error of 7e-8.

P(x) := exp(-x*x/2) *
horner([0.53070271,-0.72657601,0.71070687,-0.14224837,0.1274148,0.], 1/(1+0.23164189*x))

Note: equation constants only required 8 digits precision
I discovered this using numpy to move "/2" into the coefs.

>>> a = [1.061405429,-1.453152027,1.421413741,-0.284496736,0.254829592,0]
>>> numpy.array(a) / 2
array([ 0.53070271, -0.72657601, 0.71070687, -0.14224837, 0.1274148 , 0. ])

Turns out, 8 digits is enough to get practically the same error plots.