Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 36
Текст из файла (страница 36)
См. TotalLeast Squares в [18, 38, 100].4. Неоднородность погрешностей в исходных данных [49]. Это обобщениеприсутствует в разд. 13, 14.5. Зависимость ковариации погрешностей от неизвестного параметра(Generalized Least Squares) [106]. Это обобщение здесь отсутствует. Возможен переход к адаптивному оцениванию [73].6.
Зависимость от неизвестного параметра не только ковариации погрешностей, но и других элементов стохастической модели источника данных. Это обобщение здесь отсутствует. Возможен переход к адаптивному оцениванию [73].7. Улучшение вычислительных схем [15, 89, 93, 95, 103, 104, 107, 108, 119,123, 125, 135, 129]. Это обобщение присутствует в разд. 13, 14 и составляет их основное содержание.Вычислительные схемы МНК имеют важнейшее значение для любыхобобщений и любых приложений. Неиспользование эффективных, т. е. численно устойчивых и экономичных алгоритмов может увести исследователя всторону надуманных проблем. И наоборот, — затруднения в численной реализации могут быть с успехом преодолены, если знать и умело использовать23711 Оценивание по методу наименьших квадратовэффективные схемы вычислений.
Кроме того, как и любое знание, наука вычислений несет с собой важные, плодотворные идеи, которые полезны самипо себе. Применительно к задаче МНК, эти базовые идеи следующие: рекурсия алгоритма, скаляризация обработки векторных данных, факторизацияматриц (ковариационных и информационных), ортогонализация матриц, ивзаимная инверсность (двойственность форм) вычислительных алгоритмов.11.2Линейная задача наименьших квадратовВо многих приложениях, связанных с обработкой экспериментальныхданных, необходимо отыскивать такой вектор x ∈ Rn , линейные комбинации компонент которого, Ax, где A = A(m, n) — матрица размера (m × n),как можно более близки или, еще лучше, равны данным значениям, образующим вектор z ∈ Rm , т.
е. Ax ≈ z. Если мерой близости двух векторовсчитать квадрат евклидовой нормы разностного вектора, в данном случае,вектора v , z−Ax, то указанная задача есть линейная задача о наименьшихквадратах (см. также подразд. 10.4).Возможность сделать равным нулю вектор невязок v = z−Ax существуеттогда и только тогда, когда z ∈ R(A), где R(A) — пространство столбцовматрицы A. В этом случае имеем совместную систему уравнений Ax = z.Однако, z — вектор наблюдений, то есть экспериментальных данных и A —матрица, которую задают до того, как получат z и которую в различныхприложениях называют либо матрицей регрессоров, либо матрицей наблюдений, либо матрицей плана эксперимента. Совсем не обязательно, что условие z ∈ R(A) будет выполнено, например, из-за случайных погрешностей vво время регистрации экспериментальных данных.
Тогдаz = Ax + v,(11.1)и решение по методу наименьших квадратов (для краткости, МНКрешение) есть вектор x̄, доставляющий минимум функционалу качества (см.рис. 11.1):TJ(x) = (z − Ax) (z − Ax) =mXv(j)2 = v T v.(11.2)j=1Требуя минимум этого критерия, для искомого x̄ получаем так называемые нормальные уравнения (см. подразд. 10.4, стр. 209):AT Ax̄ = AT z.238(11.3)11.2 Линейная задача наименьших квадратоввход Aвыход zОбъектМодельx+Ax−v = z − Axmin kvk2xРис. 11.1. Линейная задача наименьших квадратовИх решение всегда существует (обе части равенства (11.3) принадлежат одному и тому же пространству R(AT ) столбцов матрицы AT ), но может бытьне единственным (если rank A < n).
В последнем случае из всех x̄ выбираютто единственное, x̄0 , которое имеет минимальную норму kx̄0k. Этот векторназывают нормальным псевдорешением. Известно (см. подразд. 10.4), чтоx̄0 = A+z,(11.4)где A+ — псевдообратная матрица к A. Как уже отмечалось, в качестве определения A+ применяют различные формулировки. Здесь для этогоиспользуем геометрический подход: A+ есть такая матрица в выражении (11.4), что для любого z ∈ Rm вектор x̄0 ∈ Rn удовлетворяет двумусловиям:Ax̄0 = p, p ∈ R(A), z − p ⊥ R(A).(11.5)x̄0 ∈ R(AT ).(11.6)Условие (11.5) требует, чтобы x̄0 отвечал совместной системе Ax̄0 = p, гдеp — проекция вектора z на R(A), а условие (11.6) требует, чтобы этот x̄0был взят из пространства R(AT ) строк матрицы A.
Условие (11.5), такимобразом, выбирает x̄0 = x̄, чтобы минимизировать функционал (11.2), аусловие (11.6) среди всех таких x̄ выбирает единственный x̄0 с минимальнойнормой.Часто матрицу A выбирают так, чтобы она имела полный столбцовыйранг, rank A = n. В этом случае m ≥ n, x̄ единственно и равно x̄0, A+ == (AT A)−1AT иx̄0 = (AT A)−1AT z.(11.8)23911 Оценивание по методу наименьших квадратовОднако иногда такое условие не выполняется, и тогда x̄0 = A+ z, где A+ —псевдообратная матрица (см. подразд.
10.4–10.5, стр. 209–212).11.3Статистическая интерпретацияПредположим, что вектор ошибок v в уравнении (11.1) образован из случайных величин с нулевым средним и известной матрицей ковариации(11.9)E {v} = 0,E vv T = Pv ,где E {·} — оператор математического ожидания (среднего) над ·, и Pv — ПО(положительно определенная) матрица. Найдем квадратно-корневое разложение Pv = SS T (например, разложение Холесского). Если теперь умножитьвектор z (11.1) на S −1, то данные z̄ = S −1z получают представлениеz̄ = Āx + v̄(11.10)с матрицей Ā = S −1A и ошибками v̄ = S −1v. Этот вектор ошибок всегдаимеет единичную ковариацию:E v̄v̄ T = E S −1vv T S −T = S −1 E vv T S −T = S −1SS T S −T = Im ,где Im — единичная матрица размера (m × m).
Вследствие этого данныеz̄ называют нормализованными экспериментальными данными. Значениепредставления (11.10) заключается в том, что оно демонстрирует, как сконструировать вектор некоррелированных между собой измерений с единичной дисперсией из вектора, элементы которого произвольно взаимно коррелированы (декоррелировать и нормализовать его). Ниже предполагаем, чтоданные z (11.1) уже декоррелированы и нормализованы, так что(11.11)E {v} = 0,E vv T = Im,где Im — единичная матрица размера (m × m). При этом из (11.3) находимAT Ax̄ = AT z = AT Ax + AT v,AT A(x̄ − x) = AT v.Отсюда, если det(AT A) 6= 0, имеемE {x̄} = x,240(AT A) E (x̄ − x)(x̄ − x)T (AT A) = AT E vv T A = AT A.(11.12)(11.13)11.4 Включение априорных статистических данныхСоотношение (11.12) выражает собой свойство несмещенности решения(оценки) x̄ относительно неизвестного (постоянного) вектора x, измеряемогов виде экспериментальных данных z, (11.1) или (11.10).
Соотношение (11.13)дает выражение для ковариации оценки x̄ в виде(11.14)Px̄ = E (x̄ − x)(x̄ − x)T = (AT A)−1.при определении x̄ по нормализованным экспериментальным данным.Обратная матрица Px̄−1 от ковариации Px̄ называется информационнойматрицей. Ее обозначение будет Λx̄ или просто Λ. При использовании нормализованных данных она равна AT A, а в более общем случае (11.9) онаравна Λ = AT Pv−1 A.11.4Включение априорных статистических данныхПредположим, что в добавление к линейной системе (11.1) мы имеемаприорную несмещенную оценку неизвестного вектора x в виде x̃ и соответствующую априорную информационную матрицу Λ̃.
Это означает, чтоE {x̃} = x и(11.15)Λ̃−1 = E (x̄ − x)(x̄ − x)T = P̃ ,где P̃ — ковариация оценки x̃. Найдем какой-нибудь квадратный корень Λ̃1/2из матрицы Λ̃, например, по одному из разложений Холесского (см. разд. 6):Λ̃ = Λ̃1/2Λ̃T /2 = R̃T R̃,где Λ̃1/2 = R̃T . Образуем вектор ṽ = (Λ̃1/2)T (x̃ − x) = R̃(x̃ − x). Этот векторимеет смысл нормализованной ошибки для априорной оценки x̃ вектора x.Действительно, его ковариация равна единичной матрице размера (n × n):E ṽṽ T = R̃ E (x̃ − x)(x̃ − x)T R̃T = Λ̃T /2Λ̃−1 Λ̃1/2 = In .Так как о векторе x, кроме экспериментальных данных z, (11.1), известнааприорная оценка x̃ с ковариацией P̃ = Λ̃−1 , эту информацию целесообразновключить в контекст задачи о наименьших квадратах, рассматривая модифицированный функционал качества J1(x) = ṽ T ṽ + v T v вместо (11.2).
Онсоединяет в себе квадрат нормы нормализованной ошибки (невязки) априорной оценкиṽ = R̃(x̃ − x) = Λ̃T /2(x̃ − x),24111 Оценивание по методу наименьших квадратовAОбъектz+v = z − Ax−МодельxAxmin(kṽk2 + kvk2 )xR̃xR̃Априорнаямодель x̃−z̃ = R̃x̃+ṽ = z̃ − R̃xРис. 11.2. Включение априорных данных в линейную задачу НКс квадратом нормы нормализованной ошибки (невязки) экспериментальныхданных v = z − Ax. Так какJ1(x) = (x̃ − x)T Λ̃(x̃ − x) + (z − Ax)T (z − Ax) == (z̃ − R̃x)T (z̃ − R̃x) + (z − Ax)T (z − Ax),(11.16)где z̃ = R̃x̃, то J1(x) может быть интерпретирован просто как критерий качества метода наименьших квадратов применительно к расширенной системе(рис.