Normal distribution on 32E and 21S

[Albert Chan]

(01-29-2021 05:46 PM)Dave Britten Wrote:  The various Application Pacs and books have typically used polynomial approximations, but I'm wondering if either of these models integrate numerically, and if so, to what accuracy. Or perhaps it's something else entirely.

Q(x) is likely implemented using rational polynomial approximation of erf(x).

see http://git.musl-libc.org/cgit/musl/tree/src/math/erf.c

If not, we still do not need to numerically integrate.
erf integral can be turn into summation (or continued fraction form, if |x| is big)

see Abramowitz & Stegun, Error Function and Fresnel Integrals, 7.1.5 , 7.1.15
see XCas implementation, translated to Python

Quote:The 32E manual only says "The HP 32E computes Q^-1 using a sophisticated iterative algorithm"

This may be Newton's method: Q^-1(p) = z       → Q(z) - p = 0

>>> from mpmath import *
>>> Q = lambda x: erf(x/sqrt(2)) / 2
>>> pdf = lambda x: exp(-x*x/2) / sqrt(2*pi)
>>>
>>> p = 0.2
>>> z = 2.6*p                   # guess for Q^-1(p)
>>> z -= (Q(z) - p)/pdf(z); print z
0.524395467595098
>>> z -= (Q(z) - p)/pdf(z); print z
0.524400512701367
>>> z -= (Q(z) - p)/pdf(z); print z
0.524400512708041
>>> print Q(z)                   # should recover p
0.2

Comment: I used Casio FX-115MS definition: Q(z) = ∫(pdf(x), x = 0.5 .. z)
Casio also named P(z) = cdf(z), R(z) = 1 - cdf(z) = erfc(z/sqrt(2)) / 2

With Casio definition: P^-1(0.7) = Q^-1(0.2) = R^-1(0.3) = 0.524400512708041

Just noticed statistical Q-function = Casio's R function. Sorry about the confusion.