Alex Herbert created STATISTICS-102:
---------------------------------------

             Summary: CauchyDistribution inverse probability loses precision in 
the tails
                 Key: STATISTICS-102
                 URL: https://issues.apache.org/jira/browse/STATISTICS-102
             Project: Commons Statistics
          Issue Type: Improvement
          Components: distribution
    Affects Versions: 1.3
            Reporter: Alex Herbert


The Cauchy distribution computes inverse probability functions in the tails 
computes using {{{}location + scale * Math.tan(Math.PI * (p - 0.5)){}}}.

For a probability close to 0 or 1 the subtraction of one half cancels the 
significant digits of {{{}p{}}}: the product {{Math.PI * (p - 0.5)}} is rounded 
to ± π/2, so the result saturates at ± 1.633 × 10{^}16{^} whatever the 
probability is.
{code:java}
final CauchyDistribution dist = CauchyDistribution.of(0.5, 2);
dist.inverseCumulativeProbability(1e-20);           // -3.266247870639074e16, 
expected -6.366197723675814e19
dist.inverseCumulativeProbability(1 - 1e-16);       //  3.957875932190438e15, 
expected  5.73416113922266e15
dist.inverseCumulativeProbability(1e-300);          // -3.266247870639074e16, 
expected -6.366197723675814e299
{code}
This is not only a saturation of the extreme tail. The accuracy decreases from 
the middle of the tail: at {{p = 1e-5}} the value is out by 5e-12 relative (at 
{{p = 1e-10}} by 5e-7); the function is not monotonic, since consecutive 
probabilities below {{1e-18}} can return the identical value; and the median is 
offset by a rounding of π (1.85e-16 relative for the scale 2.1, location 1.2 
distribution). A caller using the quantile function for a likelihood, a p-value 
or a rejection threshold in the tail gets a value that stops responding to the 
input and can decrease.

Fix: use the equivalent forms valid for the two tails, {{tan(pi (p - 1/2)) = 
-cot(pi p)}} for {{p < 1/2}} and {{cot(pi (1 - p))}} for {{{}p > 1/2{}}}, which 
need no cancellation.

 



--
This message was sent by Atlassian Jira
(v10.3.26#10030026)

Reply via email to