Модифицированные статистические алгоритмы фильтрации и прогнозирования в непрерывных стохастических системах Текст научной статьи по специальности «Математика»

Известия Института математики и информатики УдГУ

2015. Вып. 2 (46)

УДК 519.218 © К. А. Рыбаков

МОДИФИЦИРОВАННЫЕ СТАТИСТИЧЕСКИЕ АЛГОРИТМЫ ФИЛЬТРАЦИИ И ПРОГНОЗИРОВАНИЯ В НЕПРЕРЫВНЫХ СТОХАСТИЧЕСКИХ СИСТЕМАХ1

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

Ключевые слова: ветвящийся процесс, оптимальная фильтрация, прогнозирование, статистическое моделирование, стохастическая система.

§ 1. Постановка задачи

Будем рассматривать модель объекта наблюдения, описываемую стохастическим дифференциальным уравнением Ито [4]:

бХ(Ь) = / (¿,Х(Ь))б + ст(Ь,Х(Ь))бШ(Ь), Х(Ьо) = Хо, (1.1)

где Х € К" — вектор состояния, Ь € [Ьо, Т + Д(Т)| — время; /(Ь, ж): [Ьо, Т + Д(Т)| х К" — К" — вектор-функция размеров п х 1; а(Ь,ж): [Ьо,Т + Д(Т)] х К"- — К"х8 — матричная функция размеров п х з; Д(Ь): [Ьо, Т] — [0, — величина опережения по времени:

