1626435387-5d0ec7dd22f55dfcc65f3209f497fb0b (844202), страница 2
Текст из файла (страница 2)
Численно моделируем (реализуем на компьютере) достаточно большое количество n выборочныхзначений ζ1 , . . . , ζn случайной величины ζ и, используя закон большихчисел (утверждение 1.1), получаем приближение требуемойвеличины:I = Eζ ≈ Zn =Sn, где Sn = ζ1 + . . . + ζn .n(1.1)Повторение (см., например, [14]). УТВЕРЖДЕНИЕ 1.1 (закон большихчисел в форме Хинчина). Пусть ζ (1) , ζ (2) , . . . – это последовательностьнезависимых одинаково распределенных случайных величин, имеющихматематическое ожидание Eζ (i) = Iˆ и Ŝn = ζ (1) + .
. . + ζ (n) . ТогдаŜn /n → Iˆ при n → ∞.В русскоязычной литературе по методам Монте-Карло (см., например, [2, 5–13]) базовая случайная величина ζ называется оценкой величины I (в англоязычной литературе используется термин «estimator»,реже – «estimate»).Таким образом, понятие оценки в этих книгах отличается от соответствующего термина в математической статистике (см., например, [15]),где уже среднее арифметическое Zn = Sn /n из соотношения (1.1) называется оценкой величины I.В данном пособии, во избежание путаницы, случайная величина ζ,математическое ожидание которой равно интересующей нас величине I,будет называться оценивателем (это буквальный перевод англоязычного термина «estimator») или монте-карловской оценкой величины I, в то время как выборочное среднее Zn = Sn /n будет называтьсястатистической оценкой величины Eζ = I.Чаще всего оцениватель ζ из соотношения (1.1) имеет видζ = q(ξ),(d)(1.2)где ξ = ξ (1) , .
. . , ξ– это случайный вектор или отрезок последовательности случайных величин с заданными распределениями, а q(x) –функция (неслучайная) d переменных.Здесь и далее понятие случайный вектор (которое, как правило, неиспользуется в классической теории вероятностей – см., например, [14])обозначает многомерную (конкретнее, d-мерную) случайную величину,состоящую из d одномерных случайных величин (или просто – случайных величин) ξ (j) , j = 1, ..., d.8По поводу сформулированного основного алгоритма 1.1 возникаютследующие вопросы.ПРОБЛЕМА 1.1. Насколько велик класс практически значимых величин I, допускающих представление I = Eζ (см.
соотношение (1.1))?ПРОБЛЕМА 1.2. Является ли выбор оценивателя ζ единственным? Если нет, то как оптимизировать этот выбор?ПРОБЛЕМА 1.3. Какова скорость сходимости алгоритма при росте числа n выборочных значений ζ1 , . . . , ζn ?ПРОБЛЕМА 1.4. Как смоделировать выборочные значения ζ1 , ..., ζnна компьютере?В первых двух разделах данного пособия будут даны ответы на поставленные вопросы.1.5.
Использование обобщенной формулы математическогоожидания. Рассмотрим проблему 1.1. Прежде всего заметим, что дляоценивателя ζ вида (1.2) справедлива следующая обобщенная формуламатематического ожидания (см., например, [14]):ZI = Eζ = q(x)fξ (x) dx,(1.3)где fξ (x) = fξ x(1) , . . . , x(d) – это плотность распределения случайноговектора ξ.Повторение (см., например, [14]).
ОПРЕДЕЛЕНИЕ 1.1. Распределение многомернойслучайной величины (случайного вектора)ξ = ξ (1) , . . . , ξ (d) в области G ⊆ Rd является абсолютно непрерывным, если существует функция fξ (x) = fξ x(1) , . . . , x(d) такая,что для любого борелевского множества A ⊆ Rd справедливо выражение для вероятностиZP{ξ ∈ A} =fξ (x) dx.(1.4)AВ этом случае функция fξ (x) называется плотностью распределения вектора ξ.ЗАМЕЧАНИЕ 1.1. Плотность fξ (x) однозначно определена «почтивсюду» (но в случае, когда предполагается непрерывность или кусочнаянепрерывность, определение функции fξ (x) становится однозначным).
Кроме того, эта функция обладает следующими свойствами:Zfξ (x) dx = 1;(1.5)G9fξ (x) ≥ 0 для x ∈ G и fξ (x) = 0 для x ∈/ G.(1.6)Для представления плотности в дальнейшем будет использоваться записьfξ (x) = S(x), x ∈ G,где S(x) – некоторое аналитическое выражение. В этой записи даетсятолько ненулевая (положительная) составляющая плотности на областираспределения случайного вектора ξ.Также в дальнейшем мы будем достаточно часто использовать следующий одномерный аналог формулы (1.4) для случая, когда G = (a, b), аборелевское множество A представляет собой интервал (c, d):ZP{ξ ∈ (c, d)} =dfξ (x) dx(1.7)c(рис.
1.1).ЗАМЕЧАНИЕ 1.2. Из соотношения (1.7) следует, что с точностью добесконечно малых второго порядка выполнено равенствоP{ξ ∈ dx} = fξ (x) dx,где dx обозначает бесконечно малую окрестность точки x (см. далее соотношение (2.9)). Последнее соотношение иногда используется в качествеопределения вероятностной плотности распределения fξ (x).Рис. 1.1. Иллюстрация к определению плотности случайной величины ξ101.6. Приближенное вычисление интеграла методом МонтеКарло. В связи с соотношением (1.3), с математической точки зрениятеория методов Монте-Карло является специальным разделом теориичисленного интегрирования (см., например, [16]).Действительно, наиболее известным примером использования метода Монте-Карло является задача приближенного вычисленияинтегралаZI=g(x) dx; X ⊆ Rd .XВыберем плотность распределения fξ (x), x ∈ Rd , случайного вектора ξ так, что fξ (x) > 0 для g(x) 6= 0 и выборочные значения ξ 1 , .
. . , ξ nмогут быть достаточно просто смоделированы на компьютере. Тогда изсоотношений (1.1)–(1.3) вытекает следующий алгоритм.АЛГОРИТМ 1.2 (см., например, [9, 13]). Численно приближаем интеграл I следующим образом:!ZZ1 g(ξ 1 )g(ξ n )g(x)fξ (x) dx = Eζ ≈+ ... +,I=g(x) dx =n fξ (ξ 1 )fξ (ξ n )XX fξ (x)(1.8)где ζ = q(ξ) = g(ξ)/fξ (ξ) и ξ 1 , . . . , ξ n – это выборочные значения случайного вектора ξ, смоделированные на компьютере согласно выбранной плотности fξ (x).Равенства (1.8) отражают принцип построения весового оценивателя (монте-карловской оценки) интеграла I.Достаточно очевидным является то обстоятельство, что вариант выбора плотности fξ (x) (а значит, и случайных величин ξ и ζ = q(ξ)) неявляется единственным (см.
выше проблему 1.2).При выборе плотности fξ (x) в алгоритме 1.2 (и случайной величины ζ в общем алгоритме 1.1) используется следующий корректный иестественный критерий оптимальности алгоритмов.КРИТЕРИЙ 1.1. Лучшим (наиболее эффективным, экономичным)является тот алгоритм, который позволяет получить заданный уровень погрешности за меньшее компьютерное время.Из последних рассуждений следует, что сформулированная вышепроблема 1.2 связана с исследованием погрешностей численных схемтипа алгоритмов 1.1 и 1.2, то есть с проблемой 1.3 (см. далее подразделы 1.7, 1.8).11ЗАМЕЧАНИЕ 1.3.
Далеко не во всех исследованиях тех или иныхчисленных алгоритмов используется критерий 1.1. Например, в теории сеточных схем численного интегрирования (см., например, [16]) считаетсялучшим тот алгоритм, который дает более высокую степень k в множителе hk (здесь h – максимальный размер шага по всем компонентам прямоугольной сетки), входящем в верхнюю границу погрешности алгоритма;при этом не учитывается, какой ценой дается увеличение степени k.
Такойподход может быть не совсем корректным с позиций критерия 1.1.ЗАМЕЧАНИЕ 1.4. Проблема приближенного вычисления интегралаи соответствующий алгоритм 1.2 считаются «методически полезными»(но не слишком актуальными с точки зрения современных прикладныхпроблем) применениями метода Монте-Карло. Гораздо более широкометоды численного статистического моделирования применяются для вычисления линейных функционалов от решений интегральныхуравнений Фредгольма второго рода, связанных с прикладнымимарковскими процессами (см., в частности, [9, 13], а также разделы 6–8данного пособия).
В этих задачах соответствующие оцениватели ζ и численные схемы типа алгоритма 1.1 строятся для приближения бесконечныхсумм интегралов видаIh =∞X(K m f, h),(1.9)m=0Z hf y (0) k y (0) , y (1) × . . .XXi (0) (1)(m−1) (m)(m)×k y,yh ydy dy . . . dy (m−1) dy (m) ,(K m f, h) =Z...(1.10)где f (y), h(y), k(x, y) – некоторые заданные функции (см. раздел 7).1.7. Построение доверительного интервала с помощьюцентральной предельной теоремы. Займемся вопросом оценки погрешности алгоритма 1.1, которая имеет вид Sn − nI δn = |I − Zn | = n(см.
проблему 1.3). Здесь можно воспроизвести следующие простые рассуждения из математической статистики (см., например, [15]).Погрешность δn может быть представлена в виде √ Sn − nI n − nEζ . = √Dζ S√(1.11)δn = √nn Dζ n 12Рис. 1.2. Иллюстрация к «правилу трех сигма»Повторение (см., например, [14]). УТВЕРЖДЕНИЕ 1.2 (центральнаяпредельная теорема для одинаково распределенных случайных величин).Пусть ζ (1) , ζ (2) , . . .
является последовательностью независимых одинаково распределенных случайных величин, которые имеют конечноематематическое ожидание Eζ (i) = Iˆ и конечную отличную от нулядисперсию 0 < Dζ (i) = σ 2 < +∞. Рассмотрим также последовательность случайных величинΞ(n) =Ŝn − nIˆ√ , где Ŝn = ζ (1) + . . . + ζ (n) .σ nТогда Ξ(n) сходится по распределению при n → ∞ к стандартнойнормальной (гауссовской) случайной величине ξ (0,1) , имеющейплотность распределенияu2(0,1)fξ (u)e− 2= √ , −∞ < u < +∞2π(1.12)(рис. 1.2).Сходимость по распределению означает сходимость функций распре-13деления:FΞn (x) = P{Ξn < x} → Φ(x) =(0,1)Fξ (x)= P ξ (0,1) < x =Zx−∞u2e− 2 du√.2πСтандартная нормальная (гауссовская) величина ξ (0,1) известна какслучайная величина, практически не принимающая больших значений.В частности, на практикешироко применяетсятрех сигма»,√R +3 «правило(0,1)означающее, что P − 3 < ξ< 3 = −3 (1/ 2π) exp(−u2 /2) du ≈0, 997 (рис.
1.2).Из соотношения (1.11) и утверждения 1.2 получаем, что для любогомалого ε > 0 найдется константа Hε такая, что для достаточно большого n 1 может быть построен следующий доверительный интервалдля погрешности δn (см., например, [15]):√noDζP δn ≤ Hε √≈ P ξ (0,1) < Hε ≥ 1 − ε.(1.13)nНапример, для ε = 0, 003 имеем Hε ≈ 3 (что соответствует «правилутрех сигма»).1.8. Низкая скорость сходимости и «универсальность»метода Монте-Карло.
Соотношение (1.13) означает, что скоростьсходимости метода Монте-Карло по числу n выборочных значений {ζi } (см. алгоритмы 1.1 и 1.2) определяется формулой n−1/2(это ответ на вопрос сформулированной выше проблемы 1.3).Таким образом, эта скорость сходимости относительно невысока(это определяет основной недостаток метода Монте-Карло). Чтобы получить следующий десятичный знак в мантиссе величины I (т. е. уменьшить примерно в десять раз погрешность δn ), нужно численно смоделировать в сто раз больше выборочных значений {ζi }. Поэтому в практических расчетах число n берется весьма большим, а учитывая тот факт,что получение каждого ζi , как правило, является трудоемким (см., например, [9, 13], а также следующий раздел 2), можно констатировать,что методы Монте-Карло являются весьма «тяжеловесными» (затратными).Для сравнения заметим, например, что уровень погрешности выRbчисления однократного интеграла I = a g(x) dx (здесь d = 1) с гладкой подынтегральной функцией g(x) с помощью простейшего метода14центральных прямоугольниковI ≈h×nx + x Xb−aii+1; xi = a + (i − 1)h; h =g2ni=1(1.14)(это вариант интегральной суммы – рис.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
