[SciPy-dev] scipy.stats: sf for symmetric distributions.

[David Warde-Farley]

On 16-Sep-09, at 5:25 PM, Pauli Virtanen wrote:

>> sage: thepdf
>> 1/4*sqrt(2)/(1/2*x^2 + 1)^(3/2)
>> sage: n(2 * integral(thepdf, x, -infinity, -1), prec=100)
>> 0.42264973081037423549085121950
>
> Btw, this integral evaluates to
>
> 	1 - 1/sqrt(3)
>
> Now, depending on how you evaluate this in floating point, you get
> different results:
>
> 	1. - 1./sqrt(3.)        = 0.42264973081037416
> 	(sqrt(3.) - 1)/sqrt(3.) = 0.42264973081037421
> 	(3. - sqrt(3.))/3.      = 0.42264973081037427
>
> Perhaps the 0.###27 was obtained like this...

It was obtained via 1 - scipy.special(2, 1),

> It's nice if you can actually maintain this level of accuracy in
> scipy.stats in a way that works on all platforms -- doesn't work in
> scipy.special as the boost tests show.

Well, I guess we're still in trouble (depending on the special  
function) but perhaps less trouble since we're not doing an  
unnecessary subtraction.

Anyway, the ###21 result is closer than ###27 so I'm taking pleasure  
in (*extremely*) small things.

David