[SciPy-User] Rician distributions lacks sigma parameter

Warren Weckesser warren.weckesser@enthought....
Wed Apr 4 06:58:19 CDT 2012


On Wed, Apr 4, 2012 at 2:00 AM, Morten Kjeldgaard <mok@bioxray.dk> wrote:

> Hi,
>
> I am a new reader of this list, please forgive me if this issue has
> already discussed.
>
> I have been wanting to use the rician distribution (stats.rice) to
> analyze some data, but it seems that the implementation in scipy does
> not take the distribution's sigma parameter into account; rather, it
> has been set to 1.0. The wikipedia article shows the traditional
> formulation of the rician distribution [0].
>
> I understand that some distributions, e.g. stats.norm, use the scale
> parameter to define std, but this does not seem to be the case with
> stats.rice.
>
> Any ideas on how to get around this, without actually modifying
> distribution.py? I have no experience with the internals of scipy, and
> wouldn't know how to modify it correctly.
>
> Cheers,
> Morten
>
> [0] http://en.wikipedia.org/wiki/Rice_distribution
>
>

Hi Morten,

Given the parameters nu and sigma (as shown in the wikipedia article), you
use scipy.stats.rice by setting the scale=sigma and the shape parameter
b=nu/sigma.  You can use the following script to verify this:

-----
import numpy as np
from scipy.special import i0
from scipy.stats import rice
import matplotlib.pyplot as plt


def rice_pdf(x, nu, sigma):
    s2 = sigma**2
    f = x/s2 * np.exp(-(x**2 + nu**2)/(2*s2)) * i0(x*nu/s2)
    return f

x = np.linspace(0, 8, 201)

nu = 3.45
sigma = 0.35

my_pdf = rice_pdf(x, nu, sigma)
sp_pdf = rice.pdf(x, nu / sigma, scale=sigma)

np.testing.assert_allclose(my_pdf, sp_pdf)

plt.plot(x, my_pdf, label="rice_pdf")
plt.plot(x, sp_pdf, label="stats.rice.pdf")
plt.legend(loc='best')
plt.show()
-----

Warren
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20120404/17a8e88b/attachment.html 


More information about the SciPy-User mailing list