Процессы с независимыми приращениями

 

Из определения процесса с независимыми приращениями вытекает (см. [2]), что любые его конечномерные распределения вполне определяются семейством двумерных распределений { F_ts (x, y);
t, s Î T }. Это семейство, таким образом, задает процесс с независимыми приращениями { X_t, t Î T } как процесс в широком смысле. Следовательно, если для процесса с независимыми приращениями известно семейство двумерных распределений, то для его моделирования можно использовать общий подход, изложенный выше (см. 3.1.).

Для двух наиболее простых и, в то же время, наиболее важных представителей рассматриваемого класса с.п., а именно, винеровского и пуассоновского, можно указать простые алгоритмы моделирования.

5.1. Винеровский процесс. Рассмотрим винеровский процесс { w_t, t ³ 0}, который, по определению, является процессом с независимыми гауссовыми приращениями: L(w_t - w_s) = N (0, s ²(t - s)) для любых 0 £ s < t. Для моделирования этого процесса на отрезке [0, T ] делят этот отрезок на n равных частей точками 0 = t₀ < t₁ <... < t_n = T. Значения процесса w_t в этих точках получают последовательно по формулам

 

, i = 1,..., n (18)

 

где n_i, i = 1,..., n - независимые стандартные гауссовы с.в. Для построения приближенной траектории винеровского процесса следует соединить отрезками прямых соседние полученные по формулам (18) точки траектории.

Заметим, что выше (см. 3.2.) рассматривался другой способ (приближенного) моделирования винеровского процесса.

5.2. Однородный пуассоновский процесс. Этот процесс, часто называемый также простейшим потоком событий, является процессом с независимыми приращениями { p_t, t ³ 0}, начинающимся из нуля и имеющим пуассоновские приращения: L(p_t - p_s) = П (l (t - s)) для любых 0 £ s < t. Параметр
l > 0 процесса называется интенсивностью потока событий. Ясно, что этот процесс принимает целые неотрицательные значения для каждого t ³ 0. Можно показать, что значения процесса p_t возрастают в случайные моменты времени t₁ < t₂ <... < t_k <... единичными скачками, а в промежутках между соседними скачками остаются неизменными. При этом промежутки между соседними скачками являются независимыми, одинаково распределенными по показательному закону с параметром l случайными величинами. Отсюда, естественно, получается алгоритм моделирования однородного пуассоновского процесса на промежутке [ 0, T] достаточно смоделировать моменты скачков по формулам

 

, k =1, 2,... (19)

 

где t ₀ = 0, q_k, k =1,2,... - независимые с.в. с функцией распределения

 

. (20)

 

Процедура моделирования заканчивается, как только впервые t_k превысит T. Смоделированная реализация будет кусочно-постоянной функцией, возрастающей единичными скачками в моменты t_k Î T.

5.3. Неоднородный пуассоновский процесс. Так называют процесс с независимыми приращениями X_t, t ³ 0, начинающийся из нуля, приращения которого X_t - X_s для любых 0 £ s < t имеют пуассоновское распределение с параметром

 

, (21)

 

где положительная ограниченная функция l (t), t ³ 0 называется мгновенной интенсивностью (в момент t) пуассоновского потока событий. Можно показать, что, как и в однородном случае, значения процесса X_t возрастают единичными скачками в случайные моменты времени t ₁ < t ₂ <... < t _k <..., а в промежутках между соседними скачками остаются постоянными. При этом промежуток q_k+1 = t _k+1 - t _k, k = 0,1,..., t ₀ = 0 имеет закон распределения, зависящий от t _k:

 

, 0 < x < ¥. (22)

 

Как следует из (22), для того, чтобы промежутки между соседними скачками были с вероятностью 1 конечными, необходимо, чтобы мгновенная интенсивность l (t) удовлетворяла условию

 

. (23)

 

Для моделирования неоднородного процесса Пуассона с интенсивностью l(t) на промежутке [0, T], очевидно, достаточно смоделировать моменты скачков по формулам (19), где величины q _k теперь следует моделировать последовательно в соответствии с функцией распределения (22). Процедура моделирования останавливается, как только впервые t_k превзойдет величину T.

Величины q _k+1 с распределением (22) можно моделировать с помощью "универсального метода" (см. гл. 1). При этом значение q _k+1 получается как решение уравнения

 

, где L(U_k) = R_[0,1].

 

Марковские процессы

 

