Siginificant figures calculating zprob

[Sarah Wang]

Hello everyone!

I want to calculate zprob(the area under the normal curve)
with python and I managed to find a function from the
internet.

But the problem is that the function calculates the result
with only a few significant figures. If I want to get the
20th number of the result(z-prob) what should I do? I want
to get the confidence level of "6sigma" and all I get at the
moment is "1".

I remember that Python's long type has unlimited number of
significant figures as long as the memory allows. What about
floating point numbers?

BTW, I don't have enough knowledge in numerical analysis.

Thanks in advance,

Sarah Wang

p.s. Apologize me if I have made any mistake writing in the
posting -- I am not very goot at both English and mathematics.

 

[Duncan Smith]

Hello everyone!

I want to calculate zprob(the area under the normal curve)
with python and I managed to find a function from the
internet.

But the problem is that the function calculates the result
with only a few significant figures. If I want to get the
20th number of the result(z-prob) what should I do?

Why would you need this degree of precision?

I want

to get the confidence level of "6sigma" and all I get at the
moment is "1".

How about posting the code? Check out Gary Strangman's stats.py
http://www.nmr.mgh.harvard.edu/Neural_Systems_Group/gary/python.html

Duncan

 

[Josiah Carlson]

I want to calculate zprob(the area under the normal curve)

with python and I managed to find a function from the
internet.

But the problem is that the function calculates the result
with only a few significant figures. If I want to get the
20th number of the result(z-prob) what should I do? I want
to get the confidence level of "6sigma" and all I get at the
moment is "1".

I remember that Python's long type has unlimited number of
significant figures as long as the memory allows. What about
floating point numbers?

Python floats are normally IEEE 754 FP values, which generally have 53
bits of precision.

It's been a while since I did this kind of stuff, but after doing a few
tests, the summation using Simpson's rule for numerical integration
converges fairly quickly for most znorm ranges. I use 4096 divisions
beacuse it seems to be fast and accurate to at least 10 decimal places
for the ranges I've tested.

More divisions may or may not increase precision, but for most tasks, I
would imagine the code given at the end would be sufficient. There are
various float casts, and floating point constants given in the code.
With the tests that I did, it seems that some are necessary for higher
precision.

Oh, I hope this is precise enough....... print div, repr(znorm_range(0, 1, div))
.... div *= 2
....
1024 0.3413447460685452
2048 0.3413447460685422
4096 0.34134474606854304
8192 0.34134474606854265
- Josiah



from math import sqrt, e, pi

def znorm_x(x,c=1.):
return c/sqrt(2.*pi)*e**(-(x*x)/2.)

def znorm_range(low, high, divisions=4096):
#divisions needs to be even, we'll increment it if necessary
divisions += divisions%2
#necessary if either low or high are integers
low, high = float(low), float(high)
inc = (high-low)/float(divisions)
x = low+inc
t = znorm_x(low) + znorm_x(high)
c = 4.
while x < high:
t += znorm_x(x,c)
x += inc
c = 6.-c
return (high-low)*t/3./divisions

 

[Paul Rubin]

But the problem is that the function calculates the result
with only a few significant figures. If I want to get the
20th number of the result(z-prob) what should I do? I want
to get the confidence level of "6sigma" and all I get at the
moment is "1".

6 sigma is a marketing term. The number you want is on the order of
10**-9, so it describes something almost never seen in real life.

http://en.wikipedia.org/wiki/Six_Sigma

I remember that Python's long type has unlimited number of
significant figures as long as the memory allows. What about
floating point numbers?

That's for integers only. Some math libraries have a function called
"erfc" (complementary error function) which does what you want. You
can evaluate it to arbitrary precision in a system like Mathematica.
I think the number you want is 0.5*erfc(6/sqrt(2)) which according to the
Emacs calculator "calc" is 9.86587644935e-10.

 

[Sarah Wang]

Why would you need this degree of precision?

I'm doing some exploration into the statistical part of "Six Sigma"
(as someone on the thread noted) with my most powerful exploratory
tool "Python Interactive Shell".

I just wanted to see myself if six sigma is really 3.4 defects per
million opportunities but got into the limit of float type's precision.

I want

How about posting the code? Check out Gary Strangman's stats.py
http://www.nmr.mgh.harvard.edu/Neural_Systems_Group/gary/python.html

Duncan

Well, actually I used Gary's code(zprob function) and it returned 1
for 6 sigma. The precision falls short.

As Paul rightly directed me, the wikipedia page explained that I was
looking for a wrong target. I should've computed 4.5 sigma instead.
Real 6 sigma is 2 per one billion. But my curiosity drove me to
experiment that rare probability. The python code was reaching far
short but Excel's normdist function did it nicely(it's significant digits
doesn't seem very long though).

Happily with help from Josiah, I could get both what I wanted(6sigma)
and what I did not wanted(4.5sigma) at the first place.

Thank you all.

 

[Paul Rubin]

I'm doing some exploration into the statistical part of "Six Sigma"
(as someone on the thread noted) with my most powerful exploratory
tool "Python Interactive Shell".

The thing is you don't really need all that precision even still.
Math libraries usually implement the error function in two forms,
erf(x) which is basically what you call z, and erfc(x) which is 1-erf(x).
So when erf(x) is close to 1, erfc(x) is close to 0 and therefore has
an accurate floating point representation. You don't get the precision
loss of subtracting two nearly equal floats, so you don't need so much
precision to start with.