Показать pdf распределения хи-квадрат с помощью Python

Я пытаюсь восстановить PDF-файл распределения хи-квадрат с 3 степенями свободы из смоделированного образца. Вот мой код на Python:

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

norm = stats.norm(0, 1)

x1 = [x * x for x in np.random.randn(1000)]
x2 = [x * x for x in np.random.randn(1000)]
x3 = [x * x for x in np.random.randn(1000)]

f = x1 + x2 + x3

plt.hist(f, 100)
plt.show()

В результате я получил вот что.

Распределение Ци со степенью свободы 3

Конечно, это неправильно. Как показано в Википедии, PDF-файл распределения хи-квадрат с 3 степенями свободы должен сначала идти вверх от нуля, а потом идти вниз, а не то, что продолжает расти, как у меня. Что-то не так с моим кодом? Я использовал следующую формулу:

Q = x1^2 + x2^2 + x3^2

где x1, x2 и x3 - независимые стандартные нормальные случайные величины.


python pdf chi-squared
person Searene    schedule 29.11.2016    source источник
comment
Будьте осторожны со своей терминологией! Существует распределение хи, которое по определению является положительным квадратным корнем случайной величины хи-квадрат. Когда вы говорите о свободе, правильная терминология - это степени свободы. Мне кажется, что ваш код генерирует распределение хи-квадрат с 3 степенями свободы.   -  person    schedule 29.11.2016
comment
@MichaelChernick предоставленная гистограмма не показывает распределение хи-квадрат с 3 df, как заметил OP. Подход правильный, и должен возвращать правильный дистрибутив, но в коде должна быть какая-то ошибка (однако я не вижу ее в предоставленном коде - но не тестировал).   -  person Tim    schedule 29.11.2016
comment
Я думаю, что OP имеет в виду распределение $ \ chi ^ 2 $. См. Формулу внизу сообщения.   -  person utobi    schedule 29.11.2016
comment
@Searene $ x_1 $, $ x_2 $ и $ x_3 $ должны быть свежими стандартными нормалями, а не переработанными.   -  person utobi    schedule 29.11.2016
comment
@utobi Я никогда не использовал Python, но полагаю, что каждый новый вызов генератора случайных чисел должен гарантировать это. В любом случае, это проверяется OP или отвечает.   -  person Nick Cox    schedule 29.11.2016
comment
Как кажется, вопрос, в чем ошибка? это принадлежит Stack Overflow.   -  person Nick Cox    schedule 29.11.2016

Ответы (2)


arrow_upward
1
arrow_downward

Хотя я пробовал ваш код и получил тот же результат, что и вы, если вы используете свою «нормальную» переменную для генерации случайных значений, похоже, это сработает.

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

norm = stats.norm(0, 1)

x1 = norm.rvs(size=100000)**2
x2 = norm.rvs(size=100000)**2
x3 = norm.rvs(size=100000)**2

f = x1 + x2 + x3

plt.hist(f, 60, normed=True)

# Plot the theoretical density of f
x = np.arange(0, 30, .05)
plt.plot(x, stats.chi2.pdf(x, df=3), color='r', lw=2)
plt.show()

В результате я получил

Гистограмма Chi2

person Travis L    schedule 29.11.2016
comment
... это косвенно доказывает, что это какая-то ошибка в коде. - person Tim; 29.11.2016

arrow_upward
1
arrow_downward

Оператор '+' работает в списках Python иначе, чем в массивах Numpy.

f = x1 + x2 + x3

объединяет три списка в один. Однако вы хотите добавить содержимое трех списков поэлементно, что можно сделать следующим образом:

f = np.array(x1) + np.array(x2) + np.array(x3)
person cpetrich    schedule 25.09.2019
comment
Это должен быть принятый ответ - person Oliver; 29.05.2021