1626435387-5d0ec7dd22f55dfcc65f3209f497fb0b (844202), страница 26
Текст из файла (страница 26)
подраздел 3.4, статью [23] и книгу[33].Рассмотрим множество Ξ(r) дискретных случайных величин κ̃, принимающих значения из множества целых чисел Z, с фиксированнымматематическим ожиданием r = Eκ̃ (это число предполагается нецелым: r ∈/ Z) и конечной дисперсией (здесь используются обозначения изподраздела 3.4). Требуется найти величину κ, имеющую минимальнуюдисперсию среди всех κ̃ ∈ Ξ(r).УТВЕРЖДЕНИЕ 10.2 [23, 33]. Случайная величина κ принимаетдва значения K (0) = [r] (это целая часть числа r) и K (0) + 1 и имеетраспределение:P κ = K (0) + 1 = p1 = r − K (0) , P κ = K (0) = p2 = K (0) + 1 − r.(10.11)ДОКАЗАТЕЛЬСТВО.УчитываясвойстводисперсииDξ = D(ξ + c); c = const (см., например, [14]), для произвольной величины κ̃ ∈ Ξ(r) имеем1 21 2 1− r − K (0) −. (10.12)= E κ̃ − K (0) −Dκ̃ = D κ̃ − K (0) −222Вычитаемое в правой части (10.12) равно константе.
Для κ̃ = κуменьшаемое в правой части (10.12) равно1 21 211 2 (0)= K −K (0) −p1 + K (0) +1 −K (0) −p2 = .E κ̃−K (0) −2224Для любого другого κ̃ ∈ Ξ(r) с распределением P{κ̃ = z} = p(z) ;z ∈ Z имеем1 2 X 1 2 (z) 1 2 X (z)1E κ̃ − K (0) −=z − K (0) −p >p = ,2224zzтак как 21.221(z)(0)(0)(0)p=1,адляz=6K,K+1выполненоz−K−>z2P158Таким образом, для всех κ̃ ∈ Ξ(r) таких, что κ̃ 6= κ, имеем Dκ̃ > Dκ.Утверждение 10.2 доказано.Для моделирования выборочного значения κ0 случайной величиныκ с распределением (10.11) можно использовать алгоритм 10.3 или вариант формулы С. А. Роженко (10.10), имеющей следующий «компактный» вид: κ0 = [r + α0 ], где α0 ∈ U (0, 1) – стандартное случайное число(к слову, формула (10.10) была «сочинена» именно в этом виде).10.4. Случай большого числа значений.
В случае M 1 (вплотьдо M = ∞) для задания распределения (10.1) вместо конкретных значений {xi } и вероятностей {pi } используются формулы их вычисленияxi = ϕ(i); pi = ψ(i).(10.13)Для получения вероятностей часто более удобными (экономичными)являются рекуррентные формулы видаpi+1 = z(pi ), а конкретнее pi+1 = pi × r(i + 1).(10.14)При реализации алгоритма 10.1 в случае бесконечного числа значений перед вычитанием соответствующей вероятности требуется вычислить ее по одной из формул (10.13) или (10.14).АЛГОРИТМ 10.4 (см., например, [9, 13]). Моделируем стандартноеслучайное число Q := α0 ∈ U (0, 1) и полагаем m := 1 и P := p1 (илиP := ψ(1)). Производим переприсваиваниеQ := Q − P.(10.15)Если новое Q не положительно, то в качестве выборочного значенияξ0 дискретной случайной величины ξ берем ξ0 = ϕ(m) для текущего m;в противном случае полагаем m := m+1, производим пересчет вероятности P := ψ(m) (или P := z(P ), или P := P r(m)) и переприсваивание(10.15) и вновь проверяем Q на положительность и т.
д.Средние затраты алгоритма 10.4 равны!MXt=a+m pm (b + c),(10.16)i=mгде c – средние затраты на пересчет вероятности.В случае, когда пересчет вероятностей происходит по рекуррентнымформулам (10.14) и вероятность p1 задана, число t из равенства (10.16)159уменьшается на величину c. Для целочисленных случайных величин ξс распределением (10.5) при i = 1, 2, . . . формула (10.16) имеет видt1 = a + (b + c)Eη.(10.17)Существует ряд способов понизить трудоемкость (10.16), к числукоторых относится, в частности, расположение (если это возможно) вероятностей pi в порядке убывания (см. утверждение 10.1).В случае, когда пересчет вероятностей по одной из формул (10.13)или (10.14) является трудоемким (т.
е. величина c велика), можно выбрать число M0 так, чтобы сумма вероятностей RM0 = p1 +p2 +. . .+pM0была близка к единице и имелась возможность сохранить в оперативнойпамяти ЭВМ массив значений p0 , p1 , . . . , pM0 . Тогда при α0 < RM0 реализуется алгоритм 10.1 (без пересчета вероятностей), а формулы (10.13)или (10.14) будут использоваться только при α0 ≥ RM0 , т. е.
достаточноредко.Существенно снижают затраты (10.6), (10.16) алгоритмов 10.1 и 10.4для M 1 рассмотренные далее в подразделах 10.6 и 10.7 квантильный метод и метод Уолкера (алгоритмы 10.6 и 10.8). В этих алгоритмахсущественно используется специальный алгоритм моделирования равномерного дискретного распределения (алгоритм 10.5).10.5. Специальный алгоритм моделирования равномерногодискретного распределения. Реализация выборочного значения ξ0дискретной случайной величины ξ с конечным числом значений заметноупрощается, когда все значения x1 , . . . , xM равновероятны, т.
е. в соотношениях (10.1) все pi равны 1/M (такое распределение вероятностейназывается дискретным равномерным). Здесь целесообразно применить следующий специальный алгоритм (см. замечание 10.1).АЛГОРИТМ 10.5 (см., например, [9, 13]). Моделируем стандартноеслучайное число α0 ∈ U (0, 1) и полагаемµ0 = [α0 M ] + 1 = [α0 M + 1](10.18)(здесь [A] обозначает целую часть числа A) и ξ0 = xµ0 .Для обоснования корректности моделирования номера µ0 по формуле (10.18) мы должны показать, что µ = [α M + 1] является целочисленной равномерно распределенной случайной величиной, т.
е. длякаждого k = 1, 2, . . . , M выполнено равенство P{µ = k} = 1/M .Используя определение целой части числа и соотношение (9.6), получаемP{µ = k} = P{[α M + 1] = k} = P{k − 1 ≤ α M < k} =160nk − 1kk ok−11==P≤α<−=,MMMMMчто и требовалось доказать.Для достаточно большого числа значений M преимущество алгоритма 10.5 по сравнению с алгоритмом 10.1 (где требуется порядка M/2вычитаний вида (10.4)) вполне очевидно. Но наши исследования показали, что алгоритм 10.5 имеет преимущество и для малых M , в частности, он эффективней (экономичней) алгоритма 10.3 (здесь M = 2 иp1 = p2 = 1/2).10.6. Квантильный алгоритм.
Бинарный поиск. Алгоритм 10.5позволяет также построить следующую эффективную модификациюстандартного алгоритма 10.1, которая называется квантильным алгоритмом и используется, как правило, для M 1 (вплоть доM = ∞).Зададим целое число K и разобьем интервал (0, 1) на K равныхчастей [(j − 1)/K, j/K), j = 1, . . . , K.Далее построим массив целых чисел {Xj }Kj=1 такой, чтоnj − 1o,Xj = min k : Rk = p1 + p2 + . . .
+ pk ≥K(10.19)который называется массивом нижних квантилей (рис. 10.2).Этот массив задает номер k элемента массива {Ri ; i = 1, 2, . . . , M }, скоторого следует начинать поиск «вверх» (т. е., как и в алгоритмах 10.1 и 10.4, вычитать величины Rq , q = k, k+1, . . . из α0 ∈ U (0, 1) дополучения первого отрицательного значения) при (j−1)/K ≤ α0 < j/K.Рис. 10.2. Схема квантильного алгоритмаОкончательно моделирование дискретной случайной величины выглядит следующим образом.161АЛГОРИТМ 10.6 (см., например, [9, 13]).
1. Моделируем выборочноезначение α0 равномерно распределенной в интервале (0, 1) случайнойвеличины α ∈ U (0, 1).2. Вычисляем номер j0 полуинтервала [(j − 1)/K, j/K), в которыйпопадает α0 по формуле типа (10.18): j0 = [Kα0 + 1].3. Реализуем последовательный поиск «снизу вверх» начиная с RXj0 .Многочисленные тестовые вычисления показали, что при M . N0 /2(здесь N0 – максимальный суммарный объем массивов для компьютера и языка программирования, используемых в вычислениях) числоK квантилей [(j − 1)/K, j/K) следует выбирать так, что M/K ≈ z,где 2 ≤ z ≤ 3. В этом случае затраты алгоритма 10.6 практическине растут с ростом M (а у стандартного алгоритма 10.1 затраты возрастают линейно при увеличении M ).
Эффект достигается за счет использования дополнительной оперативной компьютерной памяти дляформирования массива нижних квантилей (10.19).Последнее обстоятельство особо подчеркивает преимущество алгоритма 10.6 перед другой «популярной» модификацией стандартного алгоритма 10.1, которая носит название метод деления отрезка пополам (или бинарный поиск).АЛГОРИТМ 10.7 (см., например, [9]). 1. Моделируем стандартноеслучайное число α0 ∈ U (0, 1) и полагаем s1 := 1, s2 := M + 1.2. Вычисляем целую часть k := [(s1 + s2 )/2].3.
Если Rk = p1 + . . . + pk < α0 , то полагаем s1 := k, иначе присваиваем s2 := k.4. Если s2 − s1 > 1, то идем на пункт 2, иначе полагаем m = s1 иξ0 = xm .Можно показать, что при 2r < M ≤ 2r+1 для получения выборочного значения xm требуется (r + 1) сравнений.
Таким образом, затратыt̂ алгоритма 10.7 имеют порядок log2 M (это медленнее, чем линейныйрост затрат алгоритма 10.1, но все же рост трудоемкости имеет место).Некоторое преимущество алгоритма 10.7 перед квантильным методом (алгоритм 10.6) и методом Уолкера (алгоритм 10.8) состоит в том,что он не требует использования дополнительной оперативной компьютерной памяти.Особо подчеркнем, что квантильный метод (алгоритм 10.6) работаети в случае бесконечного числа значений случайной величины ξ (т. е. дляM = ∞); здесь применяется пересчет вероятностей по формулам вида(10.13) или (10.14) и выбирается K ≈ N0 /2.162Отметим также, что можно построить эффективную реализацию алгоритма 10.6 в случае, когда число M относительно мало и1N0 > s =.mini=1,...,M piТогда можно взять число квантилей как целую часть K = [s]; приэтом в каждом квантиле [(j −1)/K, j/K) будет не более одного значенияRk и в пункте алгоритма 10.6 потребуется не более одного вычитаниявида (10.4).10.7.
Метод перераспределения вероятностей (алгоритмУолкера). С помощью алгоритма 10.5 можно сконструировать отличный от стандартного альтернативный алгоритм широкого применения – метод Уолкера [34], который применяется для конечных,но больших M , таких, что M . N0 /2, и требует проведения подготовительной работы (иногда весьма затратной).Идею алгоритма проще всего объяснить графически. Изобразим распределение (10.1) в виде диаграммы, состоящей из столбцов единичнойширины и высоты p̃i = pi (рис. 10.3).Рис.