27. Пусть двоичное представление р имеет вид (МЬЬз ... )2- Поступим далее в соответствии со следующими правилами.

81. [Инициализация.] Присвоить m <- iV <- О, j <- 1. (В этом алгоритме т обозначает количество моделируемых равномерно распределенных величин, для которых соотношение с р еще неизвестно, так как их старшие j - 1-двоичные разряды совпадают с этими же разрядами числа р, N - число моделируемых случайных величин, о которых известно, что они меньше р.)

82. [Взгляд на следующий столбец двоичных разрядов.] Генерировать случайное целое число М с биномиальным распределением (тп, 5). (Сейчас М означает число неизвестных случайных величин, таких, что их j-й разряд не совпадает с bj.) Присвоить т < т - М, если bj = 1, то присвоить N N + М.

83. [Сделано?] Если тп = О или если остающиеся двоичные разряды (.bj+ibj+2 ... )2 Р все равны нулю, алгоритм завершен. Иначе - присвоить j <- j -I-1 и возвратиться к шагу В2. I

[Когда bj = 1 для бесконечного числа j, среднее число итераций At удовлетворяет равенствам

Ао=0; " = 1 + E(fc)* для п > 1.

Положив A{z) = EAnz"/n\, получим A{z) = е - 1 + A{z)e Поэтому A{z)e~ = 1 е- + A(i2)e-/ = i:*>o(l - е-") = 1 - е" - E„>i(-)7(n!(2" - 1)) и

в обозначениях упр. 5.2.2-48.]

28. Генерируем случайную точку (у1,...,Уп) на единичной сфере и предположим, что р = УЧкУк- Генерируем независимую равномерно распределенную случайную величину и. Если "+(7 < K\/J2o-lvL то на выходе получится точка {yi/p, ,Уп/р); в остальных случаях начинаем сначала. Здесь К = min{(Eafcyfc)""V(I]ау) = 1} = Оп~\ если па„ > аи {{п + l)/(ai -t-a„))"""(oiOn/n)"-в остальных случаях.

29. Предположим, что Xn+i = 1, затем присвоим Xk <- XfjU или Xk Xk+ie~ для к = п, тг - 1, 1, где Uk - равномерно распределенная случайная величина или Yk - случайная величина с показательным распределением. [АСМ Trans. Math. Software 6 (1980), 359-364. Этот метод введен в употребление в 60-х годах Давидом Сенешолом (David Senescliol); см. работу Amer. Statistician 26,4 (October, 1972), 56-57. Альтернатива состоит в генерировании п равномерно распределенных случайных чисел и их сортировке наиболее быстрым методом. Предложенный здесь метод является особенно полезным, если только требуется несколько наибольших или наименьших Xj. Заметим, что (F~4(Ai), ..., F~(An)) соответствуют упорядоченным случайным величинам с функцией распределения F]