6.1. Общие марковские процессы. Рассмотрим случайный процесс { X _t, t Î T } со значениями в измеримом "фазовом" пространстве (E, e). Его называют марковским, если при любом фиксированном X _s = x значения процесса X _t, t > s не зависят от X _u, u < s. Условную вероятность события { X _t Î A Ì E} при условии X _s = x называют переходной функцией (вероятностью) марковского процесса и обозначают P (s, x, t, A). Переходная функция для любых s < u < t удовлетворяет уравнению

 

(24)

 

называемому уравнением Чепмена-Колмогорова. Оно выражает важное свойство марковского процесса — отсутствие последействия: если известно состояние системы в некоторый момент времени u, то вероятности перехода из этого состояния не зависят от движения системы в моменты времени, предшествующие u. Если существует минимальный элемент t ₀ множества T, то вероятности P ₀(A) = P { Î A } для всех A Î e образуют так называемое начальное распределение процесса { X _t, t Î T }. Начальное распределение и переходная функция марковского процесса вполне определяют его как процесс в широком смысле. Это значит, что при заданных начальном распределении и переходной функции марковского процесса, как отмечено в начале этой главы, есть принципиальная возможность его моделирования.

Если существует неотрицательная функция p (s, x, t, y) такая, что

 

 

для всех s, x, A s < t, то ее называют переходной плотностью марковского процесса. Марковский процесс называется однородным, если его переходная вероятность зависит не от двух временных аргументов s и t, а только от их разности t - s:

 

(25)

 

Мы рассмотрим далее моделирование некоторых подклассов марковских процессов.

6.2. Цепи Маркова. Марковский процесс { X _t, t Î T }, у которого дискретно как временное множество T, так и фазовое пространство E, называют цепью Маркова. В этом случае без потери общности можно считать, что
T = {0, 1, 2...} и E = {0, 1, 2,...}. Переходную функцию при этом достаточно определить для одноточечных множеств A _k = { k }, k = 0, 1,..., и целесообразно использовать обозначение

 

P (m, i, n, k) = p_ik (m, n). (26)

 

Определенная таким образом величина p_ik (m, n) называется переходной вероятностью цепи Маркова и представляет собой условную вероятность оказаться в k -том состоянии на n -м шаге при условии, что на m -том шаге (m < n) она находилась в состоянии i. Уравнение Чепмена-Колмогорова для переходной вероятности будет иметь вид

 

(27)

 

для любых m < r < n и i, k Î E. Если наряду с переходной вероятностью (26) задано также начальное распределение

 

p _k = P { X ₀ = k }, k = 0, 1,..., (28)

 

то будут определены и конечномерные распределения цепи Маркова: для любых n шагов с номерами m ₁ < m ₂ <... < m _n и любых состояний i_r Î E, имеем

 

(29)

Моделирование значений , ,..., цепи Маркова сводится, таким образом, к моделированию случайного вектора с заданным законом распределения вероятностей: следует сначала смоделировать с.в. в соответствии с безусловным распределением

 

, k = 0, 1,..., (30)

 

а последующие значения , 1 < r £ n, цепи - в соответствии с переходными (условными) вероятностями .

Моделирование цепи Маркова в частных случаях может упрощаться. Так для конечной цепи Маркова, у которой E = {0, 1,..., N}, суммы в правых частях формул (27) и (30) будут конечными. Соответственно, более простой будет и процедура моделирования, поскольку будут конечными как совокупности вероятностей начальных состояний, так и совокупности переходных вероятностей.

Другое упрощение следует из однородности цепи Маркова, что выражается равенством p_ik (m, m + n) = p_ik (n) для любых m ³ 0, n ³ 1. Пользуясь уравнением Чепмена-Колмогорова, легко показать, что переходные вероятности однородной цепи Маркова за n шагов вполне определяются ее переходными вероятностями p_ik (1) = p_ik, i, k Î E за один шаг. Если - матрица переходных вероятностей за n шагов, а P = - матрица переходных вероятностей за один шаг, то легко убедиться, что P ⁽ⁿ⁾ = P ⁿ, где в правой части стоит n -я степень матрицы P. Это значит, что однородная цепь Маркова как процесс в широком смысле вполне задается начальным распределением и переходными вероятностями за один шаг. Отсюда вытекает алгоритм моделирования последовательных значений цепи Маркова: X ₀, X ₁, X ₂,.... Для получения значения с.в. X ₀ моделируем дискретное распределение
{p _k, k = 0, 1,... }. Пусть получено значение i ₀. Тогда значение с.в. X ₁ моделируем в соответствии с распределением { , k = 0, 1,...}. Этот процесс продолжается до получения реализации цепи Маркова нужной "длины".

