[Scipy-tickets] [SciPy] #1893: Binomial PMF returns incorrect answer

SciPy Trac scipy-tickets@scipy....
Wed Apr 17 20:22:47 CDT 2013


#1893: Binomial PMF returns incorrect answer
-------------------------+--------------------------------------------------
 Reporter:  uscjeremy    |       Owner:  rgommers   
     Type:  defect       |      Status:  new        
 Priority:  normal       |   Milestone:  Unscheduled
Component:  scipy.stats  |     Version:  0.12.0     
 Keywords:               |  
-------------------------+--------------------------------------------------

Comment(by uscjeremy):

 Digging in a little further, it appears that this bug is because when k=n
 and p=1, one term of the PMF formula reduces to 0**0, causing an
 exception. This is a case where the math doesn't comport with the logical
 meaning of the binomial distribution, though; if the probability of each
 trial succeeding is 1, then the probability of n trials yielding n
 successes is 1 and the probability of yielding any other number of
 successes is 0.

 Mathematica and R both give 1 as the answer for binom.pmf(1,1,1) where
 scipy gives nan. It looks like R's implementation has a bunch of special
 code to deal with the edge cases where p=0 or 1 and k=0 or n, which I
 reproduce below (R v3.0.0, src/nmath/dbinom.c). Note they use 'x' where
 scipy uses 'k', and q=1-p.
 {{{
     if (p == 0) return((x == 0) ? R_D__1 : R_D__0);

     if (q == 0) return((x == n) ? R_D__1 : R_D__0);

     if (x == 0) {
         if(n == 0) return R_D__1;
         lc = (p < 0.1) ? -bd0(n,n*q) - n*p : n*log(q);
         return( R_D_exp(lc) );
     }

     if (x == n) {
         lc = (q < 0.1) ? -bd0(n,n*p) - n*q : n*log(p);
         return( R_D_exp(lc) );
     }

     if (x < 0 || x > n) return( R_D__0 );
 }}}

-- 
Ticket URL: <http://projects.scipy.org/scipy/ticket/1893#comment:1>
SciPy <http://www.scipy.org>
SciPy is open-source software for mathematics, science, and engineering.


More information about the Scipy-tickets mailing list