[SciPy-dev] scipy.special.hyp2f1 bug?

[Pauli Virtanen]

Wed, 16 Sep 2009 21:02:41 -0400, Skipper Seabold wrote:
> Since y'all are looking at the distributions vs. mpmath, it reminded me
> of this.  I haven't filed a bug report yet, but I don't think one
> existed when I found this.
> 
> I was working with the connection between the incomplete beta function
> and the hypergeometric2f1 function in scipy.special for GLM in
> statsmodels and I came across this blog post on mpmath about some corner
> cases for the hypergeometric2f1.
> 
> <http://fredrik-j.blogspot.com/2009/06/hypergeometric-2f1-incomplete-
beta.html>
> 
> Two examples of discrepancies.
[clip]

Yep, it seems that not only is the Cephes implementation of hyp2f1 buggy, 
but also the specfun.f implementation is. Oh joy.

If you want to help by writing patches adding test cases, please do! 
Possible ways to do testing are

1) Generating points and known-good values on-the-fly:

   a) mpmath (to be run only when mpmath is available, otherwise skipped)
   b) known-good expansions

2) Testing against tabulated known-good values:

   a) added directly to the unit tests
   b) tabulated to a file

Tests in 1), 2a) go in scipy/special/tests/test_basic.py.
For 2b), you can just make a text file that contains known-good values of
(a, b, c, z, hyp2f1(a,b,c,z)), one set per row.

The main task in writing the tests is in choosing representative values 
of (a,b,c,z). The hypergeometric functions have all sorts of 
cancellations and singularities, and the test points should be chosen 
with some knowledge of these -- this is the main work. Some tests can 
also be probably stolen from mpmath's test cases.

-- 
Pauli Virtanen