ISSN 2686-9667. Вестник российских университетов. Математика
Том 25, № 130 2020
© Аверина Т.А., Рыбаков К.А., 2020 DOI 10.20310/2686-9667-2020-25-130-109-122 УДК 519.676
Статистические алгоритмы фильтрации для систем со случайно изменяющейся структурой
Татьяна Александровна АВЕРИНА1'2, Константин Александрович РЫБАКОВ3
1 ФГБУН «Институт вычислительной математики и математической геофизики СО РАН» 630090, Российская Федерация, г. Новосибирск, просп. Акад. Лаврентьева, 6 2 ФГАОУ ВО «Новосибирский национальный исследовательский государственный университет» 630090, Российская Федерация, г. Новосибирск, ул. Пирогова, 1 3 ФГБОУ ВО «Московский авиационный институт (национальный исследовательский университет)» 125993, Российская Федерация, г. Москва, Волоколамское шоссе, 4
Statistical filtering algorithms for systems with random structure
Tatyana A. AVERINA12, Konstantin A. RYBAKOV3
1 Institute of Computational Mathematics and Mathematical Geophysics SB RAS 6 Lavrent'eva Ave., Novosibirsk 630090, Russian Federation 2 Novosibirsk State University 1 Pirogova St., Novosibirsk 630090, Russian Federation
3 Moscow Aviation Institute (National Research University)
4 Volokolamskoe Hwy., Moscow 125993, Russian Federation
Аннотация. Предлагаются новые алгоритмы решения задачи оптимальной фильтрации для систем со случайно изменяющейся структурой с непрерывным временем. Эта задача состоит в оценивании текущего состояния системы по результатам косвенных измерений. Математическая модель системы включает нелинейные стохастические дифференциальные уравнения, правая часть которых определяет структуру динамической системы, или режим функционирования. Правая часть может изменяться в случайные моменты времени. Число структур системы предполагается конечным, а процесс смены структуры — марковским или условно марковским. Вектор состояния такой системы состоит из двух компонент: вектора с вещественными координатами и целочисленного номера структуры. Закон изменения номера структуры определяется различными условиями: достижением непрерывной частью вектора состояния заданной поверхности в фазовом пространстве или распределением случайного промежутка времени между переключениями с заданной интенсивностью (средним числом переключений в единицу времени). Каждой упорядоченной паре режимов функционирования может отвечать свой закон перехода между ними.
Алгоритмы оценивания текущего состояния систем со случайно изменяющейся структурой относятся к типу фильтров частиц, они построены на основе метода статистического моделирования (метода Монте-Карло). Работа является продолжением исследований авторов в области статистических методов и алгоритмов анализа и фильтрации для стохастических систем с непрерывным временем.
Ключевые слова: оценивание; фильтрация; метод максимального сечения; стохастическое дифференциальное уравнение; система с переменной структурой; система со случайной структурой; статистическое моделирование; фильтр частиц
Благодарности: Работа выполнена в рамках госзадания ИВМиМГ СО РАН (проект 0315-2019-0002).
Для цитирования: Аверина Т.А., Рыбаков К.А. Статистические алгоритмы фильтрации для систем со случайно изменяющейся структурой // Вестник российских университетов. Математика. 2020. Т. 25. № 130. С. 109-122. DOI 10.20310/2686-9667-2020-25130-109-122.
Abstract. New algorithms for solving the optimal filtering problem for continuous-time systems with a random structure are proposed. This problem is to estimate the current system state vector from observations. The mathematical model of the dynamic system includes nonlinear stochastic differential equations, the right side of which defines the system structure (regime mode). The right side of these stochastic differential equations may be changed at random time moments. The structure switching process is the Markov or conditional Markov random process with a finite set of states (structure numbers). The state vector of such system consists of two components: the real vector (continuous part) and the integer structure number (discrete part). The switch condition for the structure number may be different: the achievement of a given surface by the continuous part of the state vector or the distribution of a random time period between structure switchings. Each ordered pair of structure numbers can correspond to its own switch law.
Algorithms for the estimation of the current state vector for systems with a random structure are particle filters, they are based on the statistical modeling method (Monte Carlo method). This work continues the authors' research in the field of statistical methods and algorithms for the continuous-time stochastic systems analysis and filtering.
Keywords: estimation; filtering; maximum cross section method; stochastic differential equation; system with variable structure; system with random structure; statistical modeling; particle filter
Acknowledgements: The work was carried out within the State assignment ICM&MG SB RAS (project 0315-2019-0002).
For citation: Averina T.A., Rybakov K.A. Statisticheskiye algoritmy fil'tratsii dlya sistem so sluchayno izmenyayushcheysya strukturoy [Statistical filtering algorithms for systems with random structure]. Vestnik rossiyskikh universitetov. Matematika - Russian Universities Reports. Mathematics, 2020, vol. 25, no. 130, pp. 109-122. DOI 10.20310/2686-9667-202025-130-109-122. (In Russian, Abstr. in Engl.)
Введение
В работах [1,2] были предложены статистические алгоритмы оценивания траекторий динамических систем диффузионно-скачкообразного типа, когда модель системы описывается стохастическими дифференциальными уравнениями (СДУ) с винеровской и пуассоновской составляющими. По результатам измерений, модель которых является системой диффузионного типа, т. е. описывается СДУ с винеровской составляющей, оценивалось состояние динамической системы по критерию минимума средне-квадратического отклонения ошибки оценивания в текущий (фильтрация), прошлый (сглаживание) и будущий (прогнозирование) моменты времени. Основу предложенных алгоритмов составляют методы численного решения СДУ и методы моделирования пуассоновских потоков. Для моделирования пуассоновских потоков применялись метод максимального сечения и его модификация, которые позволяют сократить количество обращений к датчику случайных чисел [2,3]. Алгоритмы, предложенные в [1,2], относятся к типу фильтров частиц [4].
В этой работе предлагается применить аналогичный подход к оцениванию траекторий стохастических динамических систем с переключением режимов функционирования (подробнее о таких системах со случайно изменяющейся структурой см. [5]). Модель таких систем описывается СДУ, правая часть которых зависит от марковского или условно-марковского случайного процесса с конечным множеством состояний. В новых алгоритмах фильтрации метод максимального сечения и его модификация применяются при моделировании композиции пуассоновских потоков для изменения режимов функционирования.
Возможные приложения разработанных алгоритмов фильтрации для систем со случайно изменяющейся структурой — задачи слежения за движущимися объектами [5,6].
Решение задачи фильтрации для систем со случайно изменяющейся структурой в самом общем случае сводится к системе СДУ с частными производными, а именно, к системе обобщенных уравнений Стратоновича для условно-апостериорных плотностей [5] (аналог уравнения Стратоновича-Кушнера для систем диффузионного типа [4]). Даже для систем диффузионного типа нахождение апостериорной плотности составляет трудную задачу, которая упрощается в линейном случае и некоторых нелинейных случаях, допускающих построение так называемых конечномерных фильтров. Для систем со случайно изменяющейся структурой нахождение условно-апостериорных плотностей более трудоемко, поэтому были разработаны различные приближенные методы, в основном связанные с аппроксимацией плотностей с помощью частичных сумм функциональных рядов [7,8], моментов и семиинвариантнов [9-11], либо с применением подхода, получившего название условно оптимальной фильтрации [12]. Такие методы относят к субоптимальным, их применение обусловлено необходимостью решения задачи фильтрации в реальном времени.
Альтернативой может служить применение метода статистического моделирования (метода Монте-Карло) и построение алгоритмов оценивания текущего состояния систем со случайно изменяющейся структурой, которые можно отнести к типу фильтров частиц. Их реализация требует использования значительных вычислительных ресурсов [13], однако при этом такие алгоритмы могут обеспечить более точный результат оценивания по сравнению с алгоритмами, основанными на субоптимальных методах. Для сложных многомерных динамических систем реализация алгоритмов фильтрации в реальном времени может оказаться невозможной, но даже в этом контексте статистические алгоритмы фильтрации могут быть использованы при сравнительном анализе с более быстрыми алгоритмами субоптимальной фильтрации для оценки точности и эффективности последних.
В публикациях описано применение фильтров частиц для систем со случайно изменяющейся структурой, однако для упрощенных математических моделей. Например, для динамических систем, изменение режима функционирования которых происходит без учета информации о непрерывной части вектора состояния, либо при рассмотрении стационарных распределений [14,15]. В этой работе рассматривается общий случай.
1. Задача оптимальной фильтрации для систем со случайно изменяющейся структурой
В задаче фильтрации предполагается, что известна математическая модель динамической системы, состояние которой измеряется. В этой работе будем полагать, что математическая модель задается векторным СДУ Ито
¿X (Ь) = / + а(0(£,Х (Ь), X (0) = Хо, (1.1)
где X € Кп - п-мерный вектор, Ь € Т = [0,Т]; /<°(-, ■): Т х Кп ^ Кп, а<^(-, ■): Т х Кп ^ КпХ5 — заданные функции для всех I € Е = {1, 2,..., Ьтах}, I — номер структуры (номер режима функционирования); Ж(Ь) — 5-мерный стандартный винеровский случайный процесс. Закон распределения начального вектора состояния Х0 задан.
Значение I в каждый момент времени — это значение случайного процесса смены структуры Ь(Ь) [5]. Как правило, случайный процесс Ь(Ь) можно задать с помощью совокупности уравнений Б1г (Ь, х) = 0, Б1г (■, ■): Т х Кп ^ К; /, г =1, 2,..., Ьтах ; I = г. Тогда при условии Ь(Ь — 0) = I и достижении вектором X поверхности (Ь,х) = 0 в момент времени Ь случайный процесс Ь(Ь) принимает значение г. Такой тип смены структуры называется сосредоточенным переходом, а система (1.1) — системой с переменной структурой.
Другой вариант — задание функций Лгг(■, ■): Т х Кп ^ К+, определяющих интенсивности смены структуры; /, г = 1, 2,..., Ьтах ; I = г. Тогда при малых ДЬ
Р(Ь(Ь + ДЬ) = г | Ь(Ь) = I, X(Ь) = х) = Л|г(Ь,х)ДЬ + о(ДЬ),
^шах
Р(Ь(Ь + ДЬ) = I | Ь(Ь) = I, X(Ь) = х) = 1 — ^ Лгг(Ь,х)ДЬ + о(ДЬ),
г=1,г=1
где Р обозначает вероятность. Такой тип смены структуры называется распределенным переходом, а система (1.1) — системой со случайной структурой.
И при сосредоточенных, и при распределенных переходах смена структуры системы, т. е. изменение правой части уравнения (1.1) происходит случайным образом. В обоих случаях будем называть такие динамические системы системами со случайно изменяющейся структурой. В общем случае можно рассматривать системы, для которых часть переходов сосредоточенная, а часть — распределенная.
Начальное условие для процесса смены структуры Ь(Ь) задается в виде ¿(0) = Ь0 € Е. Траектории Ь(Ь) — это кусочно-постоянные функции со значениями из множества Е ; X0, Ь0 и Ж(Ь) независимы.
Вектор состояния X системы со случайно изменяющейся структурой состоит из двух компонент: X — непрерывная и Ь — дискретная, т. е. информация о текущем состоянии задается в виде X = [Xт, Ь]т.
Измерения будем описывать векторным СДУ
йУ (Ь) = с<^ (Ь, X+ С<° (Ь)^(Ь), У (0) = У = 0, (1.2)
в котором У € Кт — т-мерный вектор; с<г>(■, ■): Т х Кп ^ Кт, (<г>(-): Т ^ КтХй — заданные функции для всех I € Е; V (Ь) — й -мерный стандартный винеровский
случайный процесс, не зависящий от винеровского процесса W(t). При этом предполагается, что матрица n^(t) = Z^(t)[Z^(t)]T невырождена для всех l Е L и t Е T.
Отметим, что в наиболее общем случае для систем со случайно изменяющейся структурой вводят уравнение индикатора структуры [12], либо задаются соответствующие переходные плотности [16,17], тогда вектор измерений, как и вектор состояния, состоит их двух компонент: непрерывной и дискретной. Однако в этой работе ограничимся описанием измерительной системы только уравнением (1.2).
Будем предполагать, что решение СДУ (1.1), (1.2) существует и единственно. Для этого достаточно выполнения условий на функции f ^(t,x), а® (t,x), c^(t,x), Z^(t) и Alr(t,x), которые приведены в [18]. При сосредоточенных переходах Slr(t,x) —гладкие функции на T х Rn (за исключением подмножества нулевой меры). Кроме того, E|X0|2 < то, где E обозначает математическое ожидание. Также предполагается, что существует плотность для вектора X = X(t) в каждый момент времени t Е T.
Задача оптимальной фильтрации для системы со случайно изменяющейся структурой (1.1), (1.2) заключается в нахождении оценки X вектора состояния X по доступным к текущему моменту времени t измерениям YQ = {Y(т), т Е [0,t]} в виде X = Ф(^), т. е. требуется определить функцию Ф^, ■), преобразующую измерения Yq в оценку вектора состояния в текущий момент времени t Е T.
Функцию Ф(^ ■) можно найти из условия оптимальности оценки (синтезировать оптимальный фильтр), т. е. исходя из заданного критерия качества, обеспечивающего минимум среднего риска [12,19]:
En(X,Х) ^ min,
где минимум берется по всем допустимым функциям Ф^, ■), n(x,x) — функция потерь, x = [xT,l]T.
Для систем со случайно изменяющейся структурой используют следующие функции потерь n(x,x) :
1) аддитивную квадратично-простую функцию потерь
n(x, x) = (x — x)TC(x — x) + (1 — 5l ¡),
где C — симметричная положительно определенная матрица n х n, i — символ Кронекера;
2) мультипликативную простую функцию потерь
n(x, x) = 1 — öl i i(x — x),
где £(•) — дельта-функция.
Функции потерь часто рассматривают как функции ошибки оценивания [5], т. е. разности е = x — x. В этом случае они будут выражаться через 50 — = Sl г\ Другие варианты функций потерь описаны в [19].
Для того чтобы записать решение задачи оптимальной фильтрации, требуется ввести обозначения для апостериорной плотности вектора X и вероятностей, задающих распределение номера структуры L. Обозначим через <^(t, x,l | Y?), ■, -| Y0t): T х Rn х L ^ R+, функцию, которая описывает апостериорное распределение вектора
состояния X:
1 у? =1.
1=1 ^
При фиксированном I Е Ь — это ненормированная условно-апостериорная плотность ^(1)(£,ж | У?) вектора X, а функция
Р(1)(* | У?) = | У?)*г
'В,"
Г1-3)
(1)
определяет апостериорную вероятность события ¿(¿) = / (вероятности Р(1)(0) = Р I = 1, 2,... , Ьтах, определяют распределение начального состояния Ь0 ). Далее, функ-
ция
| У?)
= у(1) (¿,х | У?) 0) Р(0(* | у?)
:1.4)
определяет условно-апостериорную плотность вектора X при ¿(¿) = /, а функция
|у?)= £ |у?)
:1.5)
1=1
является нормированной апостериорной плотностью вектора X.
Согласно [12,19] оптимальная оценка Х(£) = [XТ(£),Ь(£)]Т для аддитивной квадратично-простой функции потерь определяется соотношениями
*(*) = £ Р(1)(^IУ?^(1)(*), ад = аг^ Р(1)(^|У?)
1=1
где
XX^(^Е^(*) |£(*) = /,У?]=/ |У^ж, I = 1, 2,
а для мультипликативной простой функции потерь —
Ьт
XX(¿) = XXОД = а^шах <^(1) (¿Д(1)(£) |У?),
1еь
где
XX (1)(£) = а^шах <£(1)(;£,ж | У?), / = 1, 2,...,£п
2. Статистические алгоритмы фильтрации
Для приближенного решения задачи оптимальной фильтрации применим метод статистического моделирования. В его основе лежит построение траекторий системы со случайно изменяющейся структурой по известной математической модели, включающей уравнение (1.1) для случайного процесса X (¿) и определение процесса смены структуры Ь(£). По ансамблю траекторий (реализаций вектора состояния XX системы в различные моменты времени, которые, как правило, образуют узлы разбиения отрезка Т) с помощью статистических методов можно найти оценки X, соответствующие
разным функциям потерь. При статистической обработке результатов моделирования могут быть найдены оценки плотности вектора X, оценки условных и безусловных мо-ментных характеристик вектора X, оценка распределения для номера структуры Ь. Все это не требует вывода сложных соотношений, а базируется только на методах моделирования решений СДУ и марковских процессов с конечным множеством состояний.
Для учета имеющихся измерений У0* предлагается применить метод частиц, в котором наряду с траекториями динамической системы (1.1) и процессом смены структуры требуется моделирование траекторий весовой функции по определенному правилу. Тогда вклад каждой траектории системы со случайно изменяющейся структурой в искомую оценку вектора состояния зависит от соответствующей траектории весовой функции, которая формируется на основе модели измерений (1.2) (вес уменьшается, если траектория далека от оптимальной оценки; вес растет, если траектория близка к оптимальной оценке). Для стохастических динамических систем с непрерывным и дискретным временем определения соответствующих весовых функций отличаются [4,20]. Однако применение численных методов решения СДУ означает фактический переход от непрерывного времени к дискретному, что позволяет использовать более простое определение весовой функции и не приводить уравнения для ненормированных и нормированных условно-апостериорных плотностей вектора X [5].
Воспользуемся стохастическим методом Эйлера для построения дискретной аппроксимации решения СДУ (1.1), (1.2) при фиксированном I £ Е:
Хк+1 = Хк + Ь/,Хк) + ^а« ,Хк , (2.6)
Ук+1 = Ук + ,Хк) + ^ЬС<° (4)ДП, (2.7)
где ДШк и ДVk — в -мерный и й -мерный случайные векторы, координаты которых независимы в совокупности и имеют стандартное нормальное распределение для всех к, Ь > 0 — шаг численного интегрирования, т. е.
4+1 = tk + Ь, к = 0,1,...,Х - 1, ¿0 = 0, = Т.
Случайный вектор (Ук+1 — Ук)/Ь — с< '(¿к, Хк) имеет гауссовское распределение с нулевым математическим ожиданием и ковариационной матрицей )/Ь, т. е. с плотностью
Ь™/2 С 1
Рк ( , ) (2п)™/2^ ^>(4) Ч 2 [/ (к)] I
Эта плотность выполняет роль функции правдоподобия [4] и позволяет записать правило обновления весовых коэффициентов
^к+1 = ^к(Ук+1Ь— Ук — с<1>(1к, Хк), ь) . (2.8)
Начальное условие для весового коэффициента должно быть положительным, его конкретная величина не имеет принципиального значения, например ш0 = 1. Эти весовые коэффициенты позволяют найти апостериорные моментные характеристики вектора X. В частности, апостериорное математическое ожидание вектора X в момент
времени ¿к оценивается следующим образом (оно дает приближенное решение задачи оптимальной фильтрации при выборе аддитивной квадратично-простой функции потерь):
1 м м
Хк = —Е Хк, -к = ^ ,
к г=1 г=1
где г — номер моделируемой траектории случайного процесса Х (¿), М — число траекторий.
Аналогичным образом можно оценить апостериорные начальные и центральные мо-ментные характеристики. Например, для моментов второго порядка имеем 1 м 1 м
#к = ^ Е^кХк(Хк)т, д = ?г £<4(Хк - Хк)(хк - хХк)т,
к г=1 к г=1
и оценка Дк ковариационной матрицы вектора X в момент времени ¿к может использоваться как характеристика точности фильтрации.
Отметим, что в качестве численного метода решения СДУ (1.1) можно использовать любой подходящий метод, обеспечивающий более высокую точность по сравнению со стохастическим методом Эйлера [3,21]. Численное решение уравнения (1.2) требуется только для моделирования работы фильтра. В реальной задаче фильтрации моделирование случайного процесса У (¿) не проводится.
Для систем со случайно изменяющейся структурой наряду с моделированием траекторий случайного процесса X (¿) требуется моделирование траекторий случайного процесса ¿(¿) в соответствии с его определением [3]. При этом все приведенные выше соотношения остаются справедливыми, при их применении нужно положить I = ¿к = ¿(¿к).
При сосредоточенных переходах условие смены номера структуры с I на г — это равенство (¿,Х(¿)) = 0. С учетом дискретной аппроксимации случайного процесса X(¿) это условие можно записать как (¿к, Хк) (¿к+1,Хк+1) ^ 0 [22].
При распределенных переходах можно воспользоваться условием ординарности композиции пуасоновских потоков или построить алгоритмы фильтрации на основе метода максимального сечения и его модификации. Не будем отдельно приводить соотношения для них, поскольку эти методы в приложении к задаче фильтрации для динамических систем диффузионно-скачкообразного типа описаны в [2].
Приведем сразу алгоритмы приближенного решения задачи оптимальной фильтрации, построенные на основе алгоритмов моделирования траекторий систем со случайной структурой и методов моделирования пуассоновских потоков, выбирая аддитивную квадратично-простую функцию потерь.
Алгоритм моделирования пуассоновского потока на основе определения его интенсивности применялся в [21]. В [23,24] для моделирования пуассоновского потока применялись метод максимального сечения и его модификация. Прежде чем привести эти алгоритмы, введем необходимые обозначения.
Интенсивность композиции пуассоновских потоков смены структуры для текущего номера I задается выражением
т
Аг(¿,х)= Е А|г(¿,х).
г=1,г=1
Пусть также числа А* (/, г =1, 2,... , ¿тах ; I = г) удовлетворяют условиям:
Агг (¿,х) ^ А*, /,г = 1, 2, ...,Ьтах, I = г.
Обозначим
т
А* = X/ А*г.
г=1,г=1
Алгоритм 1 (для систем с сосредоточенными переходами)
1. Задать М —число моделируемых траекторий случайных процессов X (¿), ¿(¿) и т. е. число частиц; Ь — шаг численного интегрирования. Получить реализации векторов X0 и номеров ¿0 согласно заданным распределениям и положить к = 0, 4 = 1, * = 1, 2,..., М.
2. Оценить вектор состояния и ковариационную матрицу вектора Xk, а также номер
^к ,Ьк ,4}
структуры Ьк по реализациям ^к, ¿к,4}м=1:
1 м 1 м
^^к = ^ Е4Xk, 4 = ?Г £ 4^к — xtk)(Xk — ^^к)-, ¿к = а^тах Д<°,
к г=1 к г=1
где
О = 4 Р<1> = #{*: ¿к = I, * = 1, 2,... , М}
Ок = 2^ 4, п =-м-.
г=1
Проверить условие Т — ¿к = 0. Если оно выполнено, то завершить процесс. Положить * =1.
3. Получить реализацию вектора X в следующем узле сетки Ь = ¿к + Ь:
Xk+1 = Xk + Ь/<г>(Ьк ^к) + ^>(4 ,Xk )Д^к,
где I = ¿к, Д^к — реализация в -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Обновить весовой коэффициент
4 = 4 У (Ьк+1) — У (¿к) т(, X ) Л
4к+1 = 4кРк ( -Ь--с< (Ьк),ЬЬ
Проверить условие смены структуры: если вектор Xk+l достиг поверхности (¿, х) = 0 для некоторого номера г к моменту времени ¿к+1, т. е. (¿к,Xk) х х51г(¿к+1 ,Xk+1) ^ 0, где I = ¿к, то изменить номер структуры на г (г € Е\{/} ), т. е. положить ¿к+1 = г. Иначе положить ¿к+1 = I.
4. Если * = М, то положить ¿к+1 = ¿к + Ь, к = к + 1 и перейти к шагу 2; если * < М, то положить * = * + 1 и перейти к шагу 3.
Алгоритм 2 (для систем с распределенными переходами на основе метода максимального сечения)
1. Задать М —число моделируемых траекторий случайных процессов Х (¿), ¿(¿) и
т. е. число частиц; к — шаг численного интегрирования. Получить реализации векторов Х,5 и номеров согласно заданным распределениям и положить к = 0, = 0, = 1, г = 1, 2,... , М. Получить реализации случайных величин £г, имеющих показательное распределение с параметром А*, где I = , г = 1, 2,..., М.
2. Выполнить шаг 2 в алгоритме 1.
3. Положить I = ¿к и проверить условие + £г ^ ¿к + к. Если оно выполнено, то перейти к шагу 4, иначе положить т = ¿к, X = Хк и перейти к шагу 5.
4. Получить реализацию вектора Х в следующем узле сетки Ь = ¿к + к:
Хк+1 = Хк + к/«(¿к ,Хк) + ^«(¿к ,Хк )ДЖк,
где ДЖк — реализация 5 -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Положить ¿к+1 = I, и перейти к шагу 7.
5. Получить реализацию вектора Х в дополнительном узле сетки £ = + £г:
Х = Х + к/(т, Х + (т, Х)Д^,
где к = ¿* + £г — т, ДЖ — реализация 5 -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Положить т = т + к и получить реализацию номера структуры г с вероятностью Рг = А*/А* (г € Е\{/|). Проверить условие смены структуры: если а ^ А1г (т, Х)/А*, где а — реализация случайной величины, имеющей равномерное распределение на интервале (0,1), то изменить номер структуры на г, т. е. положить I = г. Далее положить = + £г и получить новую реализацию случайной величины £г, имеющей показательное распределение с параметром А*. Проверить условие + £ ^ ¿к + к. Если оно выполнено, то перейти к шагу 6, иначе перейти к шагу 5.
6. Получить реализацию вектора Х в следующем узле сетки ^ = ¿к + к:
Хк+1 = X + к /(0 (т, Х + (т, Х)Д Ж,
где к = ¿к + к — т, Д Ж — реализация в -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Положить ¿1+1 = 1.
7. Обновить весовой коэффициент
, = , Р<0 У(¿к+О — У(¿к) /0(, Х) Л
^к+1 = ^к Рк ( -к--^ ,Хк ),кЬ
где I = ¿к.
8. Выполнить шаг 4 в алгоритме 1.
Алгоритм 3 (для систем с распределенными переходами на основе модифицированного метода максимального сечения)
1. Задать М —число моделируемых траекторий случайных процессов X (¿), ¿(¿) и
т. е. число частиц; Ь — шаг численного интегрирования. Получить реализации векторов X0 и номеров ¿0 согласно заданным распределениям и положить к = 0, ^ = 0, 4) = 1, г = 1, 2,..., М. Получить реализации случайных величин имеющих показательное распределение с параметром А*, где I = ¿0, и положить П = Ф' = 1, * =1, 2,... , М.
2. Выполнить шаг 2 в алгоритме 1.
3. Положить I = ¿к и проверить условие + ^ ¿к + Ь. Если оно выполнено, то перейти к шагу 4, иначе положить т = ¿к, X = Xk. Если Ф' =1, то получить реализацию случайной величины а, имеющей равномерное распределение на интервале (0,1), и положить Ф' = 0. Перейти к шагу 5.
4. Выполнить шаг 4 в алгоритме 2.
5. Получить реализацию вектора X в дополнительном узле сетки £ = + :
X = X + Ь/<° (т, X) + (т, Х)Д ^,
где Ь = ¿* + £ — т и ДЖ — реализация в -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Положить т = т + Ь и получить реализацию номера структуры г с вероятностью Рг = А*/А* ( г € Е\{/} ). Положить П = П'(1 — А1г(т, X)/А*г) и проверить условие смены структуры: если 1 — а > П', то изменить номер структуры на г, т. е. положить I := г, и положить П' = Ф' =1. Далее положить = + £г и получить новую реализацию случайной величины имеющей показательное распределение с параметром А*. Проверить условие ¿* + £ ^ ¿к + Ь. Если оно выполнено, то перейти к шагу 6, иначе перейти к шагу 5.
6. Получить реализацию вектора X в следующем узле сетки ^ = ¿к + Ь:
Xk+1 = X + Ь/<° (т, X) + (т, X )Д Ж,
где Ь = ¿к+Ь—т и ДЖ —реализация в -мерного случайного вектора, координаты которого независимы и имеют стандартное нормальное распределение. Положить
¿1+1 = 1.
7. Выполнить шаг 7 в алгоритме 2.
8. Выполнить шаг 4 в алгоритме 1.
Заключение
В работе предложены новые алгоритмы решения задачи оптимальной фильтрации для систем со случайно изменяющейся структурой с непрерывным временем. В настоящее время с использованием оптимизирующих компиляторов, которые позволяют наиболее эффективно задействовать вычислительные ресурсы процессоров, можно реализовать подобные алгоритмы для сложных динамических систем, обеспечивая при этом решение задачи оптимальной фильтрации в реальном времени. Для систем диффузионного типа возможность реализации фильтров частиц в реальном времени продемонстрирована в [13]. При необходимости реализации более быстрых и не требующих обработки большого объема данных алгоритмов субоптимальной фильтрации алгоритмы, построенные на основе метода статистического моделирования, могут применяться для сравнения.
References
[1] Т. А. Аверина, К. А. Рыбаков, Применение метода максимального сечения в задачах оценивания случайных процессов, Актуальные проблемы вычислительной и прикладной математики (АПВПМ-2019, Марчуковские научные чтения), НГУ, Новосибирск, 2019. [T. A. Averina, K. A. Rybakov, Using the maximum cross section method in random processes estimation problems, Int. Conf. on Advanced Mathematics, Computations and Applications (AMCA-2019, Marchuk Scientific Readings), NSU Publ., Novosibirsk, 2019 (In Russian)].
[2] T. A. Averina, K. A. Rybakov, "Using maximum cross section method for filtering jumpdiffusion random processes", Rus. J. Numer. Anal. Math. Modelling, 35:2 (2020), 55-67.
[3] Т. А. Аверина, Статистическое моделирование решений стохастических дифференциальных уравнений и систем со случайной структурой, Изд-во СО РАН, Новосибирск, 2019. [T. A. Averina, Statisticheskoe modelirovanie reshenii stokhasticheskikh differentsial'nykh uravnenii i sistem so sluchainoi strukturoi, Siberian Branch of the Russian Academy of Sciences Publ., Novosibirsk, 2019 (In Russian)].
[4] A. Bain, D. Crisan, Fundamentals of Stochastic Filtering, Springer, New York, 2009.
[5] И. Е. Казаков, В. М Артемьев, Оптимизация динамических систем случайной структуры, Наука, М., 1980. [I.E. Kazakov, V.M Artem'ev, Optimizatsiya dinamicheskikh sistem sluchainoi struktury, Nauka, Moscow, 1980 (In Russian)].
[6] X. R. Li, V. P. Jilkov, "Survey of maneuvering target tracking. Part V: Multiple-model methods", IEEE Trans. Aerospace Electronic Syst., 41:4 (2005), 1255-1321.
[7] М. С. Ярлыков, С. М. Ярлыкова, "Оптимальные алгоритмы комплексной нелинейной обработки векторных дискретно-непрерывных сигналов", Радиотехника, 2004, №7, 18-29. M. S. Yarlykov, S. M. Yarlykova, "Optimal algorithms for complex nonlinear processing of vector discrete-continuous signals", Radioengineering, 2004, №7, 18-29 (In Russian)].
[8] М. С. Ярлыков, С. М. Ярлыкова, "Потенциальная точность синтезированных субоптимальных алгоритмов комплексной нелинейной обработки векторных дискретно-непрерывных шумоподобных сигналов", Радиотехника, 2007, №1, 46-61. [M.S. Yarlykov, S. M. Yarlykova, "The potential accuracy of synthesized suboptimal algorithms for complex nonlinear processing of vector discrete-continuous noises", Radioengineering, 2007, № 1, 46-61 (In Russian)].
[9] И. М. Косачев, Ю.Е. Кулешов, "Методология высокоточной оптимальной фильтрации случайных процессов, наблюдаемых в стохастических динамических системах со случайно изменяющейся структурой. Часть 1", Вестник Военной академии Республики Беларусь, 2016, №3(52), 57-66. [I.M. Kosachev, Yu. E. Kuleshov, "Methodology of high-precision optimal filtering of random processes observed in stochastic dynamical systems with random structure. Part 1", Military Academy of the Republic of Belarus Proceedings, 2016, № 3(52), 57-66 (In Russian)].
[10] И. М. Косачев, Ю.Е. Кулешов, "Методология высокоточной оптимальной фильтрации случайных процессов, наблюдаемых в стохастических динамических системах со случайно изменяющейся структурой. Часть 2", Вестник Военной академии Республики Беларусь, 2016, №4(53), 64-73. [I.M. Kosachev, Yu. E. Kuleshov, "Methodology of high-precision optimal filtering of random processes observed in stochastic dynamical systems with random structure. Part 2", Military Academy of the Republic of Belarus Proceedings, 2016, № 4(53), 64-73 (In Russian)].
[11] И. М. Косачев, Ю.Е. Кулешов, "Методология высокоточной оптимальной фильтрации случайных процессов, наблюдаемых в стохастических динамических системах со случайно изменяющейся структурой. Часть 3", Вестник Военной академии Республики Беларусь, 2017, №1 (54), 56-65. [I.M. Kosachev, Yu.E. Kuleshov, "Methodology of high-precision optimal filtering of random processes observed in stochastic dynamical systems with random structure. Part 3", Military Academy of the Republic of Belarus Proceedings, 2017, № 1(54), 56-65 (In Russian)].
[12] Е. А. Руденко, "Конечномерные рекуррентные алгоритмы оптимальной нелинейной логико-динамической фильтрации", Известия РАН. Теория и системы управления, 2016, №1, 43-65; англ. пер.: E. A. Rudenko, "Finite-dimensional recurrent algorithms for optimal nonlinear logical-dynamical filtering", J. Comput. Sys. Sc. Int., 55:1 (2016), 36-58.
[13] К. А. Рыбаков, А. А. Ющенко, "Непрерывные фильтры частиц и их реализация в реальном масштабе времени", Вестник ВГУ. Серия Системный анализ и информационные технологии, 2018, №3, 56-64. [K. A. Rybakov, A. A. Yushchenko, "Continuous particle filters and its real-time implementation", Proceedings of Voronezh State University. Series: Systems Analysis and Information Technologies, 2018, №3, 56-64 (In Russian)].
[14] F. Karame, "A new particle filtering approach to estimate stochastic volatility models with Markov-switching", Econometrics and Statistics, 8:C (2018), 204-230.
[15] T. Lux, "Inference for nonlinear state space models: A comparison of different methods applied to Markov-switching multifractal models", Economics Working Paper, 2018-07 (2018), 1-43.
[16] В.А. Бухалев, Распознавание, оценивание и управление в системах со случайной скачкообразной структурой, Физматлит, М., 1996. [V. A. Bukhalev, Raspoznavanie, otsenivanie i upravlenie v sistemakh so sluchainoi skachkoobraznoi strukturoi, Fizmatlit, Moscow, 1996 (In Russian)].
[17] В.А. Бухалев, А.А. Скрынников, В.А. Болдинов, Алгоритмическая помехозащи-та беспилотных летательных аппаратов, Физматлит, М., 2018. [V. A. Bukhalev, A. A. Skrynnikov, V. A. Boldinov, Algoritmicheskaya pomekhozashchita bespilotnykh letatel'nykh apparatov, Fizmatlit, Moscow, 2018 (In Russian)].
[18] M. Ghosh, A. Arapostathis, S. Marcus, "Optimal control of switching diffusions with application to flexible manufacturing systems", SIAM J. Control Optim., 31:5 (1993), 11831204.
[19] А. П. Трифонов, Ю. С. Шинаков, Совместное различение сигналов и оценка их параметров на фоне помех, Радио и связь, М., 1986. [A. P. Trifonov, Yu. S. Shinakov, Sovmestnoe razlichenie signalov i otsenka ikh parametrov na fone pomekh, Radio i svyaz', Moscow, 1986 (In Russian)].
[20] К. А. Рыбаков, "О вычислении весовых коэффициентов в непрерывном фильтре частиц", Научный вестник МГТУ ГА, 21:2 (2018), 32-39. [K.A. Rybakov, "Calculation of weight coefficients in continuous particle filter", Civil Aviation High Technologies, 21:2 (2018), 32-39 (In Russian)].
[21] Н. В. Черных, П. В. Пакшин, "Алгоритмы численного решения стохастических дифференциальных систем с переключаемой диффузией", Управление большими системами, 2012, №36, 106-143; англ. пер.: N. V. Chernykh, P. V. Pakshin, "Numerical solution algorithms for stochastic differential systems with switching diffusion", Autom. Remote Control, 74:12 (2013), 2037-2063.
[22] Т. А. Аверина, К. А. Рыбаков, "Моделирование мультиструктурных систем на многообразиях в задачах статистического анализа и фильтрации", Вестник Тамбовского университета. Серия: Естественные и технические науки, 23:122 (2018), 145-153. [T. A. Averina, K. A. Rybakov, "Modeling of multistructural systems on manifolds for statistical analysis and filtering problems", Tambov University Reports. Series: Natural and Technical Sciences, 23:122 (2018), 145-153 (In Russian)].
[23] Т. А. Аверина, "Построение и обоснование статистических алгоритмов моделирования решения систем со случайной структурой, заданной стохастическими дифференциальными уравнениями", Вестник Тамбовского университета. Серия: Естественные и технические науки, 20:5 (2015), 986-991. [T. A. Averina, "Construction and justification of statistical algorithms for simulation of switching diffusion, given by stochastic differential equations", Tambov University Reports. Series: Natural and Technical Sciences, 20:5 (2015), 986-991 (In Russian)].
[24] Т. А. Аверина, "Использование модификаций метода максимального сечения для моделирования систем со случайной структурой с распределенными переходами", Сибирский журнал вычислительной математики, 19:3 (2016), 235-247; англ. пер.: T.A. Averina, "A randomized maximum cross-section method to simulate random structure systems with distributed transitions", Numer. Anal. Appl., 9:3 (2016), 179-190.
Информация об авторах
Аверина Татьяна Александровна, кандидат физико-математических наук, старший научный сотрудник лаборатории численного анализа СДУ, Институт вычислительной математики и математической геофизики СО РАН; доцент механико-математического факультета, Новосибирский национальный исследовательский государственный университет, г. Новосибирск, Российская Федерация. E-mail: ata@osmf.sscc.ru ORCID: https://orcid.org/0000-0002-2206-289X
Рыбаков Константин Александрович,
кандидат физико-математических наук, доцент кафедры математической кибернетики. Московский авиационный институт (национальный исследовательский университет), г. Москва, Российская Федерация. E-mail: rkoffice@mail.ru ORCID: https://orcid.org/0000-0002-6839-1090
Конфликт интересов отсутствует. Для контактов:
Константин Александрович Рыбаков E-mail: rkoffice@mail.ru
Поступила в редакцию 25 марта 2020 г. Поступила после рецензирования 14 мая 2020 г. Принята к публикации 8 июня 2020 г.
Information about the authors
Tatyana A. Averina, Candidate of Physics and Mathematics, Senior Researcher, Laboratory of SDE Numerical Analysis, Institute of Computational Mathematics and Mathematical Geophysics SB RAS; Associate Professor of the Mechanics and Mathematics Department, Novosibirsk State University, Novosibirsk, Russian Federation. E-mail: ata@osmf.sscc.ru ORCID: https://orcid.org/0000-0002-2206-289X
Konstantin A. Rybakov, Candidate of Physics and Mathematics, Associate Professor of the Mathematical Cybernetics Department. Moscow Aviation Institute (National Research University), Moscow, Russian Federation. E-mail: rkoffice@mail.ru
ORCID: https://orcid.org/0000-0002-6839-1090
There is no conflict of interests.
Corresponding author: Konstantin A. Rybakov E-mail: rkoffice@mail.ru
Received 25 March 2020 Reviewed 14 May 2020 Accepted for press 8 June 2020