Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 42
Текст из файла (страница 42)
Этот факт широко известен и вошел в учебные курсы численныхметодов алгебры (разложение Холесского) — cм. выше разд. 6 и [74].13.1 Фильтрация Калмана в историческом аспектеДругая идея, также вошедшая в учебные курсы, — ортогонализация матриц (cм. выше разд. 7 и [74]) — развивалась параллельно для численногорешения плохо обусловленных и переопределенных систем по методу наименьших квадратов (МНК) [112].
Важный промежуточный результат техлет — опубликование обзора [104] и монографии Бирмана [15].Кроме исходных аэрокосмических приложений, эти и другие современныеметоды калмановской фильтрации вошли в другие области: эконометрика[93], метеорология [110], спутниковая геодезия [91], телекоммуникационныесети [121], высоконадежная спутниковая навигация [111], обработка изображений в реальном времени [124] и многое другое.Казалось бы, основные идеи исчерпаны и все полезные результаты вэтой области получены. Однако широта приложений заставила использовать новую вычислительную среду реализации этих алгоритмов — параллельные компьютеры [8], в результате чего они оказались представлены вновой, «блочной» форме с ортогонализацией [17] (эти инновации изложеныниже в разд.
14). Появилось учебное пособие, рассчитанное на широкий кругпользователей и отражающее почти все достижения в этой области к концу20 столетия [16].Указанные исследования велись за пределами СССР. По настоящее времяони остаются в России малоизвестны, применены в немногих сложныхпроектах (например, в авиационном приборостроении [80] или судостроении [78, 79]) и преподаются в немногих специальных курсах [77]. В российской научной литературе квадратно-корневые алгоритмы фильтрациисовсем не изложены [27] или освещены слабо [64].
Так, в книге Огаркова,1990 г. [64] этой теме посвящен лишь один параграф 3.7. Нехватка литературы сказывается в том, что многие отечественные специалисты, работающие в прикладных областях регрессионного моделирования или эконометрики, продолжают использовать алгоритмы, которые можно считать устаревшими, несмотря на то, что они испытывают значительные трудностив случае плохо обусловленной схемы наблюдения или при мультиколлинеарности регрессоров [26]. Даже в публикациях последних лет [33] авторыпользуются уравнениями фильтра Калмана в их классическом виде, а нев форме современных вычислительно эффективных алгоритмов. Ситуацияв России такова, что здесь лишь немногие специалисты уделяли вниманиеквадратно-корневым, а также ортогонализованным блочным алгоритмам,принадлежащим новому классу параллельно-ориентированных реализаций[64, 73, 78, 79, 86, 107, 108, 132, 135, 138].27113 Устойчивые алгоритмы фильтрацииПервая цель данного раздела — привлечь внимание студентов, магистрантов, аспирантов и специалистов к этому направлению исследований.Здесь мы ограничиваемся лишь кратким обзором алгоритмов нового класса,поскольку более подробные обзоры доступны.
Коллективная монография[73] и книга [15] содержат подробный анализ и сравнительные таблицыэффективности стандартных и квадратно-корневых алгоритмов дискретнойфильтрации. Прекрасный обзор, — начиная от зарождения идей у ГалилеоГалилея (1564–1642), через работы Карла Фридриха Гаусса (1777–1855), Норберта Винера (1894–1964) и многих других (список огромен) до РудольфаЭмиля Калмана (1930–) и последних достижений, датируемых до 2001 года,содержится в учебном пособии [16].Вторая цель данного раздела — показать прикладные возможностиквадратно-корневых алгоритмов. Для примера рассмотрено их использование в актуальной задаче обнаружения момента вхождения судна в маневр,новое решение которой содержится в [41].
Соответствующий материал с опорой на достижения в этой области [60, 90, 95, 114, 128] помещен в подразд. 13.12. Работа [95], решающая ту же задачу сопровождения маневрирующих целей, также использует квадратно-корневую реализацию фильтров.В отличие от [95], в [41] применен расширенный фильтр первого, а не второго порядка; его построение здесь включено для примера.13.2Стандартный фильтр КалманаРассмотрим источник данных, представимый в виде линейной динамической системы, возмущаемой дискретным белым шумом:x(ti+1) = Φ(ti+1, ti)x(ti) + B(ti )u(ti) + Γ(ti )w(ti),z(ti ) = H(ti)x(ti) + v(ti),(13.1)(13.2)где Φ(ti+1, ti ), B(ti), Γ(ti ), H(ti ) — известные матрицы-параметры системы;x(ti) — n-мерный вектор состояния системы, u(ti) — r-мерный вектор входного воздействия, z(ti ) — m-мерный вектор измерений, w(ti ) и v(ti) — независимые нормально распределенные векторы шумов с нулевым средним иизвестными ковариационными матрицами Q(ti ) и R(ti ) соответственно, причём Q(ti) ≥ 0, R(ti) > 0.
Начальный вектор состояния системы x0 распределен по нормальному закону с математическим ожиданием x̄0 и ковариационной матрицей P0 .При заданных матрицах-параметрах системы Φ(ti+1, ti ), B(ti ), Γ(ti ), H(ti )27213.2 Стандартный фильтр Калманаи начальных условиях x̄0, P0 решение задачи оптимального оцениваниявектора состояния x(ti) системы (13.1), (13.2) дается фильтром Калмана.Алгоритм этой рекуррентной обработки данных {z(ti ); ∀i ≥ 1} определенследующими уравнениями:xb(t−x(t+i ) = Φ(ti , ti−1 )bi−1) + B(ti−1 )u(ti−1),+−TTP (ti ) = Φ(ti, ti−1)P (ti−1)Φ (ti, ti−1) + Γ(ti−1)Q(ti)Γ (ti−1), −−T−1TK(ti) = P (ti )H (ti )[H(ti)P (ti )H (ti ) + R(ti )] ,(13.3)−−),)−K(t)H(t)P(t)=P(tP (t+iiiii),ν(ti) = z(ti ) − H(ti )bx(t−i+−b(ti ) + K(ti)ν(ti).xb(ti ) = xЗамечание 13.1.
Уравнения (13.1), (13.2) описывают широкий классмоделей систем, например, систем с управлением. Это могут быть и системыс возможными параметрическими нарушениями, — если уравнения (13.1),(13.2) рассматривать для каждого режима функционирования (с нарушением или нет) отдельно. Это могут быть и адаптивные системы с идентификацией, если матрицы, входящие в (13.1), (13.2), содержат неизвестныепараметры, которые подлежат определению в процессе функционированиясистемы для целей слежения за состоянием системы или для целей управления. В частности, управление u(ti) может отсутствовать. Его наличие илиотсутствие, так же как и наличие или отсутствие нарушений, не влияет напринципы эффективной численной реализации фильтров.
Чтобы показатьуниверсальность таких реализаций, управление не исключено из уравнений(13.1), (13.2). Наличие или отсутствие нарушений здесь не отмечено индексом режима возле обозначений матриц исключительно для упрощения записей.Общность исходной модели (13.1), (13.2) и, соответственно, уравненийфильтра Калмана (13.3), объясняется также следующим.Замечание 13.2. Если Φ(ti+1, ti ) = I, u(ti) ≡ 0 и w(ti) ≡ 0, томодель (13.1), (13.2) совпадает с той моделью, которая принята в теориирегрессионного моделирования, когда x(ti) = const — параметр, подлежащий оцениванию по наблюдениям z(ti ), i = 1, 2, . . . , N .
Алгоритм (13.3) осуществляет оптимальное оценивание этой регрессионной модели. Здесь принципиальная особенность заключается в том, что регрессионная модель ипроцесс ее оценивания — последовательные. Полная матрица «регрессоров»H = [H(ti); i = 1, 2, . . . , N ]T расщеплена на отдельные «порции» H(ti ) соответственно тем «порциям» наблюдений z(ti ), которые поступают в обработку27313 Устойчивые алгоритмы фильтрациив последовательные моменты времени i = 1, 2, . . . , N . Из теории наименьших квадратов, изложенной выше (см. подразд. 11.6), известно, что для значения текущей оптимальной оценки по экспериментальным данным безразлично, поступили эти данные в обработку целиком или порциями.
Последовательная (рекуррентная) форма алгоритма метода наименьших квадратовв форме (13.3) предпочтительнее, чем единовременная обработка всех данных «целиком» — методом решения нормальных уравнений (см. разд. 12),хотя последнее все еще (к сожалению) остается традицией в эконометрикеили при обработке данных наблюдений в астрономии или спутниковой геодезии [26].Замечание 13.3.
Еще раз отметим важное обобщение, присутствующее в модели (13.1), (13.2) и также в алгоритме фильтра (13.3). В отличие отпредыдущих разделов этой книги, здесь источник обрабатываемых данных{z(ti ); ∀i ≥ 1} — не статический, а динамический. Он содержит информацию не о неизменном по своему значению неизвестном (возможно, случайном) векторе x(ti) = const, — эта ситуация выделена в Замечании 13.2, —а об изменяющемся векторе x(ti) = var состояния динамического объекта(13.1), находящегося (в общем случае) под воздействием детерминистскогоуправления u(ti) и/или случайного возмущения w(ti).13.3Скаляризованная форма фильтра КалманаКогда матрица R(ti ) диагональная, R(ti) = diag [r1(ti ), r2(ti ), .
. . , rm(ti )],возможно произвольное расщепление системы наблюдений (вектора z(ti )) нааприорную и текущую части, — см. выше подразд. 10.8, 11.6 или [15]. Тогдапредпочтительно вместо (13.3) использовать фильтр Калмана со скалярнойобработкой измерений. Скаляризуя обработку измерений в (13.3), получаем:I. Экстраполяция:xb(t−x(t+i ) = Φ(ti , ti−1 )bi−1) + B(ti−1 )u(ti−1),+TT(13.4)P (t−i ) = Φ(ti , ti−1 )P (ti−1)Φ (ti , ti−1 ) + Γ(ti−1)Q(ti−1)Γ (ti−1 ),xb(t+P (t+0 ) = x̄0 ,0 ) = P0 .II. Обработка измерения:274А. Начальное присваивание: Pe = P (t−e=xb(t−i ), xi ).13.4 Стабилизованный фильтр Калмана–ДжозефаБ.
m-кратное повторение процедуры скалярного обновления.Для j = 1, 2, . . . , m выполнять:α := hT Peh + rj (ti ) , K := Peh/α , Pb := Pe − KhT Pe, (13.5)xb := xe + K(z − hT xe)(13.6)с экстраполяцией между повторениями: Pe := Pb , xe := xb.bВ. Завершающее присваивание: P (t+xb(t+b,i ) := P ,i ) := xTгде h — j-й столбец матрицы H (ti ); z — j-й элемент вектора z(ti ) ,j = 1, 2, . . . , m .13.4Стабилизованный фильтр Калмана–ДжозефаЭтап экстраполяции стабилизированного алгоритма совпадает с этапомэкстраполяции стандартного алгоритма Калмана.
Поэтому приведем лишьэтап обработки измерения. Джозеф [15] предложил для алгоритма (13.3)использовать общую формулу, справедливую для матрицы Pb при любом, необязательно оптимальном K: Pb = (I − KH)Pe(I − KH)T + KRK T . (Смыслприменяемых здесь обозначений матриц виден из подразд. 13.3.)Замечание 13.4.При оптимальном Kopt = PeH T (H PeH T + R)−1эта формула превращается в Pb = Pe − Kopt H Pe. При подстановке сюда фор-мулы для Kopt получаем дискретное алгебраическое уравнение Риккати(Discrete Algebraic Riccati Equation, DARE) [109].
Поэтому этот алгоритмможет также рассматриваться как процесс решения DARE.При таких вычислениях результирующая матрица сохраняет симметричность; кроме того, исчезает опасная для потери положительной определенности операция вычитания матриц, однако, сложность вычислений возрастаетпочти вдвое. Второй этап приобретает вид:II.