Известия РАН. Теория и системы управления, 2020, № 6, стр. 42-59

Введение. При оптимальном управлении логико-динамическими (многорежимными, переключательными) стохастическими системами по доступным наблюдениям их выходов возникает задача оперативного распознавания режима работы и оценивания внутрирежимного состояния этих систем. Однако применение для этого бесконечномерного абсолютно-оптимального фильтра (АОФ) [1–4] сопряжено с существенными вычислительными трудностями, а конечномерный условно-минимаксный фильтр [5] оптимален лишь параметрически и имеет небольшой фиксированный порядок (размерность вектора состояния), что препятствует достижению им значительной точности. В то же время потенциальная точность лишенных этих недостатков непараметрических конечномерных фильтров оптимальной структуры (ФОС) также ограничена, ибо они имеют или бесконечное время памяти и малый порядок, равный порядку объекта наблюдения [6, 7], или любой порядок, кратный размерности вектора измерений, но конечную память [8]. Программные реализации этих фильтров описаны в [9].

В настоящей работе предлагается способ преодоления этой ограниченности логико-динамических ФОС путем построения нового фильтра большого порядка (ФБП). Его вектор состояния предлагается формировать из нескольких векторов последних оценок, поэтому порядок ФБП кратен размерности оцениваемого вектора состояния системы наблюдения. Время памяти ФБП бесконечно, а противоречие между его точностью и сложностью регулируется выбором величины кратности. Нахождение структурных функций оптимального ФБП выполняется до начала процесса оценивания с помощью интегральных операций с условными распределениями вероятности, что можно реализовать численно методом Монте-Карло с построением выборочных приближений к искомым функциям обработкой их многомерных гистограмм.

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

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

1. Постановка задачи синтеза фильтра. Рассмотрим функционирующий в условиях неопределенности многорежимный объект. Пусть его поведение в каждый момент дискретного времени $k{\text{ = }}\overline {{\text{0,}}\,K} $ с конечным горизонтом K определяется одним из L случайных номеров его режима (структуры) ${{I}_{k}} \in \{ \overline {1,L} \} $ и случайным абсолютно непрерывным вектором его состояния ${{X}_{k}} \in {{\mathbb{R}}^{{{{n}_{x}}}}}$. При этом изменение во времени расширенного вектора состояния объекта ${{\bar {X}}_{k}} = {{[{{I}_{k}}\,\,X_{k}^{{\text{Т}}}]}^{{\text{Т}}}}$ описывается следующей системой разностных уравнений:

Здесь $\phi _{k}^{{(i)}}(x,{v}_{{}}^{\Delta }) \equiv \phi _{k}^{{}}(i,x,{v}_{{}}^{\Delta })$ – дискретнозначная скалярная функция логического блока (переключателя) объекта, $a_{k}^{{(i)}}(x,{v}) \equiv a_{k}^{{}}(i,x,{v})$ – непрерывнозначная вектор-функция его динамического блока, ${{\bar {V}}_{k}} = (V_{k}^{\vartriangle },{{V}_{k}})$ – вектор дискретных белых возмущений (порождающих шумов) этих блоков. Эти возмущения не зависят от начального состояния объекта $({{I}_{0}},\,{{X}_{0}})$, который, как и другие системы из дискретных и абсолютно непрерывных случайных величин, будем определять совместным вероятностным распределением

в виде произведения вероятности $P_{0}^{{(i)}}$ номера режима I₀ = i на условную плотность вероятности $p_{0}^{{(i)}}(x)$ вектора состояния X₀ при режиме i:

Пусть также в каждый момент времени k для уточнения состояния ${{\bar {X}}_{k}}$ системы (1.1) фиксируется случайный вектор расширенных наблюдений ${{\bar {Y}}_{k}} = {{[{{J}_{k}}\,\,Y_{k}^{{\text{Т}}}]}^{{\text{Т}}}}$, состоящий из дискретной М-значной переменной индикации режима ${{J}_{k}} \in \{ \overline {1,M} \} $ и абсолютно непрерывного вектора ${{Y}_{k}} \in {{\mathbb{R}}^{{{{n}_{y}}}}}$ измерений, причем их значения определяются по формулам

Здесь $o_{k}^{{(i)}}(x,{{w}^{\Delta }})$, $b_{k}^{{(i)}}(x,w)$ – известные функции логического и динамического блоков измерителя, ${{\bar {W}}_{k}} = (W_{k}^{\Delta },{{W}_{k}})$ – вектор дискретных белых погрешностей (шумов) измерений, не зависящий ни от ${{\bar {X}}_{0}}$, ни от ${{\bar {V}}_{k}}$.

Тогда априорное вероятностное поведение системы наблюдения (1.1), (1.2) полностью описывается двумя условными распределениями [7, 8], определяемыми по ее уравнениям и законам распределения шумов:

– распределением перехода объекта от старого состояния ${{\bar {X}}_{k}}$ к новому ${{\bar {X}}_{{k + 1}}}$:

– условным распределением измерения ${{\bar {Y}}_{k}}$ при ${{I}_{k}} = {{i}_{k}},\,\,{{X}_{k}} = {{x}_{k}}$ (функцией правдоподобия):

Здесь ${\text{M}}$ – оператор математического ожидания, выполняющий усреднение по случайным величинам, которые отмечены прописными буквами, ${{\delta }_{{i,\,l}}}$ – символ Кронекера, δ(x) – функция Дирака.

