Obsah 1.Lineární problém, notace lineární úlohy........................................................................2 2.Operace s maticemi.......................................................................................................3 3.Repertoár řešení “standardního” systému lineárních rovnic.........................................5 3.1.Gaussova eliminace...............................................................................................5 3.2.Faktorizace matic a jejich využití k řešení “standardního” systému rovnic............5 3.3.Pomocí determinantů (Cramerova metoda)...........................................................6 3.4.Pomocí inverze matice (neefektivní postup)..........................................................6 4.Přeurčená “nestandardní” lineární úloha.......................................................................7 5.Podurčená “nestandardní” lineární úloha......................................................................8 6.Derivace s vektory a maticemi.......................................................................................9 7.Lagrangeovy multiplikátory..........................................................................................12 8.Vlastní vektory a vlastní čísla......................................................................................13 9.Singular value decomposition (Singulární rozklad matice).........................................14 9.1.Up,Vp,Sp..............................................................................................................15 9.2.Pseudoinverzní matice pomocí SVD....................................................................17 9.3.SVD a vlastní čísla/vektory...................................................................................18 9.4.Zobecněná inverze...............................................................................................19 9.5.Matice rozlišení.....................................................................................................20 9.6.Aplikace SVD........................................................................................................20 10.Regularizace špatně podmíněných úloh...................................................................21 10.1.Tichonovova regularizace...................................................................................21 10.2.Hierarchie diferenčních operátorů......................................................................22 11.Inverze slabě nelineárních úloh.................................................................................23 11.1.Metoda Newton-Raphson...................................................................................24 11.2.Metoda proměnné metriky..................................................................................27 11.3.Metoda Levenberg-Marquardt............................................................................28 11.4.Metoda “Downhill simplex”.................................................................................28
1
1. Lineární problém, notace lineární úlohy Klasifikace obrácených úloh, klasifikace lineárních obrácených úloh, atd. viz. úvodní přednáška. Nadále uvažujeme diskrétní lineární problém, tj. vektor parametrů m a vektor dat d a jejich vzájemnou relaci prostřednictvím matice G. Výhodná je maticová notace: Gm = d, Ax = b ... přímá úloha m = ???, x = ??? ... obrácená úloha Gm = d nebo mTGT = dT ? Vektory budeme předpokládat jako sloupcové. Je možná i opačná varianta. Důležitá je platnost příslušného fyzikálního zákona, resp. “funkčnost” matematického popisu. Pro vyšší dimenze než 2 (tenzory nebo matice vyšších dimenzí) nemá pojem “řádka” a “sloupec” smysl. Dimenze budou značeny dim(m) = n, dim(d) = m, tj. dim(G) = [m,n]. Fyzika uvažuje vektor jako veličinu mající velikost a směr. V obrácených úlohách uvažujeme vektory jako abstraktní konstrukce seskupující paramery modelu (m) nebo jednotlivá měření (d). Jednotlivé složky vektorů spolu nemusí nijak vnitřně souviset, např. mohou mít různý fyzikální rozměr. Důvodem seskupení jednotlivých údajů do vektoru je snadná manipulace s kompaktním objektem. Nutnost pečlivé kontroly dimenzí, resp. rozlišení “řádkového” a “sloupcového” vektoru (pozn: “řádkové” vektory nechceme používat vůbec): Příklad možné chyby
[]
x1 x xc = 2 , x r =[ x1 , x3 ⋯ xk
[
a 11 a12 a 13 ⋯ a 1k a a a ⋯ a 2k x2 , x 3 , ⋯ x k ] , A= 21 22 23 ⋯ ak1 ak2 ak3 ⋯ akk
]
a) xc*(xc)T*A b) xr*(xr)T*A Oba případy a) i b) mají smysl a výsledkem jsou matice k x k, ale každá jiná. Pozor na chyby při programování! Lineární problém: Gm1 = d1, Gm2 = d2 => G(m1 + m2) = d1 + d2, Gm = d => G(c.m) = c.d, c je konstanta. Obecný lineární problém Gm = d + d0, nebo G(m + m0) = d
2
2. Operace s maticemi součet, součin, násobení skalárem, transpozice, inverze pseudoinverze obdélníkové matice Moore 1920, Penrose 1955 AA#A = A A#AA# = A# (A#A)T = (A#A) (AA#)T = (AA#) vlastní čísla a vlastní vektory Ax = λx (pravé vlastní vektory/čísla) 1. A + 0 = 0 + A = A 2. A + B = B + A 3. (A + B) + C = A + (B + C) 4. A(BC) = (AB)C 5. A(B + C) = AB + AC 6. (A + B)C = AC + BC 7. (st)A = s(tA) 8. s(AB) = (sA)B = A(sB) 9. (s+t)A = sA + tA 10. s(A+B) = sA + sB 11. (AT)T = A 12. (sA)T = s(A)T 13. (A + B)T = AT + BT 14. (AB)T = BT AT 15. (AB)-1 = B-1A-1 16. (A-1)-1 = A 17. (AT)-1 = (A-1)T 18. A=A(n,n), B = B(n,n), AB = I → A-1 = B a B-1 = A Časté chyby : následující vztahy NEPLATÍ obecně • AB = BA • AB = 0, → A = 0 | B = 0 • AB = AC & A ≠ 0 → B = C Konvence pro derivace s maticemi • d/dx(xTA) = A • d/dx(Ax) = AT • d/dx(xTAx) = Ax + Atx
3
Vlastnosti determinantů čtvercových matic 1. 2. 3. 4. 5. 6. 7.
det(In) = 1 det(AT) = det(A) det(A.B) = det(A).det(B) det(cA) = cn.det(A) det(L) = Π(L) det(U) = Π(U) det(A-1) = 1/det(A)
4
3. Repertoár řešení “standardního” systému lineárních rovnic g11 m 1+ g12 m 2+ ⋯+ g 1n mn =d 1 g21 m 1+ g22 m 2+ ⋯+ g 2n mn =d 2 g31 m 1+ g 32 m 2+ ⋯+ g 3n mn=d 3 ⋯ gm1 m1+ g m2 m2+ ⋯+ gmn mn=d m m=n , det (G)≠0 Jako “standardní” uvažujeme situaci přesně určeného systému se čtvercovou maticí soustavy G. V geofyzice převládá v některých metodách podurčený systém (např. linearizovaný krok seismické tomografie), nebo přeučený systém (např. inverze momentového tenzoru) a přesně určený sytém je spíše výjímečný. Pojem “standarní” se tak týká jen matematického aparátu. 3.1. Gaussova eliminace Při ní se provádí “ekvivalentní úpravy” rovnic tak, abychom dospěli k soustavě typu h 11 m1h12 m2 h13 m3⋯h1m mm=d 1 h 22 m2h 23 m3⋯h 2m mm =d 2 h32 m2⋯h3m mm=d 3 ⋯ h nm mm =d n z níž se neznámé parametry určují dosazováním zdola nahoru (metoda “backsubstitution”). “Ekvivalentní” úpravy jsou ekvivalentní jen pro bezchybné nebo pro statisticky konzistentní rovnice. Pokud rovnice reprezentují reálné měření, nelze libovolně sčítat/odčítat/násobit jednotlivé rovnice. Je významné, jak se při takových “ekvivalentních” úpravách šíří chyby. Pozor na to! Blíže viz “Metoda maximální věrohodnosti”. 3.2. 3.2.1.
Faktorizace matic a jejich využití k řešení “standardního” systému rovnic LU dekompozice
G = LU, L...dolní trojúhelníková matice, U … horní trojúhelníková matice Postup: Gm=d ... původní formulace, dále použijeme LU faktorizaci a zavedeme novou proměnnou y : LUm=d ,Um= y → Ly=d ...zpětným dosazováním spočteme y a z rovnice Um= y ...zpětným dosazováním požadované m
3.2.2.
QR dekompozice
G = QR, Q...ortogonální matice, R … trojúhelníková matice Postup: Gm=d ... původní formulace, využijeme ortogonality matice Q: Q−1=QT QRm=d , Rm=QT d= d̃ ...zpětným dosazováním spočteme m
5
3.2.3.
Choleskiho dekompozice
G = LLT je speciální případ LU dekopozice pro G symetrickou pozitivně definitní matici, postup řešení je proto stejný jako pro LU dekopozici. Choleskiho dekompozice je numericky velmi efektivní. Pozn. Symetrická pozitivně definitní matice je definována následovně: G=GT (symetrie) , x T G x> 0 ∀x (pozitivní definitnost)
3.2.4.
Metoda singulárního rozkladu
G = USVT (Singular Value Decomposition) 3.3.
Pomocí determinantů (Cramerova metoda)
det ([G(: ,1 :i−1) d G(:, i+ 1 :n)]) det (G) čitatel obsahuje matici G která má v i-tém sloupci vektor d mi=
Numericky velmi neefektivní, nutnost výpočtu n+1 determinantů 3.4.
Pomocí inverze matice (neefektivní postup)
Gm=d ,
m=G−1 d
6
4. Přeurčená “nestandardní” lineární úloha m > n ... přeurčená úloha, více rovnic než neznámých. Gm = d, G = G[m,n] C = C[n,m] ... “kontrakční” matice CGm = Cd (CG)-1CGm = (CG)-1Cd
=>
m = (CG)-1Cd
C = GT => m = (GTG)-1GTd (Gm - d)T(Gm - d) = min (1.Gaussova transformace, normální rovnice, metoda nejmenších čtverců) Důkaz L2 = (Gm - d)T(Gm - d) = min(m) L2 = (Gm - d)T(Gm - d) = (mTGT - dT)(Gm - d) = mTGTGm - mTGTd - dTGm + dTd dL2/dm = GTGm + (GTG)Tm - GTd - (dTG)T = 2GTGm - 2GTd = 0 GTGm = GTd m = (GTG)-1GTd
7
5. Podurčená “nestandardní” lineární úloha m < n ... podurčená úloha, méně rovnic než neznámých. Gm = d, G = G[m,n] C = C[n,m] ... “expanzní” matice m = Cz ... z je nová proměnná místo m Gm = GCz = d
=>
z = (GC)-1d
=> m = Cz = C(GC)-1d
C = GT => m = GT(GGT)-1d mTm = min (2.Gaussova transformace, metoda nejmenší normy) Důkaz mTm = min(m) Gm - d = 0 grad(mTm) = grad(Gm - d)λ grad(mTm) = 2m, grad((Gm - d)λ) = GTλ Gm - d = 0 2m = GTλ 1/2 GGTλ - d = 0 => λ = 2(GGT)-1d m = 1/2 GTλ = 1/2 GT 2 (GGT)-1d = GT(GGT)-1d
8
6. Derivace s vektory a maticemi Konvence pro derivace
[] []
x1 y1 x= x 2 , y= y 2 , ⋯ ⋯ xm yn
[ ]
dy 1 dx 1 d y dy 1 y= y ( x), D= = dx 2 dx ⋯ dy 1 dx m x : { xi } ,
y : { yi } ,
dy 2 dx 1 dy 2 dx 2 ⋯ dy 2 dx m
dy n dx 1 dy n ⋯ dx 2 ⋯ ⋯ dy n ⋯ dx m ⋯
D : { d ij }
d/dx(xTA) = A d/dx(Ax) = AT d/dx(xTAx) = Ax + ATx
• • • Důkaz
a) y = xTA, dy/dx = A D: dij =
dy j d ( akj x k ) d (a1j x 1+ a2j x 2+ ⋯+ aij x i+ ⋯+ amj x m ) = = =aij dx i dx i dxi dy D= =A dx
b) y = Ax, dy/dx = AT D: d ij=
dy j d (a jk x k ) d (a j1 x 1+ a j2 x 2+ ⋯+ a ji xi + ⋯+ a jm x m ) = = =a ji dx i dx i dx i dy T D= =A dx
9
c) y = xTAx, dy/dx = (A+AT)x D: d ij=
dy j d ( alk x k xl ) d (ali x i xl + aik x k x i) = = =ali xl + aik x k =aTil x l+ aik x k =(aTil + ail ) x l dx i dx i dx i dy D= =( AT + A) x dx (xT x) A=I → d =2 x dx (x T A x) A symetrická → d =2 Ax dx
Příklad 1 Spočtěte první a druhou derivaci skalární funkce f vektorového argumentu x
[]
x1 f =f (x)=f ( x 2 ⋯ xm
[] [
d2 f df dx 21 dx 1 d2 f df d 2 f df d df ), = , 2= ( )= dx 2 d x d x d x dx2 dx 1 dx ⋯ ⋯ df d2 f dx m dx m dx 1
d2 f dx 1 dx 2 d2 f dx 22 ⋯ d2 f dx m dx2
]
d2 f dx 1 dx m d2 f ⋯ =H dx 2 dx m ⋯ ⋯ d2 f ⋯ dx2m ⋯
H … Hessian, symetrická čtvercová matice Příklad 2 Mějme následující skalární funkci f vektorového argumentu x: 1 T T f ( x)=f (x0 )+ a ( x−x0 )+ (x−x 0) A( x−x 0) 2 Najděte x v němž je f(x) minimální. aT je nějaký vektor a A je nějaká matice (s konzistentními dimenzemi). Řešení pro 1D: f (x)=f (x 0)+ a(x −x0 )+
1 ( x−x 0)2 b 2
df =a+ b ( x−x 0)=0 dx a x=x 0− b Obecné řešení: 10
1 f ( x)=f (x0 )+ aT ( x−x0 )+ (x−x0 )T A( x−x 0) 2 df 1 T =a+ ( A+ A )( x−x0 )=0 dx 2 1 ( A+ AT )=H 2 H ( x−x0 )=−a x=x 0−H −1 a Pozn. Taylorův rozvoj skalární funkce f vektorového argumentu x je: f ( x)≈f (x0 )+ (
df T 1 T d2 f ) ( x−x0 )+ ( x−x 0) ( x−x 0)+ ⋯ dx 2 d x2
11
7. Lagrangeovy multiplikátory Uvažujme skalární funkci f vektorového argumentu x: f = f(x). K takovým funkcím patří například různé normy shody naměřených a predikovaných dat, tyto normy jsou následně předmětem optimalizace, tj. hledání min(f(x)), resp. max(f(x)). Obecnější je optimalizace s podmínkou, neboli hledání vázaných extrémů. V dalším postupu se diskutuje veličina df/dx. T
f ( x)≈f ( x 0)+ a (x− x 0) df ≈a dx Δ f =f (x+ Δ x)−f ( x)≈aT Δ x=∣aT∣.∣Δ x∣cos (ϕ) Δ f =max ≡ ϕ=0 Δ f =min ≡ ϕ=π Δ f =0 ≡ ϕ=π /2
• • •
df/dx je směr maximálního růstu funkce f izolinie funkce f jsou určeny podmínkou f(x) = const. Tečný směr k izoliniím v je určen podmínkou (df/dx)T v = 0
Obyčejná optimalizace
Jednoduchý vázaný extrém
Násobný vázaný extrém
f = f(x)
f = f(x), g = g(x)
f = f(x), g1(x), g2(x),..., gq(x)
x0: f(x0) = min
x0: f(x0) = min g(x) = 0
x0: f(x0) = min g1(x) = 0 g2(x) = 0 ... gq(x) = 0
df/dx = 0
df/dx – λdg/dx = 0 g(x) = 0
df/dx – λ1dg1/dx – λ2dg2/dx – ... – λqdgq/dx = 0 g1(x) = 0 g2(x) = 0 ... gq(x) = 0
L = f(x) + λg(x) x = [x1 x2 x3 ... λ]T dL/dx = 0
12
8. Singular value decomposition (Singulární rozklad matice) G=U∗S∗V T
=
G
T
T
U
−1
S
*
VT
*
T
U ∗U = I ,U ∗U = I ,U =U V T ∗V = I , V ∗V T = I , V −1=V T
T
G=U p∗S p∗V p
Sp G
=
Up
(Vp)T
*
*
Model Null space
Data Null space
Model NULL space Nechť Gm = d má známé řešení m = G#d. Pak m+m0 je také řešením daného problému, kde m0 je libovolný vektor z nulového prostoru:
13
G(m+m0 )=d , →
Gm0=0
n
m 0=
∑
i= p +1
,
βi v i
vektory vi jsou sloupce matice V odpovídající nulovým singulárním číslům.
Gm0=USV T m 0=U p S p V Tp m 0=
n
∑
i= p +1
βi U p S p V Tp v i =0
z podmínky ortogonality matice V , protože V Tp v i=0,i≠ p Vektory m, resp. m+m0 jsou po transformaci do datového prostoru nerozlišené. Data NULL space Při přeurčené úloze nelze část prostoru dat dosáhnout žádnou transformací vektorů z prostoru parametrů. Modelové vektory se transformují do nějaké nadroviny v prostoru dat. Doplněk zobrazitelné části prostoru dat je data_NULL_space. Data z NULL space nejsou obrazem žádného vektoru m.
8.1.
Up,Vp,Sp
U Tp ∗U p = I p , U p∗U Tp≠ I p ,U −1 p =nedefinována T T −1 V p ∗V p =I p ,V p∗V p ≠I p ,V p =nedefinována −1 det S p ≠0, S p existuje
[ ]
U =[ U p U 0 ] , U T =
[ ]
T
U U =[ U p U 0 ]∗
U T U=
U Tp U T0
[ ] U Tp U
T 0
U Tp T T T = [ U p U p U 0 U 0 ]= I ⇒ U p U p ≠ I T U0
∗[ U p U 0 ]=
[
U Tp U p U Tp U 0 T 0
T 0
U U p U U0
][ =
Ip 0 0 I m− p
14
]
⇒ U Tp U p=I p
Struktura matice S určuje typ řešené úlohy Spp: sij = 0, i ≠ j sij ≠ 0, i = j
Spp
Smn
přesně určená úloha
přeurčená úloha
podurčená úloha
15
smíšená úloha
8.2.
Pseudoinverzní matice pomocí SVD
Obecná definice pseudoinverzní matice (Penrose-Moore) G G # G=G G # G G #=G # (G G# )T =G G # (G# G)T =G # G
#
−1
T
pinv(G)=G =V p∗S p ∗U p
!!! Pozor, toto neplatí G# ≠V∗S−1∗U T a) G G # G=G T T T −1 T T G G # G=U P S P V TP . V P S−1 P U P . U P SP V P =U P S P . V P V P . S P . U P U P . SP V P = T T U P . S P S−1 P . S P V P =U P SP V P =G
b) G# G G #=G# T -1 T T -1 T G # G G# =V P S-1P U TP . U P S P V TP . V P S-1 P U P =V P S P . U P U P . SP . V P V P . SP U P = -1 T -1 T V P . S−1 P S P . SP U P =V P S P U P=G
# T
c),d)
#
(G G ) =G G # T # (G G) =G G
T T T T T (G G #)T =(U P SP V TP . V P S−1 P U P ) =(U P U P ) =U P U P T T G G #=U P S P V TP . V P S−1 P U P =U P U P
a obdobně pro (G#G)T
16
8.3.
SVD a vlastní čísla/vektory
Ax=λ x , ∣x∣=1 ( A−λ I ) x=0 x=0 ∨ det ( A−λ I )=0
G T G=V S T U T U S V T =V S T S V T G T G=V p S Tp U Tp U p S p V Tp =V p S 2p V Tp G T G V =V S T S V T V =V S T S v 1j V =v ij , v j = v 2j ⋯ v nj
[]
G T G v j =s2j v j
G G T U=U S ST u1j U=u ij , v j= u2j ⋯ u mj
[]
G G T u j=s 2j u j Ax= λx
G = USVT
Matice A čtvercová
Matice G je obecná
Vlastní čísla reálná/komplexní
Singulární čísla nezáporná
Matice vlastních vektorů X = [x1,x2,...xn]
Matice vlastních vektorů U = [u1,u2,...um] příslušející matici GGT Matice vlastních vektorů V = [v1,v2,...vn] příslušející matici GTG GTG má nenulová vlastní čísla diag(Sp)2 GGT má nenulová vlastní čísla diag(Sp)2
17
8.4.
Zobecněná inverze m=G# d T m=V p∗S−1 p ∗U p∗d
“Standardní” situace (přesně určený systém) U =U p V =V p S =S p , S regulární m=V S−1 U−1 d ... použita definice pseudoinverze −1 −1 −1 −1 −1 G=U S V , G =V S U , m=G d ...řešení je ekvivalentní užití klasické inverzní matice Přeurčený systém – metoda nejmenších čtverců U ≠U p V =V p S ≠S p , S p regulární #
T
−1
T
m=G d≡(G G) G d T T m=(V p S p U Tp U p S p V Tp )−1 V p S p U Tp d=(V p S2p V Tp )−1 V p S p U Tp d=V p S−2 p V p V p SpU p d = T # V p S−1 p U p d =G d Podurčený systém – metoda nejmenší normy U =U p V ≠V p S ≠S p , S p regulární #
T
T −1
m=G d≡G (G G ) d T m=V p S p U Tp (U p S p V Tp V p S p U Tp )−1 d=V p S p U Tp (U p S2p U Tp )−1 d=V p S p U Tp U p S−2 p Upd = T # V p S−1 p U p d =G d
18
8.5.
Matice rozlišení
Rm =G # G=V p V Tp ≠V V T , Model Resolution matrix Rd =G G #=U p U Tp ≠U U T , Data Resolution matrix Gm true=d true m+ =G # d true =G # G m true=Rm m true T Rm =G # G=V p S−1 U p S p V Tp =V p V Tp p Up
8.6.
Aplikace SVD
•
Zobecněná inverze a formálně jednotné řešení všech variant systému lineárních rovnic
•
Regularizace špatně podmíněných úloh vynulováním malých singulárních čísel (= aproximace matice)
•
Explicitní definice model NULL space a možnost vytváření ekvivalentních modelů
•
Projekce apriorního modelu do model NULL space
•
Resolution matrix
19
9. Tichonovova regularizace Základní motivace = alternativní metoda řešení podurčeného systému lineárních rovnic
G
m
=
d
a) m = G#d ... řešení existuje vždy, neboť pseudoinverze existuje vždy b) m = GT(GGT)-1d ... neexistuje vždy, GGT může být singulární
20
9.1.1.
Tichonovova regularizace nultého řádu:
Gm=d →
[ ][ ] [] Gm G d = m= λIm λI 0
... vždy přeurčené, pokud λ > 0,
řešíme převodem na normální rovnice:
[ ]
[]
λ I ] G m=[ GT λ I ] d λI 0
[ GT
−1
[ GT G+λ 2 I ] m=GT d →m=[ GT G+λ2 I ]
T
G d
řešení existuje vždy (důkaz pomocí SVD) dále. 1 lim m=0 , m≈ GT d λ λ →∞ Důkaz že (GTG + λ2I) je regulární Singulární rozklad matice G a GT G=U S V T T T T T G =V S U , S ≠S... transpozice je významná, nelze ji vynechat T ̃ G=G G+ λ2 I =V ST U T U SV T + λ 2 I =V S T S V T + λ2 I =V ST SV T + λ 2 V V T =V ( ST S+λ 2 I ) V T
Singulární rozklad matice (GTG + λ2I) ̃ U ̃ S̃ Ṽ T G= Ũ =V S̃ =ST S+ λ2 I Ṽ =U Determinant matice G
[
s21 +λ 2 0 0 2 2 0 s2 +λ 0 ̃ S= 2 0 0 s3 + λ2 ⋯ ⋯ ⋯
]
⋯ ⋯ ⋯ ⋯
̃ det ( G)=det ( Ũ ) det( S̃ )det ( Ṽ T )=1 . ∏ diag S̃ . 1=∏ ( s2i +λ 2)>0
21
Geometrie matice S a STS Původní systém podurčený
Původní systém přeurčený
S S ST ST matice soustav G G
G
22
9.1.2.
Hierarchie regularizačních operátorů
1. nultý řád (minimální hodnoty) Dává řešení s minimalizací absolutních hodnot. Má význam pro řešení formulované jako perturbace okolo referenčního modelu: 2.
první řád (“minimální sklon”, “minimální změny”) x m1
m2
m3
m4
m5
m6
m7
m8
m9
m10
m11
m12
m13
m14
m15
m16
y
Minimální variabilita řešení v horizontálním směru m2−m1=0 m3−m2=0 −1 1 0 dm m4−m3=0 ≈0 : ⇒ L= 0 −1 1 dx 0 0 −1 --(další vrstva)-⋯ ⋯ ⋯ m6 −m5=0 ⋯ Minimální variabilita řešení ve vertikálním směru m5−m1=0 m9−m5=0 dm 0 0 m13−m9 =0 ≈0: ⇒ L= −1 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ dy --(další sloupec)-m16−m12 =0 ⋯
[
[
0 0 0 0 1 0 ⋯ ⋯
][
][ ]
[
][ ]
⋯ ⋯ , GΔ m = Δd λL 0 ⋯ ⋯
]
0 1 ⋯ , GΔ m = Δd ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ λL 0
Matice L má v každém řádku jednou 1 a jednou -1, zbytek jsou 0.
23
3.
druhého řádu (“minimální křivost”, “maximální hladkost”)
y
m1
m2
m3
m4
m5
m6
m7
m8
m9
m10
m11
m12
m13
m14
m15
m16
Minimální křivost v horizontálním směru m3−2m2 +m1=0 1 −2 1 0 0 m4−3m3 +m2=0 d2 m 0 1 −2 1 0 ≈0 : --(další vrstva)-- ⇒ L= dx2 0 0 1 −2 1 m15−2m 14 +m13=0 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯
[
⋯ ⋯ , GΔm = Δd λL 0 ⋯ ⋯
Minimální křivost řešení ve vertikálním směru m9−2m 5+ m1=0 2 d m ≈0 : m13−2m9 +m5=0 ⇒ L= [ {0|1|-2} ] , G Δ m = Δ d 2 λL 0 dy --(další sloupec)-⋯
[
][ ]
Matice L má v každém řádku dvakrát 1 a jednou -2, zbytek jsou 0.
24
][
x
][ ]
9.1.3.
Matice rozlišení při aplikaci regularizace
a) klasická formulace bez regularizace Gm=d , G=USV T (SVD reprezentace) mresolv =G # d synth =G # Gmsynth mresolv ≈m synth , R=G# G b) aplikována regularizace mresolv =(GT G+ λ2 I )−1 G T d synth =(G T G+ λ2 I )−1 GT Gmsynth R=(GT G+λ 2 I )−1 GT G použijeme SVD R =(GT G+λ2 I )−1 GT G =(VST U T USV T + λ2 I )−1 VST U T USV T =(VST SV T + λ2 I )−1 VST SV T = (VST SV T +λ2 VV T )−1 VST SV T =(V (ST S+ λ2 I )V T )−1 V ST SV T = V (ST S+λ2 I )−1 V T V ST S V T = V (ST S+λ2 I )−1 ST S V T = V S̃ V T tedy
T R =V S̃ V = V
9.2.
[
s 21 0
]
0
0
⋯
s 22 s 22 +λ 2
0
⋯V
s 21+ λ2
0
0
s23 s23 + λ2
T
⋯
Volba parametru λ
L-curve RMS - level Cross-validation
25
10. Inverze slabě nelineárních úloh Původní (obecně nelineární) problém jako soustava rovnic je f 1 (m)=d 1 f 2 (m)=d 2 , ⋯ F( m)=d referenční model je m0 a příslušná referenční data d0 = F(m0). Perturbační přístup je založen na linearizaci: dF T dF F( m)≈F (m 0)+ Δ m=d 0+ J Δ m=d , J= dm dm
( )
T
( )
,
kde J je Jacobian zobrazení modelového prostoru do datového prostoru. Řešení je iterativní, při použití tlumení závisí na parametru λ
[ ][ ]
d 0 + J Δ m=d ⇒ J Δ m=d −d 0 =Δ d ⇒ J Δ m = Δ d λI 0 1. k =0, m k=m0 2. k =k +1, m k =m k + Δ m 3. STOP or GOTO 2
Lineární problém
Nelineární problém
Slabě nelineární problém
Gm = d
G(m) = d
•Jen jeden nelineární člen Taylorova rozvoje •Všechny funkce fi(m) jsou monotónní
F(m,d) = G(m) - d = 0 norm(F) = min(m)
• Newton-Raphson • Levenber-Marquardt • Proměnná metrika • Simplexová metoda • Robustní optimalizace
26
10.1. Metoda Newton-Raphson
[ ][ ]
f 1 ( m, d ) f 2 ( m, d ) F(m , d )= f 3 ( m, d ) = f 4 ( m, d ) ⋯ f n ( m, d )
0 0 0 0 ⋯ 0
Motivace: metoda tečen v 1D f(m) ≈ f(m0) + f'(m0)(m-m0) + ... Iterační cykly k = 0; mk = m0
k = k+1; mk = mk + ∆m
F(mk)
testy, revize
F(m+Δ m, d )≈F(m, d)+ ∇ FT Δ m+⋯=0 ∇ FT Δ m=−F(m , d )
[ ]
∂f1 ∂ m1 ∂f1 T ∇ F=J = ∂ m2 ∂f1 ∂mn
∂f2 ∂fm ⋯ ∂m 1 ∂ m1 ∂f2 ∂fm ⋯ ∂m 2 ∂ m2 ⋯ ∂f2 ∂fm ⋯ ∂m n ∂ mn
T
∇ F Δ m=−F (m , d) Δ m=−( ∇ FT )- F( m, d ) Δ m=−pinv (∇ FT ) F( m, d )
Ukázka polynom 1D Lokalizace – teorie + ukázka 27
Inverze s.l. Nalezení reprezentativního modelu predikujícího syntetická data konzistentní s naměřenými daty. Optimalizace
Inverze s.s. Promítání vektoru z datového prostoru do modelového prostoru.
Přiřazení normy shody dat (principiálně každému) vektoru z modelového prostoru. Vyhledání modelu s minimální normou.
G(m) = d
Datový prostor
Modelový prostor
H(d) = m
1
L p=(∑ (∣Δ xi∣)) p i
Datový prostor
Modelový prostor
{m1,(d1),D1}, {m1,(d1),D1},...,{mn,(dn),Dn}, mopt = argmin(D)
28
10.2. Metoda proměnné metriky Pozn. Vynecháváme d z argumentů funkcí. F(m) = 0 → φ(m) = FT(m)F(m) = ∑(fi(m))2 = min(m) Inverze → optimalizace m= f 1 m2 f 2 m2⋯ f n m2= ∑ f i m2 =minm i
1 ϕ( m0 + Δ m)≈ϕ(m 0 )+ ∇ ϕT (m 0 )Δ m+ Δ mT ∇ 2 ϕ(m 0) Δ m+⋯ 2
ϕ( m 0+ Δ m)≈min(Δ m)≡grad Δ m ϕ(m0 + Δ m)=0 grad Δ m ϕ(m 0 +Δ m)=
{
}
d 1 ϕ(m 0)+ ∇ ϕT (m 0)Δ m+ Δ mT ∇ 2 ϕ( m0 ) Δ m =0 d Δm 2
[]
∂ ∂ m1 ∂ ∇ = J = ∂ m2 ⋯ ∂ ∂ mn
[
∂2 ∂ m12 ∂2 ∇ 2 = H = ∂ m 2 ∂ m 1
∂2 ∂2 ⋯ ∂ m1 ∂ m2 ∂ m 1 ∂ mn ∂2 ∂ m22
⋯
∂2 ∂ m 2 ∂ mn
⋯ 2
∂ ∂ mn ∂ m1
2
∂ ⋯ ∂ mn ∂ m 2
∂2 ∂ m2n
29
]
{
}
d 1 ϕ(m 0)+ ∇ ϕT (m 0 )Δ m+ Δ mT ∇ 2 ϕ(m 0)Δ m = Δm 2 2 =∇ ϕ+ ∇ ϕ Δ m=J + H Δ m=0
grad Δ m ϕ(m 0 +Δ m)=
H Δ m=−J Δ m=−H −1 J
10.3. Metoda Levenberg-Marquardt Metoda Newton-Raphson s dynamicky nastavovaným tlumením. ( ∇ F)T Δ m=−F( m) J Δ m=−F J J Δ m=−J T F ( J T J + λ I ) Δ m=−J T F T
Δ m=−( J T J +λ I )−1 J T F
10.4. Metoda “Downhill simplex” /home/br/Documents/Lectures/Inverse/DownhilSimplex/DownhillSimplexA.pdf /home/br/Documents/Lectures/Inverse/DownhilSimplex/DownhillSimplex.pdf
30
⇚ ⇚ ⇚ ⇚ Efektivita ⇚ ⇚ ⇚ ⇚ Lineární problém
Newton-Raphson
Proměnná metrika
Simplex
Model⇰Data
Model⇨Data
Model⇨Data⇨Norma
Model⇨Data⇨Norma
Data⇨Model
Data≈⇨Model
Jakobián normy
Populace velikosti n+1
Jakobián + inverze matice m x n
Hessian normy+ inverze matice mxm
Žádné derivace
⇛ ⇛ ⇛ ⇛ Robustnost ⇛ ⇛ ⇛ ⇛ Genetický algoritmus Evoluční algoritmy Particle swarm optimization Harmony search ANN
31
