Předzpracování pohybového EEG pomocí ICA a klasifikace nezávislých komponent Lukáš Ručkay České vysoké učení technické v Praze, Fakulta elektrotechnická
[email protected] Abstrakt: Tato práce se zabývá metodou pro předzpracování EEG s cílem zvýšit klasifikační skóre pohybového EEG. Předzpracování EEG je provedeno pomocí analýzy nezávislých komponent (independent component analysis, ICA). Podle korelace mezi spektrogramy nalezených komponent a vzory vytvořenými průměrováním spektrogramů z elektrod pokrývající senzomotorickou oblast jsou vybrány komponenty související s pohybem. Tyto pohybové komponenty jsou klasifikovány systémem založeným na skrytých Markovských modelech (hiden Markov model, HMM). Systém dosahuje klasifikační skóre kolem 80% na čistém EEG a přes 95% na pohybových komponentách.
1.
Úvod
Klasifikace EEG je úloha, která nachází uplatnění ve dvou hlavních oblastech: • Klasifikace EEG je důležitá část systému rozhraní mezi mozkem a počítačem (braincomputer interface, BCI). BCI představuje alternativní metodu pro komunikaci člověka s počítačem a jeho použití je významné zejména pro lidi kteří nedokáží ovládat počítač standardní cestou (klávesnice, myš). Základní charakteristiky BCI systémů jsou shrnuty v [3] a [10]. • Druhou významnou oblast použití klasifikace EEG tvoří ověřování fyziologických hypotéz o funkci mozku a psychologie. Naším cílem v oblasti zpracování biologických signálů je zvýšení rozpoznávacího skóre v databázi EEG obsahující odlišné pohyby (proximální-distální, extenzní-flexní nebo palecmalíček databáze). Aplikace takovéhoto klasifikátoru je určena pro BCI systém. V předchozích prací [2] jsme vyvinuli klasifikátor založený na skrytých Markovských modelech (hiden Markov model, HMM) který dosahuje skóre okolo 70%. Vyšší skóre je možné dosáhnout dekompozicí EEG na pohybové a nepohybové komponenty s využitím slepé separace signálů (blind source separation, BSS), [5].
2.
Vlastnosti EEG databáze
Pro naše testy jsme si vybrali databázi EEG obsahující proximální a distální pohyby, kterou jsme získali od prof. Stančáka z III. LF UK. Podrobný popis EEG, jeho parametrů i procedury nahrávání je obsažen v [9].
Databáze obsahuje EEG záznam 7 osob provádějících v nepravidelných intervalech (12s 15s) proximální pohyb (zvednutí pravého nadpažku) a distální pohyb (rychlý ohyb pravého ukazováčku) [1], [9]. Pro záznam byl využit 64-kanálový systém, 59 kanálů bylo použito pro elektrody na hlavě podle systému 10-10. Dále byly užity 2 kanály pro EMG, 2 kanály pro EOG a 1 kanál pro ECG. Vzorkovací frekvence byla 500Hz. Zemnící elektroda byla umístěna na čele a jako referenční elektroda byla použita elektroda Cz. Pohyby byly monitorovány užitím citlivého akcelerometru, který byl přidělán na nehet prstu určeného pro distální pohyb a na vršek pravého ramene pro monitorování proximálního pohybu [9]. Do EEG záznamu byly ručně vloženy značky označující začátek pohybů podle EMG aktivity. Pohyby, které byly zasaženy očními nebo svalovými artefakty byly vyloučeny z následné analýzy. Pro další zpracování byly použity pouze pohyby se stejným tvarem. Z dat byla odstraněna stejnosměrná složka pomocí interpolace a následného odečtení polynomu sedmého řádu v intervalu jedné minuty. Dalším krokem bylo rozdělení EEG na desetisekundové úseky, vždy pět sekund před a pět sekund po pohybu. Počty realizací obou typů pohybů od všech pokusných osob se liší, avšak vždy pokusná osoba vykonala minimálně 27 pohybů.
3.
HMM klasifikátor
V této kapitole budou popsány základní vlastnosti HMM klasifikátoru, který již byl použit v předchozích pracích [1] a [2]. Prováděný výzkum je primárně zaměřen na off-line, single-trial klasifikaci EEG. Volba offline, single-trial systém zjednodušuje návrh systému ve dvou bodech: off-line klasifikace probíhá na dříve nahraných realizacích, a zpracování nemusí být prováděno v reálném čase; single-trial klasifikace neklasifikuje spojité EEG, ale EEG segmentované na úseky (jednotlivé realizace). K použití HMM pro klasifikaci EEG nás vedlo několik důvodů: • HMM poskytují inherentní podporu pro zachycení časových vazeb mezi prvky parametrizovaného signálu. Běžně používané neuronové sítě nejsou tohoto schopny. • Pokud jsou modely HMM klasifikátoru naučeny, lze je jednoduše použít i pro generování průběhů s danými statistickými vlastnostmi (generování syntetického EEG). • Struktura klasifikátoru může být navržena s ohledem na skutečnou fyziologickou podstatu časového vývoje EEG v mozku. Základní struktura HMM klasifikátoru je naznačena na obr. 1. Jako parametrizace je zvoleno lineární spektrum v intervalu 5–40Hz pro zachycení typických projevů pohybů v EEG. Těmi jsou desynchronizace (event related desynchronization, ERD) v průběhu pohybu a synchronizace (event related sznchronization, ERS) po pohybu [7], [9]. Výsledkem předchozích prací [1], [2] je optimální nastavení klasifikátoru, které zde zachováme: • parametrizace: FFT, lineární spektrum v intervalu 5–40Hz, • velikost segmentu: 500vz (1s), • posuv segmentů: 100vz (200ms), • okno: hanning, • počet učících iterací: 10.
100
50
0
−50
−100
µ ERD
TICHO
β ERS
β ERD
TICHO
−150
−200
−250 0
START
500
1000
1500
TICHO
2000
µ ERD β ERD
2500
3000
β
3500
ERS
4000
4500
TICHO
5000
KONEC
Obrázek 1: Struktura HMM klasifikátoru koresponduje s fyziologickým charakterem zpracovávaného EEG. Jednotlivé intervaly pohybového EEG korespondují se stavy modelu tak, jak je znázorněno na obrázku. Převzato z [1].
4.
Slepá separace signálů
Slepá separace signálů (Blind Source Separation, BSS) [5] je obor, který se zabývá získáváním skrytých informací obsažených v souboru dat. Analýza nezávislých komponent (Independent Component Analysis, ICA) je jedna z mnoha metod pro řešení tohoto problému. ICA aplikovaná na EEG data je velice zajímavá a perspektivní metoda získávání informací ukrytých v mozku [8]. ICA je metoda, která se snaží nalézt zdrojové signály ze znalosti jejich kombinace. Přitom předpokládáme, že zdrojové signály nejsou pozorovatelné, a že nemáme žádné informace o tom jakým způsobem jsou zdrojové signály transformovány. V oblasti slepé separace signálů se setkáme se třemi základními modely popisující transformaci: lineární model, šumový model a konvolutorní model. EEG i EMG data splňují lineární model [6], [8], který můžeme popsat takto: X = A · S, kde S je m × Nd matice vstupních signálů (nezávislé komponenty), A je mixážní matice rozměru n × m a X je n × Nd matice pozorovaných (měřených) signálů. Nd tedy značí délku jednotlivých řádků matic S a X. Počet měřených signálů a počet nezávislých komponent nemusí být shodný, avšak předpokladem pro úspěšnou separaci je n ≥ m. Zdrojové signály (nezávislé komponenty) mohou být odhadb =W c · X, kde W c je odhad matice W, W = A−1 . Problém nalezení nuty následovně: S zdrojových signálů (nezávislých komponent) přechází na problém nalezení matice W. Separace či rozklad signálů je prováděn tak, aby výsledné komponenty byly nezávislé jak jen je to možné a zároveň jejich lineární kombinace tvořila pozorované signály. To je dosaženo maximalizací funkce F , která měří nezávislost komponent.
5.
Výběr a klasifikace nezávislých komponent
Z EEG záznamu byla odstraněna stejnosměrná složka a následně bylo EEG rozloženo do 59 nezávislých komponent. Rozklad byl proveden pomocí algoritmu FastICA s následujícím nastavením: symmetric verze algoritmu, nelinearita tanh, stabilizace on, [4]. Aby byly pro všechny osoby zajištěny stejné podmínky, použili jsme pro rozklad jen prvních 27 realizací pohybů. Výběr komponent které souvisí s pohybem není snadná úloha a existuje více přístupů jak tento problém řešit. My jsme zvolili korelaci mezi spektrogramy komponent a vzory. Tyto vzory jsou tvořeny průměrováním spektrogramů z elektrod pokrývající senzomotorickou oblast. Na základě spektrální analýzy čistého EEG jsme zvolili elektrody, které jsou použity pro tvorbu vzorů. Jedná se o elektrody z centrální a parietální oblasti: 25–29, 34–38 a 43–47 (C3–C4, CP3–CP4 a P3–P4) podle systému 10-10, [9]. Tyto vybrané elektrody obsahují typické znaky pro pohybové EEG – desynchronizaci během pohybu (ERD) neboli pokles energie v daném frekvenčním pásmu a synchronizaci po pohybu (ERS) neboli opětovný nárůst energie v daném frekvenčním pásmu. Zvolené elektrody rovněž obsahují nemalé množství šumu a rušivých signálů které necharakterizují pohyb. Z tohoto důvodu je zvoleno několik elektrod a průměrování spektrogramů s cílem omezit rušivé složky ve vzorech. ERD a ERS se projevují během pohybu zejména v µ pásmu (8–13Hz) a β pásmu (13–30Hz) tak jak ukazuje obr. 2. IC 39, (mean=0.935, std=0.052), ref. matrix: 1
IC 40, (mean=0.955, std=0.030), ref. matrix: 1
40
40
1.6
35
2.4 2.2
35
2 1.4
1.2 25
1 20
30
Frequency [Hz]
Frequency [Hz]
30
1.8 1.6
25 1.4 20
1.2
0.8 15
1 15 0.8
0.6 10
10
0.6
0.4 5
−4
−3
−2
−1
0 Time [s]
a)
1
2
3
4
0.4 5
−4
−3
−2
−1
0 Time [s]
1
2
3
4
b)
Obrázek 2: Spektrogramy pohybových komponent. Svislá čára ukazuje místo pohybu. ERD je viditelné během pohybu v µ i β pásmu a ERS po pohybu v µ a β pásmu. a) distální pohyb, b) proximální pohyb. Korelace mezi jednotlivými komponentami a vzory se značně liší nejen mezi jednotlivými osobami, ale také v rámci souboru komponent pro danou osobu. Z tohoto důvodu nelze stanovit pevnou mez určující, které komponenty jsou pohybové a které nikoli. Nehledáme tedy komponenty jejichž korelace překračuje určitou mez, ale pro klasifikaci použijeme jen několik prvních komponent s největší podobností se vzorem (s největší korelací). Pro naše testy jsme použili pouze 3 páry komponent pro klasifikaci – každý pár obsahuje komponentu pro distální pohyb a komponentu pro proximální pohyb.
6.
Výsledky
V této sekci uvedeme výsledky klasifikace pohybových komponent. Pro srovnání jsou v následující kapitole rovněž uvedeny výsledky klasifikace čistého EEG. Ty považujeme za referenční a s nimi budeme porovnávat výsledky klasifikace pohybových komponent. Během procesu klasifikace je vždy množina vstupních dat rozdělena na trénovací množinu a testovací množinu a to tak, že obě množiny obsahují po 50% vstupních dat. V předchozím textu jsme zmínili, že minimální počet realizací je 27 (tj. 13 realizací pro trénovaní a 14 realizací pro testovaní). Pro dosažení spolehlivých výsledků klasifikace je však potřeba mnohem větší množina dat a proto je klasifikace pro daný soubor dat několikrát opakována. Pro každé opakování je soubor všech realizací náhodně rozdělen mezi trénovací a testovací data. Jednoduchým testem konzistentnosti výsledků jsme získali potřebný počet opakování, ten je cca 100×. Při tomto počtu opakování se liší průměrné klasifikační skóre pro různé dělení vstupních dat cca o ±2%. 6.1.
Klasifikace čistého EEG
Klasifikaci čistého EEG jsme provedli pro elektrody z oblasti senzomotorické oblasti: elektrody 16–22, 24–30, 33–39, 42–48 (FC5–FC6, C5–C6, CP5–CP6, P5–P6), [9]. Do klasifikace byly zahrnuty všechny realizace pro danou osobu (tj. od 27 do 100 realizací). Z nedostatku místa zde neuvedeme všechny výsledky, ale jen výsledky z 5 elektrod s nejlepším skóre pro každou osobu. Výsledky jsou shrnuty v tab. 1. osoba 1d 1p 2d 2p 3d 3p 4d 4p 5d 5p 6d 6p 7d 7p
85.9±07.9 72.1±11.9 87.3±05.5 85.8±07.3 68.0±15.0 71.3±12.1 86.8±09.4 96.2±06.0 69.8±09.8 71.9±10.0 90.6±06.1 89.4±06.7 84.8±06.0 62.7±11.2
(22) (22) (33) (33) (21) (21) (36) (36) (28) (28) (25) (25) (16) (16)
klas. skóre ± rozptyl (elektroda) 83.0±08.5 (21) 80.1±10.6 (30) 72.6±09.3 72.0±11.5 (21) 64.6±16.2 (30) 71.6±14.1 85.1±07.0 (42) 81.5±05.8 (48) 87.5±05.8 82.5±06.9 (42) 85.0±06.9 (48) 76.0±10.0 76.6±08.6 (24) 57.3±15.0 (28) 69.2±11.2 60.2±12.3 (24) 77.6±11.2 (28) 65.3±10.7 85.4±09.1 (35) 88.3±07.8 (37) 85.1±08.7 93.4±05.7 (35) 89.4±07.4 (37) 89.9±08.3 72.1±08.9 (29) 63.9±12.2 (20) 63.3±11.3 67.1±09.3 (29) 73.8±11.3 (20) 72.5±11.9 90.8±04.9 (29) 91.1±05.9 (24) 91.4±05.0 85.9±07.5 (29) 82.6±08.2 (24) 81.0±08.0 71.4±13.0 (34) 71.2±06.9 (21) 69.2±12.8 65.3±11.6 (34) 64.4±08.7 (21) 66.3±11.0
(24) (24) (25) (25) (18) (18) (34) (34) (18) (18) (39) (39) (35) (35)
77.3±12.5 (26) 63.2±15.8 (26) 82.4±08.2 (26) 78.9±09.4 (26) 75.0±10.7 (16) 59.4±11.2 (16) 86.7±08.9 (45) 87.2±08.5 (45) 66.5±09.5 (24) 69.1±08.7 (24) 91.2±05.0 (38) 80.9±08.3 (38) 67.7±10.0 (42) 67.8±09.7 (42)
Tabulka 1: Klasifikační skóre pro čisté EEG. Uvedeno je jen 5 elektrod s nejlepším klas. skóre. Mezi uvedenými výsledky můžeme pozorovat jistou závislost polohy elektrod s nejlepším klasifikačním skóre. Elektrody s nejlepším klasifikačním skóre totiž tvoří shluky (sousední elektrody taktéž vykazují vysoké klasifikační skóre). Jako příklad uveďme elektrody 21, 22 a 30 u osoby 1, elektrody 25, 26, 33 a 42 u osoby 2, elektrody 34, 35, 36, 37 a 45 u osoby 4, elektrody 20, 28, a 29 u osoby 5, elektrody 29, 38 a 39 u osoby 6. Tato skutečnost potvrzuje fakt, že signál nesoucí informace o pohybu se šíří od svého zdroje směrem ke skalpu a je pozorovatelný na více než jedné elektrodě. Zároveň je však na každé elektrodě přítomen šum nebo jiné rušivé signály (vše co nesouvisí s pohybem), které mají negativní dopad na klasifikaci. Od analýzy nezávislých komponent si slibujeme, že dokáže užitečný signál nesoucí informace o pohybu oddělit od ostatních.
6.2.
Klasifikace komponent
Pro klasifikaci komponent byly použity 3 páry komponent s nejvyšší korelací tak jak bylo zmíněno výše. Separační matice byly počítány ze 27 realizací avšak tyto matice byly následně použity pro rozklad celého EEG. Rozklad celého signálu pomocí separační matice získané ”jen” z 27 realizací si můžeme dovolit proto, že signál je složen z velice podobných úseků (jednotlivých realizací pohybu) a při prodloužení signálu o několik realizací zůstávají základní vlastnosti signálu stejné. Výsledky klasifikace pohybových komponent jsou shrnuty v tab. 2. osoba 1d 1p 2d 2p 3d 3p 4d 4p 5d 5p 6d 6p 7d 7p
klas. skóre ± rozptyl (pár č.) 83.0±09.4 (1) 98.5±01.9 (2) 85.1±07.3 (3) 85.1±11.0 (1) 99.0±03.5 (2) 77.2±12.4 (3) 100.0±00.0 (1) 94.7±03.5 (2) 90.3±03.4 (3) 99.6±01.5 (1) 68.3±13.3 (2) 73.8±07.1 (3) 97.5±02.7 (1) 78.9±08.5 (2) 95.9±03.4 (3) 91.9±07.1 (1) 91.0±06.6 (2) 89.3±06.1 (3) 99.1±02.0 (1) 97.9±02.3 (2) 98.3±02.3 (3) 95.9±04.5 (1) 97.7±03.1 (2) 100.0±00.0 (3) 99.2±02.3 (1) 97.2±02.7 (2) 95.0±06.6 (3) 97.7±02.3 (1) 95.4±05.3 (2) 93.7±05.3 (3) 99.4±03.1 (1) 100.0±00.0 (2) 99.7±00.8 (3) 98.5±02.3 (1) 100.0±00.0 (2) 99.7±01.0 (3) 98.1±01.7 (1) 96.0±02.5 (2) 100.0±00.0 (3) 84.0±04.5 (1) 88.0±05.7 (2) 98.7±01.3 (3)
Tabulka 2: Klasifikační skóre pro pohybové komponenty. Uvedeny jsou jen 3 páry pohybových komponent. Na první pohled je patrné z tab. 2, že klasifikace provedená na pohybových komponentách dosahuje výrazně lepších výsledků. U všech osob je možné nalézt komponenty u kterých se klasifikace blíží 100%. Na druhou stranu musíme výsledky interpretovat velice opatrně, protože jak bylo uvedeno výše v textu, je velice obtížné vybrat pohybové komponenty. Navíc komponent obsahující pohybovou aktivitu může být více a vybrat komponentu která obsahuje nejméně rušivých složek může být problém. Metoda kterou jsme zde uvedli je velice jednoduchá a dle výsledků i účinná. Nicméně, je potřeba ověřit zda vybrané komponenty skutečně obsahují pohybové vlastnosti a jsou vybrány správně. Jedna z metod pro ověření zda komponenty jsou vybrány správně může být založena na mapování komponent na skalp. Pokud se jedná o pohybové komponenty měly by mít stejné nebo alespoň velice podobné rozložení na skalpu. Touto úlohou se budeme zabývat v blízké budoucnosti.
7.
Závěr
V textu jsme popsali klasifikátor založený na HMM a poukázali na problém nedostatku dat. Tento problém je řešen tak, že klasifikace je několikrát opakována s náhodně rozdělenými daty do trénovací a testovací množiny ve stejném poměru. To že jsou data rozdělena náhodně však přináší další problém a tím je závislost výsledků klasifikace na rozdělení dat na trénovací a testovací množinu. Tento problém je odtraněn dostatečným opakováním klasifikace – testem konzistentnosti výsledků jsme tento počet odhadli na 100. Výběr komponent pro klasifikaci je řešen tak, že je spočítána korelace mezi spektrogramem dané komponenty a spektrogramem jistého vzoru. Tento vzor je získán průměrováním spektrogramů z elektrod pokrývající senzomotorickou oblast.
Klasifikační skóre na čistém EEG dosahuje hodnot 60–95% podle dané osoby. Nejvyšší klasifikační skóre (přes 90%) je dosaženo u osob 4 a 6. Naproti tomu, klasifikační skóre pro pohybové komponenty vybrané pomocí korelací spektrogramů dosahuje přes 95% pro všechny osoby. Použitím analýzy nezávislých komponent jsme tedy dosáhli významného zvýšení klasifikačního skóre. Jak již bylo uvedeno v textu je třeba ještě ověřit zda jsou komponenty vybrány správně. Pokud by totiž nebyly komponenty vybrány správně potom bysme klasifikovali různé signály nesoucí jiné informace a klasifikace by pak byla 100%.
Poděkování Tento výzkum byl podporován z grantu GAČR 102/03/H085 “Modelování biologických a řečových signálů” a z výzkumného záměru MŠMT MSM6840770012 “Transdisciplinární výzkum v biomedicínckém inženýrství 2”.
Reference [1] Šťastný, J. Analysis of states in EEG signals. PhD thesis, CTU FEE Prague, Department of circuit theory, March 2005. (In Czech). [2] Šťastný, J.; Sovka, P.; Stančák, A. EEG Signal Classification: Introduction to the Problem. Radioengineering 12, 3 (September 2003), 51–55. [3] Ebrahimi, T.; Vesin, J.-M.; Garcia, G. Brain-Computer Interface in Multimedia Communication. Signal Processing Magazine, IEEE 20, 1 (January 2003), 14–24. [4] Hyvärinen, A.; Oja, E. Independent Component Analysis - Algorithm and Application. Neural Networks 13, 4–5 (March 2000), 411–430. [5] Jung, T.-P.; Makeig, S.; Humphries, C.; Lee, T.-W.; Mckeown, M.; Iragui, V.; Sejnowski, T. Removing electroencephalographic artifacts by blind source separation. Psychophysiology 37, 2 (March 2000), 163–178. [6] Jung, T.-P.; Makeig, S.; Lee, T.-W.; McKeown, M.; Brown, G.; Bell, A.; Sejnowski, T. Independent Component Analysis of Biomedical Signals. In The 2nd International Workshop on Indeppendent Component Analysis and Signal Separation 2000, pp. 633–644. [7] Kuhlman, W. Functional Topography of The Human mu Rhythm. Electroencephalography and Clinical Neurophysiology 44 (1978), 83–93. [8] Lima, C.; Silva, C.; Tavares, A.; Oliveira, J. Blind Source Separation by Independent Component Analysis Applied to Electroencephalographic Signals. In International Workshop on Models and Analysis of Vocal Emissions for Biomedical Applications – MAVEBA, 3 2003. [9] Stančák, A.; Feigeb, B.; Lücking, C.; Kristeva-Feige, R. Oscillatory cortical activity and movement-related potentials in proximal and distal movements. Clinical Neurophysiology 111, 4 (April 2000), 636–650.
[10] Wolpaw, J.; Birbaumer, N.; McFarland, D.; Pfurtscheller, G.; Vaughan, T. Braincomputer interfaces for communication and control. Clinical Neurophysiology 113, 6 (June 2002), 767–791.
