Line´arn´ı stochastick´ y syst´em a jeho vlastnosti. Kovarianˇcn´ı funkce, v´ ykonov´a spektr´aln´ı hustota, spektr´aln´ı faktorizace, tvarovac´ı filtr ˇsumu, bˇelic´ı filtr. Kalman˚ uv filtr, formulace problemu, vlastnosti. LQG regul´ator, separaˇcn´ı princip. Metoda ”‘Loop Transfer Recovery”’. Ot´azka pokr´ yv´a partie pˇredmˇetu X35MTR. Tato ot´azka je vypracov´ana dle [1], [2], [3] a z´apisek z pˇredn´aˇsek.
Obsah 1 Line´ arn´ı stochastick´ y syst´ em a jeho vlastnosti
1
2 Kovarianˇ cn´ı funkce, v´ ykonov´ a spektr´ aln´ı hustota, spektr´ aln´ı faktorizace, tvarovac´ı filtr ˇ sumu, bˇ elic´ı filtr 2 3 Kalman˚ uv filtr, formulace problemu, vlastnosti
4
4 LQG regul´ ator, separaˇ cn´ı princip
6
5 Metoda Loop Transfer Recovery
8
1
Line´ arn´ı stochastick´ y syst´ em a jeho vlastnosti
Deterministick´ y syst´em popsan´ y x(t + 1) = Ax(x) + Bu(t) y(t) = Cx(x) + Du(t)
(1) (2)
rozˇs´ıˇr´ıme o ˇsum procesu v(t) resp. mˇeˇren´ı e(t), v z´akladn´ıch u ´vah´ach vych´az´ıme z faktu, ˇze se jedn´a o tzv. b´ıl´y ˇsum. Jeho parametry jsou n´asleduj´ıc´ı: ("
ε
v(t) e(t)
#)
("
= 0,ε
v(t1 ) e(t1 )
#"
v(t2 ) e(t2 )
#)
"
=
Q S ST S
#
δ(t1 − t2 ),
(3)
U stochastick´eho procesu n´as zaj´ım´a v´ yvoj stavu a v´ ystupu. Ten popisujeme v´ yvojem stˇredn´ı hodnoty. µx (t + 1) = Aµx (t) + Bu(t) µy (t) = Cµx (t) + Du(t) n
o
(4) (5)
kovarianˇcn´ı matici stavu lze zapsat jako Px = ε x˜x˜T , kde x˜ = x(t) − µx (t). Odhad stavu x˜ vlastnˇe ud´av´a rozd´ıl skuteˇcn´e hodnoty od stˇredn´ı hodnoty. Jinak ˇreˇceno pokud je nejlepˇs´ı situace nast´av´a pokud je hodnota kovarianˇcn´ı matice n´ızk´a ide´alnˇe nulov´a, pot´e je odhad stavu pˇresn´y, tud´ıˇz rozd´ıl stˇredn´ı hodnoty stavu a skuteˇcn´e hodnoty stavu je nulov´ y.
1
Odovozov´an´ım lze dospˇet k n´asleduj´ıc´ım rovnic´ım (rovnice pro odchylky), kter´e n´am ud´avaj´ı stochastick´e veliˇciny (deterministick´a sloˇzka zde nen´ı resp. stˇren´ı hodnota je nulov´a – nem´a cenu j´ı odhadovat jelikoˇz j´ı zn´ame) µx (t + 1) = Aµx (t) + v(t) µy (t) = Cµx (t) + e(t)
(6) (7)
D´ale kovarianˇcn´ı matice maj´ı tvar: n
o
(8)
o
(9)
o
(10)
ε x˜(t + 1)˜ xT (t + 1) = APx (t)AT + Q n
ε x˜(t + 1)˜ y T (t) = APx (t)C T + S n
ε y˜(t)˜ y T (t) = CPx (t)C T + R
(11) D´ale je moˇzn´e odvodit sdruˇzen´e kovarianˇcn´ı matice stavu a v´ ystupu. Tato odvozen´ı maj´ı z´asadn´ı vliv pro odovzen´ı Kalmanova filtru. ("
cov
x(t + 1) y(t)
#)
"
=
APx (t)AT + Q APx (t)C T + S CPx (t)AT + S T CPx (t)C T + R
#
δ(t1 − t2 ),
(12)
D´ale je moˇzn´e urˇcit rekurzivn´ı rovnice, kter´e maj´ı tvar: µx (t) = At µx (0) +
t−1 X
At−τ −1 Bu(τ )
(13)
τ =0
Px (t) = At Px (0)(AT )t +
t−1 X
Aτ Q(AT )τ
(14)
τ =0
Pˇredopkladem je ˇze matice A je stabiln´ı, pot´e ˇrady konverguj´ı k limitn´ı hodnotˇe Px = limt→∞ Px (t). Rovnice souvis´ı s Ljapunovou rovnic´ı pro ust´alen´e ˇreˇsen´ı: Px = APx (t)AT + Q s kovarianc´ı v´ ystupu Py = CPx (t)C T + R.
2
Kovarianˇ cn´ı funkce, v´ ykonov´ a spektr´ aln´ı hustota, spektr´ aln´ı faktorizace, tvarovac´ı filtr ˇ sumu, bˇ elic´ı filtr
Kovarianˇ cn´ı funkce: K popisu stochastick´ ych sign´al˚ u se obecnˇe pouˇz´ıvaj´ı tzv. momenty (obecn´e, centr´aln´ı a normovan´e). Nˇekdy se naz´ yvaj´ı t´eˇz charakteristick´e funkce. Ud´avaj´ı informaci o rozdˇelen´ı n´ahodn´e veliˇciny. Kovariance (obecn´ y moment) je stˇredn´ı hodnota souˇctu odchylek obou n´ahodn´ ych veliˇcin X, Y. Kovariance vlastnˇe poskytuje informaci o intenzitˇe vztahu mezi dvˇema veliˇcinami 2
vz´ajemn´a koerlaˇcn´ı funkce, pokud se jedn´a o jednu veliˇcinu pot´e hovoˇr´ıme o autokorelaˇcn´ı funkci. cov(X, Y ) = ε[X − ε(X)][Y − ε(Y )]
(15)
Je-li n´ahodn´ y vektor n-rozmˇern´ y, pot´e se pouˇz´ıv´a kovarianˇcn´ı matice, kde na diagon´ale jsou rozptyly jednotliv´ ych veliˇcin. Tato matice je symetrick´a!
C=
D(X1 ) C(X1 , X2 ) D(X2 , X1 ) D(X2 ) . . C(Xn , X1 ) C(Xn , X2 )
.. C(X1 , Xn ) .. C(X2 , Xn ) . . .. D(Xn )
(16)
V´ ykonov´ a spektr´ aln´ı hustota: Spektr´aln´ı hustotu lze definovat jako Fourierovu transformaci kovarianˇcn´ı funkce (posloupnosti). Jedn´a se o mˇeˇr´ıtko v´ ykonu procesu. Hodnotu v´ ykonov´e spektr´aln´ı hustoty Suu (ω) v bodˇe ω lze ch´apat jako hustotu v´ ykonu vstupn´ıho procesu v okol´ı frekvence ω. n o 1 Z ωN Suu dω = Ruu (0) = σu2 = ε u2 (t) 2ωN −ωN
(17)
D˚ ulˇeˇzit´ ym vztahem he vz´ajemn´a kovarianˇcn´ı funkce vstupu a v´ ystupu: Ryu (t) =
∞ X
g(τ )Ruu (t − τ ) = g ∗ Ruu
(18)
τ −∞
odov´ıdaj´ıc´ı konvoluci autokovarianˇcn´ı funkc´ı vstupu a imulzn´ı charakteristiky. Plat´ı vztah Ryu (t) = G ∗ Ruu → Syu (z) = G(z) ∗ Ruu (z). Vztahy pro vz´ajemnou kovarianˇcn´ı funkci a vz´ajemnou spektr´aln´ı hustotu umoˇzn ˇuj´ı identifikovat imulzn´ı charakteristiku a pˇrenos dynamick´eho syst´emu. Vz´ajemn´a kovarianˇcn´ı funkce je pˇr´ımo rovna impulzn´ı charakteristice syst´emu a vz´ajemn´a spektr´aln´ı hustota je rovna jeho frekvenˇcn´ımu pˇrenosu! Spektr´ aln´ı faktorizace: Jelikoˇz kovarianˇcn´ı funkce a tedy i spektr´aln´ı vlastnosti n´ahodn´eho procesu z´avis´ı jen na odchylk´ach od stˇredn´ı hodnoty, tak lze uvaˇzovat, ˇze deterministick´a sloˇzka je nulov´a. Tud´ıˇz, vstup je nulov´ y a za vstup lze povaˇzovat pouze b´ıl´ y ˇsum. T´ımto doch´az´ı ke zjednoduˇsen´ı. Pˇrenos mezi v´ ystupem a vstupem lze ps´at jako: G(z) = [C(zI − A)−1 I]
(19)
spektr´aln´ı hustota vstupu je d´ana kovarianˇcn´ı matic´ı "
Suu (z) =
Q S ST R
#
(20)
Dosazen´ım za vstup dost´av´ame " −1
Syu (z) = [C(zI − A) I]
Q S ST R
#
(21)
obodbnˇe lze odovdit Suy a Syy . Uvaˇzujeme-li nekorelovan´e (nez´avisl´e) b´ıl´e ˇsumy u(t), v(t) pak lze ps´at: Syy (z) = C(zI − A)−1 Q(z −1 I − A)−T C T + R (22) 3
Tvarovac´ı filtr ˇ sumu Tento filtr souvis´ı se stochastick´ ymi vlasnostmi Kalmanova filtru, uvaˇzujeme-li model ve tvaru: x(t + 1) = Ax(x) + Bu(t) + v(t) y(t) = Cx(x) + Du(t) + e(t)
(23) (24)
s nekorelovan´ ymi, b´ıl´ ymi ˇsumy. Pak s vyuˇzit´ım principu ortogonality a faktorizac´ı pozitivnˇe defintin´ı matice dost´av´ame inovaˇcn´ı posloupnost procesu y(t). Pot´e lze pˇrepsat Kalman˚ uv filtr do tvaru: xˆ(t + 1) = (Aˆ x(t|t − 1) + L(t)(t|t − 1) y(t) = C xˆ(t|t − 1) + (t|t − 1)
(25) (26)
nahraˇzen´ım chyb predikce inovacemi dostaneme xˆ(t + 1) = (Aˆ x(t|t − 1) + L(t)Γ (t)v(t) y(t) = C xˆ(t|t − 1) + Γ (t)v(t)
(27) (28)
coˇz je inovaˇcn´ı model procesu y(t) tedy model, jehoˇz vstupem je b´ıl´a posloupnost v(t) a v´ ystupem proces y(t). Tento model naz´ yv´ame jako tvarovac´ı filtr ˇsumu y(t). Inovace jsou vlastnˇe jen ”‘dekorelovan´e”’ chyby predikce coˇz je podstatn´a vˇec! Bˇ elic´ı filtr Model jehoˇz vstupem je proces y(t) a v´ ystupem je b´ıl´a posloupnost v(t) naz´ yv´ame bˇel´ıc´ı filtr pro proces y(t). Postup odvozen´ı tohoto filtru je takov´ y, ˇre z rovnic Kalmanova filtru hodnotu chyb predikce a dost´av´ame: xˆ(t + 1) = (A − L(t)C) xˆ(t|t − 1) + L(t)y(t) (t|t − 1) = −C xˆ(t|t − 1) + y(t)
(29) (30)
a po pˇrepoˇctu chyby predikce na inovace: xˆ(t + 1) = (A − L(t)C) xˆ(t|t − 1) + L(t)y(t) v(t) = −Γ−1 ˆ(t|t − 1) + −Γ−1 (t)y(t) (t)C x
3
(31) (32)
Kalman˚ uv filtr, formulace problemu, vlastnosti
Kalma˚ uv filtr odhaduje stavy syst´emu na z´akladˇe pozorov´an´ı vstup˚ u a v´ ystup˚ u. V determiistick´e formulaci se k odhadov´an´ı stav˚ u pouˇz´ıv´a pozorovatel stavu metodou um´ıst’ov´an´ı p´ol˚ u. Ve stochastick´e formulaci se k odhadu stav˚ u vyuˇz´ıvaj´ı krit´eria MS resp. LMS. Kalman˚ uv filtr je algoritmus generuj´ıc´ı posloupnost odhad˚ u stav˚ u xˆ(k|k) a koarianˇcn´ıch matic chyb odhadu P (k|k), pˇriˇcemˇz v kaˇzd´em kroku minimalizuje krit´erium JLM S = trace P (k|k). Stochastick´ y syst´ em - line´arn´ı, diskr´etn´ı, pozorovateln´ y 4
x(t + 1) = Ax(x) + Bu(t) + v(t) y(t) = Cx(x) + Du(t) + e(t)
(33) (34)
kde ˇsumy v(t), e(t) bereme jako diskr´etn´ı b´ıl´e posloupnosti s nulovou stˇredn´ı hodnotou ("
ε
v(t) e(t)
#)
=0
(35)
a se zn´amou kovarianˇcn´ı matic´ı ("
ε
v(t) e(t)
#)
("
= 0,ε
v(t1 ) e(t1 )
#"
v(t2 ) e(t2 )
#)
"
=
Q S ST S
#
δ(t1 − t2 ),
(36)
matice Q pˇredstavuje ˇsum stavu v(k) a matice R pˇredstavuje ˇsum ˇr´ızen´ı e(k). Matice S ud´av´a vz´ajemn´ y ˇsum (korelace mezi ˇsumem stavu a ˇr´ızen´ı). Pro nekorelovan´e ˇsumy je matice S = 0. Dosazen´ım xˆ(k +1|k −1), yˆ(k|k −1) do rovnice popisuj´ıc´ı syst´em dost´av´ame Ricattioho rovnici, pˇredstavuj´ıc´ı algoritmus Kalmanova filtru. h
i h
P (k+1|k−1) = M P (k|k − 1)M T + Q − M P (k|k − 1)C T
ih
CP (k|k − 1)C T + R
i−1 h
CP (k|k − 1)C T (37)
Tuto rovnici je v´ yhodn´e rozdˇelit do dvou nez´avisl´ ych krok˚ u. Datov´ y (filtraˇ cn´ı) krok: xˆ(k|k) = xˆ(k|k − 1) + L0 (k)ε(k|k − 1) P (k|k) = P (k|k − 1) − L0 (k)CP (k|k − 1) ε(k|k − 1) = y(k) − C xˆ(k|k − 1) − Du(k) h
L0 (k) = P (k|k − 1)C T CP (k|k − 1)C T + R
i−1
(38) (39) (40) (41)
ˇ Casov´ y (predikˇ cn´ı) krok: – vztahuje se pˇr´ımo k syst´emu xˆ(k + 1|k) = M xˆ(k|k) + N u(k) P (k|k) = P (k|k − 1) − L0 (k)CP (k|k − 1) P (k + 1|k) = M P (k|k)M T + Q
(42) (43) (44)
Celkov´ e zes´ılen´ı Kalmanova filtru lze popsat jako odhad stavu x(k|k − 1) ˇcasov´ ym krokem a na z´akladˇe mˇeˇren´ı v´ ystupu y(k) datov´ ym krokem upˇresnit stav x(k|k). Ladˇ en´ı KF se prov´ad´ı volbou kovarianˇcn´ıch matic Q,R. Obecnˇe se zvyˇsuje pomˇer tˇechto matic, jelikoˇz pˇri souˇcasn´em zv´ yˇsen´ı resp. sn´ıˇzen´ı obou matic se nic nedˇeje (kompenzuj´ı se). Pˇri sn´ıˇzen´ı pomˇeru Q/R ˇr´ık´am, ˇze vˇeˇr´ım stav˚ um syst´emu naopak pˇri zvyˇsov´an´ı pomˇeru vˇeˇr´ım v´ıce v´ ystupu. Vlastn´ı ˇ c´ısla matice A − L(t)C: Pˇri sniˇzov´an´ı zes´ılen´ı L se vlastn´ı ˇc´ısla pohybuj´ı od stabiln´ıch p´ol˚ u ke stabiln´ım nul´am (Chang-Letova vˇeta). Matice P : je kovarianˇcn´ı matice chyb odhadu. Jej´ım nastaven´ım v poˇca´teˇcn´ı f´azi P (0) urˇcujeme nakolik vˇeˇr´ıme modelu, mˇeˇren´ı popˇr´ıpadˇe jejich kovarianc´ım. Pokud vˇeˇr´ım modelu (stav˚ um), d´am malou hodnotu P (0), coˇz zajist´ı velk´ y roztpyl, jelikoˇz ˇsum znamen´a velk´ y rozptyl. 5
i
4
LQG regul´ ator, separaˇ cn´ı princip
Regul´ator LQG lze pouˇz´ıt na stochastick´ y syst´em popsan´ y n´asleduj´ıc´ım vztahem x(t + 1) = Ax(x) + Bu(t) + v(t), x(0) ≈ N (ˆ x(0), Px (0)
(45) (46)
Je moˇzn´e vych´azet ze dvou pˇr´ıstup˚ u ; je moˇzn´e mˇeˇrit stavy syst´emu jako zpˇetn´a vazba od stavu, nebo je rekonstruovat z mˇeˇren´eho v´ ystupu pomoc´ı Kalmanova filtru jako zpˇetnou vazbu od v´ ystupu. Zpˇ etn´ a vazba od stavu – stochastick´ y syst´ em: Vyuˇz´ıv´ame rozkladu stabu na stˇredn´ı hodnotu x = xˆ + x˜. Pot´e lze stˇredn´ı hodnotu kvadratick´e formy xT Qx vyj´adˇrit jako: n
o
n
o
n
o
n
o
n
o
E xT Qx = E (ˆ x + x˜)T Q(ˆ x + x˜) = E xˆT Qˆ x +E x˜T Q˜ x = xˆT Qˆ x+E trace˜ xx˜T Q = x˜T Q˜ x+trac (47) V´ ysledek t´eto u ´pravy ukazuje na to, ˇze stˇredn´ı hodnota kvadratick´eho krit´eria je oproti deterministick´e formulaci posunuta o vliv ˇsumu, kter´ y je d´an pr´avˇe stopou matice. Tud´ıˇz, rozep´ıˇseme-li ztr´atovou funkci je tento posun zˇrejm´ y: V ∗ (x(t), t) = xT P x +
N X
trace Qv P (k)
(48)
k=t+1
J ∗ = V ∗ (x(0), 0) = x(0)T P (0)x(0) +
N X
trace Qv P (t),
(49)
t+1
kde matice P (t) je ˇreˇsen´ı Ricattiho rovnice (totˇzn´e jako deterministick´e LQ ˇr´ızen´ı) a druh´a ˇca´st je vliv ˇsumu procesu. Prvn´ı ˇca´st rovnice pˇredstavuje tedy optim´aln´ı z´akon ˇr´ızen´ı (podobn´ y deterministick´e verzi) a druh´ y ˇclen odpov´ıd´a nepredikovan´e sloˇzce - ˇsumu v(t). Zpˇ etn´ a vazba od v´ ystupu – odhad stav˚ u: Tato ˇc´ast odpov´ıd´a formulaci Kalmanova filtru. Uvaˇzujeme, ˇze stavy syst´emu jsou nemˇeˇriteln´e. K tomu abychom data odadli, vyuˇz´ıv´am´e dostupn´ı minul´a data Dt−1 . Rekonstrukce stavu syst´emu lze vyj´adˇrit jako: xˆ(t + 1|t) = Aˆ x(t|t − 1) + Bu(t) + Lε(t|t − 1), ε(t|t − 1) = y(t) − yˆ(t|t − 1).
(50)
Stav Kalmanova filtru lze opˇet rozdˇelit na stˇredn´ı hodnotu xˆ a rozptyl x˜. Opˇet lze stˇredn´ı hodnotu kvadratick´e formy xT Qx lze vypoˇc´ıtat jako: n
o
E xT Qx = x˜T (t|t − 1)Q˜ x(t|t − 1) + trace Px (t|t − 1)Q
(51)
V n´asleduj´ıc´ım vztahu je vidˇet ˇclen Dt−1 , pˇredstavujc´ı pouˇz´ıv´an´ı dostupn´ ych dat v´ ystupu, u klasick´e zpˇetn´e vazby je zde ˇslen x(t) n
o
V ∗ (t) = min E x(t)T Q(t)x(t) + u(t)T R(t)u(t) + V ∗ (t + 1|u(t), Dt−1 ) u(t)
(52)
D´ale lze urˇcit hodnotu krit´eria a ztr´atov´e funkce jako: ∗
T
V (x(t), t) = x˜(t|t−1) P (t)˜ x(t|t−1)+
N X
T
trace LQε L P (k)+
k=t+1
N X
trace Px (k|k−1)Q(t)
k=t+1
(53) 6
kde prostˇredn´ı ˇclen rovnice pˇredstavuje vliv stochastick´e sloˇzky syst´emu (ˇsum procesu) a posledn´ı ˇclen pak vliv odhadu nam´ısto pˇr´ım´eho mˇeˇren´ı. N´asleduj´ıc´ı rovnice pˇredstavuje optim´aln´ı hodnotu krit´eria: ∗
∗
N X
T
J = V (x(0), 0) = x˜(0|−1) P (0)˜ x(0|−1)+
T
N −1 X
trace LQε L P (k)+
t+1
trace Px (t|t−1)Q(t)
t=0
(54) LQG regul´ ator - separaˇ cn´ı princip:
x(t + 1) = Ax(t) + Bu(t) + v(t) y(t) = Cx(t) + Du(t) + e(t) xˆ(t + 1|t) = Aˆ x(t|t − 1) + Bu(t) + Lε(t|t − 1) ˆ(y)(t|t − 1) = C xˆ(t|t − 1) + Du(t) u(t) = −K xˆ(t|t − 1) ε(t|t − 1) = y(t) − yˆ(t|t − 1) Uspoˇra´d´an´ı LQG regul´atoru:
Obr´azek 1: Blokov´e schema LQG regul´atoru
7
(55) (56) (57) (58) (59) (60)
h
x(t + 1) xˆ(t + 1|t)
i
"
=
A − BK BK 0 A − LC
#
h
i h
x(t) xˆ(t|t − 1) + v(t) v(t) − Le(t)
i
(61) xˆ(t|t − 1) − x(t) − xˆ(t|t − 1)
(62)
Separaˇ cn´ı princip: vlastn´ı ˇc´ısla uzavˇren´e smyˇcky: • A-BK - kvadraticky optim´aln´ı regul´ator • A-LC - Kalma˚ uv filtr
5
Metoda Loop Transfer Recovery
Pˇredpoklad je ˇze syst´em je minim´alnˇe f´azov´ y (nem´a nestabiln´ı nuly). Jelikoˇz LQG regul´ator pˇredstavuje kombinaci LQ a KF, pak se d´a pˇredpokl´adat interakce vlastnost´ı tˇechto syst´em˚ u. Rozd´ıln´e integr´aln´ı vlastnosti tˇechto syst´em˚ u maj´ı za n´asledek, ˇze se na vyˇsˇs´ıch frekvenc´ıch rozch´az´ı jejich vlastnsoti (LQ regul´ator m´a relativn´ı ˇra´d 1, KF m´a reltivn´ı ˇr´ad 2, tedy konvergence je moˇzn´a do jist´e frekvence. Tento pˇredpoklad vych´az´ı z citlivostn´ı funkce. LQ regul´ator je velmi robustn´ı, ovˇsem Kalman˚ uv filtr jiˇz tak robustn´ı nen´ı. Proto se snaˇz´ıme pˇribl´ıˇzit vlastnostem otevˇren´e smyˇcky LQ regul´atoru a to volbou parametru ρ resp. Q = Q0 + ρB T B. Pro rostouc´ı hodnotu ρ konverguje pˇrenos otevˇren´e smyˇcky LLQG (s) → LLQ (s). Tento fakt ilustruje n´asleduj´ıc´ı obr´azek:
Obr´azek 2: Princip Loop Transfer Recovery S rostouc´ı hodnotou ρ = 0, 10, 1000, 10000 doch´az´ı k ”‘vytlaˇcov´an´ı”’ kˇrivky z jednotkov´e kruˇznice a t´ım p´adem ke zlepˇsen´ı frekvenˇcn´ıch vlastnost´ıch LQG regul´atoru.
8
Reference [1] Havlena, V. (1999). Modern´ı teorie ˇr´ızen´ı – Doplˇ nkov´e skriptum. Vydavatelstv´ı ˇ CVUT, Praha. ˇ [2] Havlena, V. (2007). Modern´ı teorie ˇr´ızen´ı – Slide. CVUT, Praha. ˇ, J. Modern´ı teorie ˇr´ızen´ı [online]. Posledn´ı revize 2003-07-01 [3] Roubal, J.; Pekar [cit. 2003-07-01], hhttp://dce.felk.cvut.cz/mtr/i.
9
