Как нарисовать соответствующую кривую Белла на гистограмме? - PullRequest
1 голос
/ 14 января 2020

Пока мой код, я очень новичок в программировании и уже некоторое время пытаюсь.

Здесь я применяю преобразование Бокса-Мюллера для приближения к двум гауссовским нормальные распределения , начиная с случайной равномерной выборки. Затем я создаю гистограмму для них обоих.

Теперь я хотел бы сравнить полученные гистограммы с «реальной вещью»: стандартной кривой Белла . Как нарисовать такую ​​кривую, чтобы она соответствовала гистограммам?

import numpy as np
import matplotlib.pyplot as plt

N = 10000
z1 = np.random.uniform(0, 1.0, N)
z2 = np.random.uniform(0, 1.0, N)

R_sq = -2 * np.log(z1)
theta = 2 * np.pi * z2
z1 = np.sqrt(R_sq) * np.cos(theta)
z2 = np.sqrt(R_sq) * np.sin(theta)

fig = plt.figure()
ax = fig.add_subplot(2, 1, 1)
ax.hist(z1, bins=40, range=(-4, 4), color='red')
plt.title("Histgram")
plt.xlabel("z1")
plt.ylabel("frequency")
ax2 = fig.add_subplot(2, 1, 2)
ax2.hist(z2, bins=40, range=(-4, 4), color='blue')
plt.xlabel("z2")
plt.show()

1 Ответ

1 голос
/ 15 января 2020

Чтобы получить «оценку плотности ядра», scipy.stats.gaussian_kde вычисляет функцию для соответствия данным.

Чтобы просто нарисовать нормальную гауссову кривую, есть [scipy.stats.norm] , Вычитание среднее и деление на стандартное отклонение адаптирует положение к заданным данным.

Обе кривые будут построены так, чтобы область под кривой суммировалась в один. Чтобы настроить их в соответствии с размером гистограммы, эти кривые необходимо масштабировать по длине данных, умноженной на ширину бина. Кроме того, это масштабирование может оставаться на уровне 1, а гистограмма масштабируется путем добавления параметра hist(..., density=True).

. В демонстрационном коде данные искажены, чтобы проиллюстрировать разницу между kde и гауссовой нормалью.

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as stats

x = np.linspace(-4,4,1000)
N = 10000
z1 = np.random.randint(1, 3, N) * np.random.uniform(0, .4, N)
z2 = np.random.uniform(0, 1, N)

R_sq = -2 * np.log(z1)
theta = 2 * np.pi * z2
z1 = np.sqrt(R_sq) * np.cos(theta)
z2 = np.sqrt(R_sq) * np.sin(theta)

fig = plt.figure(figsize=(12,4))
for ind_subplot, zi, col in zip((1, 2), (z1, z2), ('crimson', 'dodgerblue')):
    ax = fig.add_subplot(1, 2, ind_subplot)
    ax.hist(zi, bins=40, range=(-4, 4), color=col, label='histogram')
    ax.set_xlabel("z"+str(ind_subplot))
    ax.set_ylabel("frequency")

    binwidth = 8 / 40
    scale_factor = len(zi) * binwidth

    gaussian_kde_zi = stats.gaussian_kde(z1)
    ax.plot(x, gaussian_kde_zi(x)*scale_factor, color='springgreen', linewidth=3, label='kde')

    std_zi = np.std(zi)
    mean_zi = np.mean(zi)
    ax.plot(x, stats.norm.pdf((x-mean_zi)/std_zi)*scale_factor, color='black', linewidth=2, label='normal')
    ax.legend()

plt.show()

resulting plot

Исходные значения для z1 и z2 очень похожи на нормальное распределение, поэтому черная линия (норма Гаусса для данных) и Зеленая линия (KDE) очень похожа друг на друга.

Текущий код сначала вычисляет реальное среднее и реальное стандартное отклонение данных. Поскольку вы хотите имитировать c идеальный гауссовский нормаль, вы должны сравнить его с кривой со средним нулем и стандартным отклонением один. Вы увидите, что они почти идентичны на сюжете.

original distribution

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...