11.7. Реализация алгоритмов

Основной, общий алгоритм Гаусса — Ньютона был представлен в форме (11.44) или (11.57). В таком виде он неудобен непосредственного применения из-за необходимости обращения на каждом шаге -матрицы . В данном разделе обсуждаются вопросы наилучшей реализации рекуррентных алгоритмов. Более подробное обсуждение изложено в книге Льюнга и Седерстрема [262, гл. 6].

Использование леммы об обращении матриц. Применяя лемму об обращении матриц (11.10) к (11.44), аналогично (11.11) получим алгоритм (для случая векторного выхода)

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

Использование факторизации. Как отмечалось в разделе 10.1, удобно представлять матрицы данных в факторизованном виде (см. (10.8)), что позволяет работать с лучше обусловленными матрицами. Для рекуррентной идентификации это означает, что представляется как произведение матриц, и заменяется алгоритмами пересчета этих матриц вместо Удобно представить

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

в которой верхняя треугольная матрица с единичными диагональными элементами, диагональная матрица. Алгоритм пересчета в (11.70) был получен Поттером [327], а устойчивые в вычислительном отношении алгоритмы для в (11.71) были разработаны Бирманом [47]. Здесь будут изложены некоторые детали соответствующих алгоритмов, непосредственно основанных на преобразовании Хаусхолдера (см. задачу и полученных Морфом и Кайлатом [298].

Шаг 1. Пусть в момент времени нижняя треугольная матрица является квадратным корнем из как в (11.70). Пусть квадратный корень из Образуем -матрицу

Шаг 2. Применяем ортогональное -преобразование

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

(Очевидно, являются нижними треугольными.)

Шаг 3. Теперь для как в имеем

Следовательно,

Проверка. Умножение матрицы (11.73) на нее же транспонированную дает

Используя тот факт, что сравнивая с убеждаемся в справедливости уравнений (11.74).

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

Решетчатые алгоритмы Соотношения (10.31) и (10.33) задают схему решетчатого фильтра, которая может быть применена при вычислении предсказаний для моделей конкретных структур (см. также задачи ). Эта схема не является рекуррентной, так как переменные на самом деле должны иметь еще индекс отражающий их зависимость от всего набора данных через (10.33), Представляется заманчивым подход, связанный с аппроксимациеи в (10.33) с помощью прошлых, уже вычисленных невязок и Тогда может вычисляться рекуррентно, и мы получаем следующую схему:

Этот алгоритм был разработан Гриффитсом [148] и Макхоулом [272] и получил название градиентный решетчатый алгоритм. Как было показано, он основан на аппроксимации и, как следствие, не реализует точное значение оценки, вычисляемой по накопленным данным. Интересно отметить, что ее можно получить посредством небольшой модификации схем пересчета Это было показано Ли, Морфом и Фридлапдером [229]. Ниже представлен полученный в результате алгоритм для случая многомерного сигнала в (10.13) (что включает приложения к динамическим системам; см. задачу Кроме того, в алгоритм введен постоянный фактор забывания

1. Начинаем при полагаем:

2. В момент запоминаем

3. В момент вычисляем для

4. Для вычисляем

5. Вычисляем для :

6. Переходим к шагу 2.

Предсказание величины основанное на задается, таким образом, выражением

Оно может быть вычислено до получения путем представления в виде