Я использую разные python для подбора функций плотности в наборе данных. Этот набор данных состоит из положительных значений времени, начиная с 1 секунды.
Я тестировал различные функции плотности из scipy.statistics
и библиотеки powerlaw
, а также свои собственные функции, используя функцию scipy.optimize
curve_fit()
.
До сих пор я получал наилучшие результаты при подборе следующей «модифицированной» степенной функции:
def funct(x, alpha, x0):
return((x+x0)**(-alpha))
Мой код выглядит следующим образом:
bins = range(1,int(s_distrib.max())+2,1)
y_data, x_data = np.histogram(s_distrib, bins=bins, density=True)
x_data = x_data[:-1]
param_bounds=([0,-np.inf],[np.inf,np.inf])
fit = opt.curve_fit(funct,
x_data,
y_data,
bounds=param_bounds) # you can pass guess for the parameters/errors
alpha,x0 = fit[0]
print(fit[0])
C = 1/integrate.quad(lambda t: funct(t,alpha,x0),1,np.inf)[0]
# Calculate fitted PDF and error with fit in distribution
pdf = [C*funct(x,alpha,x0) for x in x_data]
sse = np.sum(np.power(y_data - pdf, 2.0))
print(sse)
fig, ax = plt.subplots(figsize=(6,4))
ax.loglog(x_data, y_data, basex=10, basey=10,linestyle='None', marker='.')
ax.loglog(x_data, pdf, basex=10, basey=10,linestyle='None', marker='.')
Подгонка возвращает значение 8,48 для x0 и 1,40 для альфы. На графике журнала данные и график соответствия выглядят так:
- Мой первый вопрос технический . Почему я получаю следующее предупреждение и ошибку в
opt.curve_fit
при изменении (x + x0) на (x-x0) в функции funct
? Поскольку мои оценки для x0 равны (-inf, + inf), я ожидал, что фитинг вернет -8.48.
/ anaconda3 / lib / python3 .7 / site-packages /ipykernel_launcher.py:3: RuntimeWarning: деление на ноль, встречающееся в обратном. Это отдельно от пакета ipykernel, поэтому мы можем избегать выполнения импорта до тех пор, пока ValueError: Остатки не являются конечными в начальной точке.
- Другие мои вопросы теоретические . Является ли (x + x0) ^ (- alpha) стандартным распределением? Что представляет собой значение x0, как физически интерпретировать это значение 8.48s? Из того, что я понимаю, это означает, что мое распределение соответствует сдвинутому распределению по степенному закону? Могу ли я считать, что x0 соответствует значению xmin, классически необходимому при подгонке данных к мощности l aws?
- Что касается этого значения xmin, я понимаю, что имеет смысл рассматривать только данные, превышающие этот порог для процесс подбора, чтобы охарактеризовать хвост распределения. Однако мне интересно, каков стандартный способ характеризации полных данных с помощью распределения, которое было бы степенным законом после xmin и чем-то еще до xmin.
Это много вопросов, как и я очень незнаком с предметом, любые комментарии и ответы, даже частичные, будут очень благодарны!