Подгонка гистограммы с функцией Пуассона

Question

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

Я пытаюсь это:

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

def poisson(k, lamb):
   return (lamb**k/factorial(k))*np.exp(-lamb)

fig, ax = plt.subplots()
ax.hist(magz,bins=10,alpha=0.3)
y,x=np.histogram(magz,bins=10)

x = x + (x[1]-x[0])/2
x = np.delete(x,-1)
parameters, cov_matrix = curve_fit(poisson, x, y)
x_new = np.linspace(x[1], x[-1], 50)
ax.plot(x_new, poisson(x_new, *parameters), color='b')

результат:

Синяя линия подгонки находится внизу графика, и кажется, что она не работает.

magz значения:

magz = [ 24.638505,  20.446914,  22.10271,   21.227533,  21.761152,  18.923867,
  24.054868,  23.92457,   21.515022,  21.835458,  21.204597,  21.848573,
  24.036382,  21.126777,  21.599414,  20.044833,  24.390594,  23.772577,
  19.608918,  22.676774,  23.6312,    24.12077,   21.22321,   20.350204,
  20.548614,  22.650914,  20.561528,  24.892959,  22.49959,   22.94469,
  24.346355,  23.934491,  22.448417,  20.535562,  20.785362,  25.131568,
  24.462043,  24.173652,  19.105512,  20.641586,  19.5268,    25.0747,    23.254556,
  24.460447,  24.37759,   22.708406,  20.765025,  17.27031,   22.723192,  22.4452,
  23.366233,  23.238118,  20.72437,   22.278445,  22.231206,  20.98687,
  20.341756,  22.968032,  22.504077,  21.277349,  19.163532,  22.44034,
  22.406581,  21.999019,  23.313155,  17.945062,  23.027715,  23.640596,
  20.60174,   20.124156,  22.39343,   23.786641,  25.201003,  25.227441,
  22.537135,  23.779794,  19.416995,  23.550315,  23.225908,  22.579364,
  21.69014,   20.50741,   22.543212,  21.640962,  20.425422,  23.793624,
  20.554366,  23.522594,  21.03864,   18.282573,  17.33285,   21.015233,
  21.611837,  24.00014,   22.741971,  24.416913,  24.053067,  24.314742,
  22.400808,  19.622712,  24.911735,  20.04981,   19.128503,  23.83536,
  21.864501,  25.037927,  24.449472,  24.200244,  24.492582,  22.378057,
  25.411336,  21.727248,  20.16391,   22.045026,  20.017703,  18.534318,
  21.995092,  21.434434,  21.825042,  22.78943,   19.006103,  22.398988,
  24.198657,  25.018417,  18.829948,  23.94818,   20.813966]

Я хочу что-то подобное:

Answer 1

Я считаю, что здесь есть две отдельные проблемы.

Во-первых, при использовании curve_fit() вам действительно нужно подумать и дать начальные значения для параметров.Поскольку вы не даете начальное значение для lamb, scipy's curve_fit делает абсолютно неоправданное предположение, что значение равно 1,0.Позвольте мне прояснить: в curve_fit() есть серьезный недостаток, заключающийся в том, что он позволяет вам не указывать начальное значение и молча дает значение для вас.Этот недостаток сбивает с толку многих пользователей.

Проблема здесь в том, что начальное значение lambda = 1 не дает интенсивности при k> 10, где находятся все ваши выборки.Таким образом, при внесении небольших изменений в значение лямбды, при подгонке не наблюдается улучшения при подгонке, и поэтому нет способа исследовать значения и находить лучшее решение.Вы можете проверить это, напечатав значения, которые curve_fit() отправляет вашей функции poisson().Вы определенно хотели бы получить начальное значение около 10 или около того.Он не должен быть идеальным, но он не может быть настолько далек, что функция модели вообще не дает никакой интенсивности.

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

Короче говоря, я думаю, вы захотите что-то вроде

def poisson(k, lamb, scale):
    return scale*(lamb**k/factorial(k))*np.exp(-lamb)

, а затем вызвать curve_fit() с приличными начальными значениями для lamb и scale:

parameters, cov_matrix = curve_fit(poisson, x, y, p0=[10., 10.])

Для меня это все еще не очень подходит, но, по крайней мере, даетчто-то разумное.Надеюсь, это поможет.