Разумеется, моделирование будет еще проще, когда однородная цепь Маркова является конечной.

6.3. Марковские процессы с непрерывным временем и дискретным фазовым пространством. Рассмотрим марковский процесс { X_t, t Î T }, у которого пространство состояний E = {0, 1,...} - конечно или счетно, а
T = [0,¥). Предположим, что он однороден во времени и обозначим его переходную вероятность

 

P (s, i, s + t, k) = p _ik (t) (31)

 

для любых i, k Î E и s ³ 0, t ³ 0. Будем также предполагать, что с вероятностью 1 реализации x _t(w) рассматриваемого процесса непрерывны справа в каждой точке t ³ 0. Состояние i будем называть поглощающим, если p _ii(t) = 1 для всех t > 0. Остальные состояния будем называть непоглощающими. Пусть t - момент первого выхода из начального состояния.

Утверждение 1. Если начальное состояние X ₀ = i является непоглощающим, t - момент выхода из начального состояния, то при некотором
l_i > 0

P _i {t > t } = exp{-l_it}, (32)

где P _i обозначает условную вероятность при условии X₀ = i; при этом величины t и X_t независимы.

Доказательство этого утверждения можно найти в [7].

Замечание 1. Для поглощающих состояний j будем считать l _j = 0. Так как в этом случае для любого t P _j {t > t } = 1, то формула (32) остается справедливой.

Обозначим

 

P _i {X_t = j } = p _i (33)

 

и рассмотрим величины: t₁ - момент выхода из начального состояния, , t₂ - время, которое процесс пробудет в состоянии , q₂ = t₁ + t₂, ,.... Если , t_n-1, q_n-1 определены, то t _n - время, которое процесс находится в состоянии с момента q_n-1 до выхода из этого состояния,
q_n = q_n-1 + t _n, = .

Утверждение 2. Величины , ,... образуют однородную цепь Маркова, для которой вероятности перехода за один шаг даются равенством (33); условное распределение t₁,...,t_n при условии, что заданы совпадает с распределением n независимых показательно-распределенных с.в. с параметрами .

Утверждение 2 позволяет описать поведение рассматриваемых марковских процессов. Имеется цепь Маркова , n = 0, 1,..., называемая вложенной цепью Маркова, с вероятностями перехода за один шаг p_ij. Последовательные состояния этой цепи определяют последовательные состояния процесса с непрерывным временем. Время пребывания процесса X_t в состоянии является показательно-распределенной величиной (при условии, что задано = ) с параметром и не зависит от времени пребывания в предыдущих состояниях. Таким образом, зная вложенную цепь Маркова, можно построить процесс Маркова с непрерывным временем до момента

 

 

Процесс называется регулярным, если для всех i Î E

 

Таким образом, моделирование однородного марковского процесса с непрерывным временем и дискретным фазовым пространством сводится к моделированию вложенной однородной цепи Маркова, рассмотренному в предыдущем пункте настоящего параграфа, и моделированию последовательности независимых показательно-распределенных с.в.

6.4. Диффузионные марковские процессы. Так называется класс непрерывных марковских процессов, которые могут служить математической моделью явления диффузии - перемещения частицы под воздействием молекул среды, находящихся в хаотическом тепловом движении. Формально диффузионный марковский процесс можно определить через его переходную функцию. Пусть E - конечномерное евклидово пространство - фазовое пространство процесса X = { X_t, t ³ 0}; P (s, x, t, A) - его переходная функция. Прежде всего будем предполагать, что процесс X с вероятностью 1 имеет реализации без разрывов второго рода. Теперь сформулируем условия диффузионности: для x Î E, t ³ 0, e > 0

1)

где V _e(x) - шар радиуса e;

2) ,

где a(t, x) - измеримая функция на [0,¥)´ E со значениями в E;

3) для всех z Î E

 

,

где (·,·) - скалярное произведение в E, B (t, x) - измеримая симметричная неотрицательно определенная матричная функция на [0,¥) ´ E, действующая из E в E. Функции a (t, x) и B (t,x) называются диффузионными коэффициентами процесса: a (t,x) - вектором (коэффициентом) переноса, B(t,x) - оператором (коэффициентом) диффузии.