Замечание 1. В случае зависимости шумов ${{\bar {V}}_{k}},\,\,{{\bar {W}}_{k}}$ приходится оперировать с более общим условным распределением ${{\gamma }_{k}}({{i}_{{k + 1}}},{{x}_{{k + 1}}},{{j}_{k}},{{y}_{k}}\,{\text{|}}\,{{i}_{k}},{{x}_{k}})$ перехода от старого состояния ${{\bar {X}}_{k}}$ не только к новому состоянию ${{\bar {X}}_{{k + 1}}}$, но и к текущему наблюдению ${{\bar {Y}}_{k}}$ [7, 8]. Тогда условное распределение ${{\beta }_{k}}( \cdot )$ наблюдения ${{\bar {Y}}_{k}}$ является частным относительно ${{\gamma }_{k}}( \cdot )$, а вместо распределения перехода ${{\alpha }_{k}}( \cdot )$ нужно использовать условное распределение ${{\zeta }_{k}}({{i}_{{k + 1}}},{{x}_{{k + 1}}}\,{\text{|}}\,{{i}_{k}},{{x}_{k}},{{j}_{k}},{{y}_{k}})$ нового состояния ${{\bar {X}}_{{k + 1}}}$, уточненное еще и последним наблюдением ${{\bar {Y}}_{k}}$, которое находится как отношение ${{\gamma }_{k}}( \cdot )$ к ${{\beta }_{k}}( \cdot )$:

Требуется с целью преодоления перечисленных во Введении недостатков ФОС малого порядка (ФМП) [6, 7] и ФОС с конечной памятью (ФКП) [8] найти оценку ${{\hat {\bar {X}}}_{k}} = {{[{{\hat {I}}_{k}}\,\,\hat {X}_{k}^{{\text{Т}}}]}^{{\text{Т}}}}$ состояния ${{\bar {X}}_{k}} = {{[{{I}_{k}}\,\,X_{k}^{{\text{Т}}}]}^{{\text{Т}}}}$ системы (1.1) как некоторую функцию ${{\bar {g}}_{k}}( \cdot )$ одного последнего измерения (1.2) и не более чем $l \in \mathbb{N}$ имеющихся предыдущих оценок ${{\hat {\bar {X}}}_{{k - 1}}},{{\hat {\bar {X}}}_{{k - 2}}},\, \ldots $ Поэтому будем искать новый фильтр в виде разностного уравнения l-го порядка для вектора его оценки:

или, подробнее,

Здесь ${{\bar {g}}_{k}}( \cdot ) = {{[\begin{array}{*{20}{c}} {{{\iota }_{k}}( \cdot )}&{g_{k}^{{\text{Т}}}( \cdot )} \end{array}]}^{{\text{Т}}}}$ – искомая структурная вектор-функция фильтра, двойным верхним и нижним индексом отмечен отрезок соответствующей последовательности, например,

тогда как нижний индекс $\max (0,k - l)$ позволяет объединить в (1.5) начальный этап накопления оценок на первых l тактах работы фильтра с последующим этапом их обновления:

В оценке (1.5), в отличие от ФКП [8], старые измерения ${{\bar {Y}}_{{k - 1}}},{{\bar {Y}}_{{k - 2}}}, \ldots ,{{\bar {Y}}_{0}}$ не забываются, ибо они аккумулированы в предыдущих оценках, поэтому время памяти ФБП (1.5) как динамической системы бесконечно.

Замечание 2. Сравнивая уравнение ФБП (1.5) с уравнением двухаргументного вида

более общего трехаргументного ФМП [7] видим, что ФМП является частным случаем ФБП при исключении из аргументов функции ${{\bar {g}}_{k}}( \cdot )$ предыдущего измерения ${{\bar {Y}}_{{k - 1}}}$ и при числе учитываемых в ФБП оценок l = 1.

Условием оптимальности синтезируемой оценки (1.5) выберем минимум среднего значения аддитивной квадратично-простой функции потерь от неверного распознавания номера структуры ${{I}_{k}}$ и неточного оценивания вектора состояния X_k [10]:

где ${{C}_{k}} = C_{k}^{{\text{Т}}} > 0$ – весовая матрица квадратичной формы ошибки оценивания состояния.

2. Вектор состояния и порядок предлагаемого фильтра. Соберем учитываемые в уравнении (1.5) оценки $\hat {\bar {X}}_{{\max (0,k - l)}}^{{k - 1}}$ в блочный вектор ${{\bar {Z}}_{{k - 1}}}$, размерность которого на первых l тактах растет по мере их накопления. Осуществляя для удобства сдвиг на один такт вперед, обозначим

Процесс создания этого вектора представим как рекуррентную процедуру. Так на первых l тактах для получения вектора ${{\bar {Z}}_{k}} = \hat {\bar {X}}_{0}^{k}$ достаточно просто дополнить предыдущий вектор ${{\bar {Z}}_{{k - 1}}} = \hat {\bar {X}}_{0}^{{k - 1}}$ блоком новой оценки ${{\hat {\bar {X}}}_{k}}$, например, сверху:

На всех же последующих тактах $k \geqslant l$ дополнительно необходимо выполнить еще и обновление содержащегося в векторе ${{\bar {Z}}_{{k - 1}}}$ набора предыдущих оценок $\hat {\bar {X}}_{{k - l}}^{{k - 1}}$, удаляя из него устаревший теперь нижний блок ${{\hat {\bar {X}}}_{{k - l}}}$. Последнее запишем как операцию умножения столбца $\hat {\bar {X}}_{{k - l}}^{{k - 1}}$ на “урезающую” блочную матрицу такого удаления:

где Е, О – единичная и нулевая матрицы соответствующих порядков. Действительно,

Поэтому в дополнение к (2.2) имеем

Тогда, объединяя соотношения (2.2), (2.3), получаем для вектора (2.1) рекуррентную формулу

