Моделирование случайной величины с распределением Фишера

[MaryNak]

Добрый день.
Возникла необходимость моделирования случ. величины с распределением Фишера. При этом заранее задаются параметры $n$ и $m$. Нужен именно алгоритм/формулы.
Мои размышления: Согласно википедии, распределение случ величины равно

$$F=\frac{Y_1/d_1}{Y_2/d_2}$$

, где $Y_1$ и $Y_2$ — две независимые случайные величины, имеющие распределение хи-квадрат.
Как найти значения $Y_1$ и $Y_2$? В каком диапазоне будут случ. величины?

[Henrylee]

Исходные данные какие есть? Генератор равномерного на $[0,1]$ распределения?

[Хорхе]

Сгенерировать много нормальных величин и поскладывать квадраты.

[MaryNak]

Henrylee в сообщении #579201 писал(а):

Исходные данные какие есть? Генератор равномерного на $[0,1]$ распределения?

Исходные данные n и m любые какие-нибудь выберите.

Хорхе в сообщении #579209 писал(а):

Сгенерировать много нормальных величин

А как это сделать попроще?

[Хорхе]

Никак.

[Munuvonaza]

MaryNak в сообщении #579211 писал(а):

А как это сделать попроще?

А что сложного в приведенном методе: можно использовать таблицу нормального распределения, например здесь http://www.risktheory.ru/distr_tab_normal.htm, выбрать из нее равномерным случайным образом несколько строчек, и как было отмечено, правильно сложить квадраты
Вообще в целях ускорения моделирования, значительно проще сделать таблицу распределения Фишера, хотя ее также легко найти в открытом доступе, например здесь http://www.exponenta.ru/educat/referat/XIkonkurs/student1/index.asp
Можно вообще моделировать программно, исходя из определения самого распределения Фишера; в общем, задание тривиально

[MaryNak]

Цитата:

А что сложного в приведенном методе: можно использовать таблицу нормального распределения, например здесь http://www.risktheory.ru/distr_tab_normal.htm, выбрать из нее равномерным случайным образом несколько строчек, и как было отмечено, правильно сложить квадраты
Вообще в целях ускорения моделирования, значительно проще сделать таблицу распределения Фишера, хотя ее также легко найти в открытом доступе, например здесь http://www.exponenta.ru/educat/referat/XIkonkurs/student1/index.asp
Можно вообще моделировать программно, исходя из определения самого распределения Фишера; в общем, задание тривиально

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

А d1 и d2 в нашем случае это n и m, да?

[ewert]

MaryNak в сообщении #579289 писал(а):

Как самому сгенерировать/посчитать необходимые случайные величины?

Стандартный способ генерации стандартного нормального распределения: генерируем случайный равномерный полярный угол $\varphi$ и случайный квадрат расстояния до центра $r$ (величина $r^2$ распределена по показательному закону и потому генерируется явно, через логарифм равномерно распределённой величины). Тогда величины $x_n=r_n\cos\varphi_n$ и $y_n=r_n\sin\varphi_n$ будут распределены независимо и нормально.

[Евгений Машеров]

0. Классическая книга: Д. Кнут. Искусство программирования на ЭВМ. т.2. Получисленные алгоритмы. Гл. 3. Случайные числа.
1. Генерация F-распределённой случайной величины даётся в 3.4.1.Е(3) (стр. 144 русского изд. 1977г.) Там рекомендуется генерировать, как отношение двух случайных $\chi^2$-распределённых величин с $\nu_1$ и $\nu_2$ степенями свободы (которые являются степенями свободы для F-распределения)
2. Для генерации $\chi^2$-распределённой величины он предлагает для $\nu=2k$ генерировать k экспоненциально распределённых величин со средним 1 и брать их удвоенную сумму, для $\nu=2k+1$ генерировать k экспоненциально распределённых величин со средним 1, брать их удвоенную сумму и прибавлять к ней квадрат нормально распределённой случайной величины $Z^2$, $Z \sim N(0,1)$
3. Генерация экспоненциальной с.в. крайне проста. $X=-\ln U, где U - равномерно распределённая между 0 и 1 величина. Другие описанные Кнутом методы (случайной минимизации Марсальи, прямоугольника-клина-хвоста) для неё утратили значение, поскольку логарифм сейчас считается аппаратно и быстрее, чем по названным сложным алгоритмам, они для очень старых или очень специализированных ЭВМ.
4. Нормально распределённое число можно считать методом Бокса-Мюллера, он же полярных координат (с отбрасыванием попыток, как у Кнута, или через синус-косинус, что на современных процессорах целесообразнее)
http://en.wikipedia.org/wiki/Box%E2%80% ... _transform

[Хорхе]

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

[Евгений Машеров]

Так и то, и то - метод Бокса-Мюллера. Просто когда синус считался только через подпрограмму, придумали ускорение - ценой дополнительных попыток получать такие числа, что они равномерно распределены по кругу.
А теперь синус дёшев, и лучше первоначальный подход.

[Хорхе]

Ну где-то так. Вообще, кажется мне, что для больших значений степеней свободы самый экономный метод где-то в квантильных преобразованиях, а то и вообще использовать нормальную аппроксимацию, и кто там отличит, если вычислительная ошибка будет сравнима с ошибкой аппроксимации. Ну или нормальная аппроксимация + rejection (как это по-русски?), точно и молниеносно, и супермало величин надо сгенерировать.

А для малых значений степеней свободы все равно :)

Но вопрос был про простоту, а не про эффективность. А самый простой метод, как я уже писал, - воспользоваться определением :)

[Евгений Машеров]

Ну, для подачи оптимального совета хорошо бы конкретику. Скажем, каковы ожидаются степени свободы (отчего можно посоветовать использовать суммы или же обратную функцию) или насколько часто будет генерация (отчего можно выбрать между простотой реализации и быстродействием).
Но для решения, хоть и не оптимального, уже всё есть.

[Евгений Машеров]

Ну, если для обратной функции использовать аппроксимацию Голдстейна
http://ru.wikipedia.org/wiki/%D0%9A%D0% ... 0%B0%D1%82
то получается где-то 75-80 операций на число (включая расчёт равномерного с.ч.)
А суммирование требует $k/2$ экспоненциальных, на каждое одно с.ч. (операции 4), логарифм и сумма (условно принимаем, что все операции равноценны, для современных это куда реалистичнее, чем даже 20 лет назад), всего порядка $3k$, то есть при числе степеней свободы $k$ больше 25 есть ускорение от счёта через обратную функцию. Но практическая целесообразность ещё зависит от числа $F$- распределённых с.ч. Если монтекарлим с миллиардами реализаций - надо, если их максимум тысячи - экономии никто не заметит.

[--mS--]

(Оффтоп)

Для непонятливых объясню: автор ждёт, пока ей напишут формулу. Желательно в виде готового кода. Читать литературу или советы ей, судя по другому форуму, неинтересно.