Permalink
Browse files

Merge branch 'discrete-ppf-1802' into master.

Reviewed as #394
  • Loading branch information...
2 parents e08beea + 2bc5e74 commit 951dffca35bcd808c54997c13bbf53edc4704e7d @rgommers rgommers committed Jan 5, 2013
Showing with 30 additions and 12 deletions.
  1. +16 −12 scipy/stats/distributions.py
  2. +14 −0 scipy/stats/tests/test_distributions.py
@@ -19,10 +19,10 @@
from numpy import all, where, arange, putmask, \
ravel, take, ones, sum, shape, product, repeat, reshape, \
zeros, floor, logical_and, log, sqrt, exp, arctanh, tan, sin, arcsin, \
- arctan, tanh, ndarray, cos, cosh, sinh, newaxis, array, log1p, expm1
+ arctan, tanh, ndarray, cos, cosh, sinh, newaxis, log1p, expm1
from numpy import atleast_1d, polyval, ceil, place, extract, \
any, argsort, argmax, vectorize, r_, asarray, nan, inf, pi, isinf, \
- power, NINF, empty
+ NINF, empty
import numpy
import numpy as np
import numpy.random as mtrand
@@ -5354,10 +5354,10 @@ def _drv2_moment(self, n, *args):
return tot
def _drv2_ppfsingle(self, q, *args): # Use basic bisection algorithm
- b = self.invcdf_b
- a = self.invcdf_a
+ b = self.b
+ a = self.a
if isinf(b): # Be sure ending point is > q
- b = max(100*q,10)
+ b = int(max(100*q,10))
while 1:
if b >= self.b: qb = 1.0; break
qb = self._cdf(b,*args)
@@ -5366,7 +5366,7 @@ def _drv2_ppfsingle(self, q, *args): # Use basic bisection algorithm
else:
qb = 1.0
if isinf(a): # be sure starting point < q
- a = min(-100*q,-10)
+ a = int(min(-100*q,-10))
while 1:
if a <= self.a: qb = 0.0; break
qa = self._cdf(a,*args)
@@ -5380,7 +5380,7 @@ def _drv2_ppfsingle(self, q, *args): # Use basic bisection algorithm
return a
if (qb == q):
return b
- if b == a+1:
+ if b <= a+1:
#testcase: return wrong number at lower index
#python -c "from scipy.stats import zipf;print zipf.ppf(0.01,2)" wrong
#python -c "from scipy.stats import zipf;print zipf.ppf([0.01,0.61,0.77,0.83],2)"
@@ -5392,10 +5392,16 @@ def _drv2_ppfsingle(self, q, *args): # Use basic bisection algorithm
c = int((a+b)/2.0)
qc = self._cdf(c, *args)
if (qc < q):
- a = c
+ if a != c:
+ a = c
+ else:
+ raise RuntimeError('updating stopped, endless loop')
qa = qc
elif (qc > q):
- b = c
+ if b != c:
+ b = c
+ else:
+ raise RuntimeError('updating stopped, endless loop')
qb = qc
else:
return c
@@ -5595,8 +5601,6 @@ def __init__(self, a=0, b=inf, name=None, badvalue=None,
self.badvalue = badvalue
self.a = a
self.b = b
- self.invcdf_a = a # what's the difference to self.a, .b
- self.invcdf_b = b
self.name = name
self.moment_tol = moment_tol
self.inc = inc
@@ -6661,7 +6665,7 @@ def _rvs(self, M, n, N):
def _argcheck(self, M, n, N):
cond = rv_discrete._argcheck(self,M,n,N)
cond &= (n <= M) & (N <= M)
- self.a = N-(M-n)
+ self.a = max(N-(M-n), 0)
self.b = min(n,N)
return cond
def _logpmf(self, k, M, n, N):
@@ -938,6 +938,20 @@ def test_norm_logcdf():
finally:
np.seterr(**olderr)
+def test_hypergeom_interval_1802():
+ #these two had endless loops
+ assert_equal(stats.hypergeom.interval(.95, 187601, 43192, 757),
+ (152.0, 197.0))
+ assert_equal(stats.hypergeom.interval(.945, 187601, 43192, 757),
+ (152.0, 197.0))
+ #this was working also before
+ assert_equal(stats.hypergeom.interval(.94, 187601, 43192, 757),
+ (153.0, 196.0))
+
+ #degenerate case .a == .b
+ assert_equal(stats.hypergeom.ppf(0.02, 100, 100, 8), 8)
+ assert_equal(stats.hypergeom.ppf(1, 100, 100, 8), 8)
+
if __name__ == "__main__":
run_module_suite()

0 comments on commit 951dffc

Please sign in to comment.