Стохастический мир
Скачать 2.82 Mb.
|
T 0 x(t) t где T произвольное число. Причина неоднозначности - в нарушении условия Липшица (бесконечность производной a 0 (x) = 2/x 1/3 в x = 0). В реальном Мире, если некоторая система описывается ( 5.29 ), то она не сдвинется из начального состояния x(0) = 0, если решает уравнение итерациями. Спустя произвольный момент времени T может произой- ти внешнее воздействие (флуктуация), которое столкнјт итерационный процесс с нулевой отметки. Это и приведјт к неоднозначному решению (l C 25 ). • Кроме неоднозначности, с решениями уравнений может произойти ещј одна неприятность. Рассмотрим следующую задачу: dx dt = x 2 => x(t) = x 0 1 ? (t ? t 0 )x 0 Через конечное время t ? t 0 = 1/x 0 от начального момента решение обращается в бесконечность. Подобная ситуация называется взрывом решения. За редкими исключениями подобные модели не соответствуют реальным системам или являются лишь их грубым приближением. • С начальным условием дифференциального уравнения связана одна тонкость. Не любая функция x = f(x 0 , t 0 , t) со значением x 0 = f (x 0 , t 0 , t 0 ) удовлетворяет хоть какому-то дифференциальному уравнению первого порядка. Так, в производной по времени функции x = x 0 + sin(t ? t 0 ) ни- какими заменами и выбором a(x, t) не удастся одновременно избавить- ся и от x 0 , и от t 0 . Подставляя в уравнение решение x = f(x 0 , t 0 , t) , мы должны так его преобразовать, чтобы константы x 0 , t 0 , являющиеся внешними к уравнению начальными условиями, сократились. 146 Глава 5. • Перейдјм теперь к стохастическим дифференциальным уравнениям. Так как мы работаем со случайными процессами, различные реализа- ции траектории x(t) могут отличаться сколь угодно сильно. Говоря о единственности решения, мы подразумеваем единственность плотности вероятности P (x 0 , t 0 ? x, t) , которая влечјт за собой единственность среднего значения, волатильности, автоковариационной функции и т.д. Докажем, что для уравнения dx = a(x, t) dt + b(x, t) ?W решение будет единственным, если производные по x сноса a(x, t) и вола- тильности b(x, t) непрерывны. Непрерывность означает, что производные ограничены (нигде не обращаются в бесконечность): ?a(x, t) ?x 6 M a , ?b(x, t) ?x 6 M b Поэтому в силу формулы Лагранжа каждая из них удовлетворяет нера- венству Липшица: |a(y, t) ? a(x, t)| 6 M a |y ? x|, |b(y, t) ? b(x, t)| 6 M b |y ? x|. (5.30) Будем действовать так же, как и в детерминированном случае. Перейдјм к интегральному уравнению: x(t) = x 0 + t Z t 0 a x(s), s ds + t Z t 0 b x(s), s ?W s Пусть существуют две разные случайные функции x t = x(t) и y t = y(t) с одинаковым начальным условием x(t 0 ) = y(t 0 ) = x 0 , которые удовле- творяют этому уравнению. Найдјм дисперсию их разности: (y t ? x t ) 2 = D ? ? t Z t 0 a yx (s) ds + t Z t 0 b yx (s) ?W s ? ? 2 E , где a yx (s) = a y(s), s ? a x(s), s , b yx (s) = b y(s), s ? b x(s), s - разно- сти сноса или волатильности, вычисленные для каждого решения. Стохастические интегралы 147 Для двух n - мерных векторов {? 1 , ..., ? n } и {? 1 , ..., ? n } скалярное про- изведение всегда меньше, чем произведение их длин (косинус угла между ними меньше единицы). Поэтому справедливо неравенство Коши - Бу- няковского: (? 1 ? 1 + ... + ? n ? n ) 2 6 (? 2 1 + ... + ? 2 n ) (? 2 1 + ... + ? 2 n ). Если все ? i = 1 , имеем такой вариант этого неравенства: (? 1 + ... + ? n ) 2 6 n (? 2 1 + ... + ? 2 n ). В нашем случае n = 2, поэтому: (y t ? x t ) 2 6 2 D ? ? t Z t 0 a yx (s) ds ? ? 2 E + 2 D ? ? t Z t 0 b yx (s) ?W s ? ? 2 E Для первого интеграла (точнее, его интегральной суммы) снова приме- ним неравенство Коши-Буняковского c n = (t?t 0 )/?s . Среднее значение квадрата стохастического интеграла по ?W можно записать через сред- нее от обычного интеграла по времени ( 5.13 ), стр. 132 , поэтому: (y t ? x t ) 2 6 2(t ? t 0 ) t Z t 0 D a 2 yx (s) E ds + 2 t Z t 0 D b 2 yx (s) E ds. Воспользуемся теперь неравенствами Липшица ( 5.30 ), возведя их в квад- рат. В результате: (y t ? x t ) 2 6 M t Z t 0 (y t ? x t ) 2 ds, где M = 2(t ? t 0 )M 2 a + 2M 2 b . Среднее разности решений (y t ? x t ) 2 это обычная функция времени, поэтому, применяя лемму Гронуолла Беллмана ( 5.26 ), приходим к выводу, что (y t ? x t ) 2 = 0 Среднее является интегралом с положительной плотностью вероят- ности. Величина (y t ? x t ) 2 также положительна. Интеграл от положи- тельной функции может быть равен нулю, только когда она равна нулю. В результате мы приходим к выводу, что x(t) = y(t), и решение един- ственно. 148 Глава 5. 5.6 Метод последовательных приближений При разложении средних величин в ряд по t n (стр. 95 ) мы упоминали метод последовательных приближений как один из итерационных спосо- бов построения решения. Рассмотрим его теперь подробнее, используя стохастическое интегральное уравнение со сносом и волатильностью, не зависящими от времени: x(t) = x 0 + t Z 0 a x(? ) d? + t Z 0 b x(? ) ?W ? Идея метода состоит в выборе некоторого нулевого приближения случай- ной функции x 0 (t) , удовлетворяющего начальному условию x 0 (0) = x 0 , и получении поправок к нему по следующей схеме: x k+1 (t) = x 0 + t Z 0 a x k (? ) d? + t Z 0 b x k (? ) ?W ? В правой части стоит известная случайная функция x k (t) , найденная на предыдущей итерации. В результате интегрирований получается следую- щее приближение к решению. Заметим, что на каждой итерации текущее приближение удовлетворяет начальному условию x k (0) = x 0 . Вообще го- воря, требуется доказать, что подобная процедура при бесконечном еј применении сходится к точному решению уравнения. Мы не будем этого делать, а рассмотрим пример еј использования. В качестве нулевого приближения выберем начальное условие x 0 . То- гда постоянные величины a 0 = a(x 0 ) и b 0 = b(x 0 ) выносятся за интеграл, и первая итерация имеет вид: x 1 (t) = x 0 + a 0 t + b 0 W t (5.31) Так как W t = ? ? t , при t ? 0 мы фактически получили итерационную схему для стохастического дифференциального уравнения: x 1 (t) = x 0 + a 0 t + b 0 ? ? t, (5.32) которая активно использовалась в предыдущих главах. Понятно, что она работает тем лучше, чем меньше прошло времени t от начального момен- та t 0 = 0 . При численном решении стохастических дифференциальных уравнений выражение ( 5.32 ) часто называется схемой Эйлера. Стохастические интегралы 149 Разложим снос и волатильность в ряд Тейлора в окрестности x 0 : a(x) = a 0 + a 0 0 (x ? x 0 ) + ..., b(x) = b 0 + b 0 0 (x ? x 0 ) + ..., где a 0 0 = a 0 (x 0 ) и b 0 0 = b 0 (x 0 ) . Подставляя их и ( 5.31 ) в интегральное уравнение, для второй итерации имеем: x 2 (t) = x 1 (t) + a 0 0 a 0 t 2 2 + a 0 0 b 0 t Z 0 W ? d? + b 0 0 a 0 t Z 0 ? ?W ? + b 0 0 b 0 t Z 0 W ? ?W ? Воспользуемся формулой интегрирования по частям (см. Справочник R 54 , стр. 282 ) и известным интегралом по ?W от W ( 5.10 ), стр. 131 : t Z 0 ? ?W ? = t W t ? S t , t Z 0 W ? ?W ? = 1 2 W 2 t ? t . С их помощью перепишем второе приближение к решению: x 2 (t) = x 1 (t) + 1 2 b 0 0 b 0 (W 2 t ? t) + b 0 0 a 0 t W t + a 0 0 a 0 t 2 2 + a 0 0 b 0 ? a 0 b 0 0 S t Интеграл по времени S t от винеровской переменной через W t не выра- жается. Однако, если винеровский процесс выражен через гауссову пере- менную W t = ? ? t , то такой интеграл выражается через две независимые гауссовы переменные ?, ? ? N(0, 1), см. ( 5.4 ), стр. 125 : S t = t Z 0 W s ds = ? 3 ? + ? t 3/2 2 ? 3 Поэтому для второго приближения к решению можно записать: x 2 (t) = x 0 + b 0 ? ? t + a 0 t + 1 2 b 0 0 b 0 (? 2 ? 1) t + b 0 0 a 0 ?t 3/2 + a 0 0 b 0 ? a 0 b 0 0 ? 3 ? + ? t 3/2 2 ? 3 + a 0 0 a 0 t 2 2 (5.33) Так же, как и схема Эйлера, это соотношение работает тем лучше, чем меньше t. Однако этот ряд имеет второй порядок малости по t и явля- ется более точным. Первую строку в этом решении (точность порядка t) называют схемой Милстейна. Мы воспользуемся ею и более точным вы- ражением ( 5.33 ) в девятой главе для ускорения сходимости численного решения стохастических дифференциальных уравнений. 150 Глава 5. Глава 6 Системы уравнений Одномерные стохастические уравнения позволяют описывать только сравнительно простые системы. Даже для обычного физического осцил- лятора необходимо решать систему из двух уравнений первого порядка. Реальность в общем случае многомерна. Она дајт нам множество при- меров достаточно сложных, но исключительно интересных случайных процессов. Как и в одномерном случае, мы начнјм с дискретных процессов, обоб- щение которых на непрерывный случай приведјт нас к системе стохасти- ческих дифференциальных уравнений. Фактически, эта глава повторяет большинство результатов предыдущих глав. Для тех, кто уверенно владе- ет тензорной и матричной алгеброй, соответствующие обобщения служат лишь способом повторения уже известного материала. После вывода ос- новных многомерных уравнений будут рассмотрены решения некоторых задач. 151 152 Глава 6. 6.1 Скоррелированные блуждания • Прекрасным примером скоррелированных блужданий является фи- нансовый рынок. На нјм одновременно изменяются во времени цены ты- сяч финансовых активов. Их динамика образует n-мерный случайный процесс x(t) = {x 1 (t), ..., x n (t)} . Иногда эти активы можно рассматри- вать независимым образом, однако в большинстве својм они тесно свя- заны. Подобная связь обычно возникает за счјт внешних синхронизиру- ющих факторов новостийного или макроэкономического характера. На- пример, цены акций компаний из одного сектора экономики ежедневно изменяются достаточно синхронно. Ещј более жесткая связь возникает на рынке деривативов. В этом случае, хотя цена фьючерса или опциона испытывает случайные колебания, она, тем не менее, тесно связана с дру- гой стохастической величиной финансовым активом, лежащим в основе дериватива. В конечном счјте, финансовый рынок должен описываться единой, и достаточно большой системой стохастических дифференциаль- ных уравнений. • Рассмотрим сначала два случайных дискретных процесса, накоплен- ные изменения которых равны: x t = ?(1) + ?(2) + ... + ?(t) = ? ? t y t = ?(1) + ?(2) + ... + ?(t) = ? ? t. (6.1) Несколько тяжеловесная функциональная запись номера изменения поз- волит нам не путаться при описании множества процессов. Мы по-прежнему считаем, что каждое изменение ?(i) не зависит от предыдущего и являет- ся гауссовым. Аналогично для ?(i). Поэтому суммы равны ? t , умножен- ному на гауссову случайную величину, которую для первого процесса мы обозначили через ?, а для второго через ?. Однако теперь дополни- тельно предположим, что на каждом этапе изменения процессы скорре- лированы друг с другом: h?(i) ?(j)i = ? ? ij , h?(i) ?(j)i = h?(i) ?(j)i = ? ij , (6.2) где ? ij символ Кронекера, равный единице при i = j и нулю в про- тивном случае. Это очень похоже на фондовый рынок, на котором еже- дневные изменения цен двух акций скоррелированы h?(1)?(1)i = ?, но изменения последовательных дней, как для одной бумаги h?(1)?(2)i, так и попарно h?(1)?(2)i, практически равны нулю. Процессы ( 6.1 ) со свойствами изменений ( 6.2 ) будем называть скорре- лированными дискретными винеровскими блужданиями. Системы уравнений 153 • Аналогично одномерному процессу можно вычислить различные со- отношения между средними. Так, пусть в момент времени s мы знаем значение x s . Дајт ли эта информация возможность предсказать y t в неко- торый последующий момент t > s? Найдјм корреляционный коэффици- ент между x s и y t . Для этого разобьјм время на два интервала [1, ..., s] и [s + 1, ..., t]: x s = ? ? s y t = ? ? s + ? 0 ? t ? s. Величины ? и ? пропорциональны накопленным первым s изменениям каждого из процессов. Среднее их произведения тоже равно ?: h??i ? s ? s = ?(1) + ... + ?(s)?(1) + ... + ?(s) = s?, так как при перемножении скобок после усреднения выживут только s пар с одинаковыми индексами. Накопленные изменения ? 0 второго про- цесса не зависят от предыстории и не скоррелированы ни с ?, ни с ?. Поэтому ковариация равна: cov(y t , x s ) = hy t x s i = (? ? s + ? 0 ? t ? s) ? ? s = s?. Это означает, что в линейной регрессионной модели возникает следую- щая связь между случайными величинами: y t ? t = s? ? s ? t x s ? s + ? ? t => y t = ? x s + ?, где, как и раньше, ? ошибка предсказания модели (случайный шум нашего незнания развития процесса после момента времени s). По отдельности каждый из процессов это обычное винеровское дис- кретное блуждание с начальным значением, равным нулю. Если мы зна- ем y s в момент времени s, то прогноз его в момент t равен y t = y s + ? (см. стр. 36 ). Если нам известно только начальное значение y 0 = 0 , то лучшим прогнозом y t будет ноль. В случае, когда процесс x t не зависим от y t , знание x s нам не поможет для прогноза y t . Однако для скоррелиро- ванных процессов это не так. Если ? > 0 и x s находится в положительной области, то и y t вероятнее будет находиться там же. Хотя этот прогноз не зависит от того, совпадают ли моменты времени t = s или t s, и величина его ошибки ? увеличивается со временем t ? s. 154 Глава 6. • Чтобы при помощи компьютера смоделировать скоррелированные блуждания, нам необходимо уметь генерить случайные гауссовы числа c корреляционным коэффициентом ? между ними. Для этого проще всего взять две нескоррелированные величины ? 1 и ? 2 и вычислить следующую их линейную комбинацию (см. стр. 33 ): ? = ? 1 ? = ? ? 1 + p 1 ? ? 2 ? 2 => ? ? = ?, ? 2 = ? 2 = 1. Свойства случайных величин ? и ? легко проверить прямым вычислени- ем, используя соотношения ? i = 0 , ? 2 i = 1 и ? 1 ? 2 = 0 Ниже на каждом из рисунков приведены реализации блуждания двух скоррелированных случайных процессов: -20 0 15 0 100 =0.9 =-0.9 -6 0 12 100 В первом случае коэффициент корреляции достаточно высок ? = 0.9. Во втором он такой же по модулю, но отрицательный ? = ?0.9. • Запишем i-тый стохастический процесс в функциональном виде для индекса времени t: x i (t) = ? i (1) + ... + ? i (t) = ? i ? t. (6.3) Фактически ? i (k) можно считать не функцией, а просто способом раз- личения гауссовых чисел. Величины ? ? (k) , ? ? (k) относятся к одному мо- менту времени и равны изменению ?-го и ?-го процесса. От скоррелированных величин ? i при помощи матрицы S можно пе- рейти к нескоррелированным ? i : ? i = S ij ? j (см. § 1.6 стр. 30 ). В этом случае: x i (t) = t X k=1 ? i (k) = t X k=1 S ij ? j (k) = S ij ? j ? t = S ? t, где использованы смешанные обозначения для суммы по k в явном виде и суммы по j в результате его повторения. Случайные величины ? j (k) j -того процесса в момент времени k не зависят друг от друга. Поэтому мы, как обычно, заменяем их на итоговое суммарное случайное изменение ? j , умноженное на фактор ? t Системы уравнений 155 • В непрерывном пределе мы будем рассматривать m независимых винеровских блужданий W(t) = (W 1 (t), ..., W m (t)) , каждое из которых определяется суммой бесконечно большого числа N ? ? гауссовых изменений, каждое из которых происходит за бесконечно малое время ?t ? 0 : W ? (t) = ? ? (1) + ... + ? ? (N ) ? ?t = ? ? ? N ?t ? ? ? ? t. Изменения этих процессов ?W являются шумом для динамических пе- ременных x = (x 1 , ..., x n ) той или иной системы. Например, вектор x может представлять собой цены различных инструментов на финансо- вом рынке или координаты, задающие положение броуновской частицы в трјхмерном пространстве. Обратим внимание на то, что в общем случае количество m порожда- ющих (см. § 2.8 , стр. 72 ) процессов Винера может быть отлично от числа n интересующих нас случайных функций x. Обычно предполагается, что m 6 n. Например, возможна ситуация, когда один и тот же одномерный шум m = 1 оказывает воздействие на n переменных состояния системы. Однако чаще у каждой переменной состояния свой шум, и, следователь- но, m = n. Многомерное винеровское блуждание с произвольными сносами и во- латильностями можно записать в следующем виде: x i (t) = µ i t + ? i? W ? (t). По повторяющемуся индексу ? происходит суммирование, которое пере- мешивает при помощи матрицы ? i? независимые винеровские процессы W ? (t) . В общем случае i = 1, ..., n, а ? = 1, ..., m. Если n = m и матри- ца ? i? является диагональной, то x i (t) являются нескоррелированными. Иначе ? i? задают не только волатильности процессов, но и коэффициен- ты их корреляции. Аналогично одномерному случаю стохастическое уравнение, которому удовлетворяет этот процесс, имеет вид: dx i = µ i dt + ? i? ?W ? , или в матричных обозначениях: dx = µ dt + ? ?W, где µ = (µ 1 , ...., µ n ) вектор сноса, определяющий детерминированную составляющую изменения случайных функций x(t). 156 Глава 6. 6.2 Системы стохастических уравнений • В общем случае система стохастических уравнений записывается в виде: dx i = a i (x, t) dt + b i? (x, t) ?W ? , (6.4) где i = 1, ..., n, по повторяющемуся индексу ? = 1, ..., m предполагается суммирование, и в общем случае n 6= m. Можно опустить не только знак суммы, но и индексы, записав стохастическое уравнение в матричном виде: dx = a(x, t) dt + b · ?W, (6.5) где a векторная функция, а b матричная, размерности n x m. Вектор винеровских переменных, как и в одномерном случае, записывается через гауссовы случайные числа: ?W = {?W 1 , ..., ?W m } = {? 1 , ..., ? m } ? t = · ? t. (6.6) Мы будем считать, что h? ? ? ? i = ? ?? , а эффекты корреляции переносить на матрицу b i? . Скоррелированные величины ? 0 ? можно выразить через нескоррелированные при помощи линейного преобразования 0 = S · , поэтому стохастический член в уравнении Ито со скоррелированными винеровскими переменными b 0 · 0 ? dt эквивалентен (b 0 · S) · ? dt • Численное моделирование выполняется при помощи выбора мало- го интервала времени ?t. После этого генерится вектор m нормально распределјнных чисел = {? 1 , ..., ? m } и вычисляется набор значений процессов в следующий момент времени. Для первой итерации: x i = x 0i + a i (x 0 , t 0 ) ?t + b i? (x 0 , t 0 ) ? ? ? ?t. (6.7) Процессы x(t) = {x 1 (t), ..., x n (t)} мы всегда нумеруем, начиная с индекса 1, а нулевой индекс x 0i - это значение i-того процесса в момент времени t 0 , т.е. x 0i = x i (t 0 ) Несложно проверить, что смысл коэффициентов сноса определяется средним hx i ? x 0i i /?t = a i (x 0 , t 0 ) , а диффузия: h(x i ? x 0i ) (x j ? x 0j )i ?t = b i? (x 0 , t 0 ) b j? (x 0 , t 0 ) = (b · b T ) ij (6.8) при ?t ? 0 стремится к произведению матриц волатильности, где b T ij = b ji - операция транспонирования (перестановки) индексов. Системы уравнений 157 • Обобщим лемму Ито на n-мерный случай. Пусть F (x, t) диффе- ренцируемая функция. Разложим еј в ряд Тейлора в окрестности точки x 0 , t 0 : F (x, t) = F (x 0 , t 0 ) + ?F ?t ?t + ?F ?x i ?x i + 1 2 ? 2 F ?x i ?x j ?x i ?x j + ... (6.9) По повторяющимся индексам проводится суммирование, и все функции в правой части вычисляются в точке x 0 , t 0 . В соответствии с ( 6.7 ): ?x i = a i (x 0 , t 0 )?t + b i? (x 0 , t 0 ) ? ? ? ?t. (6.10) Изменение функции dF = F (x, t) ? F (x 0 , t 0 ) подчиняется стохастиче- скому уравнению Ито: dF = A(x 0 , t 0 ) dt + B ? (x 0 , t 0 ) ?W ? (6.11) Подставляя ( 6.10 ) в ( 6.9 ) и сохраняя члены порядка ? ?t , ?t, получаем: F ? F 0 ? ?F ?t + ?F ?x i a i + 1 2 ? 2 F ?x i ?x j b i? b j? ? ? ? ? ?t + ?F ?x i b i? ? ? ? ?t. Снос A(x 0 , t 0 ) по определению равен пределу hF ? F 0 i /?t при ?t ? 0 и находится с учјтом соотношений h? ? ? ? i = ? ?? . Для диффузии, в соответствии с ( 6.8 ), имеем: (F ? F 0 ) 2 ?t = ?F ?x i ?F ?x j b i? b j? = B ? B ? Поэтому стохастическое уравнение для скалярной функции n + 1 пере- менных F (x, t), в которую вместо аргументов x подставлены случайные процессы x(t), записывается следующим образом: dF = ?F ?t + ?F ?x i a i + 1 2 ? 2 F ?x i ?x j b i? b j? dt + ?F ?x i b i? ?W ? (6.12) Если функция F не скалярная, а векторная, то это соотношение спра- ведливо для каждой из еј компонент. Введя символ следа матрицы, равного сумме диагональных элементов Tr A = A ?? = A 11 + ... + A nn , можно записать лемму Ито в матричном виде: dF = ?F ?t + ?F ?x · a + 1 2 Tr b T · ? 2 F ?x 2 · b dt + ?F ?x · b · ?W, (6.13) где ? 2 F/?x 2 матрица вторых производных. 158 Глава 6. • Получим многомерное уравнение Фоккера-Планка. Для этого необ- ходимо повторить рассуждения из одномерной задачи. Рассмотрим слу- чайный вектор y = x(t) в момент времени t и предшествующий ему x = x(t ? ?t) в момент времени t ? ?t. Они связаны диффузным стоха- стическим процессом: y = x + a ?t + b · ? ?t, где векторная a = a i (x, t ? ?t) и матричная b = b i? (x, t ? ?t) функ- ции вычислены в момент времени t ? ?t. Предположим, что плотность вероятности случайной величины x равна P (x, t ? ?t). Распределение для гауссовой переменной нам известно. Чтобы найти распределение для величины y, необходимо вычислить среднее от произвольной функции F (y) с известными плотностями P (x, t ? ?t) и P (): hF (y)i = ? Z ?? F (y) z }| { F (x + a?t + b · ? ?t) · P (x,?) z }| { P (x, t ? ?t) P (? 1 , ..., ? m ) d n xd m ?. Разложим первый множитель в ряд по малой величине a ?t + b ? ? ?t : F (y) = F (x) + ?F ?x i (a i ?t + b i? ? ? ? ?t) + 1 2 ? 2 F ?x i ?x j (a i ?t + b i? ? ? ? ?t) (a j ?t + b j? ? ? ? ?t), где по повторяющимся индексам, как и раньше, подразумевается сумми- рование. По ?t раскладываем также P (x, t ? ?t). При интегрировании по всем ? i происходит усреднение, которое дајт h? ? i = 0 и h? ? ? ? i = ? ?? . В результате, повторяя рассуждения на стр. 107 , получаем: ?P ?t + ?(a i P ) ?x i ? 1 2 ? 2 ?x i ?x j h b i? b j? P i = 0 . (6.14) где a i = a i (x, t) , b i? = b i? (x, t) , а P = P (x 0 , t 0 ? x, t) - условная плот- ность вероятности. Если в момент времени t 0 значение x 0 известно точно, то для решения этого уравнения используется начальное условие в виде n - мерной дельта - функции Дирака, равной произведению одномерных функций по каждой координате: P (x 0 , t 0 ? x, t 0 ) = ?(x ? x 0 ) Системы уравнений 159 • Аналогично выводится уравнение для производной от среднего: d F x(t), t dt = ? Z ?? ? ?t F (x, t) P (x 0 , t 0 ? x, t) d n x. Раскрывая производную произведения и подставляя ?P/?t из уравнения Фоккера - Планка, получаем динамические уравнения для средних: d F x(t), t dt = ?F ?t + a i ?F ?x i + 1 2 b i? b j? ? 2 F ?x i ?x j (6.15) Как и лемму Ито, это соотношение можно записать в матричной форме при помощи символа следа Tr. Усреднение производится при условии, что в момент времени t вектор случайного процесса был равен x 0 = x(t 0 ) Уравнение для среднего справедливо и для векторных или тензорных функций, так как выводится независимо для каждой из компонент. Вы- бирая F = x ? и учитывая, что ?x ? /?x i = ? ?i , получаем временную дина- мику среднего в компонентной и матричной форме: ? hx ? i = ha ? (x, t)i , ? hxi = ha(x, t)i . (6.16) Только для линейных по x сносов динамика среднего значения будет совпадать с решением детерминированного уравнения. Функциональная зависимость волатильности b(x, t) при этом роли не играет. Если снос нелинеен по x, то функция hxi = x(t) будет отличаться от детерминиро- ванного решения с b = 0. Производные от произведения x µ x ? выражаются через символ Кроне- кера следующим образом: ?(x µ x ? ) ?x i = x µ ? ?i + x ? ? µi , ? 2 (x µ x ? ) ?x i ?x j = ? µj ? ?i + ? ?j ? µi Поэтому, выбирая F в тензорном виде F = x µ x ? , можно записать урав- нение для среднего от произведения случайных процессов: ? hx µ x ? i = hx µ a ? + x ? a µ + b ?? b µ? i . (6.17) В частности, для свјртки (суммирования) по индексам µ и ? имеем мат- ричное выражение для изменения квадрата ? hx 2 i = 2 hx · ai + Tr b · b T 160 Глава 6. 6.3 Стохастический осциллятор • В качестве примера решения стохастической задачи в двух измере- ниях n = m = 2 рассмотрим движение по окружности с частотой ? и по- степенным уменьшением радиуса. Детерминированная версия подобной спирали может иметь следующую зависимость координат от времени: x y t x(t) = e ??t (x 0 cos ?t ? y 0 sin ?t) y(t) = e ??t (x 0 sin ?t + y 0 cos ?t) B Начальные условия: x 0 = x(0), y 0 = y(0) B Радиус: r(t) = p x 2 + y 2 = r 0 e ??t Примеры систем с таким поведением мы рассмотрим в следующей главе. Сейчас наша цель математическое описание стохастической динамики. Для этого найдјм решение системы уравнений следующего вида: dx = (?? x ? ? y) dt + ? ?W x dy = (+? x ? ? y) dt + ? ?W y Предполагается, что шум ?W x = ? x ? t , ?W y = ? y ? t по каждой координа- те нескоррелирован. В качестве упражнения (l H 29 ) стоит записать эту же систему для скоррелированного шума. • Зависимость среднего значения от времени находим из ( 6.16 ): ?x = ?? x ? ? y ?y = +? x ? ? y. Умножим второе из уравнений на мнимую единицу ? (? 2 = ?1 ) и сло- жим их. В результате для комплексной величины z = x+?y и параметра ? = ? ? ?? получим одномерное уравнение ?z = ?? z. Оно легко инте- грируется: z(t) = e ??t z 0 = e ?? t+?? t z 0 = e ??t (cos ?t + ? sin ?t)z 0 , где z 0 = z(0) = x 0 + ?y 0 начальное условие. Приравняв действительную и мнимую части: x(t) = e ??t (x 0 cos ?t ? y 0 sin ?t) y(t) = e ??t (x 0 sin ?t + y 0 cos ?t), (6.18) получаем решение для эволюции средних значений координат в виде спирали. По каждой координате происходят затухающие периодические колебания. Параметр ? является их частотой, ? скоростью затухания. Системы уравнений 161 • Найдјм теперь, как ведјт себя среднее значение квадрата расстоя- ния от начала координат. Воспользуемся матричной записью уравнения ( 6.17 ) для среднего ? hx 2 i = 2 hx · ai + Tr b · b T . В нашем случае: x = x y , a = ??x ? ?y +?x ? ?y , b = ? 1 0 0 1 Поэтому получаем: ? x 2 = ?2? x 2 + 2? 2 => x 2 (t) = ? 2 ? + x 2 0 + y 2 0 ? ? 2 ? e ?2?t Таким образом, когда колебания затухнут (t ? ?), останется ненуле- вая дисперсия квадрата, тем большая, чем медленнее происходило зату- хание! Этот результат свидетельствует о стабилизирующей роли сильно- го трения (большие ?) при внешних стохастических толчках. Несмотря на то, что ? находится в знаменателе, особенности при ? = 0 нет. Рас- кладывая экспоненту, несложно убедиться, что при ? = 0 среднее равно x 2 = x 2 0 + y 2 0 + 2? 2 t . Это означает, что, как и в винеровском блуждании, внешние толчки со временем размывают круговую траекторию. Найти асимптотическое решение x 2 (t) , y 2 (t) , xy(t) в ситуации скоррелированно- го шума предлагается в качестве упражнения (l H 30 ). • Выразим решение задачи через гауссовы переменные. В комплекс- ных обозначениях стохастическое уравнение имеет компактный вид: dz = ?? z dt + ? ?W, где ?W = ? ? dt , ? = ? x + ?? y комплексное гауссово число, а z и ? определены выше. Перейдјм, при помощи формулы Ито, к переменной F = ze ?t . Еј динамическое уравнение не будет содержать сноса: dF = ?e ?t ?W = S(t) ?W, где S(t) = ? e ?t . Решим уравнение dF = S(t) ?W итерациями (k = 1...n): F = F 0 + X S(t k?1 )?(t k ) ? ?t. Как функция S(t), так и ? k являются комплексными величинами, поэто- му необходима определјнная осторожность по сворачиванию этой суммы в одно гауссово число. Распишем действительные и мнимые части: X [S x + ?S y ][? x + ?? y ] ? ?t = X [(S x ? x ? S y ? y | {z } |S| ? 0 x ) + ?(S x ? y + S y ? x | {z } |S| ? 0 y )] ? ?t, где |S| = q S 2 x + S 2 y модуль комплексного числа, а ? 0 x , ? 0 x новые нескор- релированные ? 0 x ? 0 y = 0 гауссовы числа. 162 Глава 6. Опуская штрихи и повторяя рассуждения одномерного случая (стр. 56 ), мы можем окончательно записать: F (t) = F (t 0 ) + ? ? t Z t 0 |S(? )| 2 d? ? ? 1/2 ?, где ? = ? x +?? y по-прежнему комплексная гауссова случайная величина. Заметим, что действительная и мнимая части выражения ? 0 = e ?? ? являются независимыми гауссовыми числами. Действительно, распишем их в явном виде: ? 0 x = ? x cos ? ? ? y sin ? ? 0 y = ? x sin ? + ? y cos ? Прямым вычислением проверяем ? 02 x = ? 02 y = 1 и ? 0 x ? 0 y = 0 . Поэто- му множители типа e ??t перед комплексным гауссовым числом можно опускать, так как e ?? ? статистически эквивалентно просто ? Проводя интегрирование для |S(?)| 2 = ? 2 e 2?? и учитывая, что z = F e ??t , z 0 = F 0 , для t 0 = 0 получаем: z = z 0 e ??t+??t + ? ? 2? p 1 ? e ?2?t ? (6.19) или в явном виде для действительной и мнимой частей: x(t) = x(t) + ? ? 2? p 1 ? e ?2?t ? x y(t) = y(t) + ? ? 2? p 1 ? e ?2?t ? y , где x(t) и y(t) средние, определяемые выражениями ( 6.18 ). В качестве упражнения стоит найти x 2 (t) , y 2 (t) , xy(t) и проверить справедливость уравнений для средних (l H 31 ). Квадрат величины |z| 2 = x 2 + y 2 является квадратом радиус-вектора, для которого уже известно среднее значение. Сделаем это ещј раз при помощи ( 6.19 ): |z t | 2 = |z 0 | 2 e ?2?t + ? 2 ? 1 ? e ?2?t . Обращаем внимание на то, что |?| 2 = h?? ? i = ? 2 x + ? 2 y = 2 , где звјз- дочка обозначает комплексное сопряжение. В решении, аналогично одномерному случаю, можно выразить z в мо- мент времени t + ? через z в момент t: z t+? = z t e ??? +??? + ? ? 2? p 1 ? e ?2?? ?, что легко позволяет вычислить, например, среднее hz t z ? t+? i (l H 32 ). Системы уравнений 163 • При больших временах t ? ? решение забывает начальные усло- вия, и средние стремятся к нулю. Распределение становится стационар- ным по каждой из координат. Однако это не означает, что периодические свойства системы исчезают. Чтобы в этом убедиться, найдјм ковариаци- онную функцию, например, по координате x. Выражая решение относи- тельно начального момента времени t, имеем: x t+? = e ??? (x t cos ?? ? y t sin ?? ) + ? ? 2? ? x p 1 ? e ?2?? Найдем cov xx (t, t + ? ) = hx t x t+? i ? hx t i hx t+? i в пределе t ? ?. Так как в этом случае hx t i = hy t i = hx t y t i = 0 , а x 2 t = ? 2 /2? , получаем, как и следовало ожидать, стационарную ковариационную функцию, завися- щую только от ? > 0: cov xx (t, t + ? ) ? cov(? ) = ? 2 2? e ??? cos ??. Она оказывается периодической функцией сдвига ?. Фурье-образ ковари- ационной функции характеризует спектральные свойства процесса (стр. 68 ): S(?) = 2 ? ? Z 0 cov(? ) cos(?? ) d? = ? 2 2? 1 ? 2 + (? + ?) 2 + 1 ? 2 + (? ? ?) 2 Таким образом, спектр имеет максимум в окрестности ? = ?. Он тем уже, чем меньше параметр затухания ?. Тем не менее, это не строго пери- одическое движение, так как типичная частота размазана и сдвинута первым слагаемым в квадратных скобках. На левом рисунке ниже представлена траектория стохастического ос- циллятора при достаточно больших временах, когда начальные условия уже забыты. Справа колебания по каждой координате: y x x y Системы, обладающие подобным поведением, мы рассмотрим в следую- щей главе, а сейчас решим многомерное линейное уравнение. 164 Глава 6. 6.4 Линейные многомерные модели • Найдјм решение линейных стохастических уравнений (по j сумма): dx i = A ij (x j ? c j ) dt + B ij ?W j Постоянный вектор c j можно убрать сдвигом x j ? x j + c j . В решении делается обратный сдвиг. Поэтому будем изучать однородное уравнение, которое запишем в матричной форме: dx = A · x dt + B · ?W, где A и B не зависящие от x и времени матрицы. Для определения среднего проще всего сразу воспользоваться соотно- шением ( 6.16 ): ?x = A · x => x = e At · x 0 , (6.20) где x 0 вектор начального значения. Если мы хотим вернуть c, то потребуются две замены: x ? x ? c и x 0 ? x 0 ? c • Монотонная зависимость от t в матричной записи решения ( 6.20 ) об- манчива. Рассмотрим стохастический осциллятор из предыдущего раз- дела: dx dy = ?? ?? ? ?? · x y + ? ?W x ?W y (6.21) В этом случае матрицу A можно разбить на сумму двух матриц: A = ?? ?? ? ?? = ? q ? ? 1, 1 = 1 0 0 1 , q = 0 ?1 1 0 Несложно проверить, что: q 2 = ?1, q 3 = ?q, q 4 = 1, q 5 = q, ... Так как матрицы 1 и q коммутируют друг с другом (q · 1 = 1 · q), экспо- нента суммы разбивается на произведение e At = e ?1?t e q ?t . Раскладывая второй множитель по степеням t и учитывая аналогичное разложение для синуса и косинуса, решение можно представить в следующем виде: Ї x Ї y = e ??t cos ?t ? sin ?t sin ?t cos ?t x 0 y 0 (6.22) Оно же выше было получено другим методом. Таким образом, монотон- ная зависимость от времени в матричных соотношениях вполне может превратиться в периодическую функцию. Системы уравнений 165 • Найдјм более практичное, чем ( 6.20 ), представление для решения линейного уравнения. Будем его искать в виде: x(t) = u e at => A · u = a u. (6.23) Постоянный вектор u является собственным вектором матрицы A, а па- раметр a еј собственным значением. Перенося (a u) в левую часть, получаем систему однородных уравнений относительно u, которая име- ет ненулевое решение, только если еј детерминант равен нулю: det(A ? a 1) = 0. Это уравнение называется характеристическим и является полиномом n -той степени по a. Обычно оно имеет n различных решений a 1 , ..., a n Часть из них может оказаться комплексными. Для каждого из них мы решаем уравнение ( 6.23 ) и находим собственные вектора u (k) . Внимание! Верхний индекс это номер собственного вектора, а не его компонента. Теперь общее решение для среднего значения вектора переменных со- стояния можно записать в следующем виде: x(t) = X k µ k u (k) e a k t , x 0 = X k µ k u (k) , (6.24) где µ k произвольные константы, выражающиеся через начальные усло- вия x 0 = x(0) . Прямой подстановкой в исходное уравнение можно про- верить справедливость этого решения. Действительная часть собствен- ных значений a k будет приводить к экспоненциально уменьшающимся (Re a k < 0 ) или увеличивающимся (Re a k > 0 ) решениям. Мнимая часть соответствует колебательным режимам. Если матрица A симметрична, то собственные вектора можно вы- брать ортогональными: u (?) · u ?(?) = ? ?? (звјздочка комплексное сопря- жение). В этом случае µ k = x 0 · u ?(k) Когда µ k выражены через x 0 , можно найти явное представление экспо- ненты от матрицы. Действительно, из ( 6.20 ), взяв производную по компо- нентам начального условия, имеем e At ?? = ?x ? /x 0? . В частности, если собственные вектора ортогональны (µ k = x 0 · u ?(k) ), то: e At ?? = X k u (k) ? u ?(k) ? e a k t (6.25) В качестве упражнения (l H 33 ) предлагается найти e At для матрицы 2x2, у которой A 12 = A 22 = 0 . Необходимо это сделать прямым разложением экспоненты при помощи собственных значений. 166 Глава 6. • Выразим теперь решение стохастической линейной системы через гауссовы переменные. Введјм новый вектор y, удовлетворяющий, по лем- ме Ито ( 6.13 ), стр. 157 , следующему уравнению: y = e ?At · x => dy = e ?At B ?W = G(t) ?W. Матрица G(t) = e ?At B зависит только от времени, поэтому решение этого уравнения легко найти при помощи итерационного метода: y µ (t) = y µ (t 0 ) + X k G µ? (t k )? ? (t k ) ? ?t = y µ (t 0 ) + g µ? ? ? Сумма независимых гауссовых чисел ? ? (t k ) снова пропорциональна гаус- совому числу, которое удобно представить в виде суммы независимых величин ? ? (второе равенство). Найдјм значения g µ? . Для этого вычис- лим среднее от y(t) ? y(t 0 )) µ (y(t) ? y(t 0 ) ? : g µ? g ?? h? ? ? ? i = X k,l G µ? (t k )G ?? (t l ) h? ? (t k )? ? (t l )i ?t. Учитывая независимость случайных величин h? ? (t k )? ? (t l )i = ? ?,? ? k,l и h? ? ? ? i = ? ?,? , а также переходя к непрерывному пределу ?t ? 0, полу- чаем (t 0 = 0 ): g µ? g ?? = X i G µ? (t i )G ?? (t i )?t = t Z 0 G µ? (? )G ?? (? ) d?, или: g(t) · g T (t) = t Z 0 e ?A? B B T e ?A T ? d?. (6.26) Напомню, что (A · B) T = B T · A T (см. стр. 304 ). Решение для y запишем в матричном виде, учитывая, что y 0 = x 0 при t = 0: y = x 0 + g(t) · Поэтому, так как x = e At y , окончательное решение системы линейных стохастических уравнений имеет вид: x(t) = Ї x(t) + S(t) · , (6.27) где S = e At g . Вектор = {? 1 , ..., ? n } представляет собой набор незави- симых случайных чисел с гауссовым распределением, имеющим нулевое среднее и единичную дисперсию, а Їx(t) среднее значение ( 6.20 ), ( 6.24 ). В качестве упражнения (l H 34 ) предлагается найти матрицу e At для двухмерного осциллятора и проверить решение ( 6.27 ). Системы уравнений 167 • Вычислим матрицу дисперсий: D ?? = h(x ? Ї x) ? (x ? Ї x) ? i = S ?i S ?j h? i ? j i = S S T ?? = h e A t g g T e A T t i ?? Учитывая ( 6.26 ), имеем: D(t) = S S T = t Z 0 e A(t?? ) B B T e A T (t?? ) d?. (6.28) Это соотношение можно (l H 38 ) сразу получить из уравнения для сред- них ( 6.17 ), стр. 159 , из которых следует матричное уравнение: ? D = A · D + D · A T + B · B T (6.29) Если существует стационарный режим, то ?D = 0 и уравнение ( 6.29 ) поз- воляет легко найти D. • Распределение для x имеет гауссовый вид, поэтому, зная матрицу дисперсий, можно записать марковскую плотность вероятности: P (x 0 , 0 ? x, t) = 1 (2?) n/2 pdet D(t) exp ? 1 2 (x ? Ї x) ? D ?1 ?? (t) (x ? Ї x) ? , где D ?1 обратная матрица дисперсий и Їx = Їx(t) средние значения ди- намических переменных. Они полностью определяют свойства процесса. В частности, характеристическая функция (l H 35 ): he ? p·x i = e ? p·Ї x? 1 2 p·D·p позволяет легко находить моменты произвольных порядков. • При помощи ( 6.27 ), ( 6.28 ) несложно (l H 39 ) найти ковариационную матрицу: cov ?? (t, t+? ) = hx ? (t)x ? (t + ? )i?hx ? (t)i hx ? (t + ? )i = D(t) e A T ? (6.30) Если в пределе t ? ? у системы существует стационарный режим, то в этом случае матрица дисперсий D становится постоянной, а ковариация зависит только от разности времјн ?. • Таким образом, алгоритм решения линейной задачи следующий: B Находим собственные значения и вектора матрицы A. B Записываем решение для средних ( 6.24 ) и выражаем µ k через x 0 B При помощи соотношения e At ?? = ?x ? /x 0? находим e At B Вычисляем матрицу дисперсий D ?? 168 Глава 6. 6.5 Многомерие помогает одномерию • Рассмотрим построение решения одномерного уравнения при помо- щи стохастического интеграла. Пусть у нас есть система уравнений n x 1 с одинаковым шумом dx ? = a ? dt + b ? ?W . Для произвольной функции F = F (t, x) в этом случае справедлива следующая версия леммы Ито ( 6.12 ), стр. 157 (суммирование по повторяющимся индексам): dF = ?F ?t + a ? ?F ?x ? + b ? b ? 2 ? 2 F ?x ? ?x ? dt + b ? ?F ?x ? ?W. Пусть x = {x, W }, где x = x(t) решение некоторого одномерного сто- хастического уравнения dx = a(x) dt + b(x) ?W , а W порождающий винеровский процесс с нулевым сносом и единичной волатильностью. В этом случае лемма Ито для функции F = F (t, x, W ) имеет вид: dF = ?F ?t + a ?F ?x + b 2 2 ? 2 F ?x 2 + b ? 2 F ?x?W + 1 2 ? 2 F ?W 2 dt+ b ?F ?x + ?F ?W ?W. Всегда подходящим выбором F : F (t, x, W ) = f (t, z) = f t, Z dx b(x) ? W можно (l H 41 ) волатильность (множителя при ?W ) сделать равной нулю. Подставляя F = f(t, z) в лемму Ито, получаем: df = ?f ?t + a(x) b(x) ? b 0 (x) 2 ?f ?z dt. Если выбором функции f удајтся добиться, чтобы множитель при dt зависел только от t и W , то это уравнение можно проинтегрировать, выразив решение в явном виде через порождающий винеровский процесс W t . Так как множитель при dt не должен зависеть от x, то частная производная по x равна нулю, и мы получаем следующее уравнение: ? 2 f ?t?z + a(x) b(x) ? b 0 (x) 2 ? 2 f ?z 2 + b(x) a(x) b(x) ? b 0 (x) 2 0 ?f ?z = 0. Это уравнение допускает разделение переменных f(t, z) = e ?t f (z) : ? + a(x) b(x) ? b 0 (x) 2 µ + b(x) a(x) b(x) ? b 0 (x) 2 0 = 0, где ?, µ некоторые константы, которые необходимо подобрать так, что- бы это соотношение обращалось в тождество. Тогда f(t, z) = e ? t+µ z Системы уравнений 169 • Рассмотрим в качестве примера логистическое уравнение (стр. 88 ): dx = x · (1 ? x) dt + p 2? x ?W. (6.31) Несложно проверить, что µ = ? ? 2? , ? = 1 ? ?, F = e (1??) t+ ? 2? W t /x , dF = e (1??) t+ ? 2? W t dt => F = 1 x 0 + t Z 0 e (1??) ? + ? 2? W ? d?, где при интегрировании учтено начальное условие F 0 = F (0) = 1/x 0 , и x 0 = x(0) . Поэтому решение ( 6.31 ) можно записать в следующем виде: x(t) = x 0 e (1??) t+ ? 2? W t ? ? 1 + x 0 t Z 0 e (1??) ? + ? 2? W ? d? ? ? ?1 (6.32) Замкнутая форма ( 6.32 ) может быть полезна при построении приближјн- ных методов, однако, к сожалению, получить с еј помощью конкретные результаты (например, среднее Їx(t)), вообще говоря, не просто. • Простое интегральное представление для решения имеет также ли- нейное по x стохастическое уравнение: dx = (? + ?x) dt + ? x ?W. В этом случае µ = ?, ? = (? 2 /2) ? ? , и для процесса y = x e ?(??? 2 /2) t??W t получаем уравнение: dy = ? e ?(??? 2 /2)t??W t dt. Поэтому решение выражается через стохастический интеграл: x(t) = e (??? 2 /2)t+?W t ? ? x 0 + ? t Z 0 e ?(??? 2 /2)s??W s ds ? ? , который позволяет вычислять средние: x(t) = x 0 D e (??? 2 /2)t+?? ? t E + ? D e (??? 2 /2)t+?W t t Z 0 e ?(??? 2 /2)s??W s ds E Первое среднее вычисляется обычным образом, а для второго необходи- мо использовать формулу ( 5.7 ), стр. 129 : x(t) = x 0 e ?t + ? t Z 0 D e (??? 2 /2)t+?(? 1 ? s+? 2 ? t?s) e ?(??? 2 /2)s??? 1 ? s E ds, или x(t) = x 0 e ?t + ? ? e ?t ? 1 . Заметим, что x(t) в данном случае проще найти при помощи динамических уравнений для средних. 170 Глава 6. • Иногда многомерные системы позволяют находить точные реше- ния одномерных стохастических уравнений. Рассмотрим блуждание в n -мерном пространстве, считая, что по каждой координате реализуется процесс Орнштейна-Уленбека с нулевым равновесным уровнем и одина- ковым притяжением к нему: dx i = ? ? 2 x i dt + ? 2 ?W i Блуждания предполагаются нескоррелированными, с одной волатильно- стью шума. Рассмотрим случайный процесс, равный квадрату радиус- вектора y(t) = x 2 1 + ... + x 2 n . Найдјм стохастическое уравнение, которому он подчиняется. В данном случае a i = ??x i /2 , b ij = ?? ij /2 . Производ- ные от y равны ?y/?x i = 2x i , ? 2 y/?x i ?x j = 2? ij , и по лемме Ито ( 6.12 ), стр. 157 имеем следующее уравнение: dy = ?? · y ? n? 2 4? dt + ? x i ?W i Сумму стохастических членов в этом уравнении можно выразить через единственную винеровскую переменную: ? i ?W i = ? i ? i ? dt = ? ? dt = ?W, где ? i = x i / ? y . Действительно, сумма гауссовых чисел снова дајт гаус- сово число, и, так как ? 2 1 + ... + ? 2 n = 1 , оно имеет единичную дисперсию. Вообще говоря, величины ? i (t) являются случайными функциями. Одна- ко на каждом шаге итерационного решения они принимают определјнное значение, но так, что сумма их квадратов всегда равна единице. Поэтому мы и переходим в уравнении к скалярной винеровской переменной ?W . В результате получается одномерное уравнение Феллера: dy = ?? · (y ? ?) dt + ? ? y ?W, с равновесным уровнем, равным ? = n? 2 /4? . Его решение выражается через известные нам случайные процессы Орнштейна-Уленбека (стр. 60 ): x i (t) = x 0i e ??t/2 + ? 2 ? ? p 1 ? e ??t ? i и n независимых гауссовых величин {? 1 , ..., ? n } Системы уравнений 171 • Решение одномерного уравнения должно зависеть от одной констан- ты начального условия y 0 = y(0) . В полученное решение входит n незави- симых констант x 0i , соответствующих начальным условиям по каждой координате. Покажем, что они, тем не менее, сворачиваются в един- ственную константу y 0 = x 2 01 + ... + x 2 0n . Для этого запишем решение в следующем виде: y(t) = n X i=1 x 0i µ(t) + 1 ? 2 s(t)? i 2 = y 0 µ 2 (t) + p 2y 0 s(t)µ(t) ? + s 2 (t) u, где µ(t) = e ??t/2 , s(t) = p?(1 ? e ??t ) , ? = ? 2 /2? и введены две новые случайные величины ? и u: ? = n X i=1 ? i ? i , u = 1 2 n X i=1 ? 2 i , n X i=1 ? 2 i = 1. Сумма квадратов весов ? i = x 0i / ? y 0 равна единице. Поэтому вели- чина ? имеет гауссово распределение с нулевым средним и единичной дисперсией, а u подчиняется ? 2 -распределению с n степенями свободы (l C 26 ). Так как обе эти величины зависят от одних и тех же гауссовых чисел ? 1 , ..., ? n , они не являются независимыми. Однако их совместная плотность вероятностей не зависит от весов ? i , и, следовательно, от констант начального условия x 0i . Действительно, найдјм производящую функцию: ?(k, p) = e k ?+ p u = n Y i=1 ? Z ? e k ? i ? i +p ?2 i 2 ? ?2 i 2 d? i ? 2? = e k2/2 1?p (1 ? p) n/2 (6.33) Введение мнимой единицы k ? ?k, p ? ?p превращает производящую функцию в характеристическую, фурье-интеграл от которой равен плот- ности вероятности P (?, u). Разложение в ряд по k и p производящей функции позволяет легко найти различные средние для случайных ве- личин ? и u. Подобное представление было получено в третьей главе при изучении процесса Феллера (стр. 82 ). Несмотря на то, что мы начали с n процес- сов Орнштейна-Уленбека, целочисленный параметр n в решении можно аналитически продолжить в область непрерывных значений n = 2?/?. Таким образом, процесс y(t) зависит от единственной константы на- чального условия y 0 = y(0) и двух случайных величин ? и u, имеющих совместное распределение ( 6.33 ). 172 Глава 6. 6.6 Некоторые точные решения ? • Рассмотрим нестационарное многомерное стохастическое уравнение: dx i = f i (t) dt + s i? (t) ?W ? (6.34) При его решении итерациями получатся ряды следующего вида: x i (t) = x i (t 0 )+ h f i (t 0 )+f i (t 1 )+... i ?t+ h s i? (t 0 )? ? (t 0 )+s i? (t 1 )? ? (t 1 )+... i ? ?t. Последний член является суммой независимых гауссовых случайных чи- сел. Поэтому решение ( 6.34 ) можно записать следующим образом: x i (t) = Ї x i (t) + S i? (t) ? ? , (6.35) где по ? по-прежнему производится суммирование, и Ї x i (t) = x i (t 0 ) + t Z t 0 f i (? ) d?, D ij = S i? S j? = t Z t 0 s i? (? )s j? (? ) d?. Явный вид матричной функции S i? (t) обычно не требуется. Нестацио- нарное гауссово блуждание полностью определяется вектором средних значений Їx i (t) и симметричной матрицей дисперсий: D = D ij = h(x i ? Ї x i )(x j ? Ї x j )i = t Z t 0 s i? (? )s j? (? ) d? = t Z t 0 s(? ) · s T (? ) d?. Через них выражается производящая функция для средних (действи- тельный аналог характеристической функции): ?(p) = he p·x i = e p·Ї x e p·S· = e p·Ї x+ 1 2 p·D·p Усреднение проводится покомпонентно для каждого гауссового числа = {? 1 , ..., ? n } при помощи формулы ( 1.11 ), стр. 16 Моменты произвольного порядка находятся взятием частных произ- водных от ?(p). Например, для: e ijkl = h(x i ? Ї x i )(x j ? Ї x j )(x k ? Ї x k )(x l ? Ї x l )i = ??(p) ?p i ?p j ?p k ?p l p=0 получаем: e ijkl = D ij D kl + D ik D jl + D il D jk Заметим, что это выражение автоматически симметрично по всем четы- рем индексам. Системы уравнений 173 • Изменения цен различных финансовых инструментов (например, ак- ций) обычно скоррелированы друг с другом. Простейшей моделью явля- ется многомерное логарифмическое блуждание. В этом случае относи- тельное изменение цены это n-мерный винеровский процесс: dx i x i = µ i dt + n X j=1 ? ij ?W j (!) Сейчас мы используем явное обозначение для суммы, а не соглашение о суммировании по повторяющимся индексам. Дисперсия относительных изменений двух акций выражается через матрицу ? ij . Действительно, для небольшого интервала времени ?t, пред- ставив ?W j = ? j ? ?t , имеем: ?x i x i ? µ i ?t ?x j x j ? µ j ?t = n X k,l=1 ? ik ? jl h? k ? l i ?t = n X k=1 ? ik ? jk ?t. Для получения решения перейдјм, как и в одномерном случае, к на- туральному логарифму от x i . Тогда по лемме Ито имеем: d ln x i = µ i ? 1 2 n X j=1 ? 2 ij ! dt + n X j=1 ? ij ?W j Решение этого уравнения с начальным условием x 0i = x i (0) , выраженное через гауссовы переменные, имеет вид: x i (t) = x 0i exp ( µ i ? 1 2 n X j=1 ? 2 ij ! t + n X j=1 ? ij ? j ? t ) Среднее значение экспоненциально изменяется со скоростью, опреде- ляемой параметром µ i : hx i (t)i = x 0i e µ i t Аналогично, среднее значение квадрата имеет вид: x 2 i (t) = x 2 0i exp ( 2µ i t + n X j=1 ? 2 ij t ) Более подробно к вопросу стохастического описания финансовых рынков мы вернјмся в восьмой главе. 174 Глава 6. • Как и в одномерном случае, некоторую систему стохастических урав- нений можно попытаться свести к простому нестационарному случаю. Для этого подберјм такую векторную функцию F = F k (x, t) , которая убирает x из уравнения (по повторяющимся индексам снова предпола- гается суммирование): dF k = ?F k ?t + ?F k ?x ? a ? + 1 2 ? 2 F k ?x i ?x j b i? b j? | {z } f k (t) dt + ?F k ?x i b i? | {z } s k? (t) ?W ? Пусть b ?1 i? обратная к b i? матрица. Тогда для функций волатильности s k? (t) можно записать: ?F k ?x i b i? = s k? (t) => ?F k ?x i = s k? (t) b ?1 ?i (6.36) Для нестационарного сноса f k (t) : f k (t) = ?F k ?t + s k? b ?1 ?? a ? ? 1 2 s k? b ?1 ?? ?b ?? ?x j b j? , (6.37) где мы подставили ( 6.36 ) и воспользовались соотношением: ?b ?1 ?x j = ?b ?1 · ?b ?x j · b ?1 которое получается дифференцированием b ?1 · b = 1 по x j Возьмјм производную выражения ( 6.36 ) по t и производную по x i от ( 6.37 ). Вычитая их, получаем условие совместности в следующем виде: ? ?t s k? (t) b ?1 ?i + s k? (t) ? ?x i b ?1 ?? a ? ? 1 2 ?b ?? ?x j b j? = 0. (6.38) Как и в одномерном случае, если при данных a i (x, t) и b ij (x, t) удајтся подобрать такие функции времени s k? (t) , что ( 6.38 ) обращается в тож- дество, то решение стохастического уравнения записывается в неявном виде: F k (x(t), t) = F k (x 0 , t 0 ) + t Z t 0 f k (? ) d? + S i? (t) ? ? , (6.39) где ? ? нормированные независимые гауссовы случайные числа, а S i? (t) S j? (t) = t Z t 0 s i? (? )s j? (? ) d?. Приведјм пример использования этого алгоритма. Системы уравнений 175 • Для системы линейных уравнений с постоянной матрицей A и зави- сящими от времени вектором c(t) и матрицей B(t): dx i = A ij x j + c j (t) dt + B ij (t) ?W j , условие совместности ( 6.38 ) и его решение имеют вид: d dt (s · B ?1 ) = ?(s · B ?1 ) · A => s(t) = s(t 0 ) · B ?1 (t 0 ) · e ?A t · B(t). При использовании алгоритма поиска точного решения нам достаточно найти частное решение условия совместности, так как фактически мы ищем наиболее простую замену F(x, t), приводящую исходное уравне- ние к нестационарному винеровскому процессу ( 6.34 ). Поэтому выберем начальное условие для матрицы s в следующем виде s(t 0 ) = B(t 0 ) и, следовательно: s(t) = e ?A t · B(t). В результате функции замены F(x, t) ( 6.36 ) и сноса f(t) ( 6.37 ) равны: F(x, t) = e ?A t · x, f (t) = e ?A t · c(t). Окончательное решение является нестационарным гауссовым процессом: x(t) = e A (t?t 0 ) x 0 + t Z t 0 e A·(t?? ) · c(? ) d? + G · , где x 0 = x(t 0 ) начальное условие. Матрица G = e A t ·S(t) удовлетворяет соотношению, определяющему матрицу дисперсии процесса: D = G · G T = e At · t Z t 0 ss T d? · e A T t = t Z t 0 e A(t?? ) · B(? )B T (? ) · e A T (t?? ) d?. Если c = 0, а B является постоянной матрицей, то эти формулы совпа- дают с результатами, полученными в предыдущем разделе. В случае, когда матрица A зависит от времени, вместо e ?At необхо- димо использовать матрицу ?(t), удовлетворяющую уравнению ??(t) = ??(t) · A(t) . Явный вид ?(t) можно выразить через A(t), однако для этого необходимы специальные обозначения, упорядочивающие матри- цы, так как интеграл от матрицы A(?) по интервалу [t 0 , ..., t] в общем случая не коммутирует с матрицей A(t) в момент времени t. 176 Глава 6. 6.7 Как решать стохастические задачи? Подведјм итоги рассмотренных выше методов анализа стохастических задач. Основными математическими инструментами описания эволюции системы, при воздействии на неј случайных факторов, являются стоха- стическое дифференциальное уравнение и уравнение Фоккера-Планка. Обычно именно стохастическое уравнение будет исходным, так как оно оперирует терминами изменения наблюдаемых переменных состояния си- стемы (координата, импульс, количество особей, цена и т.п.). Часто стоха- стическое дифференциальное уравнение является естественным обобще- нием уже известного детерминированного уравнения эволюции системы, к которому добавляется член шумового воздействия, пропорциональный ?W Функции сноса и волатильности стохастического уравнения вместе с начальными и граничными условиями полностью задают изучаемую си- стему. В качестве начальных условий может быть выбрано конкретное значение переменных состояния x 0 в момент времени t 0 , или некоторая их плотность вероятности P (x 0 ) . Знание сноса и волатильности позволя- ет легко записать уравнение Фоккера-Планка, которое обычно оказыва- ется более удобным при наличии граничных условий. Когда стохастическое описание применяется к реальной (естественной или искусственной) системе, важно понимать порядок значений парамет- ров, входящих в функции сноса и волатильности. Возможно, некоторыми членами в уравнении можно пренебречь или рассматривать их как по- правки к более простому уравнению. Аналогичным образом, необходимо понимать, является ли стохастика ведущим приближением в уравнении или небольшим возмущением детерминированного случая. Часто при помощи скейлинговых замен x ? ?x и t ? ?t и соответ- ствующего выбора констант ? и ? можно уменьшить число значимых параметров системы, сведя их к размерным величинам, характеризую- щим типичные масштабы времени и переменные состояния. Полным решением стохастической задачи является марковская плот- ность вероятности P (x 0 , t 0 ? x, t) . Она может быть задана в виде ре- шения стохастического уравнения x = f(t, ?), выраженного через ска- лярную случайную переменную ?, или в явном виде, как функция, удо- влетворяющая уравнению Фоккера - Планка. Знание плотности вероят- ности позволяет находить различные интегральные величины среднее значение, волатильность, автоковариацию и т.п. Иногда именно они ин- тересуют исследователя, а не полная плотность вероятности. Системы уравнений 177 Перечислим некоторые пријмы, позволяющие получить информацию об изучаемой стохастической системе B Если задача допускает точное решение, то иногда его можно найти непосредственно из стохастического уравнения. Для этого служит алго- ритм на стр. 57 B При помощи формулы Ито и подходящей замены функций можно попытаться свести исходное уравнение к другому, точное или приближјн- ное решение которого получить проще. В частности, всегда возможно изменить функциональную зависимость волатильности шума, изменяя при этом, естественно, и снос уравнения B Если система при длительной эволюции выходит на стационарное решение, то его удобно получать при помощи стационарного уравнения Фоккера - Планка. Для одномерных задач при этом достаточно решить обыкновенное дифференциальное уравнение, которое обычно легко ин- тегрируется. Стационарная плотность вероятности позволяет вычислить асимптотические средние значения наблюдаемых величин. B При помощи динамических уравнений для средних можно полу- чать важные соотношения между наблюдаемыми величинами. Иногда удајтся найти явное изменение во времени или стационарное решение, когда все или часть средних величин перестала изменяться. В послед- нем случае можно получить промежуточное между стационарным и ди- намическим решением. При этом одни переменные состояния являются константами, а эволюция других в этом пределе упрощается. B Выявление особых точек, в которых снос обращается в ноль, и ли- неаризация уравнений в их окрестности дајт важную качественную ин- формацию о характере решений (см. ниже). При этом необходимо про- вести анализ возможных бифуркаций системы при изменении значений еј параметров. Вообще, решение детерминированного уравнения обычно должно предшествовать анализу более сложной стохастической задачи. B Если в задаче есть граничные условия, то плотность вероятности можно представить в виде ряда по ортогональному базису собственных функций (стр. 116 ). B Естественно, многие практически интересные задачи не позволяют получить точного решения. В этом случае на помощь приходят прибли- женные или численные методы, которые мы рассмотрим в девятой главе. 178 Глава 6. • Как и в детерминированном, случае важно уметь понимать харак- тер поведения решения системы стохастических уравнений, не решая их непосредственно. Рассмотрим в качестве примера двумерную задачу с уравнениями: dx = f (x, y)dt + ? 1? (x, y) ?W ? dy = g(x, y)dt + ? 2? (x, y) ?W ? По индексу ? предполагается суммирование от единицы до двух, а компо- ненты вектора ?W ? = {?W x , ?W y } являются независимыми бесконечно малыми изменениями винеровского процесса. Предположим, что в некоторой точке x 0 , y 0 сносы обоих уравнений обращаются в ноль: f (x 0 , y 0 ) = g(x 0 , y 0 ) = 0. Такая точка называется особой. Имеет смысл выяснить поведение реше- ния в еј окрестности. Для этого функции f(x, y), g(x, y) разложим в ряд Тейлора до слагаемых первого порядка малости по отклонениям от осо- бой точки X(t) = x(t)?x 0 , Y (t) = y(t)?y 0 . Для матрицы волатильности ? i? = ? i? (x 0 , y 0 ) возьмјм нулевое приближение, вычислив еј значение в особой точке: dX = (f x X + f y Y ) dt + ? 1? ?W ? dY = (g x X + g y Y ) dt + ? 2? ?W ? , где константа f x это частная производная ?f(x 0 , y 0 )/?x 0 , вычисленная в особой точке (x 0 , y 0 ) , и аналогично f y , g x , g y . Если изучить поведе- ние решений этих уравнений, мы поймјм, как ведут себя в окрестности особой точки (при малых отклонениях от неј X, Y ) и решения общего уравнения. Подобное линейное уравнение мы рассматривали в разделе § 6.4 , стр. 164 Уравнения для среднего совпадают с детерминированными уравнени- ями. Их решение можно искать в виде Ї X(t) = Ae ?t , ЇY (t) = Be ?t , где параметр ? удовлетворяет характеристическому уравнению: det f x ? ? f y g x g y ? ? = (f x ? ?)(g y ? ?) ? f y g x = 0. Это квадратное относительно ? уравнение имеет два решения ? 1 , ? 2 Пример подобного качественного анализа системы уравнений мы рас- смотрим в седьмой главе в рамках модели Охотник - Жертва (стр. 198 ), а пока приведјм классификацию особых точек в двумерном случае. Системы уравнений 179 Возможны следующие случаи: B ? 1 < 0 , ? 2 < 0 устойчивый узел, к которому решение стремится и в котором может находиться сколь угодно долго. При небольших от- клонениях от особой точки отрицательный снос будет возвращать x, y обратно. Естественно, как положение равновесия, так и возврат к нему имеют нерегулярный, стохастический характер. В частности, стохастика может вытолкнуть решение из окрестности особой точки, уведя его в другую часть пространства переменных состояния. B ? 1 > 0 , ? 2 > 0 неустойчивый узел, который решение покидает при любом малом возмущении, которых в стохастическом мире доста- точно. Если предыдущий случай можно представлять как движение в потенциале в виде двухмерного параболоида, то в этом случае парабо- лоид перевјрнут, и решение охотно с него соскальзывает, удаляясь от особой точки. B ? 1 ? 2 < 0 комбинация двух предыдущих случаев, называемая сед- лом. В зависимости от значений параметров f x , f y , g x , g y , это седло опре- делјнным образом повјрнуто в пространстве x, y. Вдоль одного направ- ления (оси лошади) небольшие отклонения будут возвращать решение обратно к особой точке. Перпендикулярно лошади находится направ- ление наибольшей неустойчивости. B ? 1,2 = a ± i ? . Если a < 0, то это затухающий колебательный режим, называющийся фокусом. При a > 0 колебания самовозбуждаются. Если a = 0 , то в системе возникают незатухающие колебания с частотой ? (центр). Как мы видели выше, в стохастическом случае даже при a < 0 движение полностью не останавливается и происходит квазипериодиче- ское колебание вокруг особой точки. Необходимо помнить, что анализ решения в окрестности особой в ли- нейном приближении будет корректен только в случае небольшой вола- тильности. При анализе и стационарного уравнения Фоккера-Планка мы видели, что асимптотически точное решение ( 3.13 ), стр. 89 , для средних совпадает с линейным приближением только в пределе малых волатиль- ностей ?. Это же справедливо и для многомерного случая. Различные типы особых точек обладают качественно различным пове- дением решения в их окрестности. Если мы начнјм медленно изменять параметры системы, то в какой-то момент она скачком может перейти из одного типа решения в другой. Говорят, что при этом произошла би- фуркация, перестройка решения. Анализ подобных возможностей в изу- чаемых системах исключительно важен. 180 Глава 6. Глава 7 Стохастическая природа В этой главе приведены примеры природных систем, которые есте- ственным образом описываются при помощи стохастических дифферен- циальных уравнений. Эти системы охватывают широкий спектр прило- жений от физики до биологии, однако не требуют глубоких познаний в соответствующих областях. Большинство разделов не связаны друг с другом и могут быть прочитаны в любом порядке, независимо друг от друга. Первое стохастическое дифференциальное уравнение в 1908 году записал Поль Ланжевен (Paul Langevin). Именно с него начинается эта глава. 181 182 Глава 7. 7.1 Теория броуновского движения Рассмотрим сферическую частицу радиуса a. При движении со скоро- стью v в жидкости с вязкостью ? на неј действуют сила трения, про- порциональная скорости F f = ?6??av , сила тяжести F g = mg и сила Архимеда F a = ?? 0 gV = ?F g · (? 0 /?) , где ? 0 плотность воды, а ? бро- уновской частицы. Если, кроме этих сил, частица подвержена хаотиче- ским толчкам со стороны молекул воды, то систему уравнений движения можно записать в следующем виде (уравнения Ланжевена): dx = v dt dv = (?g ? v/? ) dt + ? ?W, z y |