ш«(( + Д(г)) = Т + Д(Т);

Ш(Ь) — з-мерный стандартный винеровский процесс, не зависящий от начального состояния Хо. Модель измерительной системы представим в форме

2 (Ь) = С(Ь,Х (Ь)) + С (Ь)^ (Ь), (1.2)

где 2 € К™ — вектор измерений, Ь € [Ьо, Т]; с(Ь, ж): [Ьо,Т] х К" — К™ — вектор-функция размеров т х 1; £(Ь): [Ьо,Т] — КтхЛ — матричная функция размеров т х б, удовлетворяющая условию

с (Ь)с Т(Ь)=0;

N (Ь) — б мерный стандартный гауссовский белый шум [4,8].

Задачи фильтрации и прогнозирования состоят в нахождении оценки Х(Ь + Д(Ь)) по результатам измерений 2 = {2(т), т € [Ьо, Ь)}, для случая Д(Ь) = 0 — это задача фильтрации, а для Д(Ь) > 0 — задача прогнозирования.

При формировании новых алгоритмов оптимальной фильтрации и прогнозирования, как и ранее [6], будем исходить из несмещенности оценки и минимума среднеквадратического отклонения. Тогда

X (Ь + Д(Ь)) = М[Х (Ь + Д(Ь))|^о ],

где М — знак математического ожидания.

1 Работа поддержана РФФИ (грант № 13-08-00323-а).

§ 2. Метод решения задач оптимальной фильтрации и прогнозирования

В основе предлагаемого метода к решению задач оптимальной нелинейной фильтрации и прогнозирования лежит применение уравнения Дункана-Мортенсена-Закаи, которое описывает эволюцию ненормированной апостериорной плотности вероятности, характеризующей распределение оцениваемого вектора состояния объекта наблюдения по результатам измерений на фоне помех. По сравнению с уравнением для априорной плотности вероятности, то есть уравнением Фоккера-Планка-Колмогорова, используемым для получения прогноза вектора состояния, уравнение Дункана-Мортенсена-Закаи содержит дополнительное слагаемое, а именно ненормированную плотность вероятности с коэффициентом, зависящим в общем случае от времени, вектора состояния объекта наблюдения и текущих измерений [4,7]. Если в данном коэффициенте выделить две составляющие, первая из которых содержит отрицательные значения исходного коэффициента, а вторая — положительные, то полученные компоненты будут иметь вполне определенный физический смысл — сток и источник тепла или вещества, поскольку уравнение Дункана-Мортенсена-Закаи — это линейное уравнение параболического типа, оно может описывать процессы теплопроводности и диффузии. Вероятностный смысл этих составляющих — интенсивности пуассоновских потоков поглощения и восстановления траекторий случайного процесса [1] с теми же коэффициентами сноса и диффузии, что задают исходный объект наблюдения (1.1).

В работах [5, 6] показано, что оптимальная оценка вектора состояния может быть найдена приближенно с помощью моделирования траекторий вспомогательного случайного процесса с обрывами и ветвлениями, то есть с учетом поглощения и восстановления траекторий. Каждая из ветвей, рассматриваемая отдельно, представляет собой часть траектории процесса, описываемого уравнением (1.1), результат измерений (1.2) оцениваемого вектора состояния влияет на распределение моментов времени обрывов и появления новых ветвей. Обрывы и ветвления траекторий образуют пуассоновские потоки событий. В наиболее простом случае при использовании критерия минимума среднеквадратической ошибки оценка вектора состояния объекта наблюдения может быть получена усреднением по построенному в результате моделирования ансамблю траекторий. В более сложном случае по этому же ансамблю можно оценить апостериорную функцию распределения или апостериорную плотность вероятности, а далее на основе этих характеристик получить искомую оценку вектора состояния. Без обрывов и ветвлений можно оценить априорное распределение, а с обрывами и ветвлениями — апостериорное, что можно использовать для реализации алгоритмов фильтрации при перерывах в работе измерительной системы или для задач прогнозирования.

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

где Якг(Ь) — элементы квадратной матрицы, обратной по отношению к произведению ((Ь)(т(Ь).

Функция Л(Ь, х, г) выражается через коэффициенты сноса и диффузии измерительной системы (1.2) и зависит от времени, вектора состояния вспомогательной стохастической системы с обрывами и ветвлениями траекторий, а также текущих измерений оцениваемого вектора состояния исходной системы наблюдения.

Таким образом, обрывы траекторий образуют неоднородный пуассоновский поток событий с интенсивностью Л-(Ь) = Л-(Ь,X(Ь)^(£)), определяемой формулой

а ветвления траекторий образуют другой неоднородный пуассоновский поток событий, интенсивность которого Л+(Ь) = Л+(Ь,Х(Ь)^(Ь)) задается следующим образом:

т т

Л (Ь,х,г)

-Л(Ь,х,г), Л(Ь,х,г) < 0, 0, Л(Ь,х,г) ^ 0,

Л+(Ь, х, г)

Л(Ь, х, г), Л(Ь, х, г) > 0, 0, Л(Ь, х, г) ^ 0,

следовательно, вероятность Р-(Ь, ДЬ) обрыва траектории и вероятность Р+(Ь, ДЬ) ветвления траектории на промежутке времени [Ь, Ь + ДЬ] при X (Ь) = (Ь) = г можно представить в виде

Р-(Ь, ДЬ) = Л-(Ь, х, г) ДЬ + о(ДЬ), Р+(Ь, ДЬ) = Л+(Ь,х,г)ДЬ + о(ДЬ). (2.1)

Как и в [6], чтобы получить оценку вектора состояния объекта наблюдения в текущий момент времени или его прогноз в будущий момент времени предлагается моделировать траектории вспомогательного случайного процесса с коэффициентами сноса /(Ь, х) и диффузии а(Ь, х), определяющими исходный объект наблюдения (1.1). Начальное условие для этого вспомогательного случайного процесса совпадает с заданным начальным условием X(Ьо) = Хо. На этапе фильтрации траектории могут обрываться или разветвляться, вероятности этих событий заданы выше; на этапе прогнозирования при отсутствии измерений обрывы и ветвления невозможны.

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

Хк+1 = Хк + Н!{1к,Хк) + уДь(т{1к,Хк)Шк, (2.2)

гк = фк,Хк) + ^АУк, (2.3)

4 = Ьо + Нк, к = 0,1,...,

где ДЖ^ и ДУк — случайные векторы размеров в х 1 и й х 1 соответственно, их координаты независимы и имеют стандартное нормальное распределение, Н — шаг численного интегрирования, } — равномерная сетка по времени. Для других методов численного решения стохастических дифференциальных уравнений алгоритмы формируются аналогично. В приведенных соотношениях нижний индекс к указывает на соответствие векторов Х^, Zfc, и ДУк моменту времени и не является номером координаты вектора.

Н

поднялось условие |Л(Ь, х, г)| Н ^ 1. Проверку условий обрывов или ветвлений предлагается осуществлять только в узлах сетки {Ь^}. Условие обрыва траектории на промежутке , + Н] согласно (2.1) можно записать в виде а ^ Л-(Ь&,Х&, Zfc) Н, а условие ветвления траектории на том же промежутке , + Н] — в виде а ^ Л+(^, Х&, Zfc) Н, где а — случайная величина,

