Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 40
Текст из файла (страница 40)
При R = I (что уже обсуждалось как при25811.9 Полная статистическая интерпретация рекурсии в МНКнятое неограничительное условие нормализации зашумленных наблюдений)имеемx̂ = arg min (ζ − Aξ)T (ζ − Aξ) = x̂МНК = x̂МП ,ξто есть x̂ есть МНК-решение, x̂МНК, системы (11.1), переписанной здесь ввиде ζ = Aξ + v с R = I, и одновременно, это x̂ есть оценка максимумаправдоподобия, x̂МП.Тем самым, получена полная статистическая интерпретация МНКрешения с точки зрения теории оценок:1. Если в МНК-решении учтена априорная информация, то полученноерешение совпадает с оценкой МАВ, x̂МАВ (11.55). Эту оценку называюттакже байесовской, поскольку ее вывод основан на формуле Байеса(11.40). Любой алгоритм (т. е. инструмент) вычисления оценки называют оценивателем.
Байесовский оцениватель в линейной задаче оценивания с нормальными законами распределения случайных величиндается формулами (11.55) и (11.56). Он совпадает с формулами фильтра Калмана — на этапе обновления оценок по измерениям, и он имеетнаглядную схемную интерпретацию на рис. 11.5.2. Если в МНК-решении не учтена априорная информация, то полученное решение имеет смысл оценки МП, x̂МП. Пренебрежение априорнойинформацией может быть вызвано тем, что она очень мала (недостоверна). Этот факт выражается простой записью: P̃ = ε−2I при ε → 0,как уже отмечалось на стр. 248.3. Поскольку недостоверность или скудность априорной информации выражается большим ростом диагональных элементов матрицы P̃ по сравнению с внедиагональными элементами: P̃ −1 → 0 при P̃ −1 = ε2 I, ε → 0,то P̂ −1 → AT A, что формально видно из (11.51) при R = I.
Оцениватель, который пренебрегает априорной информацией, называют фишеровским оценивателем в связи с именем R. A. Fisher — известногостатистика. В этом оценивателе−1 Tx̂МАВ → AT AA ζ = x̂МНК,(11.57)что следует из преобразований выражения (11.53) для (11.55) прииспользовании (11.51) в пределе при P̃ −1 → 0 (и с неограничительнымусловием нормализации наблюдений R = I), — если, конечно, матрица25911 Оценивание по методу наименьших квадратовAT A обратима. Если же информационная матрица AT A полученных измерений z ≡ ζ не обратима, тогда вместо (11.57) предельноесоотношение записывается в видеx̂МАВ → A+ζ + I − A+A ξ˜ = x̂МНК,(11.58)что вытекает из теоремы 10.16 (см.
подразд. 10.1) и означает сходимостьx̂МАВ к общему решению нормальных уравнений МНК.Замечание 11.6. В (11.58) первое слагаемое есть нормальное псевдорешение при минимизации kz − Axk2, где z = ζ, а второе — проекциялюбого ξ˜ на N (A).Теперь, в завершение необходимых действий по приведению (11.45)к стандартному виду плотности для нормального закона распределения,кроме (11.52), докажем, что R · P̃ = P̂ ,(11.59)AP̃ AT + R где P̂ определено формулой (11.56). Для этого введем вспомогательную матрицуTP̃P̃AP∗ =.AP̃ AT + RAP̃Выполним ее блочное U L-разложение, где U — верхняя треугольная матрица с положительно определенными блоками на диагонали, а L — нижняятреугольная матрица с единичной диагональю. ПолучаемTP̃−KAP̃P̃AI0P∗ =.KT I0AP̃ AT + RЭто позволяет найти определитель для P ∗ как произведение определителейдиагональных блоков.
С учетом (11.56), имеем ∗ P = P̂ · AP̃ AT + R .(11.60)С другой стороны, воспользуемся следующей леммой обращения симметричной положительно определенной матрицы [115].Лемма 11.2. При симметричных X1 > 0, X2 > 0 и произвольнойTX21 = X21справедливо равенство −1−1 X10X1 X12I −X1−1 X12−1 ×=TX21 X20I0X2 − X12X1−1 X12I0×.T−X12X1−1 I26011.10 Основные результатыДоказательство. Проведите самостоятельно прямой проверкой.2Беря определитель в этом выражении, получаем X1 X12 −1−1T−1 = |X1 |−1 · X2 − X12 .XX121 X T X2 12Применим этот результат к P ∗ , обозначая X1 = P̃ , X2 = AP̃ AT + R, X12 == P̃ A.
Находим ∗ P = P̃ · R .(11.61)Сопоставляя (11.60) и (11.61), получаем (11.59).11.102Основные результаты1. Критерий наименьших квадратовmin (z − Ax)T (z − Ax)x(11.62)для решения алгебраических линейных систем общего видаz = Ax + v(11.63)не является изначально статистическим критерием. Решение x̄, отвечающее этому критерию, дается нормальными уравнениямиAT A x̄ = AT z(11.64)и имеет вид суммы взаимно ортогональных векторов:x̄ = A+z + I − A+ A y(11.65)для некоторого y. Это решение единственно тогда и только тогда, когдаA+A = I, что означает, что только нулевой вектор составляет ядроN (A) матрицы A (т. е.
A имеет полный столбцовый ранг). Нормальное псевдорешение x̄0 есть единственный вектор из (11.65), имеющийминимальную норму, т. е. x̄0 = A+ z.2. Статистическая интерпретация задачи (11.62) означает следующее:(а) x есть постоянный вектор (случайный или нет), который необходимо оценить по результатам наблюдения z в виде (11.63).26111 Оценивание по методу наименьших квадратов(б) A есть матрица наблюдений, которая показывает, какие линейныекомбинации элементов вектора x включены, по условиям задачи, ввектор наблюдения z.(в) При этом v есть случайная ошибка наблюдения, ковариация которой R задана как положительно определенная матрица.
Однако,если задача (11.62), (11.63) чисто алгебраическая, R = I.(г) При статистической интерпретации, матрица AT A в (11.64) называется информационной матрицей Λ. Она обратима тогда и толькотогда, когда rank A = n (матрица A имеет полный столбцовыйранг). Обратная к Λ матрица называется ковариационной матрицей ошибки оценивания вектора x, P = Λ−1 , и она характеризуетстепень неопределенности в решении задачи оценки (в скалярномслучае она равна дисперсии).3. Статистическая задача оценивания, указанная в предыдущем пункте,приводит к таким же (алгебраически эквивалентным) результатам, какрешение алгебраической задачи (11.62), (11.63).4.
Решение задачи (11.62), (11.63), так же как и решение ее статистического эквивалента, может быть получено в рекуррентной форме, удобной с вычислительной точки зрения. Возможны две базовые эквивалентные рекуррентные формы (взаимно инверсные): информационнаясхема оценивания (неявная относительно апостериорной оценки x̂) иковариационная схема оценивания (явная относительно x̂).5. Оцениватель максимального правдоподобия, выбирает оценку x̂МП так,чтобы максимизировать функцию правдоподобия (11.42). При гауссовых независимых ошибках она является решением тех же самых нормальных уравненийAT R−1Ax̄ = AT R−1ζ ,(11.66)которые отвечают критерию взвешенных наименьших квадратовmin (z − Ax)T R−1 (z − Ax),xz = ζ.(11.67)По смыслу, x̂МП есть значение неизвестного вектора x, которое с наибольшей вероятностью обеспечивает то событие, что при случайномнаблюдении (11.63) результат окажется равным ζ (x̂МП обеспечиваетнаибольшую вероятность события z = ζ).
Эта оценка не учитываетникакой априорной информации о векторе x. С другой стороны, оценка,26211.10 Основные результатыудовлетворяющая «взвешенному» критерию (11.67) и даваемая, соответственно, решением «взвешенных» нормальных уравнений (11.66),называется марковской оценкой в связи с теоремой Гаусса–Маркова,которая утверждает, что эта оценка x̂МП , удовлетворяющая уравнениям (11.66), — эффективная в классе всех линейных несмещенныхоценок, см., например, [49, 17].6. Оцениватель максимума апостериорной вероятности выбираетоценку x̂МАВ , в отличие от x̂МП, принципиально в рекуррентной форме,т. е. учитывает априорную информацию и наилучшим образом комбинирует ее с текущими наблюдениями.
Оценка x̂МАВ вычисляется в явномвиде в алгоритме Калмана (11.54), (11.55), (11.56) (см. подразд. 11.7)или в неявном виде в информационном алгоритме (см. подразд. 11.6).При недостаточной (т. е. крайне неопределенной, скудной) априорнойинформации, оценка x̂МАВ будет такой же, как x̂МП, или, что то жесамое, x̂МНК в задаче (11.67) с решением из нормальной системы (11.66).7. Рекурсия (рекуррентная схема вычислений) в алгоритмах МНК иоптимального оценивания является краеугольным камнем теории ипрактики эффективных вычислительных алгоритмов, рассматриваемых ниже в разд.
13 и 14.Следующий разд. 12, занимающий промежуточное положение, рассматривает случай, когда рекуррентная схема вычислений в задачах МНК неиспользуется, т. е. когда нормальные уравнения решаются для всех накопленных экспериментальных данных целиком — после того, как все они накоплены и сохранены для обработки.26312Одновременное решениенормальных уравнений12.1Метод нормальных уравненийНормальные уравнения (11.3) показаны выше (см. стр. 238). Довольночасто их используют для отыскания МНК-решения x̄ [97].
Этот подход противоположен последовательным методам, берущим начало от идеи расщепления исходной переопределенной системы на априорную и текущую части(см. подразд. 11.5). Поэтому иногда его называют «одновременным решением» всех нормальных уравнений, характеризующих всю исходную переопределенную систему Ax ≈ z (см. стр. 238 и выражение (11.1)).Алгоритм использует полный векторнаблюдений z ∈ Rm и всю матрицу плана эксперимента A ∈ Rm×n , имеющуюrank(A) = n. По матрице A и вектору z алгоритм вычисляет решение x̄задачи наименьших квадратов, доставляющее min kz − Axk2.
Алгоритм:Метод Нормальных Уравнений.1.2.3.4.Вычисляет нижний треугольник матрицы Λ = AT A.Вычисляет d = AT z.Вычисляет разложение Холесского Λ = SS T .Решает сначала систему Sy = d и затем систему S T x̄ = y.Здесь в качестве разложения Холесского может быть взято либо нижнеетреугольное разложение (S = L), либо верхнее треугольное разложение (приS = U ) — см.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
