Распределение Фиттинга Пуассона на гистограмме - PullRequest
Новая
       57

Распределение Фиттинга Пуассона на гистограмме

1 голос
/ 25 марта 2019

Несмотря на подавляющее количество постов о размещении распределения Пуассона на гистограмме, следуя всем из них, ни один из них, похоже, не работает для меня.

Я ищу соответствие распределения Пуассона на этой гистограммекоторый я построил так: 

import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit
from scipy.misc import factorial

def poisson(t, rate, scale): #scale is added here so the y-axis 
# of the fit fits the height of histogram
    return (scale*(rate**t/factorial(t))*np.exp(-rate))

lifetimes = 1/np.random.poisson((1/550e-6), size=100000)

hist, bins = np.histogram(lifetimes, bins=50)
width = 0.8*(bins[1]-bins[0])
center = (bins[:-1]+bins[1:])/2
plt.bar(center, hist, align='center', width=width, label = 'Normalised data')

popt, pcov = curve_fit(poisson, center, hist, bounds=(0.001, [2000, 7000]))
plt.plot(center, poisson(center, *popt), 'r--', label='Poisson fit')
# import pdb; pdb.set_trace()
plt.legend(loc = 'best')
plt.tight_layout()

Гистограмма, которую я получаю, выглядит следующим образом:

enter image description here

Я дал предположениес масштабом 7000, чтобы масштабировать распределение до той же высоты, что и ось Y гистограммы, которую я нанес, и предположение 2000 в качестве параметра скорости, поскольку оно 2000 > 1/550e-6.Как вы можете видеть, красная пунктирная линия в каждой точке равна 0.Странно pdb.set_trace() говорит мне, что poisson(center, *popt) дает мне список значений 0. 

    126     plt.plot(center, poisson(center, *popt), 'r--', label='Poisson fit')
    127     import pdb; pdb.set_trace()
--> 128     plt.legend(loc = 'best')
    129     plt.tight_layout()
    130 


ipdb> 
ipdb> poisson(center, *popt)
array([ 0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.])

Что не имеет смысла.Я хочу подогнать распределение гистограммы по пуассону так, чтобы оно нашло наилучший коэффициент уравнения распределения Пуассона.Я подозревал, что это может быть связано с тем, что вместо этого я строю гистограмму lifetimes, которая представляет собой технически случайную выборку данных из обратного распределения Пуассона.Поэтому я попытался вычислить якобиан распределения, чтобы я мог изменить переменные, но это все равно не сработает.Я чувствую, что здесь что-то не так, что не кодирование, а математика.

1 Ответ

0 голосов
/ 25 марта 2019

Ваш расчет округляется до нуля.Со скоростью 2000 и шкалой 7000 ваша формула Пуассона уменьшается до:

7000 * 2000 ^ t / (e ^ (2000) * t!)

Используя приближение Стирлинга t!~ (2 * pi * t) ^ (1/2) * (t / e) ^ t вы получаете:

[7000 * 2000 ^ t] / [Sqrt (2 * pi * t) * e^ (2000-т) * (т ^ т)] ~ пуассон (т)

Я использовал python, чтобы получить первые пару значений пуассона (т): 

poisson(1) -> 0
poisson(2) -> 0
poisson(3) -> 0

Использование вольфрамаальфа вы обнаружите, что производная знаменателя больше, чем производная числителя для всех действительных чисел больше нуля.Следовательно, пуассон (t) приближается к нулю, когда t становится больше.

Это означает, что независимо от того, какой т, если ваш показатель равен 2000, функция Пуассона вернет 0.

Извините заформатирование.Они пока не позволяют мне публиковать TeX.

...