Diskrétní waveletová transformace, zhlazování a prahování
Dana erná1, Václav Fin¥k2
Abstrakt. Waveletová analýza m·ºe být vyuºita p°i monitoro- vání, vyhodnocování a analýze strojových proces·. Nap°íklad p°i detekci závad p°evodového systému byl získáván vibra£ní signál z ak- celerometru. Tento signál byl dále transformován pomocí diskrétní waveletové transformace a následn¥ byly pouze dominantní wavele- tové koecienty p°edány neuronové síti ke klasikaci. Tímto zp·so- bem se výrazn¥ zmen²ilo mnoºství zpracovávaných dat a p°esto bylo moºné úsp¥²n¥ detekovat v²echny chybové stavy (viz. [8]). V tomto p°ísp¥vku se budeme v¥novat diskrétním waveletovým transforma- cím, zhlazování a prahování signálu. Vhodné prahování totiº umoº-
¬uje výrazn¥ zredukovat mnoºství zpracovávaných dat p°i minimální ztrát¥ informace.
Klí£ová slova: Wavelety, waveletová transformace, zhlazování a prahování signálu, neuronové sít¥, detekce závad.
1 Wavelety úvod
Díky velmi rychlému rozvoji se v posledních letech m·ºe zdát, ºe pojem wa- velet je pouºíván v mnoha významech, které spolu zdánliv¥ v·bec nesouvisí.
P°esto je moºné najít n¥kolik sty£ných bod·. Wim Sweldens, který se za- slouºil o roz²í°ení pojmu wavelet, v jednom ze svých £lánk· uvádí následující pokus o denici:
Wavelety jsou stavební bloky, které umoº¬ují rychle dekorelovat data.
1Podporována projektem MMT £íslo LC06024.
2Podporován projektem MMT £íslo 1M06047.
Pokusme se to nyní blíºe vysv¥tlit. Ozna£íme-li X prostor funkcí nebo signál·, potom jsou wavelety stavebními bloky pro tyto funkce. Matema- ticky to znamená, ºe tvo°í bázi tohoto prostoru. Tedy libovolný prvek z pro- storu X m·ºeme vyjád°it lineární kombinací wavelet·. Ozna£íme-li wavelety ψi, pak signál f ∈ X lze napsat ve tvaru
f =X
i
ciψi
Studovaný signál je nyní reprezentován posloupností koecient· ci.
Signály okolního sv¥ta nejsou náhodné hodnoty, ale vykazují jistou ko- relaci (závislost), a to jak v oblasti prostorové, tak i v oblasti frekven£ní.
Korelace je tém¥° vºdy lokální vzorky, které jsou blízko sebe, jsou více korelované neº vzorky, které jsou daleko od sebe. Podívejme se nap°íklad na obrázek. Body, které jsou blízko sebe, budou pravd¥podobn¥ siln¥ kore- lované (budou mít t°eba velmi podobnou barvu). Naopak lze t¥ºko o£ekávat závislost mezi body leºícími na opa£ných stranách obrázku.
Wavelety mají schopnost korelaci v datech odstra¬ovat. To znamená, ºe reprezentace pomocí wavelet· je oproti p·vodní reprezentaci kompaktn¥j²í.
V praxi se tato vlastnost projevuje tak, ºe mezi waveletovými koecienty ci
je jich pouze málo významných. Toho se s úsp¥chem vyuºívá nap°íklad p°i kompresi obrázk·. Aby wavelety mohly potla£ovat korelaci v signálu, m¥ly by být zpracovávanému signálu podobné (korelované stejným zp·sobem).
V¥t²inou tedy budeme poºadovat wavelety dob°e lokalizované jak v prosto- rové, tak ve frekven£ní oblasti. Lokalizovanost v prostorové oblasti znamená, ºe funkce má (tém¥°) kompaktní nosi£. Lokalizovanost ve frekven£ní oblasti odpovídá lokalizovanosti spektra tedy lokalizovanost Fourierovy transfor- mace waveletu. Pokles spektra sm¥rem k vysokým frekvencím odpovídá hladkosti waveletu. ím je wavelet hlad²í, tím rychlej²í je pokles. Pokles spektra sm¥rem k nízkým frekvencím koresponduje s po£tem nulových mo- ment· waveletu.
Dal²í vlastnost, kterou poºadujeme od wavelet·, je moºnost rychlých výpo£t·. Jde o to, abychom byli schopni rychle p°echázet mezi p·vodní a waveletovou reprezentací dat, protoºe rychlost je ve v¥t²in¥ aplikací klí-
£ová. Za rychlé algoritmy se v¥t²inou povaºují algoritmy, jejichº sloºitost je
O(N logN ), kde N je velikost vstupních dat. Waveletové algoritmy mají ob- vykle sloºitost O(N), zatímco nap°íklad diskrétní Fourierova transformace má sloºitost O(NlogN).
2 Wavelety
Waveletem budeme rozum¥t takovou funkci ψ ∈ L2(R),ºe systém {ψj,k}j,k∈Z jejích translací a dyadických dilatací tvo°í ortonormální bázi Hilbertova pro- storu L2(R).Funkce ψj,k je pro kaºdá dv¥ celá £ísla j, k denována p°edpi- sem
ψj,k = 2j/2ψ¡
2jx − k¢ .
Písmeno k reprezentuje translaci (posun) a 2j reprezentuje dyadickou dila- taci p·vodní funkce ψ. Wavelet je obvykle konstruován tak, aby waveletový systém {ψj,k}j,k∈Z tvo°il ortonormální systém v prostoru L2(R),ale jsou i jiné moºnosti. Ortonormalita systému znamená, ºe informace obsaºená ve waveletovém koecientu není obsaºena v ºádném jiném waveletovém koe- cientu a také ºe m·ºeme kompletn¥ zrekonstruovat p·vodní signál pomocí waveletových koecient·. Nebo-li pro kaºdou funkci (signál) f ∈ L2(R)platí
f =X
j∈Z
X
k∈Z
hf, ψj,ki ψj,k
(konvergence °ad je v norm¥ prostoru L2(R)a symbol hf, ψj,kizna£í skalární sou£in v prostoru L2(R)).
Základem konstrukce waveletového systému {ψj,k}j,k∈Z zaloºeného na translaci a dyadické dilataci jediné funkce ψ ∈ L2(R) je obvykle tzv. mul- tirozklad (multiresolution analysis). Není to jediný zp·sob konstrukce or- tonormální waveletové báze, nicmén¥ kaºdý wavelet s kompaktním nosi£em lze zkonstruovat pomocí multirozkladu. Multirozkladu se zde nebudeme dále v¥novat, pouze vyuºijeme n¥které jeho vlastnosti. Hlavní je existence tzv.
²kálové funkce φ ∈ L2(R). Podobn¥ jako pro wavelet i pro ²kálovou funkci denujme op¥t pomocí translace a dilatace
φj,k = 2j/2φ¡
2jx − k¢ .
Multirozklad je konstruován tak, ºe systém funkcí {φj0,k}k∈Z reprezentuje hladkou £ást signálu (aproximaci signálu) a p°idáváním waveletových sys- tém· na r·zných úrovních rozli²ení (²kálách) {ψj0,k}k∈Z, {ψj0+1,k}k∈Z, . . . postupn¥ p°idáváme dal²í detaily zkoumaného signálu. To znamená, ºe kaºdý signál f ∈ L2(R) m·ºeme jednozna£n¥ reprezentovat kombinací ²kálových a waveletových koecient·
f =X
k∈Z
hf, φj0,ki φj0,k+ X∞ j=j0
X
k∈Z
hf, ψj,ki ψj,k, (1)
pro libovolné celé £íslo j0.
Dále z vlastností multirozkladu plyne, ºe ²kálovou funkci φ m·ºeme popsat pomocí (²kálové rovnice) lineární kombinace jejich zmen²ených a posunutých verzí takto
φ(x) =X
k
hkφ(2x − k),
kde φ(2x − k) je zmen²ená verze ²kálové funkce φ posunutá o celé £íslo k a násobená ²kálovým koecientem hk. P°edchozí rovnice nám v podstat¥
°íká, ºe ²kálovou funkci na jedné ²kále m·ºeme spo£ítat pomocí n¥kolika
²kálových funkcí na p°edchozí ²kále, coº je základem diskrétní waveletové transformace. Podobným zp·sobem mohou být ²kálové koecienty pouºity rovn¥º k výpo£tu waveletu na jedné ²kále pomocí n¥kolika ²kálových funkcí na p°edchozí ²kále prost°ednictvím tzv. waveletové rovnice
ψ(x) =X
k
(−1)kh1−kφ(2x − k).
Tato konstrukce zaji²´uje, ºe wavelety a k nim p°íslu²né ²kálové funkce jsou navzájem ortogonální. Je²t¥ poznamenejme, ºe dále budeme uvaºovat pouze wavelety s kompaktním nosi£em, které mají pouze kone£ný po£et ²kálových koecient· r·zných od nuly.
P°íklad 1. Nejjednodu²²ím p°íkladem ortonormálního waveletu je Haarova funkce. Její ²kálová i waveletová rovnice obsahují pouze dva koecienty r·zné
od nuly
φ(x) = φ(2x) + φ(2x − 1) a ψ(x) = φ(2x) − φ(2x − 1).
e²ení ²kálové rovnice je zobrazeno na obrázku 3 je to první funkce zleva a následn¥ z waveletové rovnice dostaneme Haar·v wavelet to je druhá funkce zleva na obrázku 3. Obrázek 1 prezentuje Haarovu ²kálovou bázi pro j0= 2 ze vztahu (1) na intervalu [0, 1]. Tuto bázi m·ºeme ekvivalentn¥ p°e- vést na bázi, která obsahuje Haarovu ²kálovou bázi pro j0= 1 ze vztahu (1) a k tomu je²t¥ dva wavelety na úrovni rozkladu j = 1 ze vztahu (1). Tato báze je zobrazena na obrázku 2. T°etí moºnost, jak ekvivalentn¥ reprezento- vat tuto bázi, spo£ívá v uºití Haarovy ²kálové báze pro j0 = 0ze vztahu (1) (obsahuje jedinou funkci a sice ²kálovou funkci φ), k tomu musíme p°idat wavelet ψ odpovídající úrovni rozkladu j = 0 ze vztahu (1) a op¥t dva wavelety na úrovni rozkladu j = 1 ze vztahu (1). Tato báze je zobrazena na obrázku 3. Prost°ednictvím t¥chto bází lze p°esn¥ reprodukovat funkce, které jsou po £ástech konstantní na intervalech£
0,14¤ ,£1
4,12¤ ,£1
2,34¤ ,£3
4, 1¤ .
2
0,6 1,5
1
0,4 0,5
0 0,2 0
x 1 0,8
2
0,6 1,5
1
0,4 0,5
0 0,2 0
x 1 0,8
2
0,6 1,5
1
0,4 0,5
0 0,2 0
x 1 0,8
2
0,6 1,5
1
0,4 0,5
0 0,2 0
x 1 0,8
Obrázek 1: kálové funkce φ2,0, φ2,1, φ2,2, φ2,3.
3 Diskrétní waveletová transformace
V p°edchozí £ásti jsme vid¥li, ºe wavelety a ²kálové funkce na r·zných
²kálách spolu úzce souvisí. Toho nyní vyuºijeme p°i odvozování diskrétní
0,6 0,4 0,2 0 1,4 1,2 1 0,8 0,6 0,4 0,2 0
x 1
0,8 0 0,2 0,4 0,6
1,4 1,2 1 0,8 0,6 0,4 0,2 0
x 1 0,8
0,6 0,4 0,2 0 1
0,5
0
x -0,5
-1
1
0,8 0 0,2 0,4 0,6
1
0,5
0
x -0,5
-1
1 0,8
Obrázek 2: kálové funkce φ1,0, φ1,1 a wavelety ψ1,0, ψ1,1.
2
0,6 1,5
1
0,4 0,5
0 0,2 0
x 1 0,8
1
0,6 0,5
0 0,4
-0,5
-1 0,2 0
x 1
0,8 0 0,2 0,4 0,6
1
0,5
0
x -0,5
-1
1
0,8 0 0,2 0,4 0,6
1
0,5
0
x -0,5
-1
1 0,8
Obrázek 3: kálová funkce φ = φ0,0 a wavelety ψ = ψ0,0, ψ1,0, ψ1,1.
waveletové transformace. Je-li funkce (signál) f ∈ L2(R), potom ozna£me symbolem yj,k ²kálové koecienty funkce f a symbolem xj,k waveletové ko- ecienty funkce f a tedy platí
xj,k = hf, ψj,ki = Z ∞
−∞
f (x)ψj,k(x) dx,
yj,k = hf, φj,ki = Z ∞
−∞
f (x)φj,k(x) dx.
Dosadíme-li do t¥chto vztah· ²kálovou (resp. waveletovou) rovnici a po n¥- kolika jednoduchých úpravách dostaneme:
xj,k = Z ∞
−∞
f (x)ψj,k(x) dx = Z ∞
−∞
f (x)2j/2ψ¡
2jx − k¢ dx
= Z ∞
−∞
f (x)X
l
(−1)lh1−l2j/2φ(2(2jx − k) − l) dx
= X
l
(−1)lh1−l2j/2 Z ∞
−∞
f (x) φ(2j+1x − 2k − l) dx
= X
l
(−1)lh1−l 1
√2 Z ∞
−∞
f (x) 2(j+1)/2φ(2j+1x − 2k − l) dx
= 1
√2 X
l
(−1)lh1−lyj+1,2k+l a podobn¥ pro ²kálové koecienty
yj,k = Z ∞
−∞
f (x) φj,k(x) dx = Z ∞
−∞
f (x) 2j/2φ¡
2jx − k¢ dx
= Z ∞
−∞
f (x)X
l
hl2j/2φ(2(2jx − k) − l) dx
= X
l
hl2j/2 Z ∞
−∞
f (x) φ(2j+1x − 2k − l) dx
= X
l
hl 1
√2 Z ∞
−∞
f (x) 2(j+1)/2φ(2j+1x − 2k − l) dx
= 1
√2 X
l
hlyj+1,2k+l.
Tím jsme odvodili algoritmus pro rozklad (analýzu) signálu. Ale m·ºeme jít i v opa£ném sm¥ru a postupn¥ rekonstruovat rozloºený signál. Porovnáme- li rovnici (1) pro j0 = l s rovnicí (1) pro j0 = l + 1 dostaneme následující
rovnost X
k∈Z
yl,kφl,k+X
k∈Z
xl,kψl,k=X
k∈Z
yl+1,kφl+1,k.
Na funkce φl,k a ψl,k v p°edchozí rovnosti op¥t aplikujeme waveletovou a
²kálovou rovnici a obdrºíme X
k∈Z
yl,k 1
√2 X
m
h1−mφl+1,2k+m+X
k∈Z
xl,k 1
√2 X
m
(−1)mh1−mφl+1,2k+m=
=X
k∈Z
yl+1,kφl+1,k (2)
nebo-li po p°e£íslování vnit°ních °ad (abychom mohli snáze porovnávat ko- ecienty)
X
k∈Z
yl+1,kφl+1,k =
=X
k∈Z
yl,k 1
√2 X
m
hm−2kφl+1,m+X
k∈Z
xl,k 1
√2 X
m
(−1)mh1−m−2kφl+1,m. (3)
Nakonec (vzhledem k tomu, ºe jednotlivé posunuté ²kálové funkce jsou na- vzájem lineárn¥ nezávislé) porovnáme koecienty nap°íklad u funkce φl+1,k
a obdrºíme algoritmus pro rekonstrukci signálu yl+1,k = 1
√2 X
m
hm−2kyl,m+ 1
√2 X
m
(−1)mh1−m−2kxl,m.
P°íklad 2. P°ipome¬me si je²t¥ ²kálovou a waveletovou rovnici pro Haar·v wavelet
φ(x) = φ(2x) + φ(2x − 1) a ψ(x) = φ(2x) − φ(2x − 1).
Tedy oba p°íslu²né ²kálové koecienty jsou rovny jedné a první waveletový koecient je roven jedné a druhý mínus jedné. Ostatní koecienty jsou rovny nule. Budeme nyní p°edpokládat, ºe vektor [1, 2, 3, 4] reprezentuje koeci- enty Haarovy ²kálové báze pro j0 = 2 ze vztahu (1) na intervalu [0, 1]
a budeme chtít provést waveletovou transformaci tohoto vektoru pomocí
Haarova waveletu. Postupujeme podle vzorc· odvozených v p°edchozí £ásti a dostaneme
·1 + 2
√2 ,3 + 4
√2 ,1 − 2
√2 ,3 − 4
√2
¸
=
· 3
√2, 7
√2,−1
√2,−1
√2
¸ .
Nejprve jsme spo£etli sou£et prvního a druhého prvku zadaného vektoru a sou£et jsme vyd¥lili√
2, poté jsme spo£etli sou£et t°etího a £tvrtého prvku a sou£et jsme op¥t vyd¥lili √
2. A nakonec jsme spo£etli rozdíl prvního a druhého prvku a tento rozdíl jsme vyd¥lili √
2, poté jsme spo£etli rozdíl t°etího a £tvrtého prvku a rozdíl jsme op¥t vyd¥lili√
2. Tím jsme provedli jednu úrove¬ rozkladu, ale m·ºeme pokra£ovat s první polovinou vektoru
dále ·
3 + 7 2 ,3 − 7
2 ,−1
√2,−1
√2
¸
=
·
5, −2,−1
√2,−1
√2
¸ .
Dále jiº rozkládat nelze, takºe jsme spo£etli úplný rozklad zadaného vektoru pomocí Haarova waveletu.
4 Srovnání diskrétní waveletové a Fourierovy trans- formace
Uº jsme si °ekli, ºe waveletová transformace má sloºitost O(N), zatímco Fourierova transformace má sloºitost O(NlogN). Tedy waveletová transfor- mace je rychlej²í. Dal²ím rozdílem je, ºe waveletová transformace má lokální charakter zatímco Fourierova globální. Z toho plyne, ºe waveletová transfor- mace poskytuje lep²í aproximaci nespojitostí. Podívejme se na n¥kolik ob- rázk·. Na prvním z nich je zadaný signál. Na druhých dvou obrázcích byla na signál aplikována jednak Fourierova a jednak waveletová transformace a bylo ponecháno 20 v absolutní hodnot¥ nejv¥t²ích koecient· z 512. U Fou- rierovy transformace to znamenalo, ºe byly zanedbány v²echny koecienty, které byly v absolutní hodnot¥ men²í neº 8, 1079. Zatímco u waveletové transformace to znamenalo, ºe byly zanedbány v²echny koecienty, které byly v absolutní hodnot¥ men²í neº 0, 3087. Na t°etí dvojici obrázk· bylo
0 100 200 300 400 500 600 0
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Obrázek 4: Studovaný signál
Fourier D4 wavelet
0 100 200 300 400 500 600
−0.2 0 0.2 0.4 0.6 0.8 1 1.2
0 100 200 300 400 500 600
−0.2 0 0.2 0.4 0.6 0.8 1 1.2
Obrázek 5: Ponecháno 20 v absolutní hodnot¥ nejv¥t²ích koecient· z 512.
po transformacích ponecháno 60 v absolutní hodnot¥ nejv¥t²ích koecient·
z 512. U Fourierovy transformace to znamenalo, ºe byly zanedbány v²echny koecienty, které byly v absolutní hodnot¥ men²í neº 2, 7375. Zatímco u wa- veletové transformace to znamenalo, ºe byly zanedbány v²echny koecienty, které byly v absolutní hodnot¥ men²í neº 1, 1 × 10−4.V tomto p°ípad¥ tedy sta£ilo pouºít 60 nejdominantn¥j²ích koecient· k aproximaci studovaného signálu a velikost zbylých 552 koecient· byla zanedbateln¥ malá.
Fourier D4 wavelet
0 100 200 300 400 500 600
−0.2 0 0.2 0.4 0.6 0.8 1 1.2
0 100 200 300 400 500 600
−0.2 0 0.2 0.4 0.6 0.8 1
Obrázek 6: Ponecháno 60 v absolutní hodnot¥ nejv¥t²ích koecient· z 512.
5 Zpracování signálu
Dosud jsme uvaºovali diskrétní transformaci funkce f a ukázali jsme si, jak lze tuto funkci reprezentovat pomocí waveletové báze nebo pomocí kombi- nace ²kálové a waveletové báze. Nyní budeme p°edpokládat, ºe máme dis- krétní vstupní signál kone£né délky, která je p°irozenou mocninou dvojky.
Obvykle se p°edpokládá, ºe jednotlivé hodnoty vstupního signálu odpovídají p°íslu²ným ²kálovým koecient·m. (Nebo si také m·ºeme vstupní signál p°edstavit jako vektor, na který aplikujeme lineární (waveletové) transfor- mace.) Nyní m·ºeme stejným zp·sobem, který jsme popsali v p°edchozí
£ásti, provést rozklad diskrétního signálu. Jakmile máme provedený úplný rozklad signálu (tedy rozklad aº na nejniº²í moºnou úrove¬), m·ºeme p°ed provedením rekonstrukce signálu zm¥nit libovolný koecient. M·ºeme na- p°íklad celé skupiny koecient· nahradit nulou, nebo nahradit pouze n¥- které koecienty nulou, nebo m·ºeme zm¥nit velikost n¥kterých koecient·.
Prost¥ m·ºeme s koecienty manipulovat mnoha r·znými zp·soby v závis- losti na ú£elu, kterého chceme dosáhnout. Budeme-li dále v textu mluvit o rozkladu signálu, budeme mít vºdy na mysli úplný rozklad. P°i analýze signálu p°edpokládáme, ºe velké waveletové koecienty lépe charakterizují podstatu studovaného signálu neº malé waveletové koecienty.
5.1 Volba waveletu
D·leºitým krokem je volba konkrétního waveletu, které budeme pouºívat pro rozklad signálu. Uve¤me si n¥kolik obecných doporu£ení. Víme, ºe my²- lenkou waveletového rozkladu je aproximace signálu ²kálovou funkcí a vyjá- d°ení detail· prost°ednictvím wavelet·. Pro hladké signály je tedy vhodné vybírat hladké ²kálové funkce. Pokud o£ekáváme, ºe signál bude obsahovat jemné detaily, je rozumné zvolit ²kálovou funkci s krátkým nosi£em, aby nedocházelo k rozmazání signálu. Zvlá²tní p°ípad nastává, m·ºeme-li p°ed- pokládat, ºe zpracovávaný signál je polynom stupn¥ k. Potom je nejlep²í zvolit wavelet, který má k + 1 nulových moment·, protoºe p°íslu²ná ²ká- lová funkce aproximuje polynomy do stupn¥ k p°esn¥. P°i rozkladu signálu pomocí tohoto waveletu dostaneme v²echny waveletové koecienty rovny nule. Tuto metodu m·ºeme vyuºít i s p°edpokladem, ºe p·vodní signál je polynomiální pouze na ur£itém intervalu [a, b]. Pak ov²em nulujeme pouze koecienty p°íslu²né wavelet·m, jejichº celý nosi£ je v tomto intervalu ob- saºen.
5.2 Zhlazování signálu
Zhlazování signálu je zaloºeno na skute£nosti, ºe zatímco waveletové koeci- enty signálu reprezentují jednotlivé detaily, tak ²kálové koecienty p°edsta- vují aproximaci zkoumaného signálu. Volbou waveletu a p°íslu²né ²kálové funkce lze ovlivnit dosaºitelnou hladkost £ásti signálu reprezentovanou ²ká- lovými koecienty (hladkost podstatným zp·sobem závisí na hladkosti zkou- maného signálu). Tato závislost funguje tak, ºe pokud si k analýze signálu vybereme nap°íklad wavelet se t°emi nulovými momenty, potom lze v p°í- pad¥, ºe by signál byl tvo°en konstantní, lineární nebo kvadratickou funkcí, tento signál reprodukovat p°esn¥ pouze na základ¥ ²kálových koecient·. V opa£ném p°ípad¥ budou ²kálové koecienty reprodukovat pouze polynomi- ální £ást (do stupn¥ dva) signálu. Zhlazování signálu tedy provedeme tak, ºe v úplném rozkladu ponecháme ²kálové koecienty plus waveletové koe- cienty do ur£ité ²kály a ostatní waveletové koecienty reprezentující jemné detaily nahradíme nulami (scale-dependent smoothing).
5.3 Prahování
Zmen²ování velikosti koecient· se v¥t²inou vyuºívá k odstran¥ní ²umu ze signálu. Zahrnuje zmen²ení nebo odstran¥ní vybraných waveletových koe- cient·. K dispozici je nep°eberné mnoºství metod pro výb¥r a modikaci koecient·. Asi nejpopulárn¥j²í jsou m¥kké a tvrdé prahování. Tvrdé pra- hování spo£ívá v tom, ºe v²echny waveletové koecienty, které jsou men²í neº zvolený práh λ, nahradíme nulami. To m·ºeme formáln¥ zapsat takto
xTj,k =
( 0 pro xj,k < λ, xj,k pro xj,k ≥ λ.
M¥kké prahování oproti tvrdému zohled¬uje, ºe v²echny waveletové koe- cienty obsahují jak signál tak ²um a pokou²í se odstranit ²um ze v²ech ko- ecient·. P°i m¥kkém prahování tedy v²echny koecienty, které jsou men²í neº zvolený práh λ, nahradíme nulami a ostatní posuneme sm¥rem k nule o velikost prahu
xMj,k =
( 0 pro xj,k < λ,
signum(xj,k)(|xj,k− λ|) pro xj,k≥ λ.
Nevýhodou m¥kkého prahování je, ºe p°i n¥m dochází k v¥t²í ztrát¥ energie (euklidovská norma na druhou) neº p°i tvrdém prahování.
Problém je, ºe v praxi v¥t²inou neznáme p°esný tvar signálu ani tvar zkreslujícího ²umu. Volba prahu je proto netriviální a k jeho vhodnému stanovení m·ºeme pouºít r·zná kritéria zaloºená na základ¥ signálu a (nebo)
²umu. Práh m·ºe být nap°íklad konstantní hodnota aplikovaná na v²echny koecienty, nebo jen na koecienty na n¥kterých ²kálách, nebo m·ºe být velikost prahu r·zná pro r·zné ²kály.
Velkého pokroku v oblasti potla£ení ²umu dosáhli Donoho a Johnstone ([2], [3]). Navrhli volbu prahu λ, která je zaloºena na následujícím tvrzení.
V¥ta 3. Nech´ e1, . . . , eM je posloupnost nezávislých náhodných veli£in,
které mají rozd¥lení N(0, σ2). Pak platí
M →∞lim P µ
maxi |ei| > σ√ 2 ln M
¶
→ 0.
A protoºe u ²umu p°edpokládáme práv¥ normální rozd¥lení N(0, σ2),
°íká p°edchozí v¥ta, ºe velikost ²umu v signálu bude pravd¥podobn¥ men²í neº σ√
2 ln M . Odtud vychází volba prahu λ = σ√
2 ln M . Protoºe orto- gonální waveletová transformace nemá vliv na rozd¥lení ²umu, z·stává i hodnota prahu stejná na v²ech úrovních rozkladu (viz. [10] strana 169).
Otázkou ov²em z·stává, jak ur£it rozptyl ²umu. Pokud nemáme k dispo- zici ºádné apriorní informace o moºném rozptylu, je t°eba n¥jakým zp·so- bem stanovit odhad. Praktické experimenty ukazují, ºe waveletové koeci- enty na nejjemn¥j²í úrovni (nejjemn¥j²í detaily) odpovídají tém¥° výhradn¥
²umu. Z toho d·vodu je odhad rozptylu v¥t²inou po£ítán práv¥ z t¥chto koecient·, nap°íklad jako jejich standardní odchylka. Donoho a Johnstone pouºili ve své práci [3] následující empirický odhad spo£ítali medián ab- solutních hodnot nejjemn¥j²ích waveletových koecient·, výsledek vyd¥lili hodnotou 0,6745 a tuto hodnotu pouºili jako odhad rozptylu ²umu.
Na následujících obrázcích je uveden p°íklad na odstran¥ní ²umu pomocí tvrdého prahování. Vzhledem k tomu, ºe studovaný signál je velmi hladký, je
0 100 200 300 400 500 600
−5
−4
−3
−2
−1 0 1 2 3 4 5
Obrázek 7: P·vodní signál plus bílý ²um s jednotkovým rozptylem.
vhodné pouºít wavelet, který má více nulových moment·. V tomto p°ípad¥
jsme zvolili wavelet D14, který má 7 nulových moment·. Byly uºity dv¥
volby prahu. P°i první volb¥ jsme vyuºili toho, ºe známe velikost rozptylu (jedna). Potom jsme dostali
λ1 = σ√
2 ln M =√
2 ln 512 ≈ 3, 532230068.
(512 je délka signálu.) V druhém p°ípad¥ jsme vyuºili vý²e uvedený odhad
0 100 200 300 400 500 600
−2
−1.5
−1
−0.5 0 0.5 1 1.5 2
0 100 200 300 400 500 600
−2
−1.5
−1
−0.5 0 0.5 1 1.5 2
Obrázek 8: P·vodní signál a signál po odstran¥ní ²umu pomocí waveletu D14.
rozptylu ²umu λ2 = σ√
2 ln M ≈ M ED 0, 6745
√2 ln 512 ≈ 3, 296049822.
MED je ozna£ení pro medián absolutních hodnot nejjemn¥j²ích waveleto- vých koecient·. Vidíme, ºe v tomto p°ípad¥ jsou oba odhady prahu po- dobné a oba odhady vedly k ponechání 12 nejdominantn¥j²ích koecient·.
5.4 Metody redukce dat
Technologický pokrok vytvá°í pro rmy i pro jednotlivce p°ístup k cenným informacím o výrobním procesu (o stavu automobilu, ...), které mohou být vyuºity ke zlep²ování kvality, zvy²ování efektivnosti nebo k v£asné detekci závad. S r·stem mnoºství dat rostou nároky na jejich zpracování a analýzu.
Z toho d·vodu je nezbytné vyuºívat techniky redukce dat. Jednou z moº- ností jak sladit protich·dné poºadavky na p°esnost modelu a na redukci
dat je i diskrétní waveletová transformace. Její nejjednodu²²í vyuºití spo-
£ívá v tom, ºe p°i analýze získaného signálu vyuºijeme jen n¥kolik procent (nebo dokonce jen n¥kolik desítek) nejdominantn¥j²ích waveletových koe- cient·. Dal²í moºností je vyuºít k redukci n¥kterou d°íve uvedenou (nebo podobnou) metodu na zhlazování nebo prahování. Nedávno se objevilo n¥- kolik zajímavých prací v tomto sm¥ru. První z nich je ²kálogram. kálogram m¥°í energii signálu na r·zných úrovních rozkladu (²kálách). V práci [5] byl navrºen ²kálogram s prahováním a bylo odvozeno asymptotické rozd¥lení
²kálogramu s prahováním. Dal²í moºností popsanou v £lánku [6] jsou odvo- zené metody pro odhad vhodného prahu, které nevyºadují odhad standardní odchylky.
6 Diagnostika poruch hnací soustavy
6.1 Popis problému
V¥t²ina strojních za°ízení obsahuje motor a jiné rotující £ásti. Aby p°i posu- zování stavu za°ízení nebylo t°eba toto za°ízení nejprve rozebrat, vyuºívají se k analýze stavu vnit°ních £ástí r·zné signály získané zevnit° stroje. U ro- tujících za°ízení se k diagnostice poruch obvykle vyuºívá vibra£ní signál. Na obrázku 9 je schéma hnací soustavy, p°i£emº k získávání vibra£ního signálu byl vyuºit akcelerometr.
Na získaný vibra£ní signál byla aplikována diskrétní waveletová trans- formace (s Daubechiovským waveletem D4) a cílem bylo detekovat ²est p°í- pad·:
A: p°evodovka i loºiska jsou v po°ádku,
B: loºiska jsou v po°ádku a závada £íslo 1 na p°evodovce, C: loºiska jsou v po°ádku a závada £íslo 2 na p°evodovce, D: vadná loºiska a p°evodovka je v po°ádku,
E: vadná loºiska a závada £íslo 1 na p°evodovce,
Obrázek 9: Schéma hnací soustavy F: vadná loºiska a závada £íslo 2 na p°evodovce.
P°íklady pro v²echny uvedené situace jsou uvedeny na následujících ²esti obrázcích. Z kaºdého signálu bylo vybráno deset nejdominantn¥j²ích koe-
cient· a jejich hodnoty a po°adí byly zvoleny jako vstup pro vícevrstvou neuronovou sí´.
6.2 Vícevrstvá neuronová sí´ - metoda backpropagation Vícevrstvá neuronová sí´ s algoritmem zp¥tného ²í°ení chyby (backpropa- gation) je nejznám¥j²í a nejpouºívan¥j²í model neuronové sít¥ a vyuºívá se
ve v¥t²in¥ aplikací. Sí´ má zpravidla n vstupních neuron·, m výstupních neuron· a jednu nebo dv¥ skryté vrstvy. Neurony v jednotlivých vrstvách jsou spojeny perceptrony. Aktiva£ní funkce spojitého perceptronu m·ºe být
nap°íklad sigmoida denovaná vzorcem:
yj = y(ξj) = 1 1 + e−λξj, kde
• yj jsou aktiva£ní funkce,
• xi jsou vstupní data (p°ípadn¥ aktiva£ní funkce v p°ípad¥ druhé a dal²ích vrstev),
• λje parametr strmosti sigmoidy ur£uje nelineární nár·st sigmoidy v okolí bodu nula (míru rozhodnosti neuronu); obvykle se uvaºuje λ = 1,
Obrázek 10: Schéma dvouvrstvé neuronové sít¥.
• ξj =P
iwi,jxi+ bj,
• bj jsou hledané prahové koecienty,
• wi,j jsou hledané váhové koecienty.
Dal²í moºnosti volby aktiva£ní funkce p°edstavují nap°íklad skoková funkce, lineární funkce a podobn¥.
Obrázek 11: Sigmoida (λ = 1).
Stav kaºdého neuronu se pohybuje v rozmezí od 0 do 1. Pokud se hod- nota na vstupu neuronu blíºí k ∞, blíºí se výstupní hodnota k 1. Pokud se hodnota vstupu blíºí k −∞, výstup se blíºí k 0. Spoj od neuronu i k neuronu j je ohodnocen váhou wi,j.
U£ení sít¥ metodou backpropagation spo£ívá v minimalizaci celkové chy- by pro v²echny vzory z trénovací mnoºiny, p°i£emº chyba vzoru p°edstavuje st°ední kvadratickou odchylku vektoru poºadovaných a vypo£tených hod- not, tj. matematicky
E(w) = 1 2
X
i
(y(xi) − di)2,
kde di jsou hodnoty poºadovaných výstup·. P°ed u£ením jsou váhy a pra- hové koecienty nastaveny náhodn¥. Váhy potom upravujeme nap°íklad tak, ºe od váhy v p°ede²lém kroku ode£teme ur£itý díl derivace chyby podle této váhy reprezentovaný koecientem u£ení. P°itom postupujeme od horních vrstev ke spodním - od toho byl odvozen název metoda zp¥tného
²í°ení chyby. Podrobnosti lze nalézt nap°íklad v [9] nebo v [11].
6.3 Výsledky neuronové sít¥
Pro u£ení i pro testování bylo zvoleno 48 vzork· vibra£ních signál·. Kaºdý z nich obsahoval 6 vzork· pro kaºdou situaci. Na vstupu bylo 20 neuron·
a na výstupu 6 neuron·. Neuronová sí´ byla trénována tak, aby pro u£ební vzorek produkovala na p¥ti uzlech nulu a na jednom jedni£ku. P°i£emº pro kaºdý klasikovaný p°ípad byl rezervován jeden výstupní neuron (jedni£ka na prvním výstupním neuronu potom signalizovala, ºe nastává p°ípad A - tedy ºe p°evodovka i loºiska jsou v po°ádku). Topologie skryté vrstvy byla ur£ena obvyklou metodou pokus· a omyl· (optimální bylo 14 neuron· ve skryté vrstv¥). Pr·m¥rná úsp¥²nost p°i klasikaci testovacích vzork· byla 96%. Nam¥°ené vibra£ní signály a výsledky neuronové sít¥ byly £erpány z [8].
7 Záv¥r
Diskrétní (rychlá) waveletová transformace se ukázala být velmi ú£inným nástrojem p°i analýze signál· (které nap°íklad nejsou periodické, obsahují
²um, jsou p°eru²ované, ...). Mezi nejuºite£n¥j²í vlastnosti pat°í rychlé al- goritmy, rozklad signálu na aproxima£ní £ást (reprezentovanou ²kálovými koecienty) a detaily (reprezentovanou waveletovými koecienty), p°i£emº se ukazuje, ºe v¥t²ina waveletových koecient· je malých (zanedbatelných).
Zanedbáním malých koecient· (redukcí dat) získáme kompaktn¥j²í tvar signálu, který m·ºe být vyuºit nap°íklad k úsp¥²né detekci závad hnací soustavy.
Reference
[1] P. Addison: The Illustrated Wavelet Transform Handbook, Bristol and Philadelphia, Institute of Physics Publishing, 2002.
[2] D. L. Donoho: De-noising by soft-thresholding, IEEE Transactions on Information Theory 41(3), 613-627, 1995.
[3] D. L. Donoho, I. M. Johnstone: Ideal Spatial Adaptation via Wavelet Shrinkage, Biometrika 81(3), 425-455, 1994.
[4] M. W. Frazier: An Introduction to Wavelets through Linear Algebra, Springer, 1999.
[5] M. K. Jeong, D. Chen, J. Ch. Lu: Thresholded scalogram and its appli- cations in process fault detection, Applied stochastic models in business and industry, 19(3), 231-244, 2003.
[6] M. K. Jeong, J. Ch. Lu, B. Vidakovic, D. Chen: Wavelet-Based Data Reduction Techniques for Process Fault Detection, Technometrics, 48(1), 26-40, 2006.
[7] K. Najzar: Základy teorie wavelet·, Praha, Karolinum, 2004.
[8] B. A. Paya, I. I. Esat, M. N. M. Badi: Articial Neural Network Based Fault Diagnostics of Rotating Machinery Using Wavelet Transforms as a Preprocessor, Mechanical Systems and Signal Processing, 11(5), 751- 765, 1997.
[9] J. íma, R. Neruda: Teoretické otázky neuronových sítí, Praha, MAT- FYZPRESS, 1996.
[10] B. Vidakovic: Statistical Modeling by Wavelets, New York, John Wiley
& Sons., 1999.
[11] I. Vondrák: Um¥lá inteligence a neuronové sít¥, Ostrava, VB TU, 1994.
Adresy autor·:
RNDr. Václav Fin¥k, Ph.D., Mgr. Dana erná, Technická univerzita Libe- rec, fakulta pedagogická, katedra matematiky a didaktiky matematiky, Hál- kova 6, 461 17 Liberec 1.
e-mail: vaclav.nek@tul.cz, dana.cerna@tul.cz
Tato práce byla vytvo°ena v rámci projekt· MMT £íslo LC06024 a 1M06047.