30. Генерировать случайные числа Zi = -р~ ln(7i, Z2 - Z\ - In [/2, ... до тех пор, пока не выполнится Zm+i > 1. На выходе получим {Xj,Yj) = i{Zj) для 1 < j < m, где f{{.bib2 ... b2r)2) = ((.ЬхЬг ... Ьг)2, (.br+ibr+2 • • • b2r)2). Если менее старшие двоичные разряды значительно менее случайны, чем более старший двоичный разряд, то надежнее (но медленнее) положить, что /((.6162... Ь2г)2) = {{-bibs ... b2r-i)2, (.b2b4 ... Ь2г)2)-

31. (а) Достаточно рассмотреть случай, когда к = 2, так как aiXi Н-----\-akXk = Xcos6 +

Ysind при ЛГ = Xl, cos6 = ai и У = (агХг Н-----h akXk)/ sin б. И справедливо равенство

Pr(Xcos6 + ysin6<x)= / e~~*dsdt[scose + tsine<x]

после подстановки и = s cos б + t sin , d = - s sin б + t cos б.

(b) Существуют числа q > 1 и /3 > 1, такие, что {а~* + а~)/\/2 = 1 и §/3" + /3~ = 1, поскольку числа Х„ растут экспоненциально по п согласно свойствам линейных рекуррентных соотношений.

Если отказаться от формы линейного рекуррентного соотношения, скажем, используя рекуррентное соотношение Хп = Хп-24 cos бп + Xn-ss sin бп, где бп выбрано равномерно в [О .. 2п), можно получить подходящий результат, но потребуется выполнить намного больше вычислений.

(c) Начните, пожалуй, с 2 048 нормальных случайных величин Хо, ..., Хюгз, Yo, ..., Yio23- После использования около 1/3 из них генерируйте еще 2 048 случайных величин следующим образом. Выберите целые числа а, Ь, с и d независимо в [О.. 1024), причем а и с должны быть нечетцыми. Затем присвойте

Xj Ч- X(oj>b) mod 1024 COS б + (ej+d) mod 1024 Sin б, Yj i--X(oj4-b) mod 1024 sin б + у(с+Й) mod 1024 COS б

для о < j < 1024, где cos б и sin б - случайные отношения {U - V)/({7 + v) и 2UV/{U + v"), выбранные, как в упр. 23. Можно отбросить (7 и V, кроме случаев, когда I cos б I > I и I sin б I > 5. 2 048 новых случайных величин заменят старые. Заметим, что для получения новой случайной величины необходимо выполнить лишь несколько операций.

Этот метод не расходится подобно последовательностям, содержащимся в (Ь), так как сумма квадратов iXJ + Yj) = ЦХУ + (у/)) остается постоянной со значение.м 5 « 2048, исключая незначительную ошибку округления. С другой стороны, постоянство S на самом деле является недостатком метода, поскольку сумма квадратов действительно должна иметь -распределение с 2 048 степенями свободы. Чтобы преодолеть возникшую проблему, следовало бы использовать не нормальные случайные величины Xj, а aXj, где = 5(У102з -I- 4/4095)5 является заранее вычисленным нормирующим множителем. (Величина (У102з -I- v/4095) будет приемлемо приближать требуемую случайную величину.)

Литература. С. S. Wallace, АСМ Trans, on Math. Software 22 (1996), 119-127; R. P. Brent, Lecture Notes in Сотр. Sci. 1470 (1998), 1-20.

32. (a) Преобразование (х,у) = /(X, У) является взаимно однозначным соответствием преобразованию множества {ж, у > 0} самого в себя, такого, что х + у = х + у и dx dy = dxdy. Получим

xTF = (хТу - ) ""xTF = (хТу + ) ""

(b) Это преобразование является таким соответствием (two-to-one) от двух к одному, что х ->г у = X->г у к dx dy = 2 dx dy.

(c) Достаточно рассмотреть "j-транспонированное" преобразование

х = (. . .XjJr2XjJr\Xjyj-iyj-2yj-3 . . . )2, Y -{... yj+2yj + iyjXj-iXj-2Xj.3 . . . )2

для фиксированных целых j и затем составить -транспонирования для j = О, 1, -1, 2, -2, ..., Принимая во внимание, что совместные вероятностные распределения X и У сходятся при \j\ ->• 00. Каждое j-транспонирование является взаимно однозначным с операциями х + у = х + у и dx dy - dx dy.

33. Использовать Ui как начальное значение для других генераторов случайных чисел (возможен линейный конгруэнтный генератор с другим множителем); в качестве множителя взять одну из С/г, Us,

РАЗДЕЛ 3.4.2

1. Существует {Z способов выбора п - т записей из последних N - t и {Sm\) способов выбора п - тп - 1 записей т N - t - 1 после выбора {t + 1)-й записи.

2. Переход от щага S3 к шагу S5 невозможен, если количество записей, которые осталось проверить, равно п - т.

3. Не будем путать условную и безусловную вероятности. Значение т зависит от случайного выбора первых t элементов. Если взять среднее по всем возможным выборам, которые могут появиться среди этих элементов, то можно найти, что (п - тп)/{М - t) в среднем точно равно n/N. Например, рассмотрим второй элемент. Если первый элемент отобран в выборку (это происходит с вероятностью n/N), то второй элемент выбирается с вероятностью (п - l)/(iV - 1); если первый элемент не выбирается, то второй выбирается с вероятностью n/(iV - 1). Полная вероятность выбора второго элемента равна {n/N){{n - 1)/{N - 1)) -KI - n/N)(n/{N - 1)) = n/N.

4. Из алгоритма следует, что

p(m, t+l)={l- ) р(ш, t) + """:;(m - l,t).

Требуемая формула может быть доказана индукцией по t, в частности р{п, N) = 1.

5. В обозначениях упр. 4 вероятность того, что t = к, по окончании работы алгоритма равна =р{п,к)-р(п,к-1) = (*:J)/(). Среднее равно Y.k=o4k = {N +l)n/{n + l).

6. Так же, как в упр. 5, получим Еь=о + 1)9* = {N -\- 2)(N -t- l)n/(n -I- 2); дисперсия, следовательно, равна {N -\- 1){N - п)п/{п -I- 2)(п -I- 1).

7. Предположим, есть выбор 1 < zi < жг < < Жп < Л. Пусть жо = О, Xn+i = N +1. Выбор получен с вероятностью Р = Y[i<t<NР* "Д®

{N-{t-l)-n + m)/{N -{t- 1)) для Хт <t< Xm+V,

Pt -

\{n-

m)/(A-(t-l)) ДЛЯ* = Жт+1.

Знаменатель произведения pt равен Л!, числитель содержит члены N - п, N - п- 1, 1 для тех tj, которые не равны Xj, а члены п, п-1, ..., 1 для тех tk, которые равны Хк. Следовательно, р = (N - n)\n\/N\.

Пример, п = 3, iV = 8, (Ж1,Ж2,жз) = (2,3,7); р = If 11f IIf

8. (а) p(0,fc) = (;*)/() = (D/(fc) из () выборок, если пропустить первые к записей.

(b) Присвоить X к - 1, где к является минимальным с U > Ft{X > к). Затем выполнять присвоения Л <- О, р <- iV - п, g <- iV, Д <- p/q и т. д. до тех пор, пока и < R. Присвоить Х<-Х--1,р<-р - 1, q <- q- 1, R <- Rp/q. (Этот метод хорош, когда n/N, скажем, > 1/5. Можно предположить, что n/N < 1/2, иначе лучше выбрать N - п невыбранных записей.)

(c) Pr(min(yN,...,yN-n+i) > fc) = П;;оРг(У-; > fc) = U"ZoiN-.i-k)/iN - j). (Этот метод хорош, если, скажем, п < 5.)