(0, 1)

Момент времени обрыва или ветвления траектории при выполнении условия

а < |Л(Ьк,Хк^)| Н

можно положить равным = ^ + Н, то есть взять следующий по отношению к ^ узел сетки {Ь^}, поскольку дадьнейшее усреднение по имеющейся выборке величин Х& для получения оценки вектора состояния объекта наблюдения должно проводиться именно в узлах. Отметим, что указанное условие проверяется для каждой траектории вспомогательного случайного процесса независимо при одном и том же векторе измерений Zfc оцениваемого состояния. При обрыве ветвь траектории не учитывается при дальнейшем моделировании и усреднении, при ветвлении в выборку величин Х& добавляется «копия», для которой а ^ Л+(^,Х&, Zfc) Н.

В отличие от предыдущих работ [5,6] здесь обрывы и ветвления траекторий вспомогательного случайного процесса происходят только в узлах сетки }, построенной для численного решения стохастического дифференциального уравнения (1.1) (в рассматриваемом случае это равномерная сетка), что соответствует подходу, принятому при реализации фильтров частиц [9,10]. Ранее для построения алгоритмов оптимальной фильтрации и прогнозирования при моделировании пуассоновских потоков обрывов и ветвлений траекторий вспомогательного случайного процесса применялся метод «максимального сечения», согласно которому моделировались моменты времени возможных обрывов и ветвлений с их последующим прореживанием [2,3]. Промежутки времени между такими соседними моментами времени должны иметь

показательное распределение с постоянным параметром Л*, ограничивающим сверху значения функции |Л(Ь, ж, г)|, а прореживание происходит с вероятностью

¡х^тгт

л*

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

Опишем методику моделирования траекторий вспомогательного случайного процесса с обрывами и ветвлениями и методику оценивания вектора состояния объекта наблюдения в текущий и будущий моменты времени по результатам измерений.

На первом этапе необходимо смоделировать начальные состояния Хо для траекторий вспомогательного случайного процесса с обрывами и ветвлениями. По полученной в результате моделирования выборке величин Хо заданного объема М находим среднее — оценку Х"(Ьо) вектора состояния объекта наблюдения в начальный момент времени ¿о- Используя начальные Хо

дения, с выбранным шагом численного интегрирования Н, моделируем ансамбль траекторий вплоть до момента времени Ьо + Д(Ьо) и) усредняя, получаем прогноз Х"(Ьо + Д(Ьо))-

Используя выборку величин Хо и измерение 2о, модифицируем выборку величин Х1, соответствующую моменту времени ^ = Ьо + Н, а именно элементы Х1, для которых справедливо условие обрыва а ^ |Л(Ьо,Хо, 2о)| Н и Л(Ьо,Хо,2о) < 0, удаляем из выборки, а элементы, удовлетворяющие условию ветвления а ^ |Л(Ьо,Хо,2о)| Н и Л(Ьо,Хо,2о) > 0, дублируем. Здесь

Хо Х1

а

(0, 1) Х1

этапа прогнозирования при Д(Ьо) ^ Н или смоделирована согласно разностной схеме (2.2). По модифицированной выборке величин Х1 находим среднее — оценку Х"(Ь1) вектора состояния объекта наблюдения в текущий момент времени £1. Затем, применяя приведенную выше раз-

Х1

начальную, моделируем ансамбль траекторий вплоть до будущего момента времени ^ + Д(Ь1) и, усредняя, получаем прогноз + Д(^)).

Далее повторяем описанную процедуру для следующих моментов времени Ь2,Ьз, • • • ДО конечного момента времени Т, получая таким образом оценки ) и + Д(Ь&)).

