Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 37
Текст из файла (страница 37)
11.2) z̃ṽR̃=x+,(11.17)zvAвключающей, помимо текущих экспериментальных данных z, «дополнительные» экспериментальные данные z̃, соответствующие имеющейся в наличииаприорной информации x̃, Λ̃.Обозначим через x̂ МНК-решение расширенной системы (11.17), доставляющее минимум функционалу (11.16). Из этого критерия для x̂ получаем,аналогично (11.3), нормальные уравнения(Λ̃ + AT A)x̂ = Λ̃x̃ + AT z.(11.18)Простая модификация данной в п. 11.3 статистической интерпретации приводит к выражениюΛ̂(x̂ − x) = Λ̃(x̃ − x) + AT v.242(11.19)11.5 Включение предшествующего МНК-решенияТак как E {x̃ − x} = 0 и E {v} = 0, то и E {x̂ − x} = 0, то есть x̂ есть такженесмещенная оценка: E {x̂} = x.
Если ошибкааприорнойоценки и ошибкаT= 0, то после «возвеизмерения взаимно некоррелированы, E (x̃ − x)vдения в квадрат» обеих частей (11.19) и осреднения получим ковариацию(11.20)Px̂ = E (x̂ − x)(x̂ − x)T = (Λ̃ + AT A)−1 = Λ̂−1 ,где через Λ̂ обозначена информационная матрица апостериорной оценки x̂,Λ̂ = Λ̃ + AT A.Замечание 11.1. Матрица Λ̃ не обязана быть невырожденной, хотяв (11.15) это формально предполагалось. В действительности, априорноезнание некоторых (или всех) компонент вектора x может быть исчезающемало, так что соответствующие строки и столбцы в информационной матрице Λ̃ заполняются исчезающе малыми числами или даже нулями. Приэтом соотношение (11.15) сохраняет силу в пределе, в том смысле, что в ковариационной матрице P̃ соответствующие диагональные элементы стремятсяк +∞, в то время как другие остаются ограниченными.
Заметьте, что (11.18)как условие минимума функционала (11.16) получается при произвольнойнеотрицательно определенной матрице Λ̃. Если Λ̃ = 0, то (11.18) сводитсяк (11.3) и (11.20) сводится к (11.14). Таким образом, информационная и ковариационная матрицы взаимно обратны не только формально в силу определения (11.15), но и по существу как меры достоверности/недостоверностиаприорной оценки x̃ вектора x.11.5Включение предшествующего МНК-решенияРасширенная система (11.17) показала, что априорные статистическиесведения о векторе x, поступающие в виде несмещенной оценки x̃ и ее ковариации P̃ , могут быть интерпретированы как воображаемые добавочныерезультаты z̃ = R̃x̃ некоего эксперимента.
Это наводит на мысль, что ичисто алгебраическую задачу отыскания МНК-решения системы уравненийможно решать последовательно: предварительно найти x̃ как МНК-решениечасти системы, а затем включить это x̃ в полную систему, чтобы найти ееМНК-решение x̂. Такое разделение системы на две части — априорную итекущую — можно назвать ее «расщеплением» (см. формальное решениерассматриваемого здесь вопроса о расщеплении в подразд. 10.8).24311 Оценивание по методу наименьших квадратовПусть система уравнений произвольно расщеплена на эти две подсистемы: fwR̃=x+.(11.21)zvAМНК-решение x̂ этой полной системы, доставляющее минимум функционалу J1 (x) = wT w + v T v, есть решение нормальных уранений T R̃ f TTT.(11.22)x̂ = R̃ AR̃ AzAДопустим, найдено МНК-решение x̃ для подсистемы f = R̃x + w из критерия минимума функционала J(x) = wT w. Как отмечено в (11.5)–(11.6), оноудовлетворяет двум условиям:R̃x̃ = z̃ ,z̃ ∈ R(R̃) ,f − z̃ ⊥ R(R̃) ,x̃ ∈ R(R̃T ) .Разностный вектор r = f − z̃ ортогонален пространству столбцов R(R̃) матрицы R̃ и, следовательно, лежит в левом нуль-пространстве N (R̃T ), определяемом как все векторы y, удовлетворяющие уравнению R̃T y = 0.
ПоэтомуR̃T f = R̃T (z̃ + r) = R̃T z̃ + R̃T r = R̃T z̃ .Следовательно, уравнения (11.22) совпадают с уравнениями Tz̃R̃,x̂ = R̃T ATR̃ ATzAкоторые, в свою очередь совпадают с уранениями (11.18), так как R̃T R̃ = Λ̃.Тем самым доказано, что МНК-решение x̂ данной системы (11.21) совпадаетс МНК-решением системы (11.17), отличающейся от (11.21) тем, что в неевместо f включен вектор z̃, равный проекции f на R(R̃), z̃ = R̃x̃, где x̃ —МНК-решение левой подсистемы в (11.21).11.6Рекурсия МНК в стандартной информационнойформеИнтерпретация априорных статистических данных как дополнительныхнаблюдений, или, что равносильно, учет имеющегося МНК-решения подсистемы после добавления в систему новой порции уравнений, является краеугольным камнем рекурсии для МНК.
Это дает возможность обрабатывать24411.6 Рекурсия МНК в стандартной информационной формеэкспериментальные данные по мере их поступления, то есть решать задачуо наименьших квадратах по мере поступления новых уравнений. Это равносильно также тому, что исходную большую совокупность экспериментальных данных можно «расщеплять» произвольно на порции и последовательновключать их в обработку.
Результат такой последовательной обработки, какдоказано выше (подразд. 10.8 и 11.5), будет равносилен результату обработки всей совокупности экспериментальных данных целиком.Результаты при статистической интерпретации рекурсивны потому, чтотекущие величины — оценка x̂ и ковариация Px̂ — становятся априорными икомбинируются с новыми данными, чтобы образовать обновленные оценкуи ковариацию. При этом существенно, что результаты (11.18) и (11.20) независят (теоретически) от того, какой квадратный корень R̃ в разложенииΛ̃ = R̃T R̃ использован.
Эта свобода позволяет выбирать R̃ из соображенийбольшей вычислительной точности. Кроме того, если лишь окончательнаяоценка (МНК-решение) необходима, то лучше не находитьPпромежуточныхоценок, а просто накапливать информационную матрицуATj Aj и суммуP TAj zj , и лишь в нужный момент (например, в самом конце) вычислитьрешение.Информационную форму последовательного МНК запишем, вводя матрицу (Λ | d) .I. Инициализация. Устанавливают начальные значения x0, Λ0 :d0 = Λ0 x0,Λ d := Λ0 d0 .Эти начальные значения берут из априорных данных: x0 = x̃, Λ0 = Λ̃.II. Обработка наблюдений. Вводят очередную «порцию» наблюдений z == Ax + v:(11.23)Λ d := Λ d + AT A z .В общем случае ненормализованных статистических данных z вместо (11.23) используют алгоритм:Λ d := Λ d + AT R−1 A z .(11.24)III. Выдача результата.
После последовательной обработки всех порцийнаблюдений или в нужный момент, когда Λ−1 существует, вычисляютx̂ = Λ−1d,Px̂ = Λ−1.24511 Оценивание по методу наименьших квадратов11.7Рекурсия МНК в стандартной ковариационнойформеПусть априорная и апостериорная оценки, x̃ и x̂, характеризуются невырожденными информационными матрицами Λ̃ и Λ̂, соответственно. Тогдасуществуют обратные к ним ковариационные матрицы, (11.15) и (11.20).Разрешим (11.18) относительно x̂:x̂ = (Λ̃ + AT A)−1Λ̃x̃ + (Λ̃ + AT A)−1AT z .ОбозначимL = (Λ̃ + AT A)−1Λ̃,K = (Λ̃ + AT A)−1AT(11.25)и преобразуем:x̂ = Lx̃ + Kz − KAx̃ + KAx̃ = x̃ + K(z − Ax̃) ,так как L+KA = In . Для определения K воспользуемся в (11.25) следующейважной леммой.Лемма 11.1.−1−1−1−1−1(Λ1 − Λ12Λ−1= Λ−12 Λ21 )1 + Λ1 Λ12 (Λ2 − Λ21 Λ1 Λ12 ) Λ21 Λ1 ,(11.26)где предполагается, что все матрицы согласованы по размерам и тре буемыеобращения матриц существуют.Упражнение 11.1.
Докажите лемму 11.1, рассматривая блочные матрицы−1 P1 P12Λ1 Λ12=,P21 P2Λ21 Λ2и выписывая поблочно равенства ΛP = I и P Λ = I.Применим лемму 11.1 при Λ1 = Λ̃, Λ12 = AT , Λ21 = ΛT12 = A, Λ−12 = −I.Получим(Λ̃ + AT A)−1 = Λ̃−1 − Λ̃−1AT (AΛ̃−1AT + I)−1AΛ̃−1.Обозначим: P̃ = Λ̃−1, P̂ = Λ̂−1.
Имеем из подразд. 11.5 выражение (11.23):Λ̂ = Λ̃ + AT A. Следовательно, P̂ = P̃ − P̃ AT (AP̃ AT + I)−1AP̃ . Так какK = P̂ AT , тоK = P̃ AT − P̃ AT (AP̃ AT + I)−1AP̃ AT == P̃ AT (AP̃ AT + I)−1[AP̃ AT + I − AP̃ AT ] = P̃ AT (AP̃ AT + I)−1.24611.7 Рекурсия МНК в стандартной ковариационной формеТаким образом, при статистической интерпретации (см. подразд. 11.4) получаем возможность уточнять априорную несмещенную оценку x̃ и уменьшать ее ковариацию P̃ , благодаря включению в процесс обработки вновьпоступивших результатов наблюдений z = Ax +v и применяя к ним следующий алгоритм, известный как стандартный алгоритм Калмана (рис. 11.3).AОбъектАприорнаямодельzẑ = Ax̃x̃, P̃Экстраполяциямоделиx̃ := x̂, P̃ := P̂+ν−ОбновлениемоделиKν+Апостериорнаямодельx̂, P̂Рис. 11.3.
Рекурсия МНК в стандартной ковариационной форме (схема Калмана).Этап 1 — обновление модели по наблюдениям. Этап 2 — экстраполяция модели между наблюдениями. Априорная модель дает предсказание ẑ для отклика z объекта.Разность ν имеет смысл обновляющего процесса.
Весовая матрица K, умноженнаяна ν, обеспечивает обновление априорной модели по наблюдению z, т. е. переход капостериорной моделиI. Инициализация. Начальные значения x0 , P0 :x̃ := x0 ,P̃ := P0 .(11.27)II. Обработка наблюдений. Очередная порция наблюдений z = Ax + v :K = P̃ AT (AP̃ AT + I)−1 ,P̂ = P̃ − KAP̃ ,x̂ = x̃ + K(z − Ax̃) .(11.28)24711 Оценивание по методу наименьших квадратовIII.
Экстраполяция. Распространение оценки x̂ и ее ковариации P̂ междунаблюдениями, т. е. к моменту повторения этапа II со следующей порцией наблюдений :P̃ := P̂ ,x̃ := x̂ .(11.29)Замечание 11.2. Первая и вторая формулы в (11.29) выражают правило экстраполяции только для статической МНК-задачи оценивания, т. е.для случая, когда оцениваемый вектор x не изменяется в процессе наблюдения: x = const. Они заменяются другими в случае динамической МНКзадачи оценивания — первыми двумя формулами в системе уравнений (13.3)(см. ниже стр. 273).Равным образом данный алгоритм пригоден и без статистической интерпретации (см.
подразд. 11.5), когда алгебраическая задача отыскания МНКрешения переопределенной системы решается последовательно. Такое решение может стартовать с условно «пустой» системы уравнений. Практически,это должно отвечать условию Λ̃ = 0, которое легко реализовать в информационной форме (подразд. 11.5). В ковариационной форме данное условиеможно реализовать лишь приближенно, например, полагая P0 = ε−2I, гдеε → 0. При таком выборе величина x0 практически не имеет влияния надальнейший процесс, так что она может быть взята равной нулевому вектору, x0 = 0.
После такой инициализации уравнения исходной переопределенной системы могут вводиться в этап обработки измерений последовательными порциями, и в порциях может содержаться любое, не обязательно однои то же, число уравнений, выражаемое в числе строк матрицы A.Как следует из подразд. 11.5, от указанного числа МНК-решение всейалгебраической системы уравнений не зависит. В связи с этим, с вычислительной точки зрения удобным оказывается добавление в очередной порции лишь по одному уравнению.
Тогда матрица A содержит всего однустроку, которую теперь обозначим как транспонированный вектор-столбецa, aT = (a1 , a2, . . . , an ). В этом случаеz = aT x + v,(11.30)и обработка наблюдений (11.28) принимает особенно простой вид:α = aT P̃ a+1,K = P̃ a/α,P̂ = P̃ −KaT P̃ ,x̂ = x̃+K(z −aT x̃). (11.31)Это алгоритм скалярной (последовательной) обработки.