Кажется, что-то здесь не так с stats.lognorm.fit
.
В docs упоминается альтернатива, подбирая stats.norm
журнала образцов, а затем используя exp(mu)
как масштаб. Кажется, это работает намного лучше.
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
samples = [480, 900, 1140, 1260, 1260, 1440, 1800, 1860, 1980, 2220, 2640, 2700,
2880, 3420, 3480, 3600, 3840, 4020, 4200, 4320, 4380, 4920, 5160,
5280, 6900, 7680, 9000, 10320, 10500, 10800, 15000, 21600, 25200,
39000]
samples = np.array(samples)
std, loc, scale = stats.lognorm.fit(samples) # 2.865850745357322, 479.99969879223596, 1.1400622824414484
weird_samples = stats.lognorm(std, loc=loc, scale=scale).rvs(size=2000)
mu, std = stats.norm.fit(np.log(samples)) # 8.304837454505837, 0.9720253999925554
scale = np.exp(mu) # 4043.3848507251523
loc = 0
new_samples = stats.lognorm(std, loc=loc, scale=scale).rvs(size=2000)
plt.plot(samples, 1 - np.linspace(0, 1, len(samples)), label='given samples')
plt.plot(np.sort(weird_samples), 1 - np.linspace(0, 1, len(weird_samples)), label='using stats.lognorm.fit(samples)')
plt.plot(np.sort(new_samples), 1 - np.linspace(0, 1, len(new_samples)), label='using stats.norm.fit(log(samples))')
plt.legend()
plt.show()
kdeplot Сиборна показывает следующее:
import seaborn as sns
bw = 1500
sns.kdeplot(samples, bw=bw, label='given samples')
sns.kdeplot(weird_samples, bw=bw, label='using stats.lognorm.fit(samples)')
sns.kdeplot(new_samples, bw=bw, label='using stats.norm.fit(log(samples))')
plt.xlim(-5000, 45000)
plt.show()
PS: Проблема, похоже, в том, что подгонка 3 параметров с использованием ограниченных выборок работает не очень хорошо. Вы можете заставить lognorm.fit
использовать loc=0
, что найдет гораздо более разумные параметры. Параметр loc
просто сдвигает выборки на эту величину; часто loc=0
- лучший выбор.
std, loc, scale = stats.lognorm.fit(samples, floc=0) # 0.9720253999925554, 0.0, 4043.3848507251523
Вместо того, чтобы заставлять loc
с floc
, вы также можете указать начальное предположение. Это выглядит еще лучше:
std, loc, scale = stats.lognorm.fit(samples, loc=0) # 1.0527481074345748, 203.08004314932137, 3712.4903893865644