[Date Index] [Thread Index] [Author Index]

Numerical accuracy of Hypergeometric2F1

• To: mathgroup at smc.vnet.net
• Subject: [mg55717] Numerical accuracy of Hypergeometric2F1
• From: "Christos Argyropoulos M.D." <chrisarg at fuse.net>
• Date: Mon, 4 Apr 2005 00:59:23 -0400 (EDT)
• Sender: owner-wri-mathgroup at wolfram.com

Hi,
Recently I ran into a problem in applied statistics which required 
evaluation of a specific Hypergeometric2F1 functional i.e. 
Hypergeometric2F1[k+1/2,1,3/2,x] where Element[k,Integers], k>0 , 
Element[x, Reals] and  0<=x<=1.
It appears that for "large" values of k , Mathematica returns the wrong 
answer .
If one re-writes the hypergeometric as a polynomial (using the distant 
neighbor relation  http://functions.wolfram.com/07.23.17.0007.01), get 
rid of the ratio of Pochhamer functions and then examine the numerical 
values given by the two formulations diverge.
\!\(fpol[x_, k_] :=
    Module[{}, coeff = Table[1, {k + 1}];
      For[i = 2, i \[LessEqual] k + 1, \(i++\),
        coeff[\([i]\)] =
          coeff[\([i -
                  1]\)]\ \(\(-k\) + i - 1\)\/\(i - k - 1/2\)]; \
{Hypergeometric2F1[k + 1/2, 1, 3/2,
          x], \(1\/\(2  k -
                1\)\) \(\[Sum]\+\(i = 1\)\%\(k + 1\)\((\((1 - 
x)\)\^\(-i\)\ \
coeff[\([i]\)])\)\)}]\)

In[4]:=
SetPrecision[fpol[0.7,100],16]
Out[4]=
\!\({\(-2.5109638442451095387366941924357`15.9546*^56\),
    2.0578481017803705921545468593643`15.9546*^51}\)

but

In[5]:=
SetPrecision[fpol[0.7,10],16]
Out[5]=
{57439.31727367684,57439.31727367681}

or graphically ...
ListPlot3D[Table[Apply[(#1/#2)&,fpol[x,k]],{x,0.,1.,.1},{k,1,100}]]

The problems do not seem to be caused by numerical errors secondary to 
direct evaluation of high - order polynomials in the equivalent 
expression since, using the Horner formulation for the polynomial gives 
essentially identical results

\!\(fhorn[x_, k_] :=
    Module[{}, coeff = Table[1, {k + 1}];
      For[i = 2, i \[LessEqual] k + 1, \(i++\),
        coeff[\([i]\)] =
          coeff[\([i -
                  1]\)]\ \(\(-k\) + i - 1\)\/\(i - k - 1/2\)]; \
{Hypergeometric2F1[k + 1/2, 1, 3/2, x], \(1\/\(2  k - 1\)\)
          Module[{y = 1/\((1 - x)\), res = 1}, res = 0;
            For[i = k + 1, i \[GreaterEqual] \ 1\ , \(i--\),
              res = \((res + coeff[\([i]\)])\)\ y\ ]; res]}]\)

In[6]:=
SetPrecision[fhorn[0.7,10],16]
Out[6]=
{57439.31727367684,57439.31727367682}

and

In[7]:=
SetPrecision[fhorn[0.7,100],16]
Out[7]=
\!\({\(-2.5109638442451095387366941924357`15.9546*^56\),
    2.0578481017803822228995099773782`15.9546*^51}\)

and graphically :
ListPlot3D[Table[Apply[(#1/#2)&,fhorn[x,k]],{x,0.,1.,.1},{k,1,100}]]

Has anyone else run into similar problems with the Hypergeometric2F1 for 
large values of one of the parameters? Is there a quick and dirty way of 
knowing when to mistruct the hypergeometric numerical answer returned by 
Mathematica ?

Christos Argyropoulos MD PhD
University of Cincinnati,
College of Medicine