Laboratorní úloha č.4: Elektromyogram
Úvod: Svaly jsou ke kontrakci stimulovány nervovými impulsy přicházejícími z centrální nervové soustavy (CNS) skrz míchu a motorické nervové kořeny, které obsahují jednotlivá motorická nervová vlákna. Jedno motorické nervové vlákno může inervovat několik svalových vláken. Spojení jednoho nervového vlákna a všech svalových vláken, které inervuje, se nazývá motorická jednotka (motor unit). Počet svalových vláken připadajících na jedno motorické nervové vlákno, tedy velikost motorické jednotky, se odvíjí podle umístění svalu na těle, resp. podle toho, jak jemná a precizní má být činnost daného svalu. Proces, kdy se při svalové činnosti aktivují jednotlivé motorické jednotky, nazýváme nábor motorických jednotek (motor unit recruitment). Síla stažení svalu je přitom přímo úměrná počtu zapojených motorických jednotek a dále frekvenci nervových impulsů přicházejících do každé jednotky. Amplituda elektrických impulsů přiváděných nervovými vlákny a dále vedených svalovými vlákny je velmi malá, méně než 100 µV, ale protože jsou tyto impulsy vedeny souběžně velmi mnoha vlákny, je vzniklý nasčítaný signál dostatečně velký na to, aby byl měřitelný na povrchu přiléhající pokožky pomocí páru elektrod. Elektrický signál měřitelný na povrchu pokožky a produkovaný při svalové kontrakci se nazývá elektromyografický signál (EMG) a jeho záznam se nazývá elektromyogram. Z popisovaných vlastností svalu vyplývá, že čím větší síla je od svalu požadována, tím rychleji musí do svalu přicházet impulzy z CNS a musí být zapojeno více svalových jednotek. Více svalových jednotek vytváří více myopotenciálů, tedy v součtu větší amplitudu EMG. Elektromyografie – zahrnuje skupinu elektrofyziologických metod, které registrují elektrické projevy při činnosti svalového a nervového aparátu. EMG záznam poskytuje informace umožňující diferenciální diagnostiku svalových a neuromuskulárních poruch. Cíle úlohy: 1. Analyzujte vlastní nebo jeden Vámi zvolený povrchový EMG signál při různém volním úsilí: a) vypočtěte průměrné EMG z obálek signálu získaných: s použitím filtru klouzavých průměrů (MovingAverage) s použitím integrátoru vhodně zvolte parametry filtrů, diskutujte vliv délky klouzavého okna a integrační konstanty; porovnejte filtry b) vypočtěte celkovou integraci EMG za celou dobu volního úsilí c) vypočtěte integrované EMG v intervalech po např. 200 ms d) vypočtěte integrované EMG v intervalech, kdy integrace dosáhne např. 150 mV e) zobrazte výkonová spektra jednotlivých kontrakcí: vypočtěte frekvenci mediánu výkonu (50 % výkonu) 2. Analýza únavy (stanovte parametry závislé na únavě v časové a frekvenční oblasti) vyberte dva úseky analyzovaného signálu; jeden na začátku kontrakce, druhý ke konci, kdy ale ještě nedošlo k výraznému poklesu vyvinuté síly: vypočítejte a porovnejte výkonová spektra obou úseků z analýzy celého signálu zjistěte, které konkrétní parametry poskytují informaci o únavě a zda nám více informací poskytují parametry získané z časové nebo frekvenční oblasti. Klouzavě ve zvoleném okně určete trend: směrodatná odchylka počet průchodů nulou ZCR – zerocrossing rozdíl "špička-špička" medián frekvencí 1. spektrální moment – těžiště spektra 2. spektrální moment – rozprostření spektra
1
Data potřebná k vypracování úlohy: http://sami.fel.cvut.cz/bsg/cv04/Data_lab04.zip
Pořízení biologických signálů: 1.
Snímání povrchového EMG při volním úsilí. Síla se snímá dynamometrem; EMG aktivita s použitím povrchových elektrod. Během záznamu trvajícím 10 sekund se tiskne třikrát dynamometr a to vždy za užití větší síly. Dynamometr stiskněte třikrát za sebou: nejprve slabě, poté středně a nakonec silně. dynamometr
snímací elektrody sila [kg]
60 40 20 0 0
1
2
3 t [s]
4
5
6
0
1
2
3 t [s]
4
5
6
EMG [mV]
5
-5
dynamometr
Analýza únavy. Snímání vyvinuté síly a EMG aktivity se provádí stejným způsobem jako v předcházejícím bodě. Stiskněte dynamometr při co největším volním úsilí po co nejdelší dobu. Snažte se udržet konstantní sílu stisku. dynamometr sila [kg]
40 20 0
EMG [mV]
2.
0
0
10
20
30
40
50
60
0
10
20
30 t [s]
40
50
60
5 0 -5
2
Struktura dat: 1. Snímání povrchového EMG při volním úsilí Př.: emg_3xR.txt fs= 1000 Hz 1. sloupec ... síla na dynamometru [kg] 2. sloupec ... EMG [mV] 2. Analýza únavy Př.: emg_vydrz.txt fs= 1000 Hz 1. sloupec ... síla na dynamometru [kg] 2. sloupec ... EMG [mV] Nápověda k některým úkolům: 1) Průměrné EMG (obálka): integrátor, viz 3. úloha klouzavé průměry: % např. MA-filtr 5. řádu % b=ones(5,1)/5; % a=1;
obálky je třeba časově posunout a nanormovat vhodně zvolenou konstantou Obrázek 1: Obálky EMG signálu s různým volním úsilím, porovnání obálek s úsilím (normovaně k EMG) obálky emg signálu 4 2
emg MA-filter integrator
0 -2 500
60 40
1000
1500
2000
2500
3000
3500
4000
4500
5000
5500
6000
1500
2000
2500
3000 vzorky [-]
3500
4000
4500
5000
5500
6000
síla MA-filter integrator
20 0 500
1000
Celková integrace energie kumulativní součet absolutních hodnot EMG v diskrétním čase EMG v časových intervalech: postupné sčítání absolutních amplitud EMG, v 200 ms se čítač nuluje: counter=zeros(size(emg,1)-1,1); % vytvoření čítače for i=1:size(emg,1)-1; if mod(i,0.2*fs)==0 % když je pozice na 200 ms counter(i)=0; % vynuluj čítač end counter(i+1)=counter(i)+abs(emg(i+1,2)); % čítač + |emg| end
3
EMG v intervalech po dosažení 150 mV: obdobně jako v časových intervalech, změňte jen podmínku nulování Obrázek 2: EMG, celková integrace, integrace po 200 ms, integrace do 150 mV EMG U [mV]
5 0 -5
1
2
3
4 cumsum EMG
5
6
7
1
2
3
4 integrace po 200 ms
5
6
7
1
2
3
4 integrace do 150 mV
5
6
7
1
2
3
4 čas [s]
5
6
7
4000 2000 0
1000 0
200 100 0
Spektra jednotlivých kontrakcí: ručně vyberte úseky EMG signálu jednostranné spektrum: S=abs(fft(usek_emg)); % oboustranné S=S(1:round(end/2)); % jednostranné
kumulativní součet výkonového spektra: f=linspace(0,fs/2,length(S)); % frekvence spektra med_f=f(find(cumsum(S)>=0.5*sum(S),1)); % frekvence, kde kumulativní součet výkonového spektra dosáhne 50 % Obrázek 3: Odhad výkonových spekter a mediány frekvencí úseků volního úsilí EMG 4 2 0 -2 1000
1500
2000
2500 -3
x 10
0
50
100
150
150
50
100
5500
6000
50
100
150
50 100 k 0.5 = 71Hz
150
0.04 0.02
150
0
f [Hz]
100 50 0
5000
0.06
0 cumsum PSD [%]
50
50 100 k 0.5 = 46Hz
4500
f [Hz]
100
0
4000
15 10 5
f [Hz]
0
3500
PSD
PSD
PSD
x 10 14 12 10 8 6 4 2
3000 vzorky [-]
0
50 100 k 0.5 = 59Hz
4
150
cumsum PSD [%]
500 -3
cumsum PSD [%]
[mV.s -1]
2000
100 50 0
0
2) Parametry signálů: - vyberte dva úseky analyzovaného signálu; jeden na začátku kontrakce, druhý ke konci, kdy ale ještě nedošlo k výraznému poklesu vyvinuté síly. Vypočítejte a porovnejte výkonová spektra obou úseků.
Obrázek 4: Odhad spektra dvou rozdílných úseků
U [mV]
EMG 5 0 -5 0
10
20
0
10
20
30 t [s] síla
40
50
60
30 40 t [s] PSD-pwelch
50
60
[kg]
40 20
[µW]
0
0.1 0.08 0.06 0.04 0.02 0
20
40
60
80
100
f [Hz]
a) pro celý signál EMG, kde můžeme pozorovat reálné hodnoty naměřené síly (poz: odstraňte zhruba prvních a posledních 500 vzorků signálu), nejprve segmentujte signál a poté pro každý segment spočítejte následující parametry průchody nulou ZCR: segmentujte úsek signálu, zvolte vhodné okno a překryv (winsize, noverlap) zcr=[]; index=1:winsize-noverlap:length(usek)-winsize+1; usek=usek-mean(usek); % odstranění stejnosměrné s. for i=index; segment=usek(i:i+winsize-1); polarita=sign(segment); % jen znamínka signálu (+1,-1) pruchody=diff(polarita); % změna polarity pocet_nul=sum(pruchody~=0); % součet změn = zcr zcr=[zcr,pocet_nul]; % uložení zrc daného segmentu end
rozdíl "špička-špička", myšleno rozdíl mezi maximální a minimální hodnotou amplitudy v signálu medián frekvence: f=linspace(0,fs/2,length(S)); % frekvence spektra med_f=f(find(cumsum(S)>=0.5*sum(S),1)); % frekvence, kde kumulativní součet spektra dosáhne 50 %
5
spektrální momenty: =
∑
/
× ( )
=
∑
/
× ( )
− ∑ ∑ mom1=sum(f.*S)/sum(S); % frekvence těžiště spektra mom2=sqrt(sum((f.^2).*S)/sum(S)-mom1.^2); % rozprostření spektra
b) v dalším kroku určete trend proložením přímkou a vykreslete graf viz. obrázek 5. (polynomem 1. stupně): % t_start je čas začátku úseku t=linspace(t_start,(index(end)+winsize-1)/fs,length(index)); p=polyfit(t,zcr,1); % nalezne koeficienty lineární rovnice trend_zcr=polyval(p,t); % dopočte hodnoty y ze zadaných koeficientů rovnice, dále s pomocí funkce linspace vykreslete přímku.
Obrázek 5: EMG při zkoumání svalové únavy; parametry EMG v závislosti na čase proložené přímkou.
max-min [mV] STD [mV]
60
med f [Hz]
ZCR
segmentace, parametrizace - proložení přímkou 2.2 2 1.8 1.6 1.4 1.2 10
1. mom [Hz]
20
25
30
35
40
45
50
15
20
25
30
35
40
45
50
15
20
25
30
35
40
45
50
15
20
25
30
35
40
45
50
15
20
25
30
35
40
45
50
15
20
25
30 t [s]
35
40
45
50
12 10 8 10
40 10 70 60 50 10
2. mom
15
100 80 10 150
100 10
Užitečné funkce: fft,pwelch, cumsum, linspace, sign, polyfit, polyval, diff, filter, sum, mean,median, mod, std, abs 6