[ 
https://issues.apache.org/jira/browse/NUMBERS-174?page=com.atlassian.jira.plugin.system.issuetabpanels:comment-tabpanel&focusedCommentId=17442877#comment-17442877
 ] 

Alex Herbert commented on NUMBERS-174:
--------------------------------------

{quote}To be sure: Do you always copy those explanations over to the source 
files?
{quote}
Yes.

But usually it works the other way as I start with comments in the source 
files, investigate and comment on what has been updated and why, then copy the 
findings to the Jira ticket.

In this case the easy option would be to disable this evaluation method and 
effectively set the computation back to the Boost version from 1_68_0 (pre Dec 
2018). But I wanted to find out if it could be used and when. I think there is 
a definite case for leaving the method in. But it should be used more 
conservatively.

I will raise this finding as an issue with the Boost maintainers.

I think the ported code is ready to merge, following some review. There is the 
outstanding issue of whether to merge with its own continued fraction 
implementation and then update it pending the resolution of NUMBERS-175.

I also investigated using Dfp to compute the power functions in extended 
precision. I can post up results of that too. It works very well to eliminate 
errors. I have no idea on performance implications but I think it would be 
significant. It is probably best to leave it out and think of it as a possible 
future enhancement.

 

> Update Gamma functions using the Boost implementation
> -----------------------------------------------------
>
>                 Key: NUMBERS-174
>                 URL: https://issues.apache.org/jira/browse/NUMBERS-174
>             Project: Commons Numbers
>          Issue Type: Improvement
>          Components: gamma
>    Affects Versions: 1.0
>            Reporter: Alex Herbert
>            Assignee: Alex Herbert
>            Priority: Major
>             Fix For: 1.1
>
>
> The current regularised incomplete gamma functions for P and Q compute using 
> a series representation for all input. This method is not robust to extreme 
> arguments.
> The Gamma functions are used in Commons Statistics in the Gamma and Poisson 
> distributions.
> Large values of the mean in the Poisson distribution have low
>  precision for the CDF. The distribution also has a configurable
>  threshold for convergence of the gamma function evaluation
>  (STATISTICS-38). Ideally this implementation detail should be
>  removed from the API.
> The Gamma distribution has a function switch based on a threshold to
>  compute the PDF. This threshold results in incorrect values for
>  certain parameters (STATISTICS-39).
> The function switch in the Gamma distribution is based on the
>  documentation for the Boost gamma functions [1]. However it does not
>  implement all the suggested optimisations detailed in the most recent
>  Boost documentation. This includes avoiding using the converging
>  series of the gamma function when the convergence is slow or unstable,
>  i.e. addresses the need for a configurable convergence threshold in
>  the Poisson distribution.
> One part of the evaluation of the incomplete gamma functions require
>  accuracy in the leading term:
> {noformat}
> x^a exp(-x) / gamma(a) = exp(a log(x) - x - loggamma(a){noformat}
> When x and a are large then using logs to compute this leads to
>  cancellation. Use of logs to compute this term is done in the current
>  implementation in numbers for RegularizedGamma P and Q. It is the
>  source of low precision for the CDF for the Poisson distribution when
>  the mean is large.
> I propose to port the Boost gamma functions to numbers. This will be a
>  process similar to the recent port of the error functions to the same
>  package in numbers. These functions improved both accuracy and
>  speed over the existing implementation. Once the gamma functions are
>  ported a comparison can be made between the existing and new
>  implementations. 
> [1] 
> [https://www.boost.org/doc/libs/1_77_0/libs/math/doc/html/math_toolkit/sf_gamma/igamma.html]



--
This message was sent by Atlassian Jira
(v8.20.1#820001)

Reply via email to