Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 39
Текст из файла (страница 39)
е.нужна ковариация P̃ . Имея ξ˜ и P̃ и продолжая работать лишь с нормальными распределениями, вы выражаете эту априорную информацию следующей условной (априорно принятой) плотностью распределения вероятностей вектора x:hi−1/21nT−1˜ = (2π) |P̃ |˜ P̃ (ξ − ξ)˜ .exp − (ξ − ξ)fx|x̃ (ξ|ξ)(11.39)2Здесь вы формально предполагаете, что P̃ обратима (т. е. P̃ > 0), чтобыиметь возможность записать это выражение.25311 Оценивание по методу наименьших квадратовТеперь вы проводите свое наблюдение по схеме (11.38) и получаете конкретный его результат ζ, — реализованное значение случайного вектора z.Вы хотите использовать результат ζ, чтобы выработать наилучшую оценку˜ P̃ ) сx̂ для вектора x.
Как для этого скомбинировать априорное знание (ξ,полученным результатом ζ? Ответ на этот вопрос кроется в апостериорном˜ распределения вероятностейраспределении, т. е. в плотности fx|z,x̃(ξ ζ, ξ)неизвестного вектора x при условии, что известны: (1) результат ζ для вектора наблюдений z и (2) сообщенное значение ξ˜ априорного данного вектораx̃ с ковариацией P̃ , что отражено принятием априорной плотности (11.39).Если удастся найти эту fx|z,x̃(ξ ζ, ξ̃), то по ней можно будет принять решение, поскольку именно плотность распределения для любой случайнойвеличины служит наиболее полной характеристикой.Найдем fx|z,x̃(ξ ζ, ξ̃). По формуле Байеса для условных плотностей˜ = fx,z,x̃(ξ, ζ, ξ̃) .(11.40)fx|z,x̃(ξ ζ, ξ)˜fz,x̃(ζ, ξ)Преобразуем (11.40) по известным из теории вероятностей законам, опуская, для простоты промежуточных записей, аргументы.
Получаем˜ x|x̃ (ξ ξ)˜ ξ, ξ)ff(ζf·f·fx̃z|x,x̃z|x,x̃x|x̃˜ ==.(11.41)fx|z,x̃ (ξ ζ, ξ)˜fz|x̃ · fx̃fz|x̃(ζ ξ)˜ известна как априорная и задана выражеЗдесь плотность fx|x̃ (ξ ξ)нием (11.39). Остается найти другие две. Для их нахождения опираемсяна (11.38). Имеемz = Ax + v = Aξ + v = Aξ + v ,x=ξ ˜x̃=ξx̃=ξ˜так как v от ξ˜ не зависит. Отсюда и из (11.37) получаем1−1/2fz|x,x̃(ζ ξ, ξ̃) = (2π)m|R |exp − (ζ − Aξ)T R−1(ζ − Aξ) .
(11.42)2Для другой плотности, находящейся в знаменателе (11.41), снова рассмат˜ Найдем совместнуюриваем выражение (11.38), но теперь при условии x̃ = ξ.плотность векторов x и v при этом условии:˜˜˜˜ .˜ = fx,v,x̃ (ξ, ρ, ξ̃) = fv|x,x̃ (ρ ξ, ξ)fx|x̃(ξ ξ)fx̃(ξ) = fv (ρ)·fx|x̃(ξ|ξ)fx,v|x̃ (ξ, ρ ξ)˜˜fx̃(ξ)fx̃(ξ)(11.43)25411.9 Полная статистическая интерпретация рекурсии в МНКПоследнее равенство получено вследствие принятого свойства ошибок v вэксперименте (11.38): они не зависят ни от x, ни от x̃.
В выражении (11.43)перемножаются две гауссовы плотности: (11.37) и (11.39), поэтому результат — тоже гауссова плотность. Вектор z, (11.38), образован как взвешенная сумма двух векторов, x и v, совместная плотность которых гауссова.Поэтому плотность в знаменателе (11.41) — гауссова. Чтобы ее записать,достаточно найти по (11.38) первые два условные момента. Первый момент,учитывая (11.39), равенn on o˜˜˜E z x̃ = ξ = A E x x̃ = ξ + 0 = Aξ,n o˜так как E v x̃ = ξ = 0. Находим второй момент, используя (11.39) и(11.37):ionhT ˜˜˜= AP̃ AT + R, R = ImE (z − Aξ)(z − Aξ) x̃ = ξ(последнее — в силу нормализации, см.
замечание 11.4). Теперь имеем возможность записатьh i−1/2Tm˜exp{·},(11.44)fz|x̃(ζ ξ) = (2π) AP̃ A + Im где {·} =n− 21 (ζoTT−1˜˜− Aξ) (AP̃ A + Im ) (ζ − Aξ) .Все три плотности для (11.41) найдены: (11.42), (11.39) и (11.44). Подставляя их в (11.41), получаем" #−1/2R · P̃ 1exp − (α + β − γ) ,fx|z,x̃(ξ ζ, ξ̃) = (2π)n AP̃ AT + R 2T−1α = (ζ − Aξ) R (ζ − Aξ) ,T−1˜˜ξ−ξ ,β = ξ − ξ P̃T −1 T˜˜γ = ζ − AξAP̃ A + Rζ − Aξ .(11.45)(11.46)Замечание 11.5.
Здесь и далее для удобства (прослеживания местоположения матрицы ковариаций ошибки v) вместо Im , (11.36), сохраненообщее обозначение R, как и в (11.37), которое в любой момент можно заменить на Im , если будет принято предположение о нормализации наблюдений.25511 Оценивание по методу наименьших квадратовОчевидно, (11.45), (11.46) определяют нормальную плотность распределения, однако, ее явное определение требует двух действий: (1) приведение суммы трех квадратичных форм (11.46) к одной квадратичной формеи (2) приведение отношения трех определителей к одному определителю.Проведем эти преобразования.
При первом действии активно используемлемму 11.1 (см. п. 11.7), беря ее в виде:−1−1TT−1T −1AP̃(11.47)= P̃ − P̃ A AP̃ H + RP̃ + A R Hс обозначением Λ̃ = P̃ −1 для априорной информационной матрицы. Перепишем (11.47) в нескольких эквивалентных формах. Имеем−1−1TTTT−1T −1AP̃ AT =A = P̃ A − P̃ A AP̃ A + RP̃ + A R A−1 hi−1TTTTTT= P̃ A AP̃ A + RR.AP̃ A + R − AP̃ A = P̃ A AP̃ A + RОтсюдаP̃−1P̃−1T−1+A R A−1TA R−1=ATTAP̃ A + R−1.Умножая (11.47) слева и справа на P̃ −1 , получаем−1−1−1−1T−1T−1T −1A.P̃ = P̃ − A AP̃ A + RP̃P̃ + A R A(11.48)(11.49)Умножая (11.47) слева на A и справа на AT и обозначая C = AP̃ AT + R,находим−1−1T −1AT = AP̃ AT − AP̃ AT C −1AP̃ AT =A P̃ + A R A−1−1TR.= (C − R) − (C − R)C (C − R) = R − R AP̃ A + RОтсюда−1−1T −1−1T−1T −1.A R = R − AP̃ A + RR A P̃ + A R A−1(11.50)Вычисляя квадратичную форму в (11.45) как (α + β − γ), введем промежуточные обозначения:−1 , P̂ = P̃ −1 + AT R−1AT −1a = A R ζ,(11.51)b = P̃ −1 ξ˜ ,c = a + b.25611.9 Полная статистическая интерпретация рекурсии в МНКТеперь после раскрытия скобок в (11.46), приведения подобных членов в(α + β − γ) и при подстановках (11.48), (11.49) и (11.50) получимT−1T −1TTξ − P̂ c .
(11.52)α + β − γ = ξ P̂ ξ − 2ξ c + c P̂ c = ξ − P̂ c P̂Кроме того,P̂ c ====T −1−1 ˜P̂ A R ζ + P̃ ξ =T −1 ˜T −1 ˜T −1−1 ˜P̂ A R ζ + P̃ ξ + A R Aξ − A R Aξ =ihT −1−1T −1˜˜ζ − Aξ =P̂ P̃ + A R A ξ + A Rξ˜ + P̂ AT R−1 ζ − Aξ˜ .(11.53)Из (11.48) и (11.51) видно, чтоTP̂ A R−1= K = P̃ ATTAP̃ A + R−1(11.54)есть матрица Калмана, определенная выражением (11.32) (см. п.
11.7).В свою очередь, величина (11.53), входящая в квадратичную форму (11.52)плотности распределения (11.45), есть не что иное как среднее значение вектора x, обусловленное двумя событиями: z = ζ и x̃ = ξ˜ (первое — результатом ζ измерения и второе — результатом ξ˜ априорной оценки). Обозначимэто условное среднее как x̂, тогда (11.53) перепишется в виде:˜x̂ = ξ˜ + K(ζ − Aξ).(11.55)Кроме того, из (11.52) видно, что матрица P̂ в (11.51) есть ковариация апостериорного распределения (11.45), и по лемме (11.47) она дается выражением−1TP̂ = P̃ − P̃ AP̃ A + RAP̃ = P̃ − KAP̃ .(11.56)Сравнивая (11.54), (11.55) и (11.56) с подразд. 11.7, убеждаемся в том, чтоисходя из другой — чисто статистической задачи оценивания вектора x позашумленным наблюдениям его компонент через матрицу A (см. рис.
11.4),мы получаем полное алгебраическое совпадение с рекурсией для МНК встандартной ковариационной форме, если при этом принято во внимание,˜ а в качестве искомой наилучшей оценки x̂ для векторачто z = ζ и x̃ = ξ,x мы приняли апостериорное среднее значение (11.55), что соответствуеткритерию максимума апостериорной вероятности (МАВ) (11.40), (11.41),(11.45).25711 Оценивание по методу наименьших квадратовvxAAx+текущие данныеz=ζx̂˜max fx|z,x̃ (ξ|ζ, ξ)ξзадачарешениеx̃ = ξ̃AAξ˜+−νаприорные данныеP̃AAP̃+Kx̂апостериорные данныеKKAP̃P̂−+Рис.
11.5. Статистическая задача получения оценки x̂ по критерию МАВ и ее решение(полная статистическая интерпретация МНК-решения). Между моментами получения текущих данных апостериорные данные занимают место априорных данных, ипроцесс решения повторяетсяТаким образом, при нормальных законах распределения и независимыхошибках наблюдения v МНК-решение x̂ интерпретируется как МАВ-оценкавектора x (рис.
11.5). Именно МАВ-оценка отвечает на поставленный вышевопрос (см. между (11.39) и (11.40)): как наилучшим образом скомбиниро˜ P̃ ) с текущими данными ζ.вать априорные данные (ξ,Для решения этой задачи, как видно из изложенного, пришлось вычислить апостериорную плотность (11.41) и затем найти максимизирующий ееаргумент ξ, — он оказался равен (11.55). Однако потребовались сложныевычисления, которые пока еще не коснулись отмеченного выше приведенияотношения определителей в (11.45) к одному определителю.
Из условия нормировки плотности видно, что он должен оказаться равен |P̂ |, но это ещепредстоит доказать. Отложим это доказательство на окончание этого подразд. 11.9 (см. стр. 260), а сейчас рассмотрим, что будет, если изменить критерий, то есть ради простоты максимизировать по ξ не (11.41), а только однуиз плотностей в числителе (11.41), — именно, (11.42), что означает критериймаксимального правдоподобия, МП. Очевидно, такое решение получаетсямоментально: максимум (11.42) по ξ совпадает с минимумом квадратичнойформы в показателе экспоненты.