Блог

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

Я хочу подогнать распределение Пуассона к этой гистограмме, которую я построил как таковую:

 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()
  

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

Я предположил, что масштаб равен 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 , которая технически представляет собой случайную выборку данных из инверсии распределения Пуассона. Итак, я попытался вычислить якобиан распределения, чтобы я мог внести изменения в переменные, но это все равно не сработает. Я чувствую, что мне здесь чего-то не хватает, что связано не с кодированием, а скорее с математикой.

Ваши вычисления округляются до нуля. При скорости 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-t) * (t ^ t)] ~ пуассон (t)

Я использовал python для получения первых двух значений poisson (t):

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

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

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

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