Она позволяет называть ${{\bar {Z}}_{k}}$ вектором состояния фильтра, а его размерность вначале растет и достигает величины $p = l({{n}_{x}} + 1)$.

Наконец, используя обозначение (2.1) в разностном уравнении фильтра (1.5), получаем в дополнение к уравнению состояния фильтра (2.4) еще и формулу его выхода

В результате соотношения (2.4), (2.5) являются рекуррентной версией фильтра (1.5) и доказано следующее утверждение.

Теорема 1 (о конечномерности ФБП). Предлагаемый разностный ФБП (1.5) является конечномерным, его вектор состояния ${{\bar {Z}}_{k}}$ состоит не более чем из желаемого числа l запоминаемых им оценок, так что порядок фильтра не превосходит величины $p = l({{n}_{x}} + 1)$. При этом функция ${{\bar {f}}_{k}}( \cdot )$ его уравнения состояния (2.4) фиксирована, а оптимизации по критерию (1.7) подлежат только дискретнозначная ${{\iota }_{k}}( \cdot )$ и непрерывнозначная ${{g}_{k}}( \cdot )$ функции формулы его выхода (2.5).

Замечание 3. Разностные уравнения ФБП (2.4), (2.5) отличаются от подобных уравнений ФКП [8] наполнением вектора состояния оценками $\bar {Z}_{k}^{{{\text{ФБП}}}} = \hat {\bar {X}}_{{\max (0,k - l + 1)}}^{k}$, а не измерениями $\bar {Z}_{k}^{{{\text{ФКП}}}} = \bar {Y}_{{\max (0,k - l + 1)}}^{k}$, а от уравнения ФМП (1.6), где $\bar {Z}_{k}^{{{\text{ФMП}}}} = {{\hat {\bar {X}}}_{k}}$, – произвольным количеством содержащихся в нем оценок.

3. Оптимизация функций выхода фильтра. Подставляя формулы для оценок ${{\hat {I}}_{k}},\,\,{{\hat {X}}_{k}}$ из (2.5) в критерий (1.7) и решая задачу минимизации полученного функционала, аналогично [1, 2, 7, 8, 10] можно получить такие выражения для пар искомых оптималей.

Для начальных функций имеем известные соотношения

а для последующих – аналогичные формулы

Здесь ${\rm P}_{k}^{{(i)}}({{\bar {y}}_{k}},{{\bar {z}}_{{k - 1}}}) = \Pr \left[ {{{I}_{k}} = i\,{\text{|}}{{\,}_{k}} = {{{\bar {y}}}_{k}},{{{\bar {Z}}}_{{k - 1}}} = {{{\bar {z}}}_{{k - 1}}}} \right]$ – условная вероятность i-го режима, $\rho _{k}^{{(i)}}(x\,{\text{|}}\,{{\bar {y}}_{k}},{{\bar {z}}_{{k - 1}}})$ = = ${{\rho }_{k}}(x\,{\text{|}}\,i,{{\bar {y}}_{k}},{{\bar {z}}_{{k - 1}}})$ – условная плотность вероятности состояния X_k при истинности i-го режима. Произведение двух последних образует условное распределение ${{\rho }_{k}}({{\bar {x}}_{k}}\,{\text{|}}\,{{\bar {y}}_{k}},{{\bar {z}}_{{k - 1}}})$ расширенного вектора состояния ${{\bar {X}}_{k}}$ объекта (1.1):

а сами вероятность и плотность в свою очередь выражаются через это распределение

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

При этом начальные функции (3.1) и соответствующие им оценки ${{\hat {I}}_{0}},{{\hat {X}}_{0}}$ совпадают с известными [1, 2] и определяются по распределению

Здесь и далее символом numerator в знаменателе любой дроби обозначен ее числитель.

Таким образом, формулы (3.2), (3.3) сводят задачу синтеза ФБП к определению оценивающего распределения ${{\rho }_{k}}( \cdot )$ для $k \geqslant 1$.

4. Нахождение оценивающего распределения вероятности. Используя формулу Байеса, а также марковость объекта (1.1) и безынерционность измерителя (1.2), аналогично [8] можно получить следующие соотношения. Оценивающее распределение ${{\rho }_{k}}( \cdot )$ определяется аналогично (3.4) по функции правдоподобия (1.4) и по прогнозирующему распределению ${{\pi }_{k}}( \cdot )$ состояния ${{\bar {X}}_{k}}$, которое не учитывает измерения ${{\bar {Y}}_{k}}$:

В свою очередь прогнозирующее распределение следующего такта ${{\pi }_{{k + 1}}}( \cdot )$ выражается через распределение перехода (1.3) и соответствующее условное распределение ${{\xi }_{k}}( \cdot )$:

тогда как последнее определяется только совместным распределением ${{q}_{k}}( \cdot )$ состояний ${{\bar {X}}_{k}}$, ${{\bar {Z}}_{k}}$ объекта и фильтра

Таким образом, выражения (4.1)–(4.4) задают функциональную зависимость ${{\rho }_{{k + 1}}}[{{q}_{k}}( \cdot )]$ и остается найти распределение ${{q}_{k}}({{\bar {x}}_{k}},{{\bar {z}}_{k}})$.

Для его получения заметим, что пара $({{\bar {X}}_{k}},{{\bar {Z}}_{k}})$ образует марковский вектор, так как определяется замкнутой системой из возмущаемых только белыми шумами ${{\bar {V}}_{k}}$, ${{\bar {W}}_{{k + 1}}}$ разностных уравнений состояния объекта (1.1) и фильтра (2.4):

в которых также использованы формулы для измерения (1.2) и оценки (2.5).

