user3412058
user3412058

Reputation: 325

overplot poisson distribution to histogram

I know there are previous questions related to the topic but I can't find the specifics of my question. I'm trying to overplot a Poisson distribution on a histogram. The first thing I tried was using the poisson function from the stats module in scipy:

import numpy
from scipy.stats import poisson

mu = mean(data)
n, bins, patches = pyplot.hist(data, 20, normed = 1)
pyplot.plot(bins, poisson.pmf(bins, mu), 'r-')
pyplot.show()

However as shown in the figure (in blue the histogram of my data), I get the red plot which has three weird peaks.

The histogram in blue represents my data. The red line is the Poisson function using scipy.stats

Thus I tried to write my own poisson distribution function:

def poisson(mu, x):
    from scipy.misc import factorial
    return numpy.exp(-mu) * mu**x * factorial(x)**-1

y = poisson(mu, bins)

But when I try to print it I get an array of nan's. Am I doing something wrong? Or is it that the numbers in the bins are too large?

print y    
[ nan  nan  nan  nan  nan  nan  nan  nan  nan  nan  nan  nan  nan  nan nan nan  nan  nan  nan  nan  nan]

However when printing the results from stats.poisson I get:

[3.25452236e-06   0.00000000e+00   0.00000000e+00   0.00000000e+00
 0.00000000e+00   3.63110218e-04   0.00000000e+00   0.00000000e+00
 0.00000000e+00   0.00000000e+00   5.24385396e-03   0.00000000e+00
 0.00000000e+00   0.00000000e+00   0.00000000e+00   1.06061293e-02
 0.00000000e+00   0.00000000e+00   0.00000000e+00   0.00000000e+00
 3.23183010e-03]

Upvotes: 0

Views: 4563

Answers (1)

N.G
N.G

Reputation: 367

  • for the poisson function, you should give 'int' in input instead of your 'bins', numpy.arange(1200, 1475) for example.

  • for your own poisson function, you have to be careful when you are using 'factorial', especially with large x (x>20) because it rapidly increase! I suspect to be at the origin of your nan. Also factorial of float does not exist!

try:

X = np.arange( 1200, 1450 )
plt.plot( X, poisson.pmf(X,1375), 'r-' )

Upvotes: 1

Related Questions