Т

и конечный момент времени + Д(Ь&) для нахождения прогноза согласованы с шагом интегрирования, то есть найдутся такие целые числа К1 и , что Ьо + К1Н = Т и Ьо + Н = ^ + Д(Ь&), иначе при нахождении вектора состояния в этот момент времени потребуется корректировать шаг численного интегрирования.

В следующем разделе эта методика будет представлена в виде детальных алгоритмов.

Дополнительно отметим, что если при моделировании количество обрывов существенно больше количества ветвлений траекторий, то объем выборки величин Х^ сокращается, в таком случае необходимо увеличивать ее объем искусственно. Это можно делать, например, по оценке апостериорной функции распределения на основе имеющейся выборки, генерируя новые векторы, или использовать «размножение» выборки. Если при моделировании количество ветвлений существенно больше количества обрывов траекторий, то объем выборки величин Х^ может значительно вырасти, и тогда ее следует «прореживать», например, случайным образом. На этапе фильтрации объем выборки величин Х& меняется, и для каждого значения индекса к, соответствующего времени этот объем может отличаться от начального объема М, выбран-

Хо

прогнозирования объем выборки фиксирован, он равен объему выборки на последнем шаге фильтрации.

§ 3. Алгоритмы совместного моделирования системы наблюдения, фильтрации и прогнозирования

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

§3.1. Алгоритм моделирования системы наблюдения и фильтрации

Шаг 1. Задать М — число моделируемых вспомогательных траекторий; Н — шаг численного интегрирования. Получить реализации начальных состояний Хо и Х0, где Хо — начальное состояние для основной траектории (для которой проводятся измерение и оценивание), Х0 — для вспомогательных траекторий, по которым приближенно вычисляется оптимальная оценка.

Положить к = 0 = 1 (в случае обрыва траектории с номером г при последующем моделировании = 0 к > 0), г = 1,2,..., М.

III а г 2. Положить

м

Мк = ^ (Мо = М)

г=1

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

Хк = {Х£ }г=1,...,М; ^ = 1:

Хк = ~к £

г=1,...,М; ^ = 1

Проверить условие Т — ^ = 0. Если оно выполнено, то завершить процесс, иначе — положить г = 1, ] =0 {] — количество новых ветвей на шаге к).

Далее получить реализацию оцениваемого вектора состояния в следующем узле сетки {Ьк}:

Хк+1 = Хк + Н!{1к,Хк) + ^а{1к,Хк)Ш,

получить вектор измерений:

гк = фк,хк) + ^-АУ.

В этих формулах и далее ДЖ и Д V — различные для всех к и г реализации случайных векторов размеров в х 1 и й х 1 соответственно, координаты которых независимы и имеют стандартное нормальное распределение.

Шаг 3. Проверить условие ^ = 1. Если оно выполнено, то перейти к следующему шагу, иначе — к последнему шагу.

Шаг 4. Получить реализацию вектора состояния в следующем узле сетки {Ьк}:

К+1 = хгк + Н!{1к,хгк) + ^о(1к,хгк)/хш,

а

(0, 1)

а < |Л(Ьк)|Н,

где Л(Ьк) = Л(Ьк, Х£, Zk), и если оно выполнено, то перейти к следующему шагу, иначе — к последнему шагу.

Шаг 5. Проверить следующие условия:

а) если Л(Ьк) < 0 (Л-(Ьк, Х£, 2к) > 0 обрыв траектории), то положить = 0 (траектория далее не моделируется: Е =0 г > к) и перейти к последнему шагу;

б) если Л(Ьк) > 0 (Л+(Ьк, Хк, 2к) > 0, ветвление траектории), то положить Ек+1 = 1 и j = = j + 1, продублировать реализацию вектора состояния в следующем узле сетки {^} для новой траектории (новой ветви траектории с номером г, которая далее рассматривается как другая траектория с номером М + ]):

Х М +3 = Х г Хк+1 = Хк+Р

Положить Ем+3 = 1 (Егм+3 = 0 г ^ к). 111 а г 6. Проверить следующие условия:

а) если г = М, то положить М = М + Ьк+1 = Ьк + Н, к = к + 1 и перейти к шагу 2;

