Авторегресійний процес k-го порядку з постійними коефіцієнтами визначається рівнянням регресії
|
(4.12)
|
Значення процесу (4.12) у будь-який момент часу t визначається через попередні значення та випадкове збурення et. На практиці звичайно використовують авторегресійні моделі процесів першого і другого порядку (процес Маркова і Юла-Уокера), автокореляційна функція яких є згасаючою (рис. 4.27).
Рис. 4.27. Автокореляційна функція стаціонарного процесу
Параметри процесу а1, ..., ak визначаються через коефіцієнти автокореляції. Так, для процесу Юла
де  — значення автокореляційної функції при зсувах 1 та 2.
Під час побудови рівняння авторегресії висуваються дві гіпотези. Перша — про стаціонарність процесу, друга — про те, що збудження еt є випадковим процесом у широкому розумінні слова з нормальною функцією розподілу, нульовим математичним сподіванням і дисперсією 2. На рис. 4.28, а зображено графік випадкового процесу і нормальний розподіл збурення еt (рис. 4.28, б).
Рис. 4.28. Графіки випадкового процесу (а) і нормального розподілу збурення (б)
Процес yt називається марківським, якщо для будь-яких моментів часу t1 < t2 < t3 < ... < tn умовна ймовірність значення yt залежатиме від yt-1 і не залежатиме від того, в якому стані процес знаходився в попередні моменти часу. Під час моделювання марківського стаціонарного процесу з параметрами  діють таким чином.
За початковий член ряду можна взяти будь-яке значення випадкового процесу (тому що будь-яка частина стаціонарного процесу є повноцінним представником усього процесу і має ті ж імовірнісні характеристики), наприклад  або  .
За формулою (4.12) розраховуємо значення  при заданих  без урахування збурення et і моделюємо значення нормально розподіленої випадкової величини з математичним сподіванням, що дорівнює нулю, і дисперсією 2.
Отримане значення додаємо до  і таким чином отримуємо нове значення реалізації випадкового процесу.
Повторюємо процедуру для обчислення інших значень за формулою (4.12), задаючи як початкове значення  , тобто моделюємо  .
Така методика дає змогу моделювати випадкові стаціонарні процеси з будь-якими автокореляційними та багатомірними функціями розподілів.
4.8. Статистична обробка результатів моделювання
Основою для обчислення статистичної оцінки параметра системи є реалізація випадкової величини, яка формується під час прогонів імовірнісної імітаційної моделі. Статистична оцінка також є функцією від випадкових величин, які дістають у результаті прогонів моделі, тому і згадана оцінка є випадковою величиною, закон розподілу якої залежить від закону розподілу досліджуваної випадкової величини та оцінюваного параметра. Чим більше реалізацій випадкової величини, тим точнішу статистичну оцінку параметра системи ми отримуємо. У таких умовах обробка результатів моделювання повинна провадитись лише з використанням методів і алгоритмів, які є оптимальними з погляду затрат часу та використання ресурсів комп'ютера. Під час вибору таких засобів необхідно враховувати, що всі статистичні оцінки мають бути ще й незміщеними, ефективними і спроможними.
Розглянемо деякі методи обчислення основних статистичних оцінок.
4.8.1. Оцінювання ймовірності
Оцінкою ймовірності р настання деякої події А є її частість:
де m — кількість випробувань, під час за яких випадкова подія спостерігалась; N — загальна кількість випробувань.
Для її використання зазвичай на програмному рівні організовують два лічильники, один з яких призначено для підрахунку загальної кількості випробувань N, а другий — кількості успішних випробувань m.
4.8.2. Оцінювання розподілу випадкової величини
Для оцінювання функції розподілу випадкової величини, як звичайно, будується гістограма. Під час її побудови область можливих значень випадкової величини розбивають на n діапазонів і підраховують кількість попадання значень випадкової величини в конкретний інтервал — mk (k = 1, ..., n). Оцінка ймовірності попадання випадкової величини в k-й інтервал має такий вигляд:
Цю величину називають відносною частістю. У процесі моделювання під час підрахунку значень mk кожному інтервалу ставлять у відповідність окремий лічильник.
4.8.3. Оцінювання математичного сподівання
Для оцінювання математичного сподівання випадкової величини використовується формула
де xk — значення випадкової величини, що належить k-му інтервалу;
mk — кількість попадань значень випадкової величини в інтервал;
N — загальна кількість випробувань.
У більш простому випадку для оцінювання математичного сподівання випадкової величини можна використати звичайне середнє арифметичне:
Щоб запобігти непотрібному завантаженню пам'яті, суму доцільніше підраховувати шляхом поступового накопичення.
4.8.4. Оцінювання дисперсії
Для оцінювання дисперсії випадкової величини можна використати формулу
або в простому випадку
де S2 — оцінка дисперсії випадкової величини х.
Безпосередньо використовувати ці формули в розрахунках дисперсії нераціонально, тому що зі збільшенням кількості значень xі змінюється також середнє значення, а для його обчислення потрібно запам'ятовувати всі N значень xі. Тому доцільніше використовувати таку формулу:
У цьому випадку достатньо накопичувати тільки суми двох послідовностей —  і  . Але і такий спосіб має недолік — його використання може призвести до переповнення розрядної сітки комп'ютера. Щоб запобігти цьому, потрібно змінити послідовність дій при обчисленнях, використовуючи формулу
4.8.5. Оцінювання кореляційного моменту
Для обчислення оцінки кореляційного моменту можна використовувати формулу
або, більш зручну для обчислень,
При обчисленнях за цією формулою теж доцільно змінити послідовність дій.
4.9. Визначення кількості реалізацій під час моделювання випадкових величин
Точність оцінок параметрів системи, які отримують під час обробки результатів моделювання, у першу чергу залежить від кількості випробувань N. Слід враховувати, що обсяг вибірки N завжди обмежений, тому вищезгадані оцінки матимуть різні похибки і дисперсії.
Якщо треба оцінити значення деякого параметра а за результатами моделювання xі, то за його оцінку слід брати величину  , яка є функцією від усіх значень xі. Статистична оцінка також є випадковою величиною, тому вона буде відрізнятись від а, тобто
де — точність або похибка оцінки. Імовірність того, що ця нерівність виконується, позначимо через :
|
(4.13)
|
У теорії ймовірностей — це довірчий інтервал для a, довжина якого фактично дорівнює 2, а — довірчий рівень, або надійність оцінки. Вираз (4.13) можна застосувати для визначення точності результатів статистичних випробувань.
4.9.1. Оцінювання ймовірності
Припустимо, що метою моделювання є оцінювання ймовірності настання деякої події А, яка визначає стан системи. У кожній з N реалізацій процесу настання події А є випадковою величиною , що набуває значення х1 = 1 з імовірністю р і х2 = 0 з імовірністю 1 – р. Тоді можна визначити математичне сподівання і дисперсію відповідно за формулами
|
(4.14)
|
|
(4.15)
|
Як оцінку р використовують частість настання події А. Ця оцінка є незміщеною, спроможною та ефективною. За умови, що N задано, для отримання цієї оцінки достатньо накопичувати m:
|
(4.16)
|
де хi – настання події А в реалізації і.
За формулами (4.14 – (4.16) визначимо вибіркове математичне сподівання  і дисперсію  .
Згідно з центральною граничною теоремою (у даному випадку її можна взяти у вигляді теореми Хінчина) випадкова величина m/N буде мати розподіл, близький до нормального (рис. 4.27). Тому для кожного рівня достовірності з таблиць нормального розподілу можна знайти таку величину t, при якій точність обчислюватиметься за формулою
|
(4.17)
|
Якщо = 0.05, то t = 1.96, а якщо = 0.003, то t = 3.
Підставимо у формулу (4.17) вираз дисперсії:
Звідси
|
(4.18)
|
Рис. 4.27. Функція нормального розподілу для побудови довірчого інтервалу
З формули (4.18) видно, що при р = 1 або р = 0, кількість реалізацій, які необхідно провести для підтвердження того, що подія А настає (або ні), дорівнює одиниці. Але оскільки ймовірність р заздалегідь невідома, провадять випробування (N = 50 ... 100), оцінюють частість m/N і підставляють її значення у вираз (4.18) замість р, після чого визначають остаточну кількість реалізацій. Графік залежності числа реалізацій для = 0.05 і різних значень р, якщо = 0.05, наведено на рис. 4.29.
Рис. 4.28. Залежність числа реалізацій від значень імовірності
4.9.2. Оцінювання середнього значення
Нехай випадкова величина має математичне сподівання а і дисперсію 2. У і-й реалізації вона набуває значення xi. Як оцінку математичного сподівання а використаємо середнє арифметичне:
|
(4.19)
|
Згідно з центральною граничною теоремою при великих значеннях N середнє арифметичне (4.19) буде мати нормальний розподіл з математичним сподіванням а і дисперсією  . Тоді
Звідси
|
(4.20)
|
Оскільки дисперсія 2 випадкової величини невідома, потрібно провести кілька десятків (50 ... 100) випробувань і знайти оцінку  , а потім отримане значення підставити у формулу (4.20), щоб визначити необхідну кількість реалізацій N. У цьому випадку замість нормально розподіленої величини необхідно скористатись t-розподілом Стьюдента з N – 1 ступенями вільності для визначення t. Зауважимо, що за збільшення ступенів вільності t-розподіл наближається до нормального. З практичного погляду, якщо N більше 30, користуються нормальним розподілом.
Висновки
Метод статистичних випробувань визначається як спосіб побудови і дослідження на комп'ютері моделі системи або процесу з використанням послідовностей випадкових чисел.
Метод полягає в багатократному проведенні випробувань побудованої моделі й подальшій статистичній обробці результатів моделювання з метою визначення шуканих характеристик розглядуваного процесу у вигляді оцінок його параметрів.
У методі статистичних випробувань особливе значення відіграють випадкові числа, рівномірно розподілені в інтервалі [0, 1], за допомогою яких можна отримати вибіркові значення з будь-якими статистичними властивостями.
Для генерування випадкових чисел використовуються апаратні, табличні та програмні методи.
Апаратні методи генерування випадкових чисел базуються на використанні деяких фізичних явищ і процесів - випадковий електричний сигнал перетворюють у двійковий код, який уводиться в комп'ютер за допомогою спеціальних аналого-цифрових перетворювачів.
У разі використання табличного методу випадкові числа можна зберігати на зовнішніх носіях або навіть в основній пам'яті комп'ютера.
Програмні генератори дають змогу отримувати послідовності випадкових чисел за рекурентними формулами. Більшість програмних генераторів виробляють випадкові числа, рівномірно розподілені в інтервалі [0, 1].
Серед програмних генераторів найбільш розповсюдженим є лінійний мульти- плікативний конгруентний генератор.
Основним параметром програмного генератора є повний період — кількість чисел, після якої випадкові числа починають повторюватися.
Для перевірки статистичних властивостей усіх послідовностей випадкових чисел, які будуть використовуватись під час проведення досліджень, застосовують емпіричні та теоретичні критерії.
Емпіричні критерії — це звичайні статистичні тести, в яких використовують вибіркові значення, що виробляються генератором.
Теоретичні критерії визначають характеристики послідовності за допомогою методів, які формуються на основі числових значень параметрів генератора.
Контрольні запитання та завдання
Доведіть, що незалежні випадкові величини ri та 1 – ri отримані з рівномірного розподілу в інтервалі [0, 1], мають однаковий розподіл.
Обчисліть імовірність того, що за 50 підкидань монети 10 разів випаде герб. Спочатку розв'яжіть задачу аналітично, а потім порівняйте результат з оцінкою, отриманою за допомогою методу статистичних випробувань.
Використовуючи метод статистичних випробувань, обчисліть площу круга. Зважаючи на те, що формула для визначення площі круга –  , оцініть значення константи .
Молодий пілот, виконуючи перший свій політ, не зміг здійснити вдалу посадку літака. Визначте число спроб зайти на посадку, які необхідно зробити пілоту, для того щоб політ благополучно закінчився. Ймовірність того, що за одну спробу стажист вдало завершить посадку літака, становить 0,1. Кількість палива в баках вважати необмеженою.
Скориставшись методом статистичних випробувань, знайдіть оцінку інтеграла
як площі під кривою  на відрізку 0 х 1.
Імовірність одержання заліку студентом, який не відвідував лекційні і практичні заняття, становить 0,13. Скориставшись методом статистичних випробувань, знайдіть оцінку ймовірності того, що студент одержить залік, якщо загальне число спроб здати залік не може перевищувати 3.
Знайдіть оцінку площі рівнобічного трикутника зі стороною 1 см. Визначіть рівняння прямих, які утворюють сторони трикутника і використайте ці вирази в процедурі оцінювання попадання «випадкової» точки в трикутник для методу статистичних випробувань.
Для одержання допуску до іспиту студенту необхідно отримати залік з лабораторних робіт (лабораторний курс складається з К лабораторних робіт). Припустимо, що ймовірність здавання однієї лабораторної роботи студентом становить 0,3 з однієї спроби. Оцініть кількість днів, потрібних студенту на одержання допуску до іспиту, якщо за один день він зможе здати не більше однієї лабораторної роботи (протягом семестру він не захистив жодної роботи).
Стрільцю необхідно вразити чотири мішені. Ймовірність ураження однієї мішені становить 0,14. Оцініть кількість набоїв, необхідних для того, щоб уразити всі чотири мішені.
Багдадський злодій ув'язнений у підземелля з трьома дверима [46]. Одні двері ведуть на волю, другі — у довгий тунель, а треті — у короткий. Потрапивши в один із тунелів, злодій знову опиняється в темниці. Він знову пробує вийти на волю, але не пам'ятає, в які двері входив минулого разу. Імовірність того, що злодій обере потрібні двері, дорівнює 0,3, ймовірність попадання в короткий тунель — 0,2; ймовірність попадання в довгий тунель — 0,5. Час перебування злодія в короткому тунелі — 3 хв, у довгому — 6 хв. Визначіть середній час пошуку шляху на волю. Побудуйте процедуру статистичних випробувань для розв'язування задачі.
Змоделюйте поведінку винищувача-бомбардувальника [67], який атакує об'єкт ракетами класу «повітря-земля». Кожна ракета наводиться індивідуально. Розміри об'єкта – 60150 м. Заходження на атаку відбуваються в напрямі, який збігається з напрямом довгої осі цілі, точка прицілу — геометричний центр цілі. Фактичну точку влучення для кожної ракети можна визначити горизонтальним і вертикальним відхиленнями (рис. 4.29). Для відстані, з якої запускають ракети, обидва відхилення є незалежними, нормально розподіленими щодо точки прицілювання і мають нульове середнє значення.
Рис. 4.29. Схема влучення ракети в ціль
Середньоквадратичне відхилення становить 60 м у напрямку X і 30 м — у напрямку Y. Бомбардувальник під час кожного заходження випускає шість ракет. Узявши обсяг вибірки в 10 заходжень, знайдіть оцінку середнього числа влучень для кожної атаки.
Змоделюйте випадкову двомірну дискретну величину (X, Y). Випадкова величина X може набувати значення 1, 2, 3, 4, а випадкова величина Y — 10, 20, 30, 40, у цьому разі кожній парі значень хi, yi відповідає ймовірність pij (табл. 4.5).
Таблиця 4.5. Ймовірність двомірної дискретної величини
Значення хi
|
Значення yi
|
|
10
|
20
|
30
|
40
|
1
|
0,02
|
0,05
|
0,1
|
0,05
|
2
|
0,03
|
0,1
|
0,05
|
0,03
|
3
|
0,01
|
0,15
|
0,1
|
0,02
|
4
|
0,04
|
0
|
0,15
|
0,1
|
Припустимо, що об'єм споживання води в місті є випадковою величиною, що має нормальний розподіл [67]. Знайдіть оцінку середнього споживання води в день так, щоб помилка не перевищувала 6000 л з імовірністю 0,95. Відомо, що розумна область розкиду споживання води становить 120000 л за день. Який обсяг вибірки необхідний для цього дослідження?
Використовуючи результати 20 імітаційних прогонів для оцінювання часу перебування відвідувачів у системі [46], що наведені в дужках (1,1; 2,8; 3,7; 1,9; 4,9; 1,6; 0,4; 3,8; 1,5; 3,4; 1,9; 2,1; 3,8; 1,6; 3,2; 2,9; 3,7; 2,0; 4,2; 3,3), обчисліть оцінки для вибіркового середнього, дисперсії і коефіцієнта варіації. Побудуйте гістограму, яка містить п'ять інтервалів (довжина кожного інтервалу дорівнює одиниці, нижня межа першого інтервалу дорівнює нулю).
Дано випадкові некорельовані змінні А, В, С. Змінна A має нормальний розподіл з математичним сподіванням 100 і середньоквадратичним відхиленням 20. Змінна В також розподілена нормально з математичним сподіванням 20 і середньоквадратичним відхиленням 5. Розподіл змінної С задано в табл. 4.6.
Таблиця 4.6. Функція розподілу випадкової величини
Значення С
|
10
|
20
|
30
|
40
|
Імовірність
|
0,10
|
0,25
|
0,50
|
0,15
|
Застосовуючи метод статистичних випробувань, оцініть середнє значення нової змінної D що визначається так: D = (А + В)/С. Використайте вибірку із 10 значень, яку необхідно отримати за допомогою розподілу, що демонструється на рис. 4.15.
Проаналізуйте модель наземної протиповітряної ракетної установки. Скільки необхідно виконати випробувань, кожне з яких моделює одну відсіч повітряної атаки, щоб оцінка частоти поразки цілі відрізнялася від істинної ймовірності не більше ніж на 0,02 з імовірністю не меншою за 0,85?
Ме́тод Мо́нте-Ка́рло (методы Монте-Карло, ММК) — общее название группы численных методов, основанных на получении большого числа реализаций стохастического (случайного) процесса, который формируется таким образом, чтобы его вероятностные характеристики совпадали с аналогичными величинами решаемой задачи. Используется для решения задач в различных областях физики, химии, математики, экономики, оптимизации, теории управления и др.
Случайные величины использовались для решения различных прикладных задач достаточно давно. Примером может служить способ определения числа Пи, который был предложен Бюффоном в 1777 году. Суть метода была в бросании иглы длиной  на плоскость, расчерченную параллельными прямыми, расположенными на расстоянии  друг от друга (см. Рис. 1).
Рисунок 1. Метод Бюффона
Вероятность (как видно из дальнейшего контекста, речь идёт не о вероятности, а о математическом ожидании количества пересечений за один опыт; вероятностью это становится лишь при условии, что  ) того, что отрезок пересечет прямую, связана с числом Пи:
 , где
 — расстояние от начала иглы до ближайшей к ней прямой;
 — угол иглы относительно прямых.
Этот интеграл просто взять:  (при условии, что ), поэтому подсчитав долю отрезков, пересекающих прямые, можно приближенно определить это число. При увеличении количества попыток точность получаемого результата будет увеличиваться.
В 1864 году капитан Фокс, выздоравливая после ранения, чтобы как-то занять себя, реализовал эксперимент по бросанию иглы. Результаты представлены в следующей таблице:
|
Число бросаний
|
Число пересечений
|
Длина иглы
|
Расстояние между прямыми
|
Вращение
|
Значение Пи
|
Первая попытка
|
500
|
236
|
3
|
4
|
отсутствует
|
3.1780
|
Вторая попытка
|
530
|
253
|
3
|
4
|
присутствует
|
3.1423
|
Третья попытка
|
590
|
371
|
5
|
2
|
присутствует
|
3.1416
|
Комментарии:
Вращение плоскости применялось (и как показывают результаты — успешно) для того, чтобы уменьшить систематическую ошибку.
В третьей попытке длина иглы была больше расстояния между линиями, что позволило не увеличивая числа бросаний эффективно увеличить число событий и повысить точность.
Создание математического аппарата стохастических методов началось в конце XIX века. В 1899 году лорд Релей показал, что одномерное случайное блуждание на бесконечной решётке может давать приближенное решение параболического дифференциального уравнения. Андрей Николаевич Колмогоров в 1931 году дал большой толчок развитию стохастических подходов к решению различных математических задач, поскольку он сумел доказать, что цепи Маркова связаны с некоторыми интегро-дифференциальными уравнениями. В 1933 году Иван Георгиевич Петровский показал, что случайное блуждание, образующее Марковскую цепь, асимптотически связано с решением эллиптического дифференциального уравнения в частных производных. После этих открытий стало понятно, что стохастические процессы можно описывать дифференциальными уравнениями и, соответственно, исследовать при помощи хорошо на тот момент разработанных математических методов решения этих уравнений.
Рождение метода Монте-Карло в Лос-Аламосе
Сначала Энрико Ферми в 1930-х годах в Италии, а затем Джон фон Нейман и Станислав Улам в 1940-х в Лос-Аламосе предположили, что можно использовать связь между стохастическими процессами и дифференциальными уравнениями «в обратную сторону». Они предложили использовать стохастический подход для аппроксимации многомерных интегралов в уравнениях переноса, возникших в связи с задачей о движении нейтрона в изотропной среде.
Идея была развита Уламом, который раскладывая пасьянсы во время выздоровления после болезни задался вопросом, какова вероятность того, что пасьянс сложится. Вместо того, чтобы использовать обычные для подобных задач соображения комбинаторики, Улам предположил, что можно просто поставить эксперимент большое число раз и, подсчитав число удачных исходов, оценить вероятность. Он же предложил использовать компьютеры для расчётов методом Монте-Карло.
Появление первых электронных компьютеров, которые могли с большой скоростью генерировать псевдослучайные числа, резко расширило круг задач, для решения которых стохастический подход оказался более эффективным, чем другие математические методы. После этого произошёл большой прорыв и метод Монте-Карло применялся во многих задачах, однако его использование не всегда было оправдано из-за большого количества вычислений, необходимых для получения ответа с заданной точностью.
Годом рождения метода Монте-Карло считается 1949 год, когда в свет выходит статья Метрополиса и Улама «Метод Монте-Карло». Название метода происходит от названия коммуны в княжестве Монако, широко известного своими многочисленными казино, поскольку именно рулетка является одним из самых широко известных генераторов случайных чисел. Станислав Улам пишет в своей автобиографии «Приключения математика», что название было предложено Николасом Метрополисом в честь его дяди, который был азартным игроком.
Дальнейшее развитие и современность
В 1950-х годах метод использовался для расчётов при разработке водородной бомбы. Основные заслуги в развитии метода в это время принадлежат сотрудникам лабораторий ВВС США и корпорации RAND.
В 1970-х годах в новой области математики — теории вычислительной сложности было показано, что существует класс задач, сложность (количество вычислений, необходимых для получения точного ответа) которых растёт с размерностью задачи экспоненциально. Иногда можно, пожертвовав точностью, найти алгоритм, сложность которого растёт медленнее, но есть большое количество задач, для которого этого нельзя сделать (например, задача определения объёма выпуклого тела в n-мерном евклидовом пространстве) и метод Монте-Карло является единственной возможностью для получения достаточно точного ответа за приемлемое время.
В настоящее время основные усилия исследователей направлены на создание эффективных Монте-Карло алгоритмов различных физических, химических и социальных процессов для параллельных вычислительных систем.
|