# bug #56453: [Octave-Forge] chi2inv, gaminv inaccurate for small 1st argument

 Submitted by: None Submitted on: Thu 06 Jun 2019 08:39:08 PM UTC Category: Octave Forge Package Severity: 3 - Normal Priority: 5 - Normal Item Group: Inaccurate Result Status: Confirmed Assigned to: None Originator Name: Jerry W. Lewis Originator Email: -email is unavailable- Open/Closed: Open Release: 4.4.1 Operating System: Any

You are not logged in

## Discussion

Sat 14 Sep 2019 06:15:54 PM UTC, comment #3:

I applied the patch for gammaincinv from comment #2 here (https://hg.savannah.gnu.org/hgweb/octave/rev/9a946e34553c).

The remaining two functions are part of the statistics package on Octave Forge.  I have re-categorized the bug and marked the issue as confirmed.

Rik <rik5>
Spam posted by _148912
Spam posted by _148921
Spam posted by _148912
Spam posted by _148912
Spam posted by _148912
Spam posted by _148563
Fri 07 Jun 2019 11:45:41 AM UTC, comment #2:

The attached patch seems to solve the problem for gammaincinv. Maybe a similar fix can be applied to the other two functions?

Marco Caliari <caliari>
Fri 07 Jun 2019 02:49:09 AM UTC, comment #1:

Likely the issue is in part a poor starting value from which Newton's method either diverges or converges very slowly.

It follows from Abramowitz & Stegun 26.4.6 that as p approaches zero,
chi2inv(p,df) approaches 2*[p*gamma(1+df/2)]^(2/df)

Jerry W. Lewis <jlewis>
Thu 06 Jun 2019 08:39:08 PM UTC, original submission:

The following should be analytically true
chi2inv(p,df) == gaminv(p,df/2,1)*2 == gammaincinv(p,df/2)*2

Consider chi2inv(p,2), which converges to 2*p as p approaches zero.

chi2inv(1e-15,2) returns 2.000000000000001e-15 for essentially machine accuracy;
gaminv(1e-15,1,1)*2 returns 1.900000018803639e-15 for 1 figure accuracy; and
gammaincinv(1e-15,1)*2 returns 1.998401444325283e-15 for 3 figure accuracy.

chi2inv(p,2) for p<=1e-16 returns 1.490116119384766e-09 instead of 2*p;
gaminv(p,1,1)*2 for p<=1e-16 returns 2.980232238769531e-09 == chi2inv(p,2)*2 instead of 2*p; and
gammaincinv(1e-16,1)*2 returns 2.220446049250313e-16 for almost 1 figure accuracy.

gammaincinv(p,1)*2 for p<=1e-17 returns -0 instead of 2*p.

Anonymous

## Attached Files

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:

Comment:

Attached Files
file #47051:  gammaincinv.diff added by caliari (6KiB - text/x-patch)

## Dependencies

Depends on the following items: None found

Items that depend on this one: None found

Carbon-Copy List
• -email is unavailable- added by rik5 (Posted a comment)
• -email is unavailable- added by _148921 (Posted a comment)
• -email is unavailable- added by _148912 (Posted a comment)
• -email is unavailable- added by _148563 (Posted a comment)
• -email is unavailable- added by caliari (Updated the item)
• -email is unavailable- added by jlewis (Posted a comment)
• -email is unavailable- added by None (Submitted the item)
•

Do you think this task is very important?
This task has 0 encouragements so far.

Only project members can vote.