Условие 1) вместе с отсутствием разрывов 2-го рода является условием непрерывности процесса. Условие 2) указывает на наличие движения среды, в которой находится частица: частица из точки x в момент t за время h в среднем смещается на a (t,x) h. Условие 3) показывает, что матрица B (t,x) характеризует величину среднего квадратического отклонения частицы от ее положения x в момент t. Если среда изотропна (ее свойства одинаковы по всем направлениям), то B (t,x) = b ²(t,x) I, где b ²(t,x) - числовая функция, называемая коэффициентом диффузии, I - единичная матрица.

Далее для простоты будем рассматривать процессы в одномерном фазовом пространстве E = R¹. Широкий класс диффузионных процессов с непрерывными коэффициентами образуют решения стохастических дифференциальных уравнений

 

dX _t = a (t,X_t) dt + b (t, X_t) dw _t, (34)

 

или, в интегральной форме,

 

. (35)

 

Второй интеграл в правой части последнего уравнения, называемый стохастическим интегралом по винеровсому процессу { w _t, t ³ 0} просто определяется в случае непрерывной зависимости подынтегрального процесса от s. Пусть непрерывный с вероятностью 1 случайный процесс { Y_t, t ³ 0} согласован с винеровским процессом в следующем смысле: приращения w_t - w_s при t > s не зависят от совокупности случайных величин Y_u, w _u, u £ s. Тогда существует предел в смысле сходимости по вероятности интегральных сумм

 

, 0 = s ₀ <... < s_n = t

 

при который называется стохастическим интегралом и обозначается символом

 

.

 

Решением уравнения (35) с непрерывными по совокупности переменных коэффициентами a (t,x) и b (t,x) называется согласованный с винеровским процессом процесс { X_t, t ³ 0} с непрерывными траекториями, такой, что равенство (35) справедливо для всех t > 0 с вероятностью 1.

Подробнее о стохастических интегралах и стохастических дифференциальных уравнениях см. [2].

 

Решения стохастических дифференциальных уравнений моделируются путем перехода к разностной схеме. Так, если процесс { X_t, t ³ 0} определяется уравнением (34), то простейшая разностная схема имеет вид

 

,

где k = 1, 2,... 0 = t ₀ < t ₁ <... < t_{k -1} < t_k <....

Замечание 2. Попытки формального применения метода Рунге-Кутта с целью повышения точности моделирования решения стохастического дифференциального уравнения оказываются несостоятельными, поскольку винеровский процесс не обладает никакими производными.

 

Стационарные процессы

 

7.1. Общий подход к моделированию стационарных процессов. Мы рассматриваем здесь комплекснозначные процессы X = { X_t, t Î T }, у которых параметрическое множество T либо совокупность целых чисел {0, ±1, ±2,...} (случай дискретного времени), либо множество всех действительных чисел
(-¥, ¥) (случай непрерывного времени). Если у процесса X имеются моментные функции до второго порядка включительно, причем среднее значение M X_t º m = const, а корреляционная функция

 

, (36)

 

то X называется стационарным в широком смысле. Употребляемый в дальнейшем без оговорок термин "стационарный" подразумевает стационарность в широком смысле.

В рамках теории гильбертовых случайных процессов основными характеристиками стационарного процесса являются моментные функции первого и второго порядков, а именно среднее значение M X_t º m и корреляционная функция процесса. Для простоты (и без потери общности) в дальнейшем считаем, что m = 0. Центральным результатом теории стационарных процессов является теорема А. Я. Хинчина о спектральном представлении процесса и его корреляционной функции.

Теорема 2. Каждый стационарный процесс X допускает спектральное представление

(37)

в виде стохастического интеграла по комплекснозначному процессу с ортогональными приращениями Z. При этом корреляционная функция R(t) процесса X допускает представление

, t Î T (38)

в виде интеграла Стилтьеса по неубывающей непрерывной слева функции

ограниченной вариации F(w). Интегрирование в (37) и (38) производится в пределах -p £ w £ p в случае дискретного времени t, и в пределах -¥ < w < ¥ - в случае непрерывного t.

Замечание 1. Функция F (w) в представлении (38) называется спектральной функцией процесса X. Она однозначно может быть восстановлена по корреляционной функции R (t). Особенно просто это делается в случае, когда F (w) абсолютно непрерывна. Тогда ее производная f (w) = dF (w)/dw называется спектральной плотностью процесса. В этом случае "формулы обращения" имеют вид

 

, w Î [-p,p] (39)

в случае дискретного времени и

 

, w Î (-¥,¥) (40)

 

в случае непрерывного времени.

