[SciPy-User] gamma.pdf numerically unstable?

[Skipper Seabold]

On Mon, Aug 6, 2012 at 9:29 AM, Tom Furmston <tfurmston at googlemail.com> wrote:
> Hi,
>
> I'm trying to use the gamma distribution in the SciPy library, but it
> appears to be numerically unstable.
>
> First, I've imported the library with
>
> from scipy.stats import gamma
>
> I want to evaluate the gamma distribution are x = 1.6, where the gamma
> distribution has the parameters shape = 142, scale = 0.0098. When I
> type
>
> gamma.pdf(1.6, 142, loc=0, scale=0.0098)
>
> I get a value of inf.
>
> According to matlab the value should be around 0.7030.
>
> Am I using this distribution correctly? Is this distribution known to
> suffer from numerical issues?
>

What SciPy version?

[~/]
[1]: from scipy.stats import gamma

[~/]
[2]: gamma.pdf(1.6, 142, loc=0, scale=.0098)
[2]: 0.69705038205764014

[~/]
[3]: np.exp(gamma.logpdf(1.6, 142, loc=0, scale=.0098))
[3]: 0.69705038205763981

[~/]
[4]: import scipy

[~/]
[5]: scipy.__version__
[5]: '0.12.0.dev-8647010'

In R

> dgamma(1.6, 142, scale=.0098)
[1] 0.6970504


Skipper