Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 33
Текст из файла (страница 33)
. . ≥ λr > 0, а остальные λi = 0,i = r + 1, r + 2, .√. . , n. Это возможно, так как rank(AT A) =√rank A = r.Вычислим µi = λi для i = 1, 2, . . . , r. Заметим, что µi = λi = 0 дляi = r + 1, r + 2, . . . , n. Для µi , i = 1, 2, . . . , r снова запишем:AT Axi = λi xixTi AT Axi= 1,µ2i=⇒xTi (AT A)xi = λi = µ2i ,xTi AT Axi·= 1,µiµiAxi 6= 0.Обозначим yi = Axi/µi , i = 1, 2, . . . , r. ТогдаxTi AT AxjxTi λj xjxTi AT AxjT·===yi yj =µiµjµi µjµi µj1, i = j = 1, 2, .
. . , r,0, i =6 j.Следовательно, {y1 , y2, . . . , yr } — ортонормированы в Rm , так как yi ∈∈ Rm . Все yi — линейные комбинации столбцов матрицы A: yi ∈ R(A).Они могут быть дополнены в Rm до полного ортонормированногобазиса:{y1, y2, . . . , yr , yr+1, . . . , ym }.Таким образом, имеем:•216в Rm — ортонормированный базис из векторов {y1, y2 , . . . , ym },10.5 Отыскание псевдообратной матрицы•в Rn — ортонормированную систему векторов {x1, x2, .
. . , xn}.Соберем эти векторы в две матрицы:Q1 = [y1, y2 , . . . , ym ],Q2 = [x1, x2, . . . , xn].a) Имеем для i = 1, 2, . . . , r:yi µi = Axi,µi > 0,Axi 6= 0.(10.11)б) Имеем для i = r + 1, r + 2, . . . , m:kAxik2 = λi = 0=⇒Axi = 0.Поэтому в записиyi µi = Axi следует считать µi = 0, i = r + 1, r + 2, . . . , m. (10.12)Это означает, что (10.11) и (10.12) равносильны записи:D 0= A [x1, x2, . . . , xn][y1, y2 , . .
. , ym ]|{z} 0 0{z}|| {z }Q1Q2Σгде D = diag (µ1 , µ2, . . . , µr ), то естьQ1 Σ = AQ2 ,или Q1 ΣQT2 = A ,или Σ = QT1 AQ2 .Утверждение 1 теоремы доказано.2. Умножение на ортогональную матрицу не изменяет евклидовой нормывектора, поэтому имеем цепочку равенств:kAx−bk2 = kQ1ΣQT2 x−bk2 = kΣQT2 x−QT1 bk2 = kΣy−QT1 bk2 = kΣy−ck2.Введем новый векторy = QT2 x = Q−12 x.Для него kyk = kxk, так как QT2 — ортогональная матрица.Из приведенной цепочки равенств получаем следующий вывод: отыскание нормального псевдорешения системы Ax = b эквивалентно отысканию нормального псевдорешения системы Σy = c, где c = QT1 b. Имеемэти нормальные псевдорешения: ȳ0 и x̄0. Они равны:ȳ0 = Σ+ c,x̄0 = A+b.21710 Теоретические основыВ силу связи x и y имеемȳ0 = QT2 x̄0=⇒x̄0 = Q2 ȳ0 .Поэтомуx̄0 = Q2 Σ+c = Q2Σ+QT1 b = A+ b .Но b — любой вектор в Rm .
ОтсюдаA+ = Q2Σ+QT1 .234Пример 10.4. Дана матрица [13] A =. 3Тогда AT A = 3 4= 25 . Поэтому единственный собственный4вектор матрицы AT A есть: x1 = [1]. Это и есть матрица Q2 , Q2 = [1].Найдем λ1 :(AT A)x1 = λ1 x1=⇒Найдем µ1 :µ1 =Найдем вектор y1 :[25] [1] = λ1 [1]p1Ax1=y1 =µ15λ1 =34√=⇒λ1 = 25.25 = 5.1=3/54/5.y1 — есть первый столбец матрицы Q1 = [y1, y2]. Второй столбец y2 матрицы2Q1должен быть выбран так, чтобы быть в R ортогональным вектору y1 =3/5=.4/5y12В общем виде y2 =. Из условия ортогональности y1T y2 = 0 имеемy2234y12 + y22 = 0 .55Кроме того, y2 должен быть ортонормирован, то есть y2T y2 = 1.
Отсюдаследует22y12+ y22= 1.21810.5 Отыскание псевдообратной матрицыИз первого условия имеем y12 = −y224/3 . Подставляя значение y12 во второеусловие, выводим y22:222y22(4/3)2 + y22= 1 =⇒ y22[1 + (4/3)2] = 1 =⇒1192=⇒ y22===.1 + (4/3)21 + 16/9 25Можно выбрать одно из двух значений:(a) y22 = 3/5 либо(б) y22 = −3/5.Отсюда, соответственно,(a) y12 = (−3/5)(4/3) = −4/5 , либо(б) y12 = −(−3/5)(4/3) = 4/5 .Следовательно, возможны два варианта:Вариант 1Q1 =3/5 −4/54/5 3/5В варианте 1 имеемΣ = QT1 AQ2 =Вариант 2,3/5 4/5−4/5 3/5Q1 =341=3/5 4/54/5 −3/59/5 + 16/5−12/5 + 12/5.=50.В варианте 2 имеем 9/5 + 16/553/5 4/53 T1 ==.Σ = Q1 AQ2 =12/5 − 12/504/5 −3/54 5В любом варианте QT1 AQ2 = Σ, где Σ =, что и должно быть.0++Теперьнайдем(из примеров 10.1, 10.2 и 10.3): Σ+ = Σ и A .
Имеем= 1/5 | 0 . Поэтому ищем A+ = Q2Σ+QT1 .В варианте 1 3/5 4/5++ TA = Q2 Σ Q1 = 11/5 0= 3/25 4/25 .−4/5 3/521910 Теоретические основыВ варианте 2+A = Q2 Σ+QT1=Здесь также имеем+1A A=1/5 03/25 4/253/5 4/54/5 −3/534=1=3/25 4/25 .= I.Это произошло в силу того, что rank A = r = n = 1. Как мы знаем, в этомслучае−1 A+ = (AT A)−1AT = 253 4 = 3/25 4/25 .Вывод из примера 10.4. Пользоваться сингулярным разложениемA = Q1 ΣQT2 , чтобы находить A+, непросто. Для этого нужно:T1. Найти собственные векторы x1, x2, . . .
, xn матрицы A A.Тогда Q2 = x1 x2 . . . xn .2. Найти собственные значения матрицы AT A: λ1 > λ2 > . . . > λn ≥ 0.3. Оценить ранг r матрицы A, т. е. разграничить {λ1 > λ2 > . √.. >> λr > 0} и {λr+1 = λr+2 = . .
. = λn = 0} и образовать µi = λi ,D = diag (µ1 , µ2, . . . , µr ).√4. Найти yi = Axi/µi, где µi = λi , i = 1, 2, . . . , r.5. Доопределить каким-либо образом систему {yi } до ортонормированного базиса в Rm . 6. Образовать Q1 = y1 y2 . . . ym .В качестве альтернативного решения для нахождения A+ рассмотрим наследующем примере так называемое «усеченное» L̄Ū -разложение матрицы.Пример 10.5.
Дана1A=23матрица [13]25 ,rank A = 2 = n.7Выполним полное LU -разложение и затем «усеченное» L̄Ū -разложение: 1 21 21112 0 1 = 2 1 A=2 5=2 1= L̄Ū ,0 13 73 1 13 1 | {z }0 0|{z} | {z } | {z }ŪL220UL̄10.5 Отыскание псевдообратной матрицыrank L̄ = rank Ū = rank A = 2.В общем случае: пусть rank A = r, A = LU , тогда:1. L — всегда квадратная с единичной диагональю, т.
е. L = L(m, m) иdet L = 1,2. U — имеет (m − r) нулевых нижних строк.Отбросим последние (m − r) строк в U и последние (m − r) столбцов в L (вэтом и заключается переход к «усеченному» разложению), тогда получим:A = L̄Ū ,L̄ = L̄(m, r),Ū = Ū (r, n),rank L̄ = r = rank U = rank A.Свойства L̄ и Ū (так как A = L̄Ū и rank A = rank L̄ = r):1◦ Ū :R(Ū T ) = R(AT ) — совпадение пространств строк,2◦ L̄:R(L̄) = R(A) — совпадение пространств столбцов.Теорема 10.10 ( [13]).A+ = Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T .Доказательство.
Пусть b — произвольный вектор в Rm . Рассмотримвекторy = [ Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T ] b .Вектор y ∈ R(U T ) = R(AT ) . Умножим y слева на A = L̄Ū :Ay = L̄Ū Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T b = L̄(L̄T L̄)−1L̄T b .По определению, L̄(L̄T L̄)−1L̄T , есть матрица проектирования на R(L̄), т. е.на R(A), поэтомуA Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T b = p ,где p — проекция b на R(A) .Таким образом, вектор y удовлетворяет следующим условиям:1. y ∈ (AT ) .2. Ay = p, где p — проекция любого вектора b на R(A) .22110 Теоретические основы3. y = Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T b .Отсюда Ū T (Ū Ū T )−1(L̄T L̄)−1L̄T = A+.2Приведем заключительные примеры для лучшего понимания структурыи свойств псевдообратной матрицы A+ .Пример 10.6.
Если A = A(1, 1), то0, если A = 0,A+ =1/A, если A =6 0.Пример 10.7. Если A = diag (λ1, λ2 , . . . , λn ), то0,если λj = 0,++A+ = diag (λ+где λ+1 , λ2 , . . . , λn ),j =1/λj , если λj 6= 0.Пример 10.8. Матрицы4 0A1 =0 0иA2 =4 00 10−10почти одинаковы. Однако, их псевдообратные матрицы1/4 01/4 0++и A2 =A1 =0 10100 0сильно отличаются.Эти примеры, имеющие общее объяснение в теореме о сингулярном разложении, говорят о том, что операция псевдообращения — разрывная.
Действительно, функция1/λ, если λ 6= 0 ,λ+ =0, если λ = 0имеет разрыв в точке λ = 0.Отсутствие непрерывности операции псевдообращения приводит к серьезным вычислительным трудностям и, возможно, к большим вычислительнымошибкам при решении задачи МНК, особенно, если r < n или r < m.К счастью, на практике задача МНК обычно возникает при m > n = r,т. е. с матрицами A = A(m, n) максимального столбцового ранга для переопределенных систем Ax = z. В этом случае A+ = (AT A)−1AT и МНКрешение x̄ единственно и равно x̄0.22210.6 Основные теоремы по МНК и псевдоинверсии10.6Основные теоремы по МНК и псевдоинверсииТеорема 10.11 ( [1]).Пусть z ∈ Rm и A = A(m, n) .
Тогда:(1) найдется x̄0 ∈ Rn единственный, если это вектор с минимальной нормой, минимизирующий kz − Axk2 ,(2) вектор x̄0 является единственным вектором из R(AT ), удовлетворяющим уравнению Ax̄ = ẑ, где ẑ = p — проекция вектора z на R(A) .Теорема 10.12 ( [1]).Среди всех векторов x̄, минимизирующих2kz−Axk , вектор x̄0 , имеющий минимальную норму, является единственнымвектором видаx̄0 = AT y,удовлетворяющим уравнениюAT Ax̄0 = AT z.Иными словами, вектор x̄0 может быть получен с помощью любого вектора y0 , удовлетворяющего уравнениюAT AAT y = AT z,по формуле x̄0 = AT y.Теорема 10.13 ( [1]).
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
