AOP_Tom2 (1021737), страница 39
Текст из файла (страница 39)
Марсалья. Здесь исполь- зуется функция Е(х) =его(х/~/2) = )/ — / е '~~Э, х > О, (12) яв которая является функцией распределения абсолюшного значения нормальной случайной величины. Затем Х вычисляется в соответствии с распределением (12). 0.8 0.7 О.б 0.5 0.4 о.з 0.2 ол 0.0 о з 4 Рис. 9. Площадь оод графиком плотности распределения разделена на 31 часть. Площадь каждой части равна среднему числу вычислений случайной величины с такой плотностью. Припишем случайный знак ее значению, и это сделает ее действительно нормальной случайной величиной.
Метод прямоугольника-клина-хвоста основан на важных общих технических приемах, которые будут рассмотрены ниже по мере построения алгоритма. Первая клк>чевая идея — рассматривать г'(х) как смесь нескольких других функций, т. е. записать Е(Х) = Р>т>(Х) + Рзг2(Х) + '' + Рптп(Х), (13) ГДЕ Г>, Гз,..., Ä— ПОДХОДЯЩИЕ РаСПРЕДЕЛЕНИЯ И Р>, Рз,..., Рп — НЕОтРИЦатЕЛЬНЫЕ вероятности, сумма которых раина 1.
Если генерировать случайную переменную Х, выбирая распределение г с нероятностью р„то легко видеть, что Л точно будет иметь г'-распределение. С некоторыми распределениями гз(х), пожалуй, трудно иметь дело, даже труднее, чем с Г, но мы обычно уст'раиваем так, что вероят>>ости р, в этол> случае очень малы.
Большинство распределений Р;(х) будут довольно хорошо устроены, поскольку они будут простой модификацией равномерного распределения. Изложение завершается чрезвычайно эффективной программой, поскольку среднее время счета этой программы очень малб. Рассматриваемый здесь метод легче понять, если работать с производными распределений, а не с самими распределениями.
Пусть П*) = ~" (~), Ь( ) = й" (*) будут плотност ми распределений. Тогда равенство (13) можно записать как 2 (х) = Р>2 > (х) + Р222(х) + ' ' + Раз п(х). (14) Каждая Д> (х) есть > О, и общая площадь под графиком Д(х) равна 1. Поэтому существует подходящий графический метод отображения зависимости (14): площадь под 2(х) разделена на и частей н части, соответствующие ~,(х), имеют площадь р,. На рис. 9 иллюстрируется интересующий нас случай при 1(х) = вп(х) »>2>>хе ' 72; площадь нод этой кривой разделена на п = 31 часть. Существует 15 прямоугольников, представляющих р> у> (х), ..., р>21>з(х), 15 клинообразных частей, представляющих р>аЛе(х), . °, рзоУзо(х), и оставшаяся часть рз> 12>(х) — "хвост", т.
е. график 1(х) при х ) 3. о в+Л Рмс. 10. Плотность распределения, для которой алгоритм 1 может использоваться при генерированни случайных чисел. Прямоугольные части /~(х), ..., /ш(т) представляют равномерное распределение. Например, /з(я) представляет случайную равномерно распределенную величину, лежащую меислу 1 и -'. Высота р„/,(х) равна /(//5); следовательно, площадь у'-го прямоугольника равна —,'уае р = -/(у/5) = ~/ — е 1 У"е для 1 < у < 15. (15) Чтобы генерировать распределение, соответствующее таким прямоугольным частям, просто вычислим Х вЂ” 1г;+д (16) где У равномерно и Я принимает значение (у — 1)/5 с вероятностью р, и не зависит от У. Так как р~ + + ры — —,9183, можно просто использовать равномерно распределенные случайные величины, подобные этим приблизительно в 92% случаев. В оставшихся 8% случаев будем обычно генерировать одно из клинообразных распределений Рш, ..., Гзе.
Типичный пример, который показывает, что необходимо делать, представлен на рис. 10. Когда х < 1, часть кривой вогнутая, а когда х > 1, она выпуклая, но в каждом случае часть кривой достаточно близка к прямой линии и, как показано, может быть заключена между двумя параллельными линиями. Чтобы перебрать эти клиновидные распределения, будем использовать другой общий технический прием: метод ошбрпкоеки фон Неймана получения сложной плотности из другой плотности, которая ее "заключает". Описанный выше метод полярных координат является простым примером такого подхода; шаги Р1-РЗ получают случайную точку внутри единичного круга, генерируя ее в большем круге, отбраковывая ее и начиная снова, егчи точка была вне круга.
Общий метод отбраконки является даже более сильным, чем этот. Пусть нужно генерировать случайную величину Х с плотностью / и пусть д — другая плотность распределения, такая,что /(1) < сд(С) для всех 1, где с — константа. Тогда генерируем случайную величину Х с плотностью д, а также независимую равномерно распределенную случайную величину Г Если У ) /(Х)/сд(Х), отбрасываем Х и начинаем снова с другими Х и Г Когда — +р Рис. 11. Область "принятия гипотезы" алгоритма 1 .
условие ьГ < /(Х)/сд(Л ) в конце концов выполняется, на выходе Х будет иметь требуемую плотность /. (ГДокаэагпельспэао. Х < х произойдет с вероятностью р(х) = [ (д(1) г11 /(1)/сд(1))+др(х), где величина д = ( (д(Е) 41 (1 — /(г)/сд(1))) = 1 — 1/с равна вероятности отбраковки; следовательно, р(х) = [* ДС)г111 Техника отбраковки наиболее эффективна, когда с мвлб, так как должно быть в среднем с итераций, прежде чем значение будет принято (см. упр. 6).
В одних случаях /(х)/сд(х) всегда равно 0 или 1 и нет необходимости генерировать П. В других случаях, если /(х)/сд(х) трудно вычислить, следует постараться 'втиснуть" его между двумя более простыми граничными функциями (18) г(х) < /(х)/сд(х) < в(х) и точное значение /(х) /сд(х) не нужно вычислять, если не выполняется неравенство г(х) ( П < в(х). Следующий алгоритм разрешает проблему клина, совершенствуя метод отбраковки.
Алгоритм Е (Плогпносгпп, блиэкээе к лпнейнььи). Этот алгоритм можно использовать для генерирования случайной величины Л с любым распределением, плотность ,/(х) которого удовлетворяет следующим условиям (см. рис, 10) /(х) =О для х ( в и для х > в+ Ь; (19) а — Ь(х — в)/И < /(х) < 6 — 6(х — в)/И для в < х ( в+ И.
Ь1. [Получить ьг < Ц Генерировать две независимые случайныс величины П и Г, равномерно распределенные между 0 и 1. Если 1Г > 1', заменить П е+ Г, Е2. [Простой случай?) Если 1х < а/6, перейти к шагу Е4. ЕЗ. [Попытаемся снова?] Если 1' > (? + (1/6)/(в+ Ибэ), возвратиться к шагу 1.1. (Если а/6 близко к 1, данный шаг алгоритма нужен не очень часто). 1,4. [Вычисление Х.[ Присвоить Х е — в+ ИП. 1 Когда достигнут шаг Е4, точка ((г, 1г) является случайной точкой на площади, заштрихованной на рис. 11, а именно — 0 < П < 1х < Б+(1/6)/(в+ИВ ).
условия (19) гарантируют, что а, 1 — < П + - /(в+ ИГГ) < 1. Ь Ь Рис. 12. Алгоритм прямоугольника-клина-хвоста для генерирования нормальных случай- ных величин. Сейчас вероятность того, что Х < в + Ьх для 0 < х < 1, .— это площадь, лежащая слева от вертикальной линии (7 = т (см, рис, 11) и деленная нв общую площадь, т.
е. *1 7ге1 -1(в+ йи) ди/ / -1(в+ Ьи)ди = / 1(о) до. , ь / ,/, Ь Следовательно, Х имеет нужное распределение. С подходящими константами авз Ь и в„алгоритм Е можно применять к клинообразным плотностям ~1+~в (см. рис. 9) для 1 < 1 < 15. Последнее распределение Кз~ нуждается в обработке приблизительно один раз из 370; оно используется тогда, когда получаем результат Х ) 3. В упр.
11 показано, что стандартная схема отбраковки может использоваться для этого "хвоста". Рассмотрим процедуру в целом. Алгоритм М (Метод прямоугольника-хлино-хвоста длл нормильно распределенных случайных величин). Для этого алгоритма (рис. 12) используем вспомогательные таблицы (Ро,...,Рв~), (Юм...,1~~в), (1о,тэ~), (ло,,Уз~): %, ~л) Фьв 11зо), (Я~в,,Езо), построенные так, как в упр. 10; примеры выбираются из табл. 1. Предполагается, что используется бинарный компьютер (подобные процедуры могут быть построены для десятичных машин).
М1.(Получить Ц Генерируем равномерно распределенное случайное число У = ( ЬоЬеЬт . Ь~)т. (Здесь Ьд — двоичные разряды в бинарном представлении К Для приемлемой точности г должно быть по крайней мере равно 24.) Присвоить ф в — Ьо. (Позже ф будет использовано для определения знака результата.) М2. (Прямоугольник?) Присвоить с' +- (ЬсЬзЬзЬ4Ьв)с (двоичное число определяется старшими двоичными разрядами Р) и присвоить 1 в — (.ЬвЬг ..
Ьс)з (дробная часть определяется оставшимися двоичными разрядами). Если У > Р,, присвоить Х в — 1' + з'е и перейти к шагу М9. Иначе, если у < 15 (т. е. Ьс — — О), присвоить Х в — 8, + 7Ц1 и перейти к шагу М9. (Это использование метода псевдонимов Уолкера (%а)йег) (3).) М3. (Клин или хвост?] (Сейчас 15 < с < 31 и каждое определенное значение 1 появляется с вероятностью рл) Если у = 31, перейти к шагу М7.
М4. (Получить сс < Ц Генерировать две новые независимые равномерно распределенные случайные величины П и 1', если П > 1; заменить (? в+ Ъ'. (Сейчас выполняется алгоритм 1.) Присвоить Х в — 51 гв + -',(?. М5. (Простой случай?) Если T < Р, перейти к шагу М9. вз -х' М6. (Попытаемся снова?) Если К > (? + Е.(есзс- в х 1сз — 1), вернуться к шагу М4; иначе — перейти к шагу М9. (Вероятность, что этот шаг нужно будет выполнять, мала.) М7. (Получить мажорирующую случайную величину.) Генерировать две новые независимые равномерно распределенные случайные величины П и 1'и присвоить Х '9 — ЬУ М8.
(Отбраковать?) Если ПХ > 3, возвратяться к шагу М7. (Это будет случаться приблизительно в одном из двенадцати случаев после достижения шага М8.) М9. [Присоединить знак.) Если ф = 1, присвоить Х (- -Х. в Этот алгоритм является очень миленьким примером тесного переплетения математической теории с изобретательностью программирования — прекрасная иллюстрация искусства программирования! На выполнение шагов М1, М2 и М9 уходит почти все время; остальные шаги осуществляются намного быстрее. Впервые метод прямоугольннка-клина-хвоста опубликовали Дж.