Замечание 2. Спектральная функция F (w) стационарного процесса X служит структурной функцией процесса с ортогональными приращениями Z, участвующего в представлении (37) процесса X:

 

(41)

 

для любых точек w₁ < w₂ спектра процесса X.

Спектральное представление (37) стационарного процесса через процесс более простой структуры (процесс с ортогональными приращениями) можно положить в основу приближенного моделирования процесса X. Для этого разобьем интервал, в котором 0 < F (w) < F(+¥), на интервалы D_k и возьмем в D_k произвольную точку w _k. Тогда, если F_k - приращение F (w) на интервале D_k, то можно записать:

 

 

и, следовательно

 

, (42)

 

где { Z_k } - ортогональные приращения процесса Z на интервалах { D_k }. Таким образом, формула (42) указывает алгоритм моделирования значений процесса X. Для работы этого алгоритма необходимо уметь моделировать процесс с ортогональными приращениями Z со структурной функцией F (w). Моделирование приращений Z_k в формуле (42) возможно, в частности, когда Z является центрированным процессом с независимыми приращениями с заданным семейством двумерных распределений.

Изложенный выше общий подход к моделированию стационарных процессов далеко не всегда удается реализовать на практике. Поэтому ниже мы рассмотрим частные случаи стационарных процессов, для которых задача моделирования решается в силу специфической структуры этих процессов.

7.2. Моделирование процесса с дробно-рациональной спектральной плотностью. Будем рассматривать одномерный стационарный гауссов процесс X_t с непрерывным временем, обладающий спектральной плотностью вида

 

(43)

 

где P (z) и Q(z) - многочлены от z с действительными коэффициентами:

 

, , a ₀ = 1,

 

причем m < n и нули Q(z) лежат в левой полуплоскости. Представление f (w) в виде (43) означает, что процесс X_t можно получить в виде линейной комбинации процесса Y_t и m - 1 его производных:

 

, (44)

 

где Y_t является решением линейного стохастического дифференциального уравнения с постоянными коэффициентами

 

, (45)

 

где w _t, t ³ 0 - винеровский процесс. Если ввести вектор (столбец)
V _t = , то уравнение (45) запишется в векторной форме

 

d V _t = AV _t dt + d W _t, (46)

 

где

 

,

 

Заметим, что стационарное решение уравнения (45) имеет спектральную плотность

 

. (47)

 

Если вектор V _t будет найден, то значение процесса X_t находится по формуле (44) как линейная комбинация компонент вектора V _t.

Ниже излагается метод моделирования процесса X_t, предложенный Франклиным. Согласно уравнению (46), если известно значение процесса в момент t, то значение через промежуток времени D _k определяется следующим образом:

 

, (48)

 

где r (t) - гауссов случайный вектор,

 

.

 

Этот вектор некоррелирован со случайным вектором

 

.

 

Чтобы смоделировать r (t), нужно найти его корреляционную матрицу

 

R _r = = ,

 

где T - знак транспонирования, а

 

 

Можно показать, что R _r удовлетворяет соотношению

если i + j нечетно в противном случае

,

 

откуда можно найти симметрическую матрицу R _r. Следовательно, можем смоделировать гауссов вектор r (t) и получить очередное значение V _{t+ D t}, а с ним, по (44), и значение процесса X.

Чтобы смоделировать начальное значение процесса V ₀, достаточно знать корреляционную матрицу в нуле:

 

.

 

Можно показать, что компоненты этой матрицы находятся по формуле

 

 

где k _l можно определить из соотношений

 

 

и a ₀ = 1. После нахождения матрицы K (0,0), гауссов вектор V ₀ может быть смоделирован, как описано в гл.II, и по формуле (44) получается X ₀. Последующие значения процесса V _t, а следовательно, и процесса X_t, будут получаться по формулам (48) и (44) соответственно.

7.3. Случайные последовательности АРСС. Простота структуры случайных последовательностей авторегрессии и скользящего среднего (АРСС) и возможность использования их для аппроксимации широкого класса случайных последовательностей определяют как практический, так и теоретический интерес к ним. Моделирование этих последовательностей позволяет решать самые разнообразные прикладные задачи, связанные с изучением реальных процессов в науке и технике. Случайная последовательность (с. п.) АРСС с p вещественными параметрами авторегрессии b ₁,..., b_p и q вещественными параметрами скользящего среднего a ₀ = 1, a ₁,..., a _q описывается уравнением

 

, (49)

 

