Skip to content

Commit

Permalink
ENH: allow b=0 for rice distribution
Browse files Browse the repository at this point in the history
  • Loading branch information
ev-br committed Sep 15, 2013
1 parent f2677c3 commit b7e0dec
Show file tree
Hide file tree
Showing 2 changed files with 24 additions and 0 deletions.
3 changes: 3 additions & 0 deletions scipy/stats/distributions.py
Original file line number Diff line number Diff line change
Expand Up @@ -5573,6 +5573,9 @@ class rice_gen(rv_continuous):
%(example)s
"""
def _argcheck(self, b):
return b >= 0

def _pdf(self, x, b):
return x * exp(-(x-b)*(x-b)/2.0) * special.i0e(x*b)

Expand Down
21 changes: 21 additions & 0 deletions scipy/stats/tests/test_continuous_extra.py
Original file line number Diff line number Diff line change
Expand Up @@ -132,5 +132,26 @@ def test_rdist_cdf_gh1285():
values, decimal=5)


def test_rice_zero_b():
# rice distribution should work with b=0, cf gh-2164
x = [0.2, 1., 5.]
npt.assert_(np.isfinite(stats.rice.pdf(x, b=0.)).all())
npt.assert_(np.isfinite(stats.rice.logpdf(x, b=0.)).all())
npt.assert_(np.isfinite(stats.rice.cdf(x, b=0.)).all())
npt.assert_(np.isfinite(stats.rice.logcdf(x, b=0.)).all())

q = [0.1, 0.1, 0.5, 0.9]
npt.assert_(np.isfinite(stats.rice.ppf(q, b=0.)).all())

mvsk = stats.rice.stats(0, moments='mvsk')
npt.assert_(np.isfinite(mvsk).all())

# furthermore, pdf is continuous as b\to 0
# rice.pdf(x, b\to 0) = x exp(-x^2/2) + O(b^2)
# see e.g. Abramovich & Stegun 9.6.7 & 9.6.10
b = 1e-8
np.allclose(stats.rice.pdf(x, 0), stats.rice.pdf(x, b),
atol = b, rtol=0)

if __name__ == "__main__":
npt.run_module_suite()

0 comments on commit b7e0dec

Please sign in to comment.