Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 34
Текст из файла (страница 34)
Для всякой матрицы A = A(m, n) существуетпсевдообратная матрица A+ = A+ (n, m), такая что для произвольного вектора z ∈ Rmx̄0 = A+ zявляется вектором с минимальной нормой среди всех векторов x̄, минимизирующих kz − Axk2.Теорема 10.14 ( [1]).
Для любой матрицы A = A(m, n) псевдообратная матрица A+ обладает свойствами:(1) A+ = (AT A)+AT ,(6) (AAT )+ = (AT )+A+,(2) (AT )+ = (A+)T ,(7) AA+ A = A,(3) A+ = AT (AAT )+ ,(8) A+ AA+ = A+ ,(4) (A+)+ = A,(9) AA+ = (AA+)T ,(5) (AT A)+ = A+ (AT )+,(10) A+ A = (A+A)T .22310 Теоретические основыОсновное правило обращения произведения матриц, (BA)−1 = A−1B −1, невыполняется для псевдообратных матриц, то есть(BA)+ 6= A+B + .Теорема 10.15 ( [1]).Для любой симметрической матрицы A == A(n, n) с действительными элементами предельная матрицаPA = lim(A + δI)−1 A = lim A(A + δI)−1δ→0δ→0существует.
Она является матрицей проектирования на R(A) = R(AT ). Этоозначает, что для любого вектора z ∈ Rn векторẑ = PA zявляется проекцией z на R(A) = R(AT ) — на пространство строк или столбцов данной симметрической матрицы.Теорема 10.16 ( [1]).(1) Для произвольных z ∈ Rm и A = A(m, n) вектор x̄ минимизирует критерий kz − Axk2 тогда и только тогда, когда x̄ имеет вид:x̄ = A+ z + (I − A+ A)yдля некоторого y ∈ Rn .(2) Вектор x̄, минимизирующий kz − Axk2, является единственным тогдаи только тогда, когда A+ A = I. Последнее условие выполняется тогдаи только тогда, когда только нулевой вектор составляет ядро (нульпространство) матрицы A, т.
е. при r = n.(3) Уравнение Ax = z имеет решение тогда и только тогда, когда A+ Az == z. Это выполняется в том и только том случае, когда z ∈ R(A).Вектор x является решением уравнения Ax = z тогда и только тогда,когда он задается в видеx = A+ z + (I − A+ A)yдля произвольного y ∈ Rn . Это решение единственно (и тогда оно равноA+z), если и только еслиAA+z = z224иA+A = I.10.7 Вычисление матриц проектирования10.7Вычисление матриц проектированияДана матрица A. Возьмем B = AT .
ИмеемR(AT ) — пространство строк матрицы A,R(AT )— пространство столбцов матрицы AT .R(B)Определение 10.8. Матрица проектирования P на R(B) есть такаяматрица, которая обладает свойстом:(b − P b) ⊥ R(B),где P b = b̂ — проекция вектора b на R(B), (b − P b) = b̃ — перпендикуляр кR(B).Запишемb̃ = b − P b = (I − P )b,где I − P — матрица проектирования на R⊥ (B) = N (B T ) = N (A).Из определения псевдообратной матрицы (см. определение 10.5) следует,что матрица проектирования P вектора b на R(A) в общем случае такова:P b = p,P b = Ax̄0 = AA+b=⇒P = AA+.Если взять B = AT и проектировать вектор b на R(B), то мы должнывзять матрицу проектирования в виде P = BB + . Но B + = (AT )+ = (A+)T ,следовательно P = AT (A+)T .Если A = A(m, n), то P имеет размеры: (n × m)(n × m)T = (n × n).
ЕслиB = B(n, m), то B + = B + (m, n), тогда P имеет размеры: (n × m)(m × n) == (n × n).Выводы:1. Если ищут матрицу вида P = I −AT (AAT )−1A, то это есть матрицапроектирования любого вектора на ядро, т. е. на нуль-пространствоN (A) матрицы A. Но в таком виде ее можно определить, толькоесли (AAT )−1 существует, т. е. если rank A = m (A имеет полныйстрочный ранг).2. Если rank A = r < m, что возможно иногда при m ≤ n, то (AAT )−1не существует.
В этом случае для этой же матрицы P справедливо22510 Теоретические основынаиболее общее выражение, а именно: P = I − AT (A+)T , где такжеможно иметь в виду, что всегда (A+)T = (AT )+. Это означает, чтонужно уметь отыскивать A+. Для этого есть различные, не оченьпростые, вычислительные методы (см. книгу [1], а также книгу [13]).Пример 10.9. Определение проектора PA = AA+ . Дано:1 −1 0 1A = A1 = 2 −2 −1 3 = A(n, m) = A(3, 4),−1 1 −2 1n = 3, m = 4, r = rank(A) = 2 < n = 3.Делаем LU -разложение:1 −1 0 11 0 0 −1 1 .A= 2 10 0 0 0−1 2 1{z}|{z}|LUНаходим усеченное LU -разложение:11−101A = L̄Ū = 2 1 .0 0 −1 1−1 2 |{z}| {z }ŪL̄Отсюда следует:10 14 −221 −10 −14 22 .A+ = Ū + L̄+ = Ū T (Ū Ū T )−1 (L̄T L̄)−1L̄T ={z} | {z } 150 |5 −8 −41 Ū +L̄+5 22 19Найдем матрицы проектирования как на пространство строк матрицы A,т. е. на R(AT ), так и на пространство столбцов матрицы A, т. е.
на R(A).1. Матрица проектирования на R(AT ):1 2 −110−1055 −1 −2 1 1 14 −14 −8 22 =PA = AT (A+)T = 0 −1 −2 150−22 22 −41 191 3 122610.7 Вычисление матриц проектирования60 −60 30 302 −2 1 11 −60 60 −30 −30 = 1 −2 2 −1 −1=150 30 −30 90 −60 5 1 −1 3 −230 −30 −60 901 −1 −2 33 −2 1 11 −2 3 −1 −1 .I − PA = I − AT (A+)T = 5 1 −1 2 −2 1 −1 −2 2PA проектирует на R(AT ) ,Здесь:I − PA проектирует на N (A) ..Найдем базис пространства N (A). Он образован следующими двумявекторами (так как n − r = 2): 1−11 0b1 = иb=20 1 .01Любой вектор в N (A) задается в виде: 1−11 + y2 0 , где y1 , y2 − числа.y1 0 1012.
В то же время матрица проектирования на R(A):10 14 −221 −1 0 11 + −10 −14 22 =PA = AA =2 −2 −1 3 5 −8 −41 150 −1 1 −2 15 22 1925 50 −255 10 −51 1 =50 130 10 =10 26 2 .15030−25 10 145−5 2 2925 −10 51 +I − PA = I − AA =−104 −2 .305 −2 122710 Теоретические основыЗдесь:PA проектирует на R(A) ,I − PA проектирует на N (AT ) = R(A)⊥ .Найдем базис пространства R(A):1v1 = 2 ,−10v2 = −1 .−2Эти два вектора оказались взаимно ортогональны:0v1T v2 = 1 2 −1 −1 = −2 + 2 = 0.−2Найдем базис пространства N (AT ).
Это пространство определяется каксовокупность векторов y, таких что: y T A = 0. Или иначе: y T v1 = 0 иy T v2 = 0. Отсюда найдем:5y = −2 , если выбрать y3 = 1.1Действительно,15 −2 1 2 = 0 и−105 −2 1 −1 = 0 .−2Базис пространства N (AT ) состоит из одного вектора. Возьмем в качестве базисного вектора v3 = 21 y, т. е. вектор2.5v3 = −1 .0.5В качестве произвольного вектора для проектирования возьмем (рис. 10.3) 61z = 0 = 6 0 .0022810.7 Вычисление матриц проектированияe321N (AT )v3y234z1v2R(A)1234v1 = ẑe2перпендикуляр на R(A)56e1Рис.
10.3. Проектирование вектора на линейное подпространствоНайдем его проекцию на R(A): 5 10 −5111 10 26 2 6 0 = 2 = v1 (совпало с v1).ẑ = PA z =30−5 2 290−1Найдем его проекцию на N (AT ): 156z̃ = z − ẑ = 0 − 2 = −2 = y (совпало с y).0−11Пример 10.10.A = A2 = n = 4,1 2 −1−1 −2 1 = A(n, m) = A(4, 3),0 −1 −2 1 3 1m = 3, r = rank(A) = 2 < n.По сравнению с примером 10.9 данная матрица совпадает с транспонированной матрицей примера 10.9.
Используем свойство (A+)T = (AT )+ для22910 Теоретические основынахождения псевдообратной матрицы для нашего примера:10 −105 51 A+ =14 −14 −8 22 .150−22 22 −41 19Отсюда, найдем матрицу проектирования на пространство строк матрицыA, т. е. на R(AT ):10 14 −221 −1 0 11 T+ T −10 −14 22 =PA = A (A ) = 2 −2 −1 3 5 −8 −41 150 −1 1 −2 15 22 1925 50 −255 10 −51 1 =50 130 10 =10 26 2 .15030−25 10 145−5 2 29Следовательно,25 −10 51 −104 −2 .I − PA = I − AT (A+)T =305 −2 1Здесь:PA проектирует на пространство строк матрицы A, R(AT ) ,I − PA проектирует на нуль-пространство матрицы A, N (A) .Геометрическая интерпретация.
Для примера 10.9 построить картинкуневозможно (так как там b ∈ R4), а для примера 10.10 — можно (в этом случае b ∈ R3 ). Существует связь матриц двух последних примеров, а именно:A1 = AT2 . Соответственно, имеемR(AT1 ) = R(A2 )N (A1) = N (AT2 )R(AT2 ) = R(A1)N (A2) = N (AT1 )Из примера 10.9 имеем матрицу P1 , которая проектирует на N (A1), и матрицу I −P1 , которая проектирует на R(AT1 ). Из примера 10.10 имеем матрицуP2 , которая проектирует на N (AT1 ), и матрицу I − P2 , которая проектируетна R(A1) (см. рис.
10.2 на стр. 213).23010.8 Рекурсия в задаче МНК10.8Рекурсия в задаче МНКПостановка вопроса. Даны матрица A = A(m, n) и вектор z ∈ Rm .Требуется найти единственный вектор x̄0 с минимальной нормой, минимизирующий kz − Axk2.Можно ли искать его последовательно?1◦ Известно, что Ax̄0 = ẑ, где ẑ ∈ R(A), z − ẑ ⊥ R(A), x̄0 ∈ R(AT ).Нормальное псевдорешение x̄0 несовместной системы Ax = z удовлетворяет равенствуx̄0 = A+z ,где A+ — псевдообратная матрица.A1A22◦ Произвольно расщепим матрицу A на блоки, A =и, соответz1ственно, вектор z на подвекторы, z =, так чтобы A1 = A1 (k, n),z2A2 = A2(s, n), k + s = m, z1 ∈ Rk , z2 ∈ Rs .3◦ На первом шаге найдем нормальное псевдорешение только первой системы Ax1 = z1 .
Оно равно x̃0 = A+1 z1 . Проекция вектора z1 на R(A1 ):ẑ1 = A1x̃0 = (A1A+1 )z1 .4◦ Рассмотрим системуA1ẑ1x=,A2z2которая отличается от исходной системыA1z1x=A2z2(10.13)(10.14)тем, что вместо z1 используется ẑ1 . Для нее нормальное псевдорешениеобозначим x̂0 . Оно равно+ A1ẑ1x̂0 =.A2z2Вопрос: верно ли равенство x̂0 = x̄0? Напомним, что+ A1z1x̄0 =.A2z223110 Теоретические основыРешение вопроса:Для системы (10.13) запишем нормальные уравнения: TAẑ11A1 AT2x = AT1 AT2,A2z2 TA1 A1 + AT2 A2 x = AT1 ẑ1 + AT2 z2 .(10.15)Затем для системы (10.14) запишем нормальные уравнения: TAz11A1 AT2x = AT1 AT2,A2z2 TA1 A1 + AT2 A2 x = AT1 z1 + AT2 z2 .(10.16)Сравним в правой части AT1 ẑ1 и AT1 z1 . ИмеемAT1 ẑ1 = AT1 (A1A+1 )z1 .Более простое доказательство:AT (AA+) = AT (AA+)T = (AA+A)T = AT .Лемма 10.2. Докажем, чтоAT = AT (AA+) .Доказательство.
Имеем (см. стр. 220):A = L̄Ū ,AT = Ū T L̄T ,A+ = Ū +L̄+ ,AA+ = L̄Ū Ū + L̄+,AT (AA+) = Ū T L̄T · L̄Ū · Ū +L̄+ = Ū T L̄T = AT .Лемма доказана. Поэтому AT1 ẑ1 = AT1 z1 .Докажем то же самое другим способом:AT1 (ẑ1) = AT1 (z − z̃) = AT1 z − AT1 z̃,где z̃ ∈ N (AT1 ),т. е. AT1 z̃ = 0 .Поэтому AT1 ẑ1 = AT1 z1 . Таким образом, правые части уравнений (10.15)и (10.16) совпадают. Поэтому решения уравнений (10.15) и (10.16) одинаковы, то естьx̂0 = x̄0 .23210.9 Основные свойства симметрических / эрмитовых матриц10.9Основные свойства симметрических / эрмитовыхматрицСвойства симметрических матриц, A = AT , в теории МНК служат основой многих результатов.