Поэтому для совместного распределения ${{q}_{k}}( \cdot )$ существует рекуррентная формула, которую и запишем далее в виде двух цепочек соотношений, соответствующих этапу накопления оценок в векторе состояния фильтра и этапу их дальнейшего обновления. Учтем также, что известность распределения ${{q}_{{k - 1}}}( \cdot )$ позволяет по соотношениям (4.1)–(4.4) найти ${{\rho }_{k}}( \cdot )$, а последнее по формулам (3.2) определяет и очередные структурные функции ${{\iota }_{k}}( \cdot )$, ${{g}_{k}}( \cdot )$.

Сначала найдем начальное совместное распределение ${{q}_{0}}({{\bar {x}}_{0}},{{\bar {z}}_{0}})$ случайной величины ${{\bar {X}}_{0}}$ и ее оценки ${{\bar {Z}}_{0}} = {{\hat {\bar {X}}}_{0}} = {{\bar {g}}_{0}}({{\bar {Y}}_{0}})$, где функция ${{\bar {g}}_{0}}({{\bar {y}}_{0}}) = {{[\begin{array}{*{20}{c}} {{{\iota }_{0}}({{{\bar {y}}}_{0}})}&{g_{0}^{{\text{Т}}}({{{\bar {y}}}_{0}})} \end{array}]}^{{\text{Т}}}}$ получена по (3.1), (3.4). Используя известные свойства согласованности и умножения распределений, а также краткие обозначения повторяющихся интегрально-суммирующих операций

и произведения символа Кронекера на функцию Дирака последовательно вычисляем или окончательно

Чтобы выразить следующее распределение ${{q}_{1}}({{\bar {x}}_{1}},{{\bar {z}}_{1}})$ через найденное ${{q}_{0}}({{\bar {x}}_{0}},{{\bar {z}}_{0}})$, учтем, что ${{\bar {Z}}_{1}} = \hat {\bar {X}}_{0}^{1}$, и также используем известную связь ${{\hat {\bar {X}}}_{1}} = {{\bar {g}}_{1}}({{\bar {Y}}_{1}},{{\hat {\bar {X}}}_{0}})$. Действуя подобным образом, получим

В результате имеем подобное (4.6) соотношение

в котором распределение ${{p}_{1}}({{\bar {x}}_{1}},{{\hat {\bar {x}}}_{0}})$ легко представим через ${{q}_{0}}({{\bar {x}}_{0}},{{\bar {z}}_{0}})$ как

Повторяя эти преобразования для любого такта $k = \overline {1,l - 1} $, получаем, что для этапа накопления оценок распределения ${{q}_{1}}( \cdot )$, …, ${{q}_{{l - 1}}}( \cdot )$ находятся по рекуррентной цепочке из двух формул:

Здесь ${{\bar {g}}_{k}}( \cdot )$ – уже найденная по предыдущему распределению ${{q}_{{k - 1}}}( \cdot )$ и формулам (3.2), (3.3), (4.1)–(4.4) структурная функция фильтра, а начальное условие ${{q}_{0}}( \cdot )$ для этой процедуры вычисляется по (4.6).

Затем для этапа обновления оценок подобным же образом можно найти и все последующие распределения ${{q}_{l}}( \cdot )$, ${{q}_{{l + 1}}}( \cdot )$ и т.д. В основном это приводит к изменению в цепочке (4.7) только второго соотношения:

Здесь дополнительная операция “интегрирования” типа (4.5) по ${{\hat {\bar {x}}}_{{k - l}}}$ соответствует удалению из вектора ${{\bar {Z}}_{{k - 1}}} = \hat {\bar {X}}_{{k - l}}^{{k - 1}}$ его устаревшего блока ${{\hat {\bar {X}}}_{{k - l}}}$ и отражает процесс последовательной потери информации о самых старых оценках.

Итог разд. 3 и 4 представим в виде следующего утверждения.

Теорема 2 (о синтезе оптимального ФБП). Оптимальные в смысле минимума среднего риска (1.7) структурные функции предлагаемого фильтра (1.5) вычисляются по формулам (3.2) рекуррентно вместе с определяющими их оценивающим ${{\rho }_{k}}( \cdot )$ и совместным ${{q}_{k}}( \cdot )$ распределениями вероятности по интегральным соотношениям (4.1)–(4.8).

5. Алгоритм синтеза фильтра. Полученные выше соотношения задают последовательность действий по получению структурных функций фильтра (1.5). Они выполняются заранее, до появления результатов измерений, а потому практически могут быть осуществлены скорее всего численно, но на мощном вычислителе или достаточно долго. Последнее выгодно отличает этот фильтр, как и все другие виды ФОС, от классического АОФ.

Действительно, прежде всего, найдя по формулам (3.1), (3.4) начальную вектор-функцию ${{\bar {g}}_{0}}({{\bar {y}}_{0}})$, следует по (4.6) вычислить стартовое совместное распределение ${{q}_{0}}( \cdot )$. Последнее позволяет по (4.1)–(4.4) получить первое оценивающее распределение ${{\rho }_{1}}( \cdot )$, с помощью которого затем по (3.2), (3.3) находится следующая вектор-функция ${{\bar {g}}_{1}}( \cdot )$.

Далее таким же образом для любого такта $k \geqslant 1$ следует чередовать такие операции:

– нахождение совместного распределения ${{q}_{k}}( \cdot )$ по предыдущему распределению ${{q}_{{k - 1}}}( \cdot )$ и функции ${{\bar {g}}_{{k - 1}}}( \cdot )$ с помощью соотношений (4.7) для тактов $k = \overline {1,l - 1} $ накопления оценок или с помощью формул (4.8) для последующих тактов $k \geqslant l$ обновления их набора,

