Я реализую Оценщик максимального правдоподобия для дискретных данных счета с целью подбора кривой, реализуя функцию curve_fit, используемую в качестве исходного параметра угадывания для функции минимизации в scipy.Я определил и попробовал эти методы для нескольких дистрибутивов, но для простоты включу только один - дистрибутив logseries.
В этот момент я также попробовал следующие методы из методов statsmodels: 1. statsmodels.discrete.discrete_model.fit 2. statsmodels.discrete.count_model.fit 3. statsmodels.base.model.GenericLikelihoodModel
Большинство подгонок кривой имеют тенденцию встречаться с ошибками переполнения или с нулями и нулями внутри.Я подробно опишу эти ошибки в другом посте
#Import a few packages
import numpy as np
from scipy.optimize import curve_fit
from scipy.optimize import minimize
from scipy import stats
from numpy import log
import numpy as np
import matplotlib.pyplot as plt
#Given data
x=np.arange(1, 28, 1)
y=np.array([18899, 10427, 6280, 4281, 2736, 1835, 1158, 746, 467, 328, 201, 129, 65, 69, 39, 21, 15, 10, 3, 3, 1, 1, 1, 1, 1, 1, 1])
#Define a custom distribution
def Logser(x, p):
return (-p**x)/(x*log(1-p))
#Doing a least squares curve fit
def lsqfit(x, y):
cf_result = curve_fit(Logser, x, y, p0=0.7, bounds=(0.5,1), method='trf')
return cf_result
param_guess=lsqfit(x,y)[0][0]
print(param_guess)
#Doing a custom MLE definition, minimized using the scipy minimize function
def MLERegression(param_guess):
yhat = Logser(x, param_guess) # predictions based on a parameter value
sd=1 #initially guessed for fitting a normal distribution error around the regressed curve
# next, we flip the Bayesian question
# compute PDF of observed values normally distributed around mean (yhat)
# with a standard deviation of sd
negLL = -np.sum( stats.norm.logpdf(y, loc=yhat, scale=sd) ) #log of the probability density function
return negLL
results = minimize(MLERegression, param_guess, method='L-BFGS-B', bounds=(0.5,1.0), options={'disp': True})
final_param=results['x']
print(final_param)
Я ограничил оптимизатор выводом результатов, аналогичных ожидаемым (значение параметра около 0,8 или 0,9). В противном случае алгоритм выдает ноль