Téma: Řešení soustav lineárních algebraických rovnic Zpracoval Doc. RNDr. Zdeněk Hlaváč, CSc
Mějme soustavu lineárních algebraických rovnic tvaru Ax = b .
(1)
Matice soustavy A je typu (m, n), vektor pravých stran b je dimenze m, vektor řešení x je dimenze n.
1
Frobeniova věta a její důsledky
Základním tvrzením pro řešení soustavy (1) je tzv. Frobeniova věta: Řešení soustavy (1) existuje právě když platí h(A) = h([A , b]) . Poznámka: Symbolem h(A) je míněna hodnost matice A. Matice [A , b] je tzv. rozšířená matice soustavy. Vzniká z původní matice soustavy připojením vektoru pravých stran na místo (n + 1)−tého sloupce. Přiřaďme k soustavě (1) tzv. homogenní soustavu s nulovým vektorem pravých stran, tedy soustavu Ax = O .
(2)
Pro řešení této soustavy platí následující základní tvrzení: Nechť h(A) = h(≤ n). Potom množina řešení soustavy (2), na kterou pohlížíme jako na podmnožinu Euklidovského prostoru Rn , je podprostor tohoto prostoru dimenze n − h. Poznámka: V předchozím tvrzení vzpomínaný podprostor Rn−h je tzv. nulový prostor matice A a číslo η = n − h, jakožto dimenze tohoto prostoru, se nazývá nulitou této matice. Jestliže tedy h(A) = n (tedy matice soustavy má alespoň tolik řádek jako sloupců a je navíc regulární), je nulita matice rovna nule a jejím nulovým prostorem je Euklidovský prostor dimenze nula skládající se pouze z nulového vektoru. Za této situace soustava (2) má pouze triviální řešení x = O. Jak již bylo řečeno, může tento případ nastat pouze pro m ≥ n (pro m > n se jedná o tzv. přeurčenou soustavu rovnic). Jestliže platí m < n (tzv. nedourčená soustava rovnic) je vždy nulita matice soustavy kladná a množina řešení je vždy Euklidovský prostor kladné dimenze (tedy nekonečná množina). Pro η > 0 zřejmě existuje η bázových vektorů nulového prostoru matice A. Označme (j) je x(1) , . . . , x(η) , kde x(j) = [xi ] ; j = 1, . . . , η , i = 1, . . . , n . Taková řešení soustavy (2) nazýváme prvky fundamentálního systému řešení této soustavy. Každé řešení diskutované soustavy lze potom vyjádřit ve tvaru lineární kombinace prvků fundamentálního systému jako x=
η X
αi x(i) .
i=1
Toto vyjádření lze rovněž přepsat do tvaru
x = Xα , 1
(3)
kde X je tzv. fundamentální matice typu (n, η), do níž jsou po sloupcích uloženy jednotlivé prvky fundamentálního systému, tedy X = [x(1) , . . . , x(η) ] . Vektor α dimenze η je vektor koeficientů v lineární kombinaci. Bázové vektory x(i) lze volit tak, že vybereme h lineárně nezávislých sloupců matice A (jež jistě existují, protože h(A) = h) a na i−té místo ze zbylých míst, kterých je právě η, píšeme jedničku a na ostatní místa nuly. Provedeme-li výběr báze tímto způsobem, je fundamentální matice X ortonormální. Platí pro ni tedy X T X = Eη .
(4)
Protože soustava (1) je lineární, lze její obecné řešení vyjádřit ve tvaru x = x0 + xp , kde x0 je libovolné řešení homogenní soustavy (2) a xp je jediné, tzv. partikulární, řešení soustavy (1), pokud ale existuje. Podle (3) pak pro obecné řešení soustavy (1) platí x = Xα + xp ,
(5)
kde X lze volit tak, že platí (4) a α je libovolný vektor dimenze η. Diskutujme nyní tvrzení Frobeniovy věty s ohledem na typ matice soustavy a výsledek (5). Mohou nastat dva základní případy podle druhu ”obdélníkovosti” matice A a uvnitř nich vždy dva podpřípady podle její regularity. 1. Případ m ≤ n (pro m < n nedourčená soustava) (a) h = h(A) = m = min(m, n) - tedy regulární matice. Přidáním sloupce pravých stran se hodnost matice zřejmě nemůže změnit. Podle Frobeniovy věty pak řešení (1) existuje a má tvar (5). Řešení je η−parametrické množství, kde η je nulita matice soustavy, rovna n − m. Speciálně pro čtvercovou matici, kdy je m = n a tedy η = 0, je toto řešení jediné. (b) h = h(A) < m = min(m, n) - tedy singulární matice. Přidáním sloupce pravých stran hodnost matice může zůstat stejná (jestliže vektor pravých stran je lineární kombinací nezávislých sloupců matice soustavy) nebo může o jedničku vzrůst (v ostatních případech). V prvním případě podle Frobeniovy věty řešení existuje a podle (5) je jich (n − h)−parametrické množství. Je totiž n − h > 0 vzhledem k platnosti nerovnic h < m ≤ n. Nekonečné množství řešení existuje i v případě čtvercové matice soustavy. Ve druhém případě podle Frobeniovy věty řešení neexistuje. 2. Případ m > n (přeurčená soustava) (a) h = h(A) = n = min(m, n) - tedy regulární matice. Přidáním sloupce pravých stran hodnost matice opět, v případě sloupce pravých stran ve tvaru lineární kombinace (nyní všech) sloupců matice soustavy, může zůstat stejná nebo vzrůst o jedničku. V prvním případě podle Frobeniovy věty řešení existuje a protože η = n − h = 0, je podle (5) jediné. Ve druhém případě podle Frobeniovy věty řešení neexistuje. (b) h = h(A) < n = min(m, n) - tedy singulární matice. pak situace je zcela analogická jako v případě singulární matice nedourčené soustavy rovnic. 2
Z předchozího rozboru je patrno, že v některých (přesně specifikovaných) případech řešení soustavy (1) existuje právě jedno, v některých případech neexistuje řešení žádné, v některých případech jich existuje nekonečně mnoho. Otázkou je, zda v případě neexistence řešení nelze najít jakési (v jistém smyslu zobecněné) pseudořešení. Dále má smysl se ptát, zda v případě existence nekonečně mnoha řešení nelze z nich vybrat jediné, které by bylo nějakou vlastností významné. Ukazuje se, že na obě otázky je kladná odpověď. Bližším upřesněním těchto otázek se budeme zabývat ve zbytku tohoto tématu.
2
Pseudoinverzní matice k regulární matici
Definice: Nechť A je matice typu (m, n). Matice AL typu (n, m) se nazývá zleva pseudoinverzní k A, jestliže platí AL A = E n .
(6)
Matice AP typu (n, m) se nazývá zprava pseudoinverzní k A, jestliže platí AAP = E m .
(7)
Pro takto definované matice platí následující tvrzení: Nechť B je libovolná matice typu (n, m) taková, že BA je regulární matice (čtvercová, řádu n). Potom matice AL = (BA)−1 B
(8)
je k A zleva pseudoinverzní. Nechť dále matice B je libovolná matice typu (m, n), že AB je regulární matice (čtvercová řádu m). Potom matice AP = B(AB)−1
(9)
je k matici A zprava pseudoinverzní. Ověření tohoto tvrzení je zřejmé. Matice (8) totiž splňuje podmínku (6) a matice (9) splňuje podmínku (7). Pro existenci pseudoinverzních matic platí následující důležité tvrzení: Nutnou a postačující podmínkou existence zleva pseudoinverzní matice je podmínka h = h(A) = n = = min(m, n). Nutnou a postačující podmínkou existence zprava pseudoinverzní matice je podmínka h = h(A) = m = min(m, n). Poznámka: První případ předchozího tvrzení je případ regulární matice nedourčené soustavy, druhý případ je případ regulární matice přeurčené soustavy. Oběma případům společný je případ regulární čtvercové matice. Je-li splněna podmínka regularity, pak zřejmě pro obdélníkovou matici s větším počtem řádků existuje nekonečně mnoho zleva pseudoinverzních matic a žádná matice zprava pseudoinverzní. Každá ze zleva pseudoinverzních matic se získá prostřednictvím vztahu (8) pro vhodně volenou matici B. K regulární čtvercové matici existuje právě jedna zleva pseudoinverzní matice. Tato matice je pro libovolnou (vhodnou) matici B v (8) rovna inverzní matici A−1 . Jestliže totiž má být BA regulární, musí (věta o součinu determinantů) být i B regulární. V tom případě je AL = (BA)−1 B = A−1 B −1 B = A−1 . Poznámka: Předchozí úpravu nelze provést pro obdélníkové (byť regulární) matice. Inverzní matice k jednotlivým obdélníkovým činitelům totiž neexistuje. 3
Příklad: Matice A = [1; 2; 3]T je regulární matice typu (3, 1). Existuje k ní tedy nekonečně mnoho zleva pseudoinverzních matic a žádná matice zprava pseudoinverzní. Například pro matici B = [1; 0; 0] dostaneme BA = [1] (matice typu (1, 1)). Potom (BA)−1 = [1] a tedy podle (8) AL = B = [1; 0; 0]. Pro matici Bh= [1; 1;i 1] dostaneme stejným postupem jinou zleva pseudoinverzní matici AL = 61 B = 16 ; 61 ; 61 . Pro matici B = [4; 5; 6] dostaneme h
i
1 5 3 další zleva pseudoinverzní matici AL = 32 B = 81 ; 32 ; 16 atd. Analogicky pro regulární obdélníkovou matici s větším počtem sloupců existuje nekonečně mnoho zprava pseudoinverzních matic a žádná matice zleva pseudoinverzní. Každá ze zprava pseudoinverzních matic se získá aplikací výrazu (9) pro vhodně zvolenou matici B. K regulární čtvercové matici existuje právě jedna zprava pseudoinverzní matice. Tato matice je pro libovolnou (vhodnou) matici B v definici (9) rovna inverzní matici A−1 . Platí totiž
AP = B(AB)−1 = BB −1 A−1 = A−1 . Ze Sylvestrovy věty o hodnosti součinu matic vyplývá, že hodnost matice A je stejná jako hodnost matice AAT i hodnost matice AT A. Odtud plyne, že pro regulární matici (na jejím typu nezáleží) lze v definicích (8), (9) za matici B dosazovat matici AT . Získáme tím právě jednu významnou zleva i zprava pseudoinverzní matici (pokud samozřejmě vůbec nějaká existuje). Označme tuto speciální zleva pseudoinverzní matici A− (= AL pro B = AT ; m ≥ n). Podobně speciální zprava pseudoinverzní matici označme A+ (= AP pro B = AT ; m ≤ n).
3
Pseudoinverzní matice a soustavy rovnic
V dalším ukážeme souvislost takto definovaných pseudoinverzních matic s řešením soustavy lineárních algebraických rovnic. Mějme soustavu (1) s regulární maticí soustavy. Případ čtvercové matice je nezajímavý, protože jediné řešení takové soustavy je x = A−1 b. Ostatní případy rozdělíme do dvou skupin podle typu matice soustavy. 1. Nedourčená soustava (m < n). Víme, že v takovém případě existuje nekonečně mnoho ((n − m)−parametrické množství) řešení. Zároveň víme, že k matici A existuje nekonečně mnoho zprava pseudoinverzních matic. Má potom smysl psát x = AP b .
(10)
Každé takto vyjádřené x je řešením soustavy (1). Ověření se provede snadno dosazením (10) do (1) při použití definice (9). Máme totiž b = Ax = AAP b = AB(AB)−1 b = E n b = b . Poznámka: Dokonce lze dokázat (podrobnosti vypustíme), že každé řešení z výše zmíněného (n − m)−parametrického množství řešení lze vyjádřit ve tvaru (10) pro vhodně volenou matici B. Speciální tvar řešení dostaneme, klademe-li AP = A+ . Ukazuje se totiž, že řešení xp = A+ b má ze všech řešení x = AP b nejmenší Euklidovu normu. Pro ověření tohoto tvrzení si představíme libovolné řešení soustavy ve tvaru (5) a partikulárním řešením A+ b. Je tedy 4
x = Xα + A+ b ,
(11)
kde X je ortonormální fundamentální matice a α vektor konstant. Nyní stačí dokázat, že kvadrát Euklidovy normy vektoru (11) nabývá svého minima pro α = O. Užitím základních vlastností Euklidovy normy, distributivního zákona pro matice a věty o transpozici součinu matic postupně dostáváme
||x||2 = xT x = (Xα + A+ b)T (Xα + A+ b) = [αT X T + bT (A+ )T ](Xα + A+ b) =
= αT X T Xα + bT (A+ )T Xα + αT X T A+ b + bT (A+ )T A+ b .
(12)
Sloupce fundamentální matice xi však jsou řešeními soustavy (2). Proto je Axi = O m ; i = 1, . . . , η . Zapíšeme-li tuto rovnici v maticovém tvaru pro všechna i najednou, dostaneme AX = O m,η . Použitím vztahu pro transpozici součinu matic a definice matice A+ postupně máme (A+ )T X = [AT (AAT )−1 ]T X = (AAT )−T AX = = (AAT )−1 AX = (AAT )−1 O m,η = O m,η .
(13)
Transpozicí tohoto vztahu vzniká [(A+ )T X]T = X T A+ = O η,m .
(14)
Protože druhý a třetí sčítanec v poslední rovnosti výrazu (12) jsou podle (13) a (14) bilineární formy (provedené na vektory b a α) s nulovými maticemi, jsou nulové. Proto (12) dává ||x||2 = αT X T Xα + bT (A+ )T A+ b . První sčítanec tohoto výrazu je podle (4) skalární součin αT α. Proto je 2
||x|| =
η X
αi2 + bT (A+ )T A+ b .
i=1
T
+ T
+
Ke vztahu b (A ) A b se tudíž přidává nezáporné číslo. V závislosti na vektoru α je proto ||x||2 minimální právě když α = O η . Poněvadž norma vektoru je nezáporné číslo, plyne z minima jejího kvadrátu minimum i její první mocniny. Tvrzení je tím ověřeno. Poznámka: V místech, kde je potřeba rozlišení, znamená O a nulový vektor dimenze a, zatímco O a,b nulovou matici typu (a, b). 5
2. Přeurčená soustava (m > n) Víme, že v takovém případě neexistuje žádné řešení (pro vektor pravých stran, jenž není lineární kombinací sloupců matice A) nebo existuje právě jediné řešení. Přesto lze formálně napsat ˜ = A− b = (AT A)−1 AT b . x
(15)
˜ , když řešení neexistuje? Ukazuje se, že se Co v takovém případě vyjadřuje vektor x jedná o vektor, který na celém prostoru Rn minimalizuje Euklidovu normu reziduálního vektoru tvaru Ay − b. Platí tedy ||A˜ x − b|| = ||AA− b − b|| = minn ||Ay − b|| . y ∈R
(16)
Poznámka: Matice AA− není obecně jednotkovou maticí, protože je pseudoinverzní zleva a zde je jí násobeno zprava. Za účelem ověření tvrzení (16) budeme metodami matematické analýzy extremalizovat funkci f (y) = ||b − Ay||2 . Postupně zřejmě platí f (y) = (b − Ay)T (b − Ay) = (bT − y T AT )(b − Ay) = = bT b − y T AT b − bT Ay + y T AT Ay .
(17)
Protože všechny sčítance v této rovnici jsou skaláry, transponováním se nemění. Protože však (y T AT b)T = bT Ay , je možno (17) přepsat na tvar f (y) = bT b − 2y T AT b + y T AT Ay .
(18)
Tento vztah je zřejmě kvadratickou funkcí proměnné y. Sčítanci tohoto výrazu tvoří po řadě absolutní, lineární a kvadratický člen této funkce. Pro její gradient potom platí gradf (y) = −2AT b + 2AT Ay .
(19)
Nutnou podmínkou extrému vzpomínané funkce je podmínka gradf (y) = O, která dává rovnici AT b = AT Ay . Protože matice AT A je regulární, dostáváme pro stacionární bod funkce f (y) vztah y = (AT A)−1 AT b . 6
Podle definice matice A− odtud ˜. y = A− b = x Dokažme ještě, že nalezený stacionární bod je skutečně minimem. Hessova matice funkce f (y) je podle (19) zřejmě konstantní maticí tvaru H f (y) = 2AT A . Matice AT A je, stejně jako matice AAT , pro libovolnou matici A positivně semidefinitní. Pro regulární matici A je dokonce positivně definitní. Plyne to např. ze vztahu xT AT Ax = (Ax)T Ax = ||Ax||2 ≥ 0 , takže kvadratická funkce vlevo, a tím i zkoumaná matice, je positivně definitní vzhledem k základní vlastnosti normy. Je-li x = O, je i Ax = O a ||Ax||2 = 0, takže xT AT Ax = 0. Jestliže naopak ||Ax||2 = 0, je zřejmě z vlastnosti normy i Ax = O a vzhledem k regularitě matice A i x = O . Tím je ověřena positivní definitnost Hessovy matice. Stacionární bod je proto lokálním minimem a tvrzení je ověřeno. Poznámka: Výrazem (15) bylo tedy získáno jakési ”pseudořešení”, které, v případě, že skutečné řešení neexistuje, minimalizuje Euklidovu normu reziduálního vektoru. V tomto smyslu se jedná o jakési nejlepší přiblížení se ke splnění rovnice (1) Jestliže skutečné řešení přesto existuje (tedy když vektor pravých stran soustavy je lineární kombinací sloupců matice soustavy), vše, co bylo výše řečeno, zůstává v platnosti. Pouze minimální hodnota Euklidovy normy reziduálního vektoru je nulová a popsané ”pseudořešení” je řešením. Levá a pravá pseudoinverzní matice k matici A mají ještě následující vlastnosti: 1. Pravá pseudoinverzní matice AP k matici A příslušející k (vhodné) matici B v (9) je transponovanou levou pseudoinverzní maticí (AT )L k matici AT příslušející k matici B T v definici (8) a naopak. Platí tedy
B(AB)−1 = [(B T AT )−1 B T ]T , (BA)−1 B = [B T (AT B T )−1 ]T , pro každou takovou matici B, pro níž (BA)−1 resp. (AB)−1 existuje k regulární matici A. Ověření je okamžitě patrno užitím vztahu pro transpozici součinu matic. Speciálně pro matice A− a A+ platí A− = [(AT )+ ]T ; A+ = [(AT )− ]T .
(20)
2. Pro každou zprava pseudoinverzní matici AP k matici A (pokud existuje) platí vztahy AAP A = A ; AP AAP = AP . 7
(21)
Podobně pro každou zleva pseudoinverzní matici AL k matici A platí vztahy AAL A = A ; AL AAL = AL .
(22)
Ověření této vlastnosti je zřejmé. Podle definice levé a pravé pseudoinverze totiž platí AAP = E resp. AL A = E. Odtud vhodným uzávorkováním součinů na levých stranách v (21) a (22) tvrzení bezprostředně plyne. Speciálně platí tyto vztahy i pro matice A− a A+ , pokud tyto matice existují. 3. Matice A+ A i AA− jsou symetrické. Po dosazení z jejich definice totiž platí
(A+ A)T = [AT (AAT )−1 A]T = AT (AAT )−1 A = A+ A , (AA− )T = [A(AT A)−1 AT ]T = A(AT A)−1 AT = AA− .
(23)
Poznámka: Obecně matice AP A ani matice AAL symetrická být nemusí, pokud není symetrická matice B v definici (8) a (9).
4
Obecná pseudoinverzní matice
Definujme nyní pseudoinverzní matici ještě obecněji. Tato definice bude mít platnost i pro singulární (obdélníkové) matice. Nechť A je matice typu (m, n) hodnosti h ≤ min(m, n). Matice X, která je řešením soustavy rovnic AXA = A ; XAX = X ; (XA)T = XA ; (AX)T = AX ,
(24)
se nazývá pseudoinverzní matice k matici A. Existence pseudoinverzní matice ve smyslu výše popsané definice se dokazuje přes jiné velice hluboké tvrzení teorie matic, a sice přes větu o singulárním rozkladu matice (singular value decomposition - SVD). Formulujeme ji jako následující tvrzení: Nechť A je matice typu (m, n) hodnosti h ≤ min(m, n). Potom existují ortogonální čtvercové matice U řádu m a V řádu n a diagonální matice S řádu h s kladnými prvky na diagonále tak, že je A = U ΣV
T
=U
"
S O m−h,h
O h,n−h O m−h,n−h
#
VT .
(25)
Poznámky: 1. Pokud některá z nulových matic ve vztahu (25) má některý typ nulový, vynechává se. 2. Ukazuje se, že diagonální prvky matice S, tzv. singulární čísla matice A jsou nenulovými vlastními čísly matice AT A. Tato matice má podle Sylvestrovy věty o hodnosti součinu matic hodnost h. Má tedy (n − h)−násobné nulové vlastní číslo. Matice V souvisí s vlastními vektory této matice a matice U souvisí s vlastními vektory matice AAT uloženými po sloupcích. Pokud jsou uvedené matice jednoduché struktury, jsou to přímo matice vlastních vektorů výše uvedených matic. Blíže se těmito otázkami zabývat nebudeme. 8
Užijeme nyní singulárního rozkladu pro vyjádření řešení X soustavy rovnic (24). Pro praktické použití bývají matematické knihovny procedur pro použití na číslicové výpočetní technice vybaveny procedurou pro singulární rozklad matice obvykle pod identifikátorem SVD. Základem je následující tvrzení: Nechť (25) je singulární rozklad matice A. Potom matice X se singulárním rozkladem X=V
"
diag
1 , . . . , s1h s1
O n−h,h
#
O h,m−h O n−h,m−h
UT ,
(26)
kde S = diag(s1 , . . . , sh ) splňuje rovnice (24) a tedy je její pseudoinverzní maticí. Ověření tohoto tvrzení provedeme dosazením (25) a (26) do jednotlivých rovnic v (24). Po příslušných úpravách vždy obdržíme identity. Pro první rovnici (24) máme AXA = U
·
"
diag
1 , . . . , s1h s1
O n−h,h
"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
O h,m−h O n−h,m−h
#
T
U U
"
#
V TV ·
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#
VT .
(27)
Protože matice U i V jsou ortogonální (ortonormální), platí U T U = E m a V T V = E n . Poněvadž platí "
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#"
diag
1 , . . . , s1h s1
O n−h,h
O h,m−h O n−h,m−h
#
=
"
Eh O m−h,h
#
O h,m−h O m−h,m−h
a dále "
Eh
O h,m−h O m−h,m−h
O m−h,h
#"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#
=
"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#
,
upravíme pravou stranu výrazu (27) na tvar pravé strany výrazu (25). Poněvadž se rovnají též levé strany uvedených rovnic, dostáváme odtud první rovnici (24). Zcela analogicky bychom ověřili rovněž druhou rovnici tohoto výrazu. Pro levou stranu třetí rovnice (25), podle vztahu pro transpozici součinu matic, máme výraz
T
(XA) = V
"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#T
T
U U
"
diag
1 , . . . , s1h s1
O n−h,h
O h,m−h O n−h,m−h
Protože pro ortogonální matici U platí U T U = E m a " "
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
diag
1 , . . . , s1h s1
O n−h,h
O h,m−h O n−h,m−h
#T #T
=
"
diag(s1 , . . . , sh ) O h,m−h O n−h,h O n−h,m−h
=
"
diag
je 9
1 , . . . , s1h s1
O m−h,h
O h,n−h O m−h,n−h
#
,
#
,
#T
VT .
T
(XA) = V
"
diag(s1 , . . . , sh ) O h,m−h O n−h,h O n−h,m−h
#"
diag
1 , . . . , s1h s1
O m−h,h
O h,n−h O m−h,n−h
#
V T . (28)
Pro pravou stranu třetí rovnice v (24) vychází
XA = V
"
diag
1 , . . . , s1h s1
O n−h,h
O h,m−h O n−h,m−h
#
T
U U
"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
#
VT .
Protože je U T U = E m , dostáváme odtud pravou stranu rovnice (28). Tím je ověřena třetí rovnice (24). Zcela analogicky bychom ověřili též rovnici čtvrtou. Poznámka: Blokovým roznásobením prostředních činitelů ve výrazu (28) dostaneme T
(XA) = XA = V
"
Eh O h,n−h O n−h,h O n−h,n−h
#
VT .
V předchozím tvrzení byla ověřena existence pseudoinverzní matice definované vztahy (24). Její jednoznačnost je obsahem následujícího tvrzení: Pseudoinverzní matice X řešící soustavu rovnic (24) je jediná. Ověření tohoto tvrzení provedeme sporem. Nechť tedy existují dvě matice X 1 a X 2 splňující soustavu (24). Z první rovnice dostaneme při vhodném uzávorkování A = AX 1 A = A(X 1 A) , odkud transpozicí máme AT = (X 1 A)T AT . S ohledem na třetí rovnici (24) odtud AT = X 1 AAT . Podobně bychom pro řešení X 2 dostali vztah AT = X 2 AAT . Odečtením posledních dvou rovnic vzniká O = (X 1 − X 2 )AAT . Dále roznásobením X 1 AAT = X 2 AAT . Přenásobme tento vztah zprava maticí X T1 . Dostaneme X 1 AAT X T1 = X 2 AAT X T1 . Aplikací vztahu pro transpozici součinu matic odtud dostáváme X 1 A(X 1 A)T = X 2 A(X 1 A)T . Podle třetí rovnice (24) však (pozor, každá strana má jiné uzávorkování) 10
X 1 AX 1 A = X 2 AX 1 A . Podle první rovnice (24) však dále X 1A = X 2A .
(29)
Upravme nyní rozdíl X 1 − X 2 . Podle druhé rovnice (24) je X 1 − X 2 = X 1 AX 1 − X 2 AX 2 . Podle čtvrté rovnice (24) však dále X 1 − X 2 = X 1 (AX 1 )T − X 2 (AX 2 )T . Užitím vztahu pro transpozici součinu a vztahu pro roznásobení matic je X 1 − X 2 = X 1 X T1 AT − X 2 X T2 AT = (X 1 X T1 − X 2 X T2 )AT . Přenásobme levou a pravou stranu této soustavy zprava maticí (X 1 − X 2 )T . Dostaneme (X 1 − X 2 )(X 1 − X 2 )T = (X 1 X T1 − X 2 X T2 )AT (X 1 − X 2 )T . Užitím vztahu pro transpozici součinu matic odtud dostaneme (X 1 − X 2 )(X 1 − X 2 )T = (X 1 X T1 − X 2 X T2 )[(X 1 − X 2 )A]T . Podle (29) však lomená závorka na pravé straně je nulovou maticí. Odtud plyne (X 1 − X 2 )(X 1 − X 2 )T = O . Levá strana této rovnice znamená, že X 1 − X 2 je ortogonální matice (protože mimodiagonální prvky této matice jsou nulové) s nulovým kvadrátem Euklidovy normy všech řádků (protože i diagonální prvky této matice jsou nulové. Každý řádek matice X 1 − X 2 má tedy nulovou Euklidovu normu. Je tedy sám nulový. Proto je X 1 − X 2 nulová matice a tudíž X 1 = X 2 . Jednoznačnost řešení soustavy maticových rovnic (24) je tím ověřena. V předchozích dvou tvrzeních byla dokázána existence a jednoznačnost řešení soustavy maticových rovnic (24). Jedinou existující matici X splňující tyto rovnice nazýváme (obecnou) pseudoinverzní maticí k matici A a značíme ji AM . Poznámka: Maticové rovnice (24) poprvé formuloval Moore. Existenci a jednoznačnost jejich řešení poprvé dokázal Penrose. Proto se takto definované pseudoinverzní matici říká pseudoinverzní matice Mooreova-Penroseova typu. Výše bylo dokázáno, že pokud existuje matice A− , má vlastnosti matice X v (24). První dvě vlastnosti jsou ve vztahu (22) pro AL = A− , třetí vlastnost je triviální, protože A− A = E n a čtvrtá vlastnost je dokázána ve druhé rovnici (23). Rovněž bylo dokázáno, že pokud existuje matice A+ , má rovněž vlastnosti matice X ve vztahu (24). První dvě vlastnosti jsou opět dokázány v (21) pro AP = A+ , čtvrtá vlastnost je triviální, poněvadž platí AA+ = E m a třetí vlastnost je dokázána v první rovnici (23). Odtud a z věty o existenci matic A+ a A− plyne, že pro regulární matici A je AM = A+ = AT (AAT )−1 pro m ≤ n , AM = A− = (AT A)−1 AT pro m ≥ n . 11
(30)
Pseudoinverzní matici Moorova-Penroseova typu lze pro regulární matici A počítat s výhodou podle (30) bez použití singulárního rozkladu. V případě singulární matice A je nutno pro výpočet matice AM použít procedury singulárního rozkladu. Toto platí bez ohledu na typ matice A. Jestliže matice A je typu (m, n) hodnosti h ≤ min(m, n), má matice AM následující jednoduché vlastnosti: 1. AM je typu (n, m) hodnosti h; 2. (AM )T = (AT )M a je to matice typu (m, n) hodnosti h; 3. (AM )M = A. Ověření těchto vlastností provedeme přes singulární rozklad matice A. Jestliže A=U
"
diag(s1 , . . . , sh ) O h,n−h O m−h,h O m−h,n−h
"
diag
#
VT,
(31)
pak je A
M
=V
1 , . . . , s1h s1
O n−h,h
O h,m−h O n−h,m−h
#
UT .
Protože U je řádu m a V je řádu n, plyne odtud ihned první tvrzení. Užitím vztahu pro transpozici součinu matic zřejmě dostaneme
(AM )T = U
"
diag
= U
"
diag
O n−h,h
O h,m−h O n−h,m−h
#T
O h,n−h O m−h,n−h
#
1 , . . . , s1h s1 1 , . . . , s1h s1
O m−h,h
VT =
(32)
VT.
Transpozicí (31) dostaneme T
A =V
"
diag (s1 , . . . , sh ) O h,m−h O n−h,h O n−h,m−h
#
UT .
Jestliže na tento vztah provedeme konstrukci pseudoinverze matice, dostáváme T M
(A )
=U
"
diag
1 , . . . , s1h s1
O m−h,h
O h,n−h O m−h,n−h
#
VT,
což je vztah totožný se (32). Druhé tvrzení je tím ověřeno. Třetí tvrzení je zřejmé z definice (obecné) pseudoinverzní matice.
5
Obecná pseudoinverze a soustavy rovnic
Pro soustavu rovnic (1) napišme nyní formálně vztah x0 = AM b . Pro takto určený vektor x0 platí následující tvrzení: 12
(33)
Pro libovolný vektor x ∈ Rn , x 6= x0 je ||Ax0 − b|| ≤ ||Ax − b|| .
(34)
Jestliže dále existuje vektor x 6= x0 takový, že v rovnici (34) platí znaménko rovnosti, pak pro všechna taková x je (35)
||x0 || ≤ ||x|| ,
kde || . . . || označuje Euklidovu normu vektoru. Abychom mohli toto tvrzení kvalifikovaně ověřit, odvoďme důležité pomocné tvrzení: (AAM )2 = AAM ; (AM A)2 = (AM A) .
(36)
Ověření tohoto tvrzení je velmi jednoduché. Je totiž (AAM )2 = AAM AAM = A(AM AAM ) = AAM podle druhé rovnice (24). Dále platí (AM A)2 = AM AAM A = AM (AAM A) = AM A podle první rovnice v (24). Pomocné tvrzení je tím ověřeno. Ze vztahu (36) plyne, že matice AAM , chápána jako lineární operátor z Rm do Rm , je projektor (nebudeme se blíže zabývat otázkou, na který podprostor prostoru Rm tento projektor promítá, není to pro náš účel důležité). Protože pro libovolný vektor x ∈ Rm je x = E m x = AAM x + (E m − AAM )x , je zřejmě matice E m − AAM projektorem na ortogonální doplněk podprostoru, na který promítá operátor AAM v prostoru Rm . Pro libovolné vektory y ∈ Rm a z ∈ Rm jsou tedy vektory u = AAM y a v = (E m − AAm )z na sebe kolmé. Jejich skalární součin (u, v) je tedy nulový. Platí proto
(37)
(38)
AAM y, (E m − AAM )z = 0 .
Zcela analogicky bychom dokázali podobné tvrzení pro operátor AM A, jenž promítá na podprostor prostoru Rn . Pro libovolná y ∈ Rn a z ∈ Rn tak je
AM Ay, (E n − AM A)z = 0 .
Nyní máme připravenu situaci pro ověření vlastnosti (34) pro vektor x0 v (33). Vezmeme libovolný vektor x ∈ Rn a upravíme reziduální vektor Ax − b. Přičtením a odečtením vektoru AAM b uvážením první vlastnosti pseudoinverzní matice v (24) postupně dostaneme Ax − b = Ax − AAM b + AAM b − b = AAM Ax − AAM b + AAM b − E m b = = AAM (Ax − b) + (AAM − E m )b ,
(39)
kde jsme sdružili první dva a poslední dva sčítance levé strany poslední rovnosti a použili distributivní zákon pro matice. Upravíme nyní vztah pro kvadrát Euklidovy normy 13
reziduálního vektoru za použití vztahu (39) a základní vlastnosti skalárního součinu ||z||2 = (z, z), platné pro libovolný vektor z libovolného prostoru se skalárním součinem. Je potom ||Ax − b||2 = (Ax − b, Ax − b) =
= AAM (Ax − b) − (E m − AAm )b, AAM (Ax − b) − (E m − AAm )b . Protože skalární součin je lineární v obou činitelích, lze poslední výraz dále upravit, čímž získáváme
||Ax − b||2 = AAM (Ax − b), AAM (Ax − b) − AAM (Ax − b), (E m − AAM )b −
− (E m − AAM )b, AAM (Ax − b) + (E m − AAM )b, (E m − AAM )b . Vzhledem ke komutativitě skalárního součinu (nad reálnými čísly) a ke vztahu souvislosti skalárního součinu s normou máme dále
||Ax−b||2 = ||AAM (Ax−b)||2 +||(E m −AAM )b||2 −2 AAM (Ax − b), (E m − AAM )b . Podle (37) v němž položíme y = Ax − b a z = b je třetí sčítanec posledního výrazu nulový. Dostáváme odtud konečný tvar ||Ax − b||2 = ||AAM (Ax − b)||2 + ||(E m − AAM )b||2 .
(40)
Odtud plyne zanedbáním (nezáporného) prvního sčítance ||Ax − b||2 ≥ ||(E m − AAM )b||2 . Roznásobením matic na pravá straně této nerovnice a zavedením vektoru x0 ze (33) odtud dostáváme ||Ax − b||2 ≥ ||Ax0 − b||2 . Tím je (34) ověřeno, protože norma je nezáporné číslo a tudíž lze beze změny znaménka nerovnice přejít od kvadrátů k prvním mocninám. Pro ověření další vlastnosti vektoru x0 se ptejme, kdy v (34) platí znamení rovnosti. Zřejmě právě tehdy, když zanedbaný sčítanec v (40) je nulový. Protože norma je nulová právě když příslušný vektor je nulový, je rovnost v (34) ekvivalentní výrazu AAM (Ax − b) = O . Roznásobením levé strany vzniká AAM Ax = AAM b . Podle první rovnice v (24) odtud Ax = AAM b . Přenásobme tuto rovnost zleva maticí AM . S ohledem na druhou rovnici v (24) odtud máme 14
AM Ax = AM b .
(41)
Pro případ rovnosti v (34), neboli pro případ platnosti (41) nyní upravme výraz pro libovolný vektor x ∈ Rn a posléze pro jeho normu. Je totiž x = E n x + AM b − AM b . Podle (41) dále x = E n x + AM b − AM Ax = (E n − AM A)x + AM b . Jestliže podle druhé rovnice v (24) píšeme AM b = AM AAM b, dostáváme x = (E n − AM A)x + AM AAM b .
(42)
Položíme-li nyní v (38) z = x a y = AM b vidíme, že skalární součin sčítanců na pravé straně předchozí rovnice je nulový. Je tedy
(E n − AM A)x, AM AAM b = 0 .
(43)
Pro kvadrát normy libovolného vektoru x s ohledem na linearitu skalárního součinu zřejmě platí podle (42)
||x||2 = (x, x) = (E n − AM A)x + AM AAM b, (E n − AM A)x + AM AAM b =
= (E n − AM A)x, (E n − AM A)x + 2 (E n − AM A)x, AM AAM b + +(AM AAM b, AM AAM b) .
S ohledem na (43) při využití souvislosti skalárního součinu s normou odtud dostáváme konečný vztah ||x||2 = ||E n − AM A)x||2 + ||AM AAM b||2 . Zanedbáním nezáporného prvního sčítance odtud při použití druhé rovnice v (24) dostaneme ||x||2 ≥ ||AM AAM b||2 = ||AM b||2 . Podle definice vektoru x0 v (33) odtud ||x||2 ≥ ||x0 ||2 . Protože norma je nezáporné číslo, plyne odtud vztah (35). Obě základní vlastnosti vektoru x0 jsou tím ověřeny. Vektor x0 v (33) se nazývá zobecněným řešením soustavy (1). Tento vektor je tedy buď ostře nejlepším ”řešením” soustavy (1) ve smyslu minima Euklidovy normy reziduálního vektoru Ax − b, nebo existuje více nejlepších ”řešení” popsaného typu. V takovém případě má vektor x0 mezi nimi minimální Euklidovu normu. Jestliže soustava (1) má klasické řešení (někdy říkáme též, že je konzistentní), pak toto řešení x zřejmě splňuje podmínku ||Ax − b|| = 0. Pokud je jediné, platí x = x0 a pro všechny ostatní vektory y ∈ Rn je 15
||Ay − b|| > 0 . Pokud je klasických řešení více (η−parametrické množství), je x0 jedním z nich a sice to, které má minimální (svoji) Euklidovu normu. Poznámka: Vztah (16) je slabší analogií vztahu (34). Vztah (16) byl však dokazován za předpokladu regularity matice A. Kdybychom chtěli důkaz (34) provést stejným způsobem, narazíme na problém řešitelnosti soustavy rovnic tvaru AT Ay = AT b , vyplývající z nulovosti gradientu funkce, která definuje kvadrát Euklidovy normy reziduálního vektoru. Matice AT A totiž obecně regulární není, takže exaktní řešení předchozí rovnice nemusí existovat. Navíc postačující podmínka extrému popisované funkce dává singulární (semidefinitní) Hessovu matici. Pro ověření jejího minima bychom museli použít její vyšší derivace.
16