– преобразование совместного распределения ${{q}_{k}}( \cdot )$ в оценивающее ${{\rho }_{{k + 1}}}( \cdot )$ снова по (4.1)–(4.4) и получение новой вектор-функции ${{\bar {g}}_{{k + 1}}}( \cdot )$ опять по (3.2), (3.3).

Выборочное приближение к оптимали ${{\bar {g}}_{k}}( \cdot )$ можно найти методом Монте-Карло, осуществляя потактовое статистическое моделирование уравнений объекта (1.1), измерителя (1.2) и фильтра (2.4), (2.5). Его целью является получение достаточно больших пакетов реализаций случайных величин $\{ {{\bar {X}}_{k}}\} $, $\{ {{\bar {Y}}_{k}}\} $, $\{ {{\bar {Z}}_{{k - 1}}}\} $ и построение по ним многомерных гистограмм искомых функций (3.2), так как они определяются условными вероятностями и плотностями (3.3). Одновременно с этим выполняется и анализ точности фильтра. Краткая иллюстрация этой процедуры приведена ниже:

Однако оба указанных способа довольно сложны, поэтому далее рассмотрим два численно-аналитических ковариационных приближения к ФБП. Они позволяют находить вид субоптимальных функций оценивания аналитически, а их параметры – численно, но уже с более простым получением выборочных значений только вероятностей дискретных случайных величин и соответствующих им двух первых условных моментов (средних и ковариаций) непрерывных случайных величин.

6. Гауссовское приближение к фильтру. Воспользуемся известным методом нормальной аппроксимации плотностей вероятности [11].

6.1. Замена распределений вероятностями и плотностями. Конкретизируем соотношения (4.1)–(4.8), заменив в них вероятностные распределения систем из смешанных случайных величин произведениями соответствующих вероятностей дискретных величин и условных плотностей вероятности непрерывных величин.

Так, условное распределение (4.1) из числителя формулы Байеса (4.2) представим в виде произведения

где условная вероятность $\Omega _{k}^{{({{i}_{k}},{{j}_{k}})}}({{\bar {z}}_{{k - 1}}}) = \Pr \left[ {{{I}_{k}} = {{i}_{k}},{{J}_{k}} = {{j}_{k}}\,{\text{|}}\,{{{\bar {Z}}}_{{k - 1}}} = {{{\bar {z}}}_{{k - 1}}}} \right]$ есть интеграл от распределения ${{\omega }_{k}}( \cdot )$:

Для подобного представления совместного распределения ${{q}_{k}}({{\bar {x}}_{k}},{{\bar {z}}_{k}}) = {{q}_{k}}({{i}_{k}},{{x}_{k}},{{\bar {z}}_{k}})$ выделим в векторе состояния фильтра (2.1) его дискретный ${{U}_{k}} \in {{\mathbb{N}}^{l}}$ и непрерывный ${{Z}_{k}} \in {{\mathbb{R}}^{{l{{n}_{x}}}}}$ блоки:

Тогда распределения из (4.4) принимают вид, аналогичный (6.1):

Кроме того у сомножителей оценивающего распределения (3.3) в верхний индекс выделим еще и целочисленное условие j_k:

В результате подставляя (6.1) в (4.2) и используя (6.4), легко получаем запись формулы Байеса в виде двух отдельных выражений для оценивающей вероятности ${\rm P}_{k}^{{({{i}_{k}}\,|{{j}_{k}})}}( \cdot )$ и для оценивающей условной плотности $\rho _{k}^{{({{i}_{k}}\,|{{j}_{k}})}}( \cdot )$:

В свою очередь подобное (4.2) отношение (4.4) позволяет таким же образом выразить и аналогичные сомножители условного распределения ${{\xi }_{k}}( \cdot )$ через сомножители совместного распределения ${{q}_{k}}( \cdot )$:

Кроме того, из произведения (4.1) находим такие представления используемых в (6.5) функций через прогнозирующее распределение ${{\pi }_{k}}( \cdot )$:

Здесь ${\rm B}_{k}^{{(j\,|\,i)}}(x) = \Pr \left[ {{{J}_{k}} = j\,{\text{|}}\,{{I}_{k}} = i,{{X}_{k}} = x} \right]$ – известная из (1.4) условная вероятность индикации

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

6.2. Аппроксимация условной плотности и функций оценивания. С целью получения аналитических выражений для функций оценивания ${{\iota }_{k}}( \cdot )$, ${{g}_{k}}( \cdot )$ выполним гауссовскую аппроксимацию условной плотности вероятности из числителя оценивающей плотности (6.5):

Здесь $N(\theta \,{\text{||}}\,m;D)$ – плотность нормального закона распределения случайного вектора Θ с параметрами $m = {\text{M}}\,[\Theta ]$ и $D = {\text{cov}}\,[\Theta ,\Theta ]$. В качестве условия близости этих плотностей выберем стандартные равенства как их математических ожиданий:

так и ковариаций:

Видно, что моменты $\lambda _{k}^{{(ij)}}( \cdot )$, $\Psi _{k}^{{(ij)}}( \cdot )$ имеют смысл характеристик условного прогноза вектора состояния X_k объекта (1.1) по предыдущему состоянию фильтра ${{\bar {Z}}_{{k - 1}}}$, моменты $\mu _{k}^{{(ij)}}( \cdot ),\Phi _{k}^{{(ij)}}( \cdot )$ – аналогичные характеристики прогноза вектора измерения Y_k, тогда как взаимная ковариация $\Delta _{k}^{{(ij)}}( \cdot )$ – характеристика их совместного прогноза.