б) если г < М, то положить г = г + 1 и перейти к шагу 3.

§ 3.2. Алгоритм моделирования системы наблюдения и прогнозирования

МН

Хо Хог Хо

Хог

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

Положить к = 0 Ег = 1 (в случае обрыва траектории с номером г при последующем моделировании Ек = 0 к > 0), г = 1,2, • • •, М.

Шаг 2. Положить к = к (запомнить текущее время), г = 1 j = 0 (] — количество новых ветвей на шаге к), 2 = 1 (прогнозирование).

Ш а г 3. Проверить условие + Д(ЬК) — Ьк > 0. Если оно выполнено, то перейти к шагу 4. Иначе положить

м

Мк = еК (Мо = М)

г=1

и найти прогноз Хк как выборочное среднее реализаций Хк = {Хк}г=1;...;м;=1:

= ^ Е

г=1,...,М; ^Д = 1

Проверить условие Т — = 0. Если оно выполнено, то завершить процесс, иначе — положить к = к, г = 15 = 0 (фильтрация).

Далее, получить реализацию оцениваемого вектора состояния в следующем узле сетки {Ьк}

Хк+1 = Хк + 1г/(и,Хк) + ^Да(и,Хк)А]¥,

получить вектор измерений:

гк = фк,хк) + ^-АУ.

В этих формулах и далее ДШ и ДУ — различные для всех к и г реализации случайных векторов размеров § х 1 и б х 1 соответственно, координаты которых независимы и имеют стандартное нормальное распределение.

Шаг 4. Проверить условие = 1. Если оно выполнено, то перейти к следующему шагу, иначе — к последнему шагу.

Шаг 5. Получить реализацию вектора состояния в следующем узле сетки {Ьк}

Хгк+1 = Х1 + Н!(1к,Х1) + у/На{Ьк,Х{)Ш.

Если 5 = 1, то перейти к последнему шагу, иначе — получить реализацию а случайной

(0, 1)

а < |Л(Ьк)|Н,

где A(ifc) = A(ifc, Zfc), и если оно выполнено, то перейти к следующему шагу, иначе — к последнему шагу.

Шаг 6. Проверить следующие условия:

а) если A(tk) < 0 (A-(tk, Xk, Zk) > 0 обрыв траектории), то положить = 0 (траектория далее не моделируется: F? = 0 r > k) и перейти к последнему шагу;

б) если A(tk) > 0 (A+ (tk,Xk, Zk) > 0 ветвление траектории), то положить Fk+1 = 1 и j = = j + 1, продублировать реализацию вектора состояния в следующем узле сетки {tk} для новой траектории (новой ветви траектории с номером г, которая далее рассматривается как другая траектория с номером M + j):

X M+j = X Xk+1 = Xk+1-

Положить FM+j = 1 (FrM+j = 0 r ^ k). Ill а г 7. Проверить следующие условия:

а) если г = M, то положить M = M + j, tk+1 = tk + h k = k + 1 и перейти к шагу 3 при 5 = 1 или к шагу 2 при 5 = 0;

г < M г = г + 1

Список литературы

1. Казаков II.Е., Артемьев В.М., Бухалев В.А. Анализ систем случайной структуры. М.: Физматлит, 1993. 272 с.

2. Михайлов Г. А., Аверина Т. А. Алгоритм «максимального сечения» в методе Монте-Карло // ДАН. 2009. Т. 428. № 2. С. 163-165.

3. Михайлов Г.А., Рогазинский С.В. Модифицированный метод «мажорантной частоты» для численного моделирования обобщенного экспоненциального распределения // ДАН. 2012. Т. 444. № 1. С. 28-30.

4. Пантелеев А.В., Руденко Е.А., Бортаковский А.С. Нелинейные системы управления: описание, анализ и синтез. М.: Вузовская книга, 2008. 312 с.

5. Рыбаков К.А. Сведение задачи нелинейной фильтрации к задаче анализа стохастических систем с обрывами и ветвлениями траекторий // Дифференциальные уравнения и процессы управления. 2012. № 3. С. 91-110.

