On Sun, Jun 7, 2009 at 6:05 AM, David Cournapeau <da...@ar.media.kyoto-u.ac.jp> wrote: > (Please do not send twice your message to numpy/scipy ML, thank you) > > wierob wrote: >> Hi, >> >> int64 and float seem to work for the stderr calculation. Now, the >> calculation of the p-value causes an underflow. >> >> File "C:\Python26\lib\site-packages\scipy\stats\distributions.py", line >> 2829,in _cdf >> return special.stdtr(df, x) >> FloatingPointError: underflow encountered in stdtr >> >> It seems that the error occurs in scipy.special.stdtr(df, x) if df = >> array([13412]) and x = array([61.88071696]). >> > > The error seems to happen in the cephes library (the library we use for > many functions in scipy.special, including this one). I don't know > exactly what triggers it (your df is quite high, though, and the cdf is > already quite past the 1-eps at your point). >
The result is correct, and we don't get the exception at the default warning/error level, so this looks ok to me. >>> stats.t._cdf(61.88071696, 13412) 1.0 >>> stats.norm._cdf(61.88071696) 1.0 >>> stats.norm._sf(61.88071696) 0.0 >>> stats.norm._ppf(1-1e-16) 8.2095361516013874 Using a double precision library won't help it find out whether the answer is (1 - 1-e30) or (1 - 1-e80) I made comments on a similar case yesterday where the overflow actually shows up in the result without any warning or exception. http://projects.scipy.org/scipy/ticket/962 But I only know about floating point precision what I learned in these news groups, and playing with limiting cases in stats.distributions. using error level to raise doesn't seem useful if you want for example to work with include floating point inf >>> np.seterr(all="raise") {'over': 'ignore', 'divide': 'ignore', 'invalid': 'ignore', 'under': 'ignore'} >>> x = np.array([0,1]) >>> x array([0, 1]) >>> x/x[0] Traceback (most recent call last): File "<pyshell#140>", line 1, in <module> x/x[0] FloatingPointError: divide by zero encountered in divide >>> x = np.array([0.,1.]) >>> x/x[0] Traceback (most recent call last): File "<pyshell#143>", line 1, in <module> x/x[0] FloatingPointError: divide by zero encountered in divide >>> np.seterr(all="ignore") {'over': 'raise', 'divide': 'raise', 'invalid': 'raise', 'under': 'raise'} >>> x/x[0] array([ NaN, Inf]) Josef _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