Вследствие аппроксимации (6.9) приближение к распределению (4.1) почти, с точностью до множителя, гауссово

так как полностью определяется не только моментами (6.10)–(6.13), но и условной вероятностью (6.2), которая уже является характеристикой совместного прогноза значений двух дискретных случайных величин – номера режима I_k и выхода индикатора J_k. Тогда в результате применения теоремы о нормальной корреляции [12] оказывается справедливым следующее свойство.

Лемма 1 [6, 7]. Пусть в любой момент времени $k \geqslant 1$ существуют конечные условные моменты (6.10)–(6.13) непрерывных случайных величин X_k, Y_k, а их совместно-условная плотность вероятности $\omega \,_{k}^{{(ij)}}( \cdot )$ почти гауссовская (6.9). Тогда оценивающие вероятность и плотность (6.5) имеют следующие аналитические приближения:

причем параметры последней определяются по формулам в которых $ \oplus $ – символ псевдообращения матрицы по Муру–Пенроузу.

В итоге получаем соответствующие приближения к функциям оценивания (3.2).

Лемма 2 [6, 7]. В условиях леммы 1 оптимальные функции оценивания фильтра (3.2) выражаются через моменты (6.10)–(6.13) и вероятность (6.2) по следующим приближенным формулам с функциями (6.15), (6.17):

6.3. Связи между параметрами аппроксимирующей плотности. В случае приближения (6.14) прогнозирующее распределение (4.3), являясь частным по отношению к ${{\omega }_{k}}( \cdot )$:

становится полигауссовским с параметрами в виде трех характеристик прогноза состояния:

Подставляя же в характеристики прогноза измерения (6.11), (6.13) содержащую это распределение дробь (6.7), несложно получить зависимости параметров коррекции $\mu _{k}^{{(ij)}}( \cdot ),$ $\Phi _{k}^{{(ij)}}( \cdot ),$ $\Delta _{k}^{{(ij)}}( \cdot )$ от параметров прогноза $\lambda _{k}^{{(ij)}}( \cdot )$, $\Psi _{k}^{{(ij)}}( \cdot )$, $\Omega _{k}^{{(ij)}}( \cdot )$.

Для компактной записи подобных зависимостей будем использовать такое обозначение интегральной операции усреднения некоторой нелинейности $f( \cdot )$ с гауссовским весом $N( \cdot )$:

Ее результат $\varphi ( \cdot )$ известен как характеристика статистической линеаризации нелинейности $f( \cdot )$.

Лемма 3 [6, 7]. В условиях леммы 1 функциональные характеристики прогноза измерения (6.11), (6.13) выражаются через характеристики прогноза состояния (6.10), (6.12) и вероятность прогноза режима-индикации (6.2) по алгебраическим формулам

с тремя функциями гауссовской коррекции $h_{k}^{{( \cdot )}}( \cdot )$, $F_{k}^{{( \cdot )}}( \cdot )$, $G_{k}^{{( \cdot )}}( \cdot )$. Две первые, основные из них, определяются как характеристики статистической линеаризации (6.20): двух функций $\nu _{k}^{{( \cdot )}}( \cdot )$, ${\rm H}_{k}^{{( \cdot )}}( \cdot )$, полученных условным (при X_k = x) усреднением по шумам $W_{k}^{\Delta },{{W}_{k}}$ индикатора-измерителя (1.2) его функций $o_{k}^{{( \cdot )}}( \cdot )$, $b_{k}^{{( \cdot )}}( \cdot )$: тогда как третья находится с помощью дифференцирования функции $h_{k}^{{( \cdot )}}( \cdot )$ по ее первому аргументу

Следствие 1 [6]. В случае независимости белых шумов индикатора $W_{k}^{\Delta }$ и измерителя W_k условные средние (6.23) распадаются на произведения условной вероятности индикации (6.8) и условных средних только непрерывнозначной функции измерителя (1.2):

В оставшихся пока интегральными формулах для прогнозных характеристик как режима-индикации $\Omega _{k}^{{(ij)}}( \cdot )$, так и состояния $\lambda _{k}^{{(ij)}}( \cdot )$, $\Psi _{k}^{{(ij)}}( \cdot )$ можно выделить внутренние интегралы по переменной y, которые в соответствии с (4.1), (6.1) явно выражаются через прогнозирующее распределение ${{\pi }_{k}}( \cdot )$:

Поэтому из (6.2), (6.10), (6.12) получаем

В свою очередь ${{\pi }_{{k + 1}}}( \cdot )$ выражается через условное распределение ${{\xi }_{k}}( \cdot )$ по (4.3) или, используя (6.3), через его сомножители

Подставляя это в (6.25), для значений этих трех параметров на следующем такте получим интегральные выражения:

в которых $l = {{i}_{{k + 1}}},f = {{j}_{{k + 1}}}$, а новые подынтегральные функции имеют вид

Используя здесь распределение перехода (1.3), эти условные средние можно непосредственно выразить через вероятность индикации (6.8) и две функции объекта наблюдения (1.1):

В результате доказано следующее утверждение.

Лемма 4. Если моменты (6.10), (6.12) существуют, то значения прогнозов режима-индикации $\Omega _{k}^{{(ij)}}( \cdot )$ и состояния $\lambda _{k}^{{(ij)}}( \cdot )$, $\Psi _{k}^{{(ij)}}( \cdot )$ на следующем такте k + 1 формулами (6.26) выражаются через условную вероятность $\Xi _{k}^{{(i\,|\,u)}}(z)$ и результаты усреднения с весом $\xi _{k}^{{(i\,|\,u)}}(x\,{\text{|}}\,z)$ условных средних индикации и состояния объекта (6.27).