где Z = { Z_t, t = 0, ±1, ±2,... } - последовательность некоррелированных, одинаково распределенных с.в. с M Z_t = 0 D Z_t = . Последовательность Z часто называют входной последовательностью или шумом. В дальнейшем с. п., задаваемую уравнением (49), будем обозначать АРСС(p,q). Существует обширная литература по теории с.п. АРСС. Необходимые для целей данного пособия факты теории процессов АРСС можно найти в [9]. В частности, известно, что с.п. АРСС(p,q) будет стационарной при условии, что все корни характеристического уравнения

 

 

лежат внутри единичного круга | z | < 1. Условие стационарности в дальнейшем предполагается выполненным.

Умножая (49) на X_t-k и переходя к математическим ожиданиям, получаем

 

, (50)

 

где R_XZ (k) = M Z_tX_t-k. Так как X_t-k зависит только от членов входной последовательности до момента t - k, то очевидно, что R_XZ (k) = 0 при k>0 и R_XZ (k) ¹ 0 для k £ 0. Значения R_XZ (k), k = 0, -1, -2,... находятся последовательно следующим образом. Умножим на Z_t обе части равенства (49) и перейдем к математическим ожиданиям в полученном равенстве. В результате получаем R_XZ (0) = . Затем обе части (49) умножим на Z_t-1 и перейдем к математическим ожиданиям. Получим R_XZ (-1) = (a ₁ + b ₁) . Продолжая этот процесс, найдем необходимое число корреляционных моментов R_XZ (k), k £ 0.

Из (50) следует, что

 

, для k ³ q + 1 (51)

 

Для решения разностных уравнений (51) (для больших k) в качестве начальных необходимы p значений, например, R_X (q), R_X (q - 1),..., R_X (q - p + 1). Дисперсию = R_X (0) с.п. АРСС(p, q) вместе с R_X (1),..., R_X (p) получим, решая систему уравнений, получающихся из (50) при k = 0, 1, 2,..., p. Спектральная плотность с.п. АРСС(p, q) имеет вид

 

, w Î [-p, p]. (52)

 

При моделировании с.п. АРСС(p,q) предполагаются известными (заданными) векторы параметров a и b, а также предполагается известным распределение величин Z_k входной последовательности Z. Кроме того, величины входного случайного шума Z_k обычно считаются независимыми, а не только некоррелированными.

Для того, чтобы с помощью уравнения (49) можно было сразу начать моделирование стационарной последовательности, необходимы предварительные данные. Действительно, если на первом шаге (t = 1) мы хотим получить значение X₁, являющееся элементом стационарной последовательности X, мы должны в формулу (49) подставить q + 1 значений входной последовательности: Z_{-q+ 1},..., Z ₀, Z₁ и p предшествующих X₁ значений выходной последовательности: X_{-p+ 1},..., X _-1, X ₀. Следовательно, для начала процесса моделирования последовательности X нужно смоделировать случайный вектор

 

(53)

 

Это практически невозможно, если "шум" Z является негауссовым, поскольку невозможно найти распределение вектора U.

Моделирование с.п. АРСС(p,q) в случае гауссова шума начинается с моделирования вектора (53). Поскольку на уравнение (49) можно смотреть как на линейное преобразование "шума" Z в выходной процесс X, то в случае гауссова шума, как известно, линейное преобразование не выводит за пределы гауссовой системы. Таким образом, вектор U будет гауссовым с нулевым средним и корреляционной матрицей R _U, имеющей вид

 

 

где R _Z - корреляционная матрица вектора Z ^q+1= (Z_{-q+ 1},..., Z ₁) - это диагональная матрица, у которой диагональные элементы одинаковы и равны , R _X - корреляционная матрица вектора X ^p = (X_{-p+ 1},..., X ₀), элементы которой находятся из уравнений (50) и (51), а матрица R _XZ размерности
(q + 1) ´ p имеет своими элементами корреляционные моменты R_XZ (k) между координатами векторов Z ^{q+ 1} и X ^p. Алгоритм вычисления корреляционных моментов R _XZ(k) описан выше. Итак, корреляционную матрицу R _U можно считать известной. Используя алгоритм моделирования гауссовского вектора (см. гл. II), смоделируем вектор U. Подставляя его в уравнение (49) при
(t = 1), найдем значение X₁. На следующем шаге моделируется только очередное значение шума Z ₂; затем в правую часть (49) при (t = 2) подставляется вектор

 

 

и в результате имеем значение X ₂ выходной последовательности. Процесс продолжается до получения "отрезка" стационарной с.п. АРСС(p,q) заданной длины.