[Rd] pgamma discontinuity (PR#7307)

Morten Welinder terra at gnome.org
Mon Oct 25 18:04:08 CEST 2004


A little code study, formula study and experimentation reveals that the
situation is mostly fixable:

1. Get rid of the explicit alpha limit.  (A form of it is implicit in
   (2) and (3) below.)

2. Use the series formula when

       (x < alph + 10 && x < 0.99 * (alph + 1000))

   This guarantees that the sum converges reasonably fast.  (For extremely
   large x and alpha that will take about -53/log2(0.99) iterations for
   53 significant bits, i.e., about 3700 iterations.)

3. Use the continued fraction formula when

       (alph < x && alph - 100 < 0.99 * x)

   Aka, you don't want to use the formula either near the critical point
   where alpha/x ~ 1 unless the numbers are small.

4a. Go to a library and figure out how Temme does it for alpha near x,
    both large.  In this case the 0.99 from above could probably be
    lowered a lot for faster convergence.

or

4b. Use the pnorm approximation.  It seems to do a whole lot better for
    alpha near x than it did for the 10:1 case I quoted.

Comments, please.

Morten



More information about the R-devel mailing list