Несложно заметить, что условные средние (6.27) упрощаются в следующих распространенных частных случаях.

Следствие 2. Если в уравнениях объекта (1.1) порождающие возмущения $V_{k}^{\Delta }$, V_k независимы, то в (6.27) выделяется сомножитель в виде условной вероятности переключения

и они принимают более простой вид:

Следствие 3. Если индикация режима ${{J}_{k}}$ происходит без учета состояния объекта ${{X}_{k}}$, так что первое уравнение измерителя (1.2) имеет вид ${{J}_{k}} = o_{k}^{{({{I}_{k}})}}(W_{k}^{\Delta })$, то вероятность индикации (6.8) не зависит от состояния объекта: ${\rm B}_{k}^{{(j\,|\,i)}} = M[{{\delta }_{{j,o_{k}^{{(i)}}(W_{k}^{\Delta })}}}]$. Тогда из (6.28) получаем простые произведения:

которые содержат условные средние только непрерывнозначной функции объекта.

Также из (6.27) легко видеть, что наибольшего упрощения эти функции достигают в случае отсутствия в измерении (1.2) переменной индикации J_k, когда ${\rm B}_{k}^{{(j\,|\,i)}}(x) = {{\delta }_{{j,1}}}$.

Следствие 4. Если в наблюдении (1.2) отсутствует переменная индикации J_k, то функции (6.27) не зависят от целочисленной переменной  f и выражаются только через вероятность переключения и непрерывнозначную функцию объекта наблюдения:

6.4. Аппроксимация безусловной плотности и параметры фильтра. Так как условное распределение ${{\xi }_{k}}( \cdot )$ определяется, согласно (4.4), безусловным распределением ${{q}_{k}}( \cdot )$, то с целью получения алгебраических зависимостей и для характеристик прогноза $\Omega _{k}^{{(ij)}}( \cdot )$, $\lambda _{k}^{{(ij)}}( \cdot )$, $\Psi _{k}^{{(ij)}}( \cdot )$ аппроксимируем гауссовской, кроме условной плотности (6.9), еще и безусловную плотность

В качестве условия близости и этих двух плотностей тоже выберем аналогичные (6.10)–(6.13) равенства всех их математических ожиданий и ковариаций:

Тогда само распределение ${{q}_{k}}( \cdot )$, как и условное распределение (6.14), почти гауссово:

и снова по теореме о нормальной корреляции аналогично (6.15)–(6.17) легко получаем следующие связи.

Лемма 5 [6, 7]. Если справедлива гауссовская аппроксимация безусловной плотности вероятности (6.30), то условные вероятность и плотность (6.6) имеют следующие приближения:

При этом математическое ожидание $\vartheta _{k}^{{(i\,|\,u)}}(z)$ приближения к плотности $\xi _{k}^{{(i\,|\,u)}}( \cdot )$ является линейной функцией ее условия

а ее ковариация $T_{k}^{{(i\,|\,u)}}$ от этого условия не зависит и вместе с двумя параметрами из (6.33) определяется по формулам

Тогда, используя в интегральных формулах (6.26) леммы 4 гауссовость условной плотности $\xi _{k}^{{( \cdot )}}( \cdot )$, получим следующий результат.

Лемма 6 [7]. Если имеет место гауссовское приближение (6.30), то вероятность (6.2) прогноза номера структуры объекта ${{I}_{{k + 1}}}$ и его индикации ${{J}_{{k + 1}}}$, а также характеристики (6.10), (6.12) условного прогноза следующего состояния ${{X}_{{k + 1}}}$ выражаются через функции (6.32), (6.33) предыдущих измерения Y_k и состояния фильтра ${{\bar {Z}}_{k}}$ по алгебраическим формулам:

с тремя функциям гауссовского прогноза: определяемыми по условным средним (6.27) с помощью интегральной операции (6.20).

Отметим, что в построении зависимостей (6.35) по функциям (6.36) участвуют две тройки числовых параметров гауссовских приближений (6.32). Это вероятность дискретных компонентов векторов состояния объекта и фильтра $Q_{k}^{{(i\,u)}} = \Pr \left[ {{{I}_{k}} = i,{{U}_{k}} = u} \right]$ вместе с двумя условными моментами $m_{k}^{{z\,|\,i\,u}}$, $D_{k}^{{z\,|\,i\,u}}$ непрерывной части вектора состояния фильтра из (6.31), а также определяемая по (6.34) с учетом всех моментов (6.31) тройка параметров $\Gamma _{k}^{{(i\,|\,u)}},$ $\kappa _{k}^{{(i\,|\,u)}}$, $T_{k}^{{(i\,|\,u)}}$. Они находятся методом Монте-Карло в процессе последовательного статистического моделирования приведенных ниже уравнений предлагаемого фильтра так же, как и гауссовского ФМП [7].

6.5. Уравнения гауссовского фильтра. Приведенные выше леммы позволяют сформулировать следующий результат.

Теорема 3 (об уравнениях гауссовского ФБП). Если имеют место гауссовские аппроксимации плотностей вероятности (6.9), (6.30), то процедура вычисления оценок состояния многорежимного объекта (1.1) по наблюдениям (1.2) субоптимальным ФБП состоит из рекуррентной цепочки следующих шагов. После нахождения классических начальных оценок ${{\hat {I}}_{0}}$, ${{\hat {X}}_{0}}$ и задания начальных условий $\Omega _{0}^{{(ij)}},$ $\lambda _{0}^{{(ij)}},$ $\mu _{0}^{{(ij)}},$ $\Psi _{0}^{{(ij)}},$ $\Delta _{0}^{{(ij)}},$ $\Phi _{0}^{{(ij)}}$ $i = \overline {1,L} $, $j = \overline {1,M} $, на каждом такте $k \geqslant 0$ последовательно осуществляются следующие операции.

