Создание логнормальных выборок, соответствующих данным, из которых он был создан - PullRequest
2 голосов
/ 16 июня 2020

Я пытаюсь создать новый сэмпл на основе некоторых других сэмплов, но я что-то делаю / понимаю что-то не так. У меня есть 34 образца, которые, как я полагаю, распределены относительно логнормально. На основе этих образцов я хочу сгенерировать 2000 новых образцов. Вот код, который я запускаю:

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]
plt.plot(samples, 1 - np.linspace(0, 1, len(samples)))
std, loc, scale = stats.lognorm.fit(samples)
new_samples = stats.lognorm(std, loc=loc, scale=scale).rvs(size=2000)

a = plt.hist(new_samples, bins=range(100, 40000, 200),
             weights=np.ones(len(new_samples)) / len(new_samples))
plt.show()

Вот график, и, как вы можете видеть, действительно мало образцов выше 1000, хотя образец содержал довольно много выше 1000.

plotPlot2

Как мне лучше всего сгенерировать образец, который лучше представляет ожидаемые значения?

1 Ответ

3 голосов
/ 16 июня 2020

Кажется, что-то здесь не так с 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()

resulting plot

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()

kdeplots

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

fitting with initial guess

...