6. Рыбаков К.А. Алгоритмы прогнозирования состояний в стохастических дифференциальных системах на основе моделирования специального ветвящегося процесса / / Дифференциальные уравнения и процессы управления. 2015. № 1. С. 25-38.

7. Синицын 11.11. Фильтры Калмана и Пугачева. М.: Логос, 2007. 776 с.

8. Тихонов В.П., Кульман Н.К. Нелинейная фильтрация и квазикогерентный прием сигналов. М.: Советское радио, 1975. 704 с.

9. Bain A., Crisan D. Fundamentals of stochastic filtering. New York: Springer, 2009. XIII+390 p.

10. Candy J.V. Bayesian signal processing: classical, modern and particle filtering methods. New York: John Wiley & Sons, 2009. 472 p.

Поступила в редакцию 23.09.2015

Рыбаков Константин Александрович, к. ф.-м.н., доцент, кафедра математической кибернетики, Московский авиационный институт, 125993, Россия, г. Москва, Волоколамское шоссе, 4. E-mail: [email protected]

K. A. Rybakov

Modified statistical algorithms for filtering and extrapolation in continuous-time stochastic

systems

Keywords: branching process, optimal filtering, extrapolation, statistical modeling, stochastic system. MSC: 60G35, 93E11

To solve the optimal filtering and extrapolation problems for continuous-time stochastic systems, which are described by stochastic differential equations, we propose new algorithms based on modeling the paths of the special random process with terminating and branching. In previous papers the "maximal section algorithm" has been used, but in this paper we use the definition only for inhomogeneous Poisson flows modeling. Therefore, it is not necessary to build a separate grid for each path of the random process with terminating and branching in the numerical solution of stochastic differential equations. The proposed algorithms are easier, they are preferred for real-time optimal filtering and extrapolation.

REFERENCES

1. Kazakov I.E., Artem'ev V.M., Bukhalev V.A. Analiz sistem sluchainoi struktury (Analysis of systems with random structure), Moscow: Fizmatlit, 1993, 272 p.

2. Mikhailov G.A., Averina T.A. The maximal section algorithm in the Monte Carlo method, Doklady Mathematics, 2009, vol. 80, no. 2, pp. 671-673. DOI: 10.1134/S1064562409050111

3. Mikhailov G.A., Rogazinskii S.V. The modified majorant frequency method for numerical simulation of the generalized exponential distribution, Doklady Mathematics, 2012, vol. 85, no. 3, pp. 325-327.

DOI: 10.1134/S1064562412030064

4. Panteleev A.V., Rudenko E.A., Bortakovskii A.S. Nelineinye sistemy upravleniya: opisanie, analiz i sintez (Nonlinear control systems: description, analysis, and synthesis), Moscow: Vuzovskaya kniga, 2008, 312 p.

5. Rybakov K.A. Reducing the nonlinear filtering problem to the analysis of stochastic systems with terminating and branching paths, Differentsial'nye uravneniya i protsessy upravleniya, 2012, no. 3, pp. 91-110 (in Russian).

6. Rybakov K.A. Extrapolation algorithms for stochastic differential systems based on modeling special branching process, Differentsial'nye uravneniya i protsessy upravleniya, 2015, no. 1, pp. 25-38 (in Russian).

7. Sinitsyn I.N. Fil'try Kalmana i Pugacheva (Kalman and Pugachev filters), Moscow: Logos, 2007, 776 p.

8. Tikhonov V.I., Kul'man N.K. Nelineinaya fil'tratsiya i kvazikogerentnyi priem signalov (Nonlinear filtering and quasi-coherent reception of signals), Moscow: Sovetskoe radio, 1975, 704 p.

9. Bain A., Crisan D. Fundamentals of stochastic filtering, New York: Springer, 2009, XIII+390 p.

10. Candy J.V. Bayesian signal processing: classical, modern and particle filtering methods, New York: John Wiley & Sons, 2009, 472 p.

Received 23.09.2015

Rybakov Konstantin Aleksandrovich, Candidate of Physics and Mathematics, Associate Professor, Department of Mathematical Cybernetics, Moscow Aviation Institute (National Research University), Volokolam-skoe shosse, 4, Moscow, 125993, Russia. E-mail: [email protected]