>>>>> "PS" == Petr Savicky <savi...@cs.cas.cz>
>>>>>     on Sat, 23 May 2009 18:22:26 +0200 writes:

    PS> On Wed, May 20, 2009 at 11:10:11PM +0200, wolfgang.re...@gmail.com 
wrote:
    PS> ...
    >> Strange behavior of qbinom:
    >> 
    >> > qbinom(0.01, 5016279, 1e-07)
    >> [1] 0
    >> > qbinom(0.01, 5016279, 2e-07)
    >> [1] 16
    >> > qbinom(0.01, 5016279, 3e-07)
    >> [1] 16
    >> > qbinom(0.01, 5016279, 4e-07)
    >> [1] 16
    >> > qbinom(0.01, 5016279, 5e-07)
    >> [1] 0
    >> 

    PS> There is a bug in function do_search() in file src/nmath/qbinom.c. This
    PS> function contains a cycle
    PS> for(;;) {
    PS>   if(y == 0 ||
    PS>   (*z = pbinom(y - incr, n, pr, /*l._t.*/TRUE, /*log_p*/FALSE)) < p)
    PS>   return y;
    PS>   y = fmax2(0, y - incr);
    PS> }
    PS> When this cycle stops, *z contains pbinom(y - incr, ...), but is used
    PS> as if it is pbinom(y, ...) for the resulting y.

    PS> In the example qbinom(p=0.01, size=5016279, prob=4e-07), we get at 
    PS> some step y=15 as a result of a search left with incr=50, so we have
    PS> *z=pbinom(-35, ...)=0. Hence, y=15 is treated as too low and is 
increased
    PS> to 16. Since 16 is detected to be sufficient, the search stops with 
y=16,
    PS> which is wrong.

    [............]

Thanks to Wolfgang and Petr,

Petr's analysis of the problem seems right on spot to me,
and I'm currently testing the patch (and will also add some
regression tests along Petr's example)
which will make it into R-patched (to be R 2.9.1 in a while) and
R-devel.

Gratefully,
Martin Maechler, ETH Zurich

______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel

Reply via email to