Регрессия распределения по формуле

#python #pandas #dataframe #numpy #модель_данных

Есть DataFrame, значения индекса которого соответствуют x:

df = pd.Series({0: 712, 1: 254, 2: 96, 3: 85, 4: 74, 5: 65, 6: 66, 7: 60, 8: 62,
9: 63, 10: 58, 11: 58, 12: 47, 13: 38, 14: 39, 15: 35, 16: 33, 17: 29, 18: 25, 19:
21, 20: 22, 21: 23, 22: 16, 23: 18, 24: 13, 25: 13, 26: 11, 27: 14, 28: 12, 29: 10,
30: 6, 31: 7, 32: 6, 33: 3, 34: 2, 35: 5, 36: 5, 37: 2, 38: 2, 39: 3, 40: 1, 41: 2,
42: 3, 43: 2, 44: 1, 45: 1, 48: 1, 49: 1, 50: 0, 51: 0, 52: 0, 54: 1, 56: 0, 57: 0,
63: 0, 208: 0}).to_frame('y')


Есть информация, что он соответствует некоторой сложной формуле с 3-мя неизвестными:
n, z и b:

y ~ (n ** x / factorial(x) * exp(-n))*b + (z ** (x) / factorial(x) * exp(-z))*(1-b)


Каким образом можно узнать эти 3-и неизвестные? 
    

Ответы

Ответ 1

Для этого можно использовать функцию scipy.optimize.curve_fit.

Пример:

import pandas as pd
import numpy as np
from scipy.optimize import curve_fit
from scipy.special import factorial
import matplotlib.pyplot as plt
import matplotlib

matplotlib.style.use('ggplot')

def func(x, n, z, b):
    return (n ** x / factorial(x) * np.exp(-n))*b + (z ** (x) / factorial(x) * np.exp(-z))*(1-b)

df = pd.Series({0: 712, 1: 254, 2: 96, 3: 85, 4: 74, 5: 65, 6: 66, 7: 60, 8: 62,
9: 63, 10: 58, 11: 58, 12: 47, 13: 38, 14: 39, 15: 35, 16: 33, 17: 29, 18: 25, 19:
21, 20: 22, 21: 23, 22: 16, 23: 18, 24: 13, 25: 13, 26: 11, 27: 14, 28: 12, 29: 10,
30: 6, 31: 7, 32: 6, 33: 3, 34: 2, 35: 5, 36: 5, 37: 2, 38: 2, 39: 3, 40: 1, 41: 2,
42: 3, 43: 2, 44: 1, 45: 1, 48: 1, 49: 1, 50: 0, 51: 0, 52: 0, 54: 1, 56: 0, 57: 0,
63: 0, 208: 0}) \
       .to_frame('y')

popt, pcov = curve_fit(func, df.index, df['y'], maxfev=10**6)

print('n :\t{0}\nz :\t{1}\nb :\t{2}'.format(*tuple(popt)))


plt.scatter(df.index, df['y'], s=20, color='orange')
plt.plot(df.index, func(df.index, *popt))


Вывод:

n :     -1.999873816372032
z :     -2.000426132777202
b :     0.655564183300677




PS Судя по графику, "информация" о предполагаемой функции не очень достоверна.

Простая функция:

y = 1/x


подошла бы лучше.