1. Регистрация новых наблюдений ${{J}_{k}}$, ${{Y}_{k}}$ и формирование вектора состояния фильтра:

Здесь ${{C}_{u}} = \left[ {\begin{array}{*{20}{c}} {{{E}_{{l - 1}}}}&{{{O}_{{(l - 1) \times 1}}}} \end{array}} \right]$, ${{C}_{z}} = \left[ {\begin{array}{*{20}{c}} {{{E}_{{(l - 1){{n}_{x}}}}}}&{{{O}_{{(l - 1){{n}_{x}} \times {{n}_{x}}}}}} \end{array}} \right]$ – матрицы удаления из векторов ${{U}_{{k - 1}}}$ и ${{Z}_{{k - 1}}}$ их последних блоков.

2. Конкретизация условных прогнозов состояния и измерения по индикации J_k:

и подобная конкретизация параметров фильтра по вектору ${{U}_{k}}$ имеющихся оценок режимов:

3. Вычисление реализаций параметров гауссовского приближения к оценивающему распределению ${{\rho }_{k}}( \cdot )$ в соответствии с (6.15), (6.17) (лемма 1):

4. Определение оценок состояния объекта по (2.5) (теорема 1) и (6.18) (лемма 2):

5. Нахождение модификаций $\Xi _{k}^{{(i)}},$ $\vartheta _{k}^{{(i)}}$, соответственно, вероятности $Q_{k}^{{(i)}}$ по (6.32) и вектора ${{Z}_{k}}$ по (6.33) (лемма 5):

6. Вычисление характеристик условных прогнозов следующего состояния объекта ${{\bar {X}}_{{k + 1}}}$ и выхода индикатора ${{J}_{{k + 1}}}$ по (6.35) с использованием функций прогноза (6.36) (леммы 4, 6):

7. Нахождение значений характеристик условных прогнозов следующего измерения ${{Y}_{{k + 1}}}$ по (6.21) с использованием функций коррекции (6.22) (лемма 3):

8. Если дискретное время $k < K$, то положить $k: = k + 1$ и вернуться к п. 1, иначе завершить вычисления.

7. Линеаризованное приближение к фильтру. К сожалению, построение гауссовского фильтра требует выполнения довольно сложных интегральных операций нахождения как условных средних измерителя (6.23) и всей системы наблюдения (6.27), так и определения характеристик статистической линеаризации (6.20) этих нелинейностей. Поэтому рассмотрим и другое известное приближение, существенно более простое, а потому и более популярное.

Ограничимся частными случаями (следствия 1–3), когда, согласно (6.24), (6.28), (6.29), функции коррекции (6.22) и прогноза (6.36) находятся гауссовским усреднением более простых условных средних:

Они образуются как произведения постоянной вероятности индикации ${\rm B}_{k}^{{(j\,|\,i)}}$ соответственно на условные средние непрерывнозначной функции измерителя (1.2):

на переменную вероятность переключения режима ${\rm A}_{k}^{{(l\,|\,i)}}(x)$ и на произведения последней на условные средние непрерывнозначной функции объекта (1.1):

Пусть три следующие функции дифференцируемы по своим аргументам и каждую из них с достаточной точностью можно линеаризовать по Тейлору:

Здесь ${{x}^{o}} = x - m$ – малое отклонение, $m_{k}^{\delta }$ – среднее значение случайной величины ${{\Delta }_{k}}$, штрихами и нижними индексами обозначены матрицы частных производных по соответствующим переменным:

а дугами сверху – сечения функций при фиксированных значениях второго аргумента:

Тогда аналогично [6] получим такое утверждение.

Следствие 5 (о структурных функциях линеаризованного ФБП). Пусть имеют место следствия 1–3 и справедливы приближения (7.1). Тогда две основные функции гауссовской коррекции (6.22) аппроксимируются следующими выражениями:

а три функции гауссовского прогноза (6.36) имеют такие приближения:

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

Заключение. Предложен метод синтеза быстрого логико-динамического ФБП, запоминающего несколько последних оценок вектора состояния дискретного многорежимного объекта. Этот фильтр обладает наивысшей точностью в своем классе, а его структурные функции находятся заранее, после чего сам фильтр может быть реализован в реальном времени на дешевом вычислителе. Получены формулы для нахождения его оптимальных функций оценивания, найдена цепочка интегральных формул для соответствующих распределений вероятности, представлены способы вычисления структурных функций фильтра как по интегральным соотношениям, так и методом Монте-Карло. Однако последний численный алгоритм громоздок технически из-за необходимости получения гистограмм искомых функций.

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

Отметим, что в частном случае однорежимного объекта, когда в уравнениях (1.1), (1.2) отсутствуют целочисленные переменные I_k, J_k, предложенный логико-динамический фильтр вырождается в обычный динамический вроде [13].

В последующих статьях на эту тему автор планирует продемонстрировать эффективность полученных ковариационных приближений к ФБП в их сравнении с подобными приближениями не только к классическому АОФ, но и с другими версиями ФОС, как малого порядка [6, 7], так и с конечной памятью [8]. Ранее подобные аналитические сравнения сложности самих алгоритмов дискретных АОФ и ФОС уже приводились в [13], а численные сопоставления их точности были выполнены как на модельном логико-динамическом примере [6], так и на рассмотренной в [13, 14] задаче автономного оценивания траектории однорежимного спускаемого аппарата.