quantiles of the hypergeometric distribution (PR#502)

Torsten Hothorn <hothorn@ci.tuwien.ac.at> writes:

 term = exp(lfastchoose(NR, xr) + lfastchoose(NB, xb) - lfastchoose(N, n));

The expression inside the exp() function gives -1415.411856 and exp()
returns 0. So term and sum in 

 while(sum < p && xr < xend) {
        xr++;
        NB++;
        term *= (NR / xr) * (xb / NB);
        sum += term;
        xb--;
        NR--;
    }

are 0 all the time and xend = 600 is returned. 

R-Core: How to fix this?

Ick! That is of course only to be expected with recurrence relations
like this applied to large numbers. Since the term formula (as far as
I remember!) applies to all the point probabilities, I suppose that
one could start the recurrence from a "known positive" value and work
both forwards and backwards from there. Since you still want to do the
sums from the low end (although that is perhaps not all that
important?) you might need an intermediate array to hold the values.
One idea could be to start from the normal approximation of the
desired quantile. Probably a bit more than an afternoons work, but
hardly several weeks...