Ошибка в распределении хи-квадрат

Я пытаюсь сгенерировать случайные величины хи-квадрат по следующему алгоритму:

введите описание изображения здесь

где a (i) - независимые стандартные нормальные случайные величины с m четными и нечетными соответственно.

Википедия дает следующее определение:

Код, который я написал:

dch = double(1000)
t = double(1)
for(i in 1:1000) {
  for(j in 1:m) {
    x = runif(1, 0, 1)
    t = t + x*x
  }
  dch[i] = t
}

но я получаю неправильный график плотности. Итак, где находятся ошибки и как их исправить?


r distribution chi-squared
person 8txra    schedule 19.10.2016    source источник
comment
Ну, вы не определяете m. Вы используете однородные случайные величины, а не обычные случайные величины. Это те ошибки, которые я вижу.   -  person Gregor Thomas    schedule 19.10.2016
comment
Кроме того, я не могу очень хорошо прочитать ваше изображение PNG, но продукт, кажется, переходит от 2 до m / 2, тогда как вы переходите от 1 до m. А журнал продукта можно переписать как сумму логов. Вы подводите итоги, а логов нет ...   -  person Gregor Thomas    schedule 19.10.2016
comment
Конечно, обычно вы не можете вести журнал обычного RV, поэтому, возможно, ваш код (единообразная случайная величина) правильный, и ваше описание (стандартная нормальная случайная величина) является проблемой.   -  person Gregor Thomas    schedule 19.10.2016
comment
Я сделал согласно вики. потому что я не понимаю, как это могло быть от 2 до m / 2, если m = 2 (что является моим случаем). поэтому я меняю m на 2, но это все еще неверно. даже больше :(   -  person 8txra    schedule 19.10.2016
comment
Если у вас есть две степени свободы, вики говорится, что пусть X будет случайным образом однородным на (0, 1), тогда -2 * log (X) ~ ChiSq (2). Так что, может быть, тебе стоит просто сделать x = runif(1000); dch = -2 * log(x). Сравнивая hist(dch, 20) с hist(rchisq(1000), 20), это выглядит правильно.   -  person Gregor Thomas    schedule 19.10.2016

Ответы (1)


arrow_upward
2
arrow_downward

Как предположил Грегор в комментариях, вы неправильно интерпретируете входные данные алгоритма. Один из способов получить хи-квадрат с m степенями свободы - это сложить m независимых квадратов стандартных нормалей, но это не единственное известное нам соотношение распределения. Оказывается, что хи-квадрат (2) совпадает с экспоненциальным распределением со средним значением 2, а экспоненты легко генерировать с помощью выборка с обратным преобразованием, также известная как инверсия. Итак, в принципе, если m четное, вы хотите сгенерировать экспоненты m / 2 (2) и просуммировать их. Если m нечетное, сделайте то же самое, но добавьте еще один стандартный нормальный квадрат.

Все это означает, что простая реализация потребовала бы выполнения логарифмических вычислений m / 2 для генерации экспонент. Оказывается, вы можете применить свойство суперпозиции экспонент, так что у вас будет только сделать одну оценку журнала. Поскольку трансцендентные функции требуют больших вычислительных ресурсов, это повышает эффективность алгоритма.

Когда пыль уляжется - z во второй строке вашего алгоритма является стандартным нормальным, но a являются однородными (0,1), а не нормалями.

person pjs    schedule 19.10.2016