[SciPy-dev] scipy.stats: sf for symmetric distributions.
Wed Sep 16 19:02:38 CDT 2009
On 16-Sep-09, at 6:40 PM, firstname.lastname@example.org wrote:
> Since you have multi-precision available with sage, could you check
> the precision of the normal cdf in the lower tail?
> In http://projects.scipy.org/scipy/ticket/962 there was a discussion
> about using the truncated normal in the tails. In the last comment
> there, I mentioned that it would be possible to use the lower tail if
> the precision of scipy.special in this case is high enough. But I
> wasn't able to verify how good the normal cdf really is.
> These numbers are in the range 1e-20 to 1e-80, so here precision
> really matters.
Sadly, Sage works out the integral involving the erf() function, which
is only implemented for double-precision.
I did however do it with numerical_integral which provides an error
bound (not a very good one, it turns out):
sage: numerical_integral(1/sqrt(2*pi) * exp(-x^2/2), -infinity, -10.2)
In : norm.cdf(-10.2)
sage: numerical_integral(1/sqrt(2*pi) * exp(-x^2/2), -infinity, -20)
In : norm.cdf(-20)
sage: numerical_integral(1/sqrt(2*pi) * exp(-x^2/2), -infinity, -10)
In : norm.cdf(-10)
We seem to be getting about ten digits of agreement, but that could be
sage: numerical_integral(1/sqrt(2*pi) * exp(-x^2/2), -infinity, -30)
In : norm.cdf(-30)
(gotta love that bound)
Maybe GSL has an arbitrary precision erf taht we could compare against?
More information about the Scipy-dev