bugGNU Octave - Bugs: bug #35827, The hygepdf and hygecdf functions...


bug #35827: The hygepdf and hygecdf functions are numerically unstable for large populations

Submitted by:  None
Submitted on:  Tue 13 Mar 2012 01:19:38 PM UTC  
Category:  Libraries Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Postponed Assigned to:  None
Originator Name:  Joachim Arts Originator Email:  -email is unavailable-
Open/Closed:  Open Release:  dev
Operating System:  Any

Add a New Comment(Rich Markup)

You are not logged in

Please log in, so followups can be emailed to you.


Sun 20 Nov 2016 02:51:26 PM UTC, comment #5:

These functions should be part of the Octave-Forge statistics package which will happen at the next release.

Rik <rik5>
Project Administrator
Sat 19 Nov 2016 08:33:37 PM UTC, comment #4:

This issue is still present in Octave 4.2.0.

Hartmut <hardy>
Wed 14 Mar 2012 04:09:53 PM UTC, comment #3:


That's the expression I used as well, but perhaps I mistyped something. At any rate, I think this is worth adopting as it is definitely more robust as a numerical algorithm.

Rik <rik5>
Project Administrator
Wed 14 Mar 2012 09:37:59 AM UTC, comment #2:


the correct answer is 7.60282189600926e-03 with Octave too.

but is it safe to always use this expression? I think so, the difference should be O(eps).


Marco Caliari <caliari>
Project Member
Tue 13 Mar 2012 11:22:33 PM UTC, comment #1:

Just to be clear, I get the following when running the sample code

If I replace the cdf with pdf then Octave produces NaN.

What does Matlab return for

I'm guessing it returns Inf which indicates that they are not using the straightforward implementation of two binomial coefficients divided by a third.

The algorithm for hygepdf is pretty straightforward with two binomial coefficients divided by a third (http://en.wikipedia.org/wiki/Hypergeometric_distribution). Octave implements this directly which you can see by looking at the code in hygepdf.m.

The actual binomial coefficient in bincoeff.m is the exponential of some gammaln functions. Thus, it would be pretty easy to pull the definition from bincoeff into hygepdf. You could then make use of the rules for multiplying and dividing exponentials to add or subtract terms BEFORE applying the final exponential.

I just did a quick test with the example you provided and it does not produce NaN. Instead the answer is 0.0043171. This is a lot better but at this point I'm not sure which is more accurate: Matlab or Octave? How did you determine that the correct answer was 0.0076?

Rik <rik5>
Project Administrator
Tue 13 Mar 2012 01:19:38 PM UTC, original submission:

The hypergeometric pdf and cdf return NaN for large populations.
For example hygecdf(10,10000,200,1000) returns NaN, when the correct result, also given in MATLAB is 0.007602821895.
Fixing by returning the binomial approximation for large populations is dangerous here because for the example above, the binomial approximation is binocdf(10,1000,200/10000)=0.010237. Thus this approximation can be quite poor already.



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

Attach Files:

Attached Files
file #34233:  numHygeCdf.m added by ajf (203B - text/x-objcsrc - Wrapped the expression given in comment 2 into a ready-to-use numerical hygecdf function file)


   bug dependencies.

Items that depend on this one: None found


Carbon-Copy List
  • -email is unavailable- added by hardy (Posted a comment)
  • -email is unavailable- added by ajf (Updated the item)
  • -email is unavailable- added by mtmiller (Updated the item)
  • -email is unavailable- added by caliari (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by None (Submitted the item)

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.




    Follow 5 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2015-06-16 ajf Attached File- => Added numHygeCdf.m, #34233
    2015-06-04 mtmiller StatusNone => Postponed
        Release3.6.1 => dev
        Operating SystemMicrosoft Windows => Any
    2012-03-13 rik5 Dependencies- => Depends on bugs #34362

    Back to the top

    Powered by Savane 3.3