Pracownia Sygnałów Biologicznych/Zajecia 5 6

Pracownia Sygnałów Bioelektrycznych/EMG

Pomiar EMG

Wstęp

Filmik ilustrujący działanie mięśni

Sygnały elektro-fizjologiczne pochodzące z mięśni nazywa się elektromiogramem (EMG). Elektromiografia jest jednym z podstawowych badań w rozpoznawaniu chorób mięśni i nerwów obwodowych. Ma ona również wiele zastosowań naukowych. Amplituda sygnału EMG wynosi od około kilkudziesięciu μV do 10 mV, zaś pasmo sygnału obejmuje zakres częstości od 2 do 5000 Hz, przy czym największa energia sygnału znajduje się w przedziale od 50 do 150 Hz.

Istnieją dwa sposoby pomiaru sygnałów EMG — badanie igłowe i powierzchniowe. W badaniu igłowym EMG, elektroda igłowa lub igła z dwoma elektrodami wbijana jest w mięsień lub w nerw ruchowy. Następnie obserwuje się aktywność elektryczną mięśni w spoczynku i podczas wysiłku.

Przykłady elektromiogramów. Panel górny — pacjent zdrowy. Panel środkowy — pacjent ze zmianami w nerwach obwodowych (neuropatia). Panel dolny — pacjent ze zmianami w mięśniach (miopatia).

Badanie powierzchniowe EMG wykonuje się z użyciem elektrod samoprzylepnych, umieszczonych na powierzchni skóry, zwykle elektrody bipolarne są rozmieszczone na linii równoległej do włókien. Ocenie podlegają mięśnie położone powierzchownie lub grupy mięśni. Obydwie metody mają swoje wady i zalety. Metoda „igłowa” umożliwia rejestrację sygnału EMG z wybranego mięśnia, podczas gdy metoda powierzchniowa rejestruje zbiorczą aktywność wielu jednostek ruchowych. Jednakże, w przeciwieństwie do metody powierzchniowej, metoda igłowa jest badaniem inwazyjnym i czasem bolesnym, które wykonywane jest w ośrodku klinicznym.

Na zajęciach przeprowadzimy powierzchniowy pomiar EMG (w literaturze często takie badanie oznacza się skrótem sEMG, s od ang. surface — powierzchnia). W tym celu umieścimy elektrody na skórze, nad mięśniami, których aktywność chcemy zbadać.

Filmik o zastosowaniach klinicznych EMG

Lektura uzupełniająca: Reaz MBI, Hussain MS and Mohd-Yasin F. Techniques of EMG signal analysis: detection, processing, classification and applications. Biol. Proced. Online 2006; 8(1): 11-35. https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1455479/pdf/bpo_v8_p11_m115.pdf

Źródła błędu

Rejestracja sygnału EMG, podobnie jak pomiar innych sygnały bioelektrycznych podlega pewnym zakłóceniom. Są to przede wszystkim

• artefakty ruchowe,
• artefakty związane z obecnością zewnętrznych pól elektromagnetycznych.

Artefakty ruchowe

Artefakty ruchowe

Aktywacja mięśnia powoduje jego skrócenie i przemieszczenie względem elektrody oraz skóry. Skutkuje to trzema rodzajami zakłóceń sygnału EMG:

1. Zmiana amplitudy Elektroda oddala się od aktywnego mięśnia, co obniża rejestrowany sygnał.

2. Zmiana potencjału DC Przemieszczenie wpływa na potencjał elektrochemiczny na granicy skóra–elektrolit.

3. Artefakt rozciągnięcia skóry Potencjał przeznaskórkowy (PTP) w spoczynku wynosi około 30 mV. Podczas rozciągania

Artefakty sieciowe

Artefakty sieciowe stanowią poważny problem w przypadku rejestracji sygnału EMG. Zwykle bowiem nie obserwujemy tylko zakłócenia o jednej częstości, równej częstości zmian napięcia w gniazdku zasilającym (np. w Polsce jest to 50 Hz), lecz również wyższe harmoniczne tej częstości (w Polsce będą to 100 Hz, 150 Hz, 200 Hz itd.). Jak można zauważyć, częstości 50 Hz, 100 Hz i 150 Hz znajdują się w paśmie, w którym sygnał EMG ma największą energię. Stosowanie filtrów pasmowo zaporowych w takim przypadku nie jest wskazane, bowiem filtry jak wiemy nie tłumią ściśle określonych częstości tylko pasma o pewnej szerokości (np. 45-55 Hz, 95-105 Hz, itd). W efekcie znaczna część interesującego nas pasma sygnału EMG zostałaby odrzucona.

Aby zminimalizować przenikanie od rejestrowanego sygnału EMG sygnału sieciowego o częstości 50 Hz, należy:

• zadbać o niską impedancję pomiędzy elektrodą a skórą,
• mierzyć sygnały różnicowe.

Artefakty sieciowo-ruchowe

Podczas ruchu mięśnia ciało może ulec przemieszczeniu, co spowoduje ruch kabla w przestrzeni pomiędzy elektrodą a wzmacniaczem. W przestrzeni tej istnieje pole elektromagnetyczne wywołane zasilaniem sieci elektrycznej. Ruch kabla w polu elektromagnetycznym może powodować zaburzenia sygnału o częstości 50 Hz. Ponadto układ elektrody-kable-wzmacniacz, tworzy pewne ramki, które w trakcie ruchu zmieniają kształt a także powierzchnię. Zgodnie z prawem indukcji Faradaya zmiana strumienia magnetycznego powoduje powstawanie siły elektromotorycznej, która również może zakłócać pomiar.

• Artefakty ruchowe można w większości wyeliminować przez zastosowanie filtru górnoprzepustowego, którego częstość odcięcia ustawia się w granicy od 10 do 20 Hz.
• Aby wyeliminować artefakty związane z ruchem kabla w polu elektromagnetycznym można zastosować tzw. elektrody aktywne. W elektrodach tych (miniaturowy) wzmacniacz znajduje się na elektrodzie. Wzmacnianie sygnału na elektrodzie zwiększa względną czułość układu na sygnał mierzony na elektrodzie w stosunku do zaburzenia związanego z ruchem kabla.
• Innym rozwiązaniem tego problemu jest zastosowanie kabli ekranowanych. Kable TMSI, których używamy na Pracowni, są kablami ekranowanymi. W kablach tych, pomiędzy dwiema warstwami izolatora, istnieje dodatkowa osłona z przewodnika podłączona do wzmacniacza. Dzięki temu zewnętrzne pole elektromagnetyczne nie przenika do środka kabla. Dodatkowo aby wyeliminować indukowanie się ładunków w wyniku tarcia pomiędzy izolatorami a osłoną, jest ona częściowo pokryta warstwą węgla.

Ćwiczenia

Ćwiczenie I: Badanie zależności sygnału EMG od obciążenia

Pomiar

• Umieść elektrodę GND na obojczyku,
• Umieść dwie elektrody do rejestracji sygnału EMG na mięśniu dwugłowym ramienia (popularnie zwanym bicepsem). Kable tych elektrod połącz z unipolarnymi wejściami wzmacniacza numer 1 i 2.
• Ustaw częstość próbkowania sygnału na 2048

Obserwacje wstępne

W tej części proszę wykonać wskazane poniżej obserwacje i zarejestrować fragmenty sygnałów, które mogłyby ilustrować ciekawe obserwacje, aby można je było wykorzystać w końcowej prezentacji.

• W SVAROGU:
  • badać będziemy mięsień przedramienia
  • przetestuj sygnały rejestrowane w kilku różnych lokalizacjach względem mięśnia
  • przyjrzyj się przebiegowi i widmu sygnałów dla różnych wartości siły napięcia mięśnia - czy uda się zaobserwować "salwy" odpowiadające aktywacji pojedynczej jednostki ruchowej?
  • dobierz odpowiednio filtry górnoprzepustowe
  • zaobserwuj czy i kiedy w sygnale pojawiają się artefakty
  • porównaj sygnały rejestrowane przez elektrody bipolarne i przez elektrody monopolarne zmontowane bipolarnie

Rejestracja 1: napinanie i rozluźnianie ręki

Zarejestruj sygnał, w którym naprzemiennie ściskasz 1 s i rozluźniasz pięść 2 s (dla badania mięśnia przedramienia). Powtórz 30 razy.

Rejestracja 2: zwiększanie obciążenia ręki

Rejestrujemy dwa stany: spoczynek, narastajace obciążenie:

• zarejestruj około 1 min. sygnału spoczynkowego
• w dalszej części będziemy stopniowo zwiększać obciążenie: co około 20 s kolejne zwiększenie obciążenia - 5 lub 6 poziomów

Rejestracja 3: stałe obciążenie ręki

Rejestrujemy dwa stany: spoczynek, narastajace obciążenie:

• zarejestruj około 20 s. sygnału spoczynkowego
• w dalszej części dokładamy obciążenie (środkowy poziom z poprzedniego ćwiczenia) i rejestrujemy 1 min. sygnału

Analiza

• Dalsze analizy robimy w pythonie (notebook).
  • wykreśl sygnał
  • sygnały filtrujemy: 30 Hz high-pass, 500 Hz low-pass and 50 Hz notch,
  • Przedstaw przebieg średniej mocy sygnałów. Aby to zrobić podnieś próbki do kwadratu i uśrednij za pomocą średniej biegnącej o długości 0.1 s (filtrowanie oknem prostokątnym). Porównaj uzyskane przebiegi z przebiegiem sygnałów wejściowych.
  • Zaproponuj algorytm detektora, wykrywającego ruch ręki.

• Zbadaj zależność parametrów sygnału EMG od obciązenia mięśnia:
  • Predstaw wykres średnią amplitudę (odchylenie standardowe) od obciążenia
  • Przedstaw widma dla kolejnych obciążeń (na wspólnym wykresie). Do estymacji wykorzystaj metodę Welcha, dobierając parametry tak, aby rozdzielczość częstotliwościowa była 1 Hz.

• Zbadaj zależność parametrów sygnału EMG od zmęczenia mięśnia.Badanie efektu zmęczenia mięśnia w EMG
  • Zmęczenie mięśnia w elektromiografii bada się głównie poprzez analizę zmian sygnału podczas długotrwałego skurczu izometrycznego. Kluczowe metody:
    • Analiza widmowa (spektralna): Podczas zmęczenia obserwuje się przesunięcie widma mocy w kierunku niskich częstotliwości, co opisują dwa wskaźniki:
      • Częstotliwość medianowa (MNF/MDF) – wartość częstotliwości, poniżej której i powyżej której leży po 50% całkowitej mocy widma. Podczas zmęczenia systematycznie spada.
      • Częstotliwość średnia (Mean Power Frequency, MPF) – ważona średnia częstotliwości widma. Również maleje wraz ze zmęczeniem.
  • Spadek tych parametrów wynika z wolniejszego przewodzenia potencjałów wzdłuż włókien mięśniowych (obniżenie prędkości przewodzenia) oraz akumulacji metabolitów (mleczany, jony H⁺, K⁺).
• Przedstaw ewolucję czasową środka ciężkości widma częstości

Ćwiczenie III: zapoznanie się z sygnałami rejestrowanymi przez inwazyjne EMG

Proszę wczytać i przyjrzeć się sygnałom zdrowy.bin, miopatia.bin, neuropatia.bin. Sygnały są zapisane jako dtype='float64' Częstość próbkowania 4000Hz, amplitudy zapisane są w mV. Dane pochodzą z bazy Physionet: https://physionet.org/content/emgdb/1.0.0/

Proszę przeczytać informację o tych danych. W raporcie proszę napisać istotne informacje o tych sygnałach i wybrać charakterystyczne fragmenty typowe dla danego stanu klinicznego.

Ćwiczenie V: Wykorzystanie pomiaru EMG do sterowania on-line

Przykładowy fragment kodu example.py umożliwiający odbieranie sygnału on line w pythonie przedstawiony jest poniżej.

Proszę dodać fragment analizujący sygnał on-line i wykrywający moment napięcia mięśnia. Po wykryciu w najprostszej wersji niech w terminalu pojawia się komunikat o wykryciu kliknięcia. W wersji max proszę zrobić wizualizację napięcia mięśnia albo podpiąć ten sygnał do sterowania jakimś prostym interfejsem.

Wersja działająca z pythonem z pakietu Svarog-Lab

#!/opt/braintech-svarog-lab-python/bin/python3

# powinno się odpalać pythonem z pakietu Svarog-Lab
# lub python w którym są zainstalowane sterowniki z
# https://braintech.pl/pliki/svarog/svarog-streamer-src/svarog-streamer-src-latest.zip
# zaletą jest dostęp wprost do sterownika wzmacniacza
# /opt/braintech-svarog-lab-python/bin/python3


# Alternatywne wzmacniacze
# from braintech.drivers.perun32.amplifier import Perun32Amplifier as Amplifier # duży wzmacniacz na 32 kanały na USB
# from braintech.drivers.perun8.amplifiers import PerunCppAmplifier as Amplifier # headset na 8 kanałów
from braintech.drivers.tmsi.amplifiers import TmsiCppAmplifier as Amplifier # wzmacniacze TMSI

import numpy as np

# Szukamy wzmacniaczy
amps = Amplifier.get_available_amplifiers()
if len(amps) < 1:
    raise Exception("Amplifier not connected")
amp = Amplifier(amps[0])

amp.sampling_rate = 1024 # dla TMSI
# amp.sampling_rate = 500 # dla Perun8, Perun32



def samples_to_microvolts(samples):  # amplifier podaje próbki w integerach, wprost z ADC
    return samples * gains + offsets


amp.start_sampling()
gains = np.array(amp.current_description.channel_gains)
offsets = np.array(amp.current_description.channel_offsets)
while True:
    # pobieramy 16 próbek
    # Proponuje używać IPython dla eksperymentowania
    packet = amp.get_samples(16)
    print(samples_to_microvolts(packet.samples))
    print(packet.ts[0])
    print(packet.samples.shape, amp.current_description.channel_names)

Wersja dla protokołu LSL

# gdy nie ma możliwości używać pythona wbudowanego w svarog-lab możliwe jest używanie streamowania do standardu LSL

# włączenie streamowania do LSL:
# w terminalu:
# svarog_streamer -l 
# wypisze listę wzmacniaczy
# szukamy ID odpowiedniego wzmacniacza np:
# * Perun-8 Headset
#       id: "Perun8 1"

# odpalamy stream LSL danego wzmacniacza:
# svarog_streamer -a "Perun8 1" -n "nazwa_streamu"
# nazwa streamu jest ważna, ponieważ streamy są widoczne w sieci LAN
# zaleta streamu jest też taka, że można podglądać go w Svarogu jednocześnie z naszym skryptem python
# po odpaleniu streamu można odpalać poniższy skrypt dowolnym pythonem z zainstalowanym numpy oraz pylsl

from pylsl import StreamInlet, resolve_stream
import time

nazwa_streamu = "nazwa_streamu" # należy odpowiednio zmienić na nazwę użytą w svarog_streamer -n

# znajdujemy streamy
print("szukamy streamy LSL")
streams = resolve_stream('type', 'EEG')

selected_stream = None
# wybieramy nasz
for stream in streams:
	if stream.name() in nazwa_streamu:
		selected_stream = stream
if selected_stream is None:
	print("Nie znalesiono streamu", nazwa_streamu, "w liście", [i.name() for i in streams])
	exit()

# używamy streama
inlet = StreamInlet(selected_stream)

while True:
	# pobieramy próbki (w mikrowoltach)
    sample, timestamp = inlet.pull_chunk(timeout=1.0, max_samples=10)
    print(sample, timestamp, time.monotonic())

Dodatek

Do działania on-line przydatne może być filtrowanie sygnału w sposób biegnący. Najlepiej zastosować do tego funkcję lfilter z ustalonymi warunkami początkowymi. Używa się tego w następujący sposób:

import scipy.signal as ss
import numpy as np
import matplotlib.pyplot as plt

Fs = 256
T = 1
t = np.arange(0,T,1/Fs)
f0 = 10
f1 = 17
f2 = 23
x = (np.sin(2*np.pi*f0*t) +
     np.sin(2*np.pi*f1*t ) +
     np.cos(2*np.pi*f2*t))
xn = x + np.random.randn(len(t)) * 0.08

b, a = ss.butter(3, 11/(Fs/2))

zi = ss.lfilter_zi(b, a)
z, _ = ss.lfilter(b, a, xn, zi=zi*xn[0])
plt.plot(xn)
plt.plot(z)
plt.show()

Filtrowanie online

Symulacja zastosowania tego sposobu on-line; filtrujemy za każdą iteracją pętli to co przychodzi ze wzmacniacza [s]:

zi = ss.lfilter_zi(b, a)
y = np.zeros(xn.shape)
for  ind, s  in enumerate(xn):
  y_tmp, zi = ss.lfilter(b, a, [s], zi=zi)
  y[ind]=y_tmp[-1]

  
plt.plot(y)
plt.plot(z)
plt.plot(xn)
plt.show()

Dodatek 2

Dla zabawy detektor napięcia mięśni mógłby generować kliknięcia myszki, aby np. zagrać w:

https://dino-chrome.com/

W tym celu trzeba doinstalować bibliotekę pynput:

 pip3 install --user pynput

a potem zastosować kod ze strony: https://pynput.readthedocs.io/en/latest/keyboard.html

from pynput.keyboard import Key, Controller

keyboard = Controller()

def spacja():
# Press and release space
	keyboard.press(Key.space)
	keyboard.release(Key.space)


from pynput.mouse import Button, Controller
mouse = Controller()

def click():
#lewy klik myszy
	mouse.press(Button.left)
	time.sleep(500)
	mouse.release(Button.left)

wersja dla systemu dostępnego w Labie 4.59 na Ubuntu 18

from obci_cpp_amplifiers.amplifiers import TmsiCppAmplifier
import numpy as np
amps = TmsiCppAmplifier.get_available_amplifiers('usb')
amp = TmsiCppAmplifier(amps[0])

amp.sampling_rate = 512

amp.start_sampling()
gains = np.array(amp.current_description.channel_gains)
offsets = np.array(amp.current_description.channel_offsets)

def samples_to_microvolts(samples):  # z jednostek wzmacniacza do mikrowoltów
    return samples * gains + offsets
    
while True:
    # 16 próbek w pakiecie, nieodebrane próbki się bufurują i można odebrać je później
    packet = amp.get_samples(16)
    print(samples_to_microvolts(packet.samples))
    print(packet.ts[0])
    print(packet.samples.shape, amp.current_description.channel_names)

Aby wykonać go w terminalu należy uruchomić polecenie, Uwaga, aby zadziałał trzeba wyłączyć SVAROGa. Proszę przetestować czy po podłączeniu wzmacniacza i uruchomieniu tego skryptu pojawiają się w terminalu wartości próbek.

/opt/braintech/bin/python3 example.py

Moduł keystroke: 

# -*- coding: utf-8 -*-
#!/usr/bin/env python
# Author:
#     Mateusz Kruszyński <mateusz.kruszynski@gmail.com>
#

from Xlib import X, display, Xutil, XK
import Xlib
special_X_keysyms = {
    ' ' : "space",
    '\t' : "Tab",
    '\n' : "Return",  # for some reason this needs to be cr, not lf
    '\r' : "Return",
    '\e' : "Escape",
    '!' : "exclam",
    '#' : "numbersign",
    '%' : "percent",
    '$' : "dollar",
    '&' : "ampersand",
    '"' : "quotedbl",
    '\'' : "apostrophe",
    '(' : "parenleft",
    ')' : "parenright",
    '*' : "asterisk",
    '=' : "equal",
    '+' : "plus",
    ',' : "comma",
    '-' : "minus",
    '.' : "period",
    '/' : "slash",
    ':' : "colon",
    ';' : "semicolon",
    '<' : "less",
    '>' : "greater",
    '?' : "question",
    '@' : "at",
    '[' : "bracketleft",
    ']' : "bracketright",
    '\\' : "backslash",
    '^' : "asciicircum",
    '_' : "underscore",
    '`' : "grave",
    '{' : "braceleft",
    '|' : "bar",
    '}' : "braceright",
    '~' : "asciitilde"
    }

display = display.Display()
window = display.screen().root

def wait(p_keys_list):
    """Block the whole keyboard!!! And wait until some key from p_keys_list
    is pressed. By now p_keys_list is a list of strings, so use single
    ascii symbols.
    There is a way out of this hell - hit 'Escape'.
    The function returns hit button`s string representation
    Eg. for p_keys_list == ['1','2','3'] the function will hand untill
    1,2 or 3 key is preseed or Escape is pressed."""
    ds = display
    window.grab_keyboard(1, X.GrabModeAsync, X.GrabModeAsync, X.CurrentTime)
    while True:
        ev = ds.next_event()
        if ev.type == X.KeyPress:
            keysym = ds.keycode_to_keysym(ev._data['detail'], 0)
            keystr = XK.keysym_to_string(keysym)
            print("Got keysym/keystr: "+str(keysym)+ ' / '+str(keystr))
            if keystr in p_keys_list:
                ds.ungrab_keyboard(X.CurrentTime)
                ds.flush()
                return keystr
            elif str(keysym) in p_keys_list:
                ds.ungrab_keyboard(X.CurrentTime)
                ds.flush()
                return keysym
            elif keysym == 65307:
                ds.ungrab_keyboard(X.CurrentTime)
                ds.flush()
                return 'Escape'

def char_to_keysym(ch) :
    keysym = Xlib.XK.string_to_keysym(ch)
    if keysym == 0 :
        # Unfortunately, although this works to get the correct keysym
        # i.e. keysym for '#' is returned as "numbersign"
        # the subsequent display.keysym_to_keycode("numbersign") is 0.
        keysym = Xlib.XK.string_to_keysym(special_X_keysyms[ch])
    return keysym

def keysym_to_keycode(keysym):
    keycode = display.keysym_to_keycode(keysym)
    shift_mask = 0
    return keycode, shift_mask

def char_to_keycode(ch) :
    keysym = char_to_keysym(ch)
    keycode = display.keysym_to_keycode(keysym)
    if keycode == 0 :
        print "Sorry, can't map", ch

    if (is_shifted(ch)) :
        shift_mask = Xlib.X.ShiftMask
    else :
        shift_mask = 0

    return keycode, shift_mask


def is_shifted(ch) :
    if ch.isupper() :
        return True
    if "~!@#$%^&*()_+{}|:\"<>?".find(ch) >= 0 :
        return True
    return False


def send_string(str) :
    """I am not working. I dont know why:("""
    for ch in str :
        #print "sending", ch, "=",
display.keysym_to_keycode(Xlib.XK.string_to_keysym(ch))
        keycode, shift_mask = char_to_keycode(ch)
        event = Xlib.protocol.event.KeyPress(
            time = int(time.time()),
            root = display.screen().root,
            window = window,
            same_screen = 0, child = Xlib.X.NONE,
            root_x = 0, root_y = 0, event_x = 0, event_y = 0,
            state = shift_mask,
            detail = keycode
            )
        window.send_event(event, propagate = True)
        event = Xlib.protocol.event.KeyRelease(
            time = int(time.time()),
            root = display.screen().root,
            window = window,
            same_screen = 0, child = Xlib.X.NONE,
            root_x = 0, root_y = 0, event_x = 0, event_y = 0,
            state = shift_mask,
            detail = keycode
            )
        window.send_event(event, propagate = True)


if __name__ == "__main__":
    import sys,time
    print(wait(sys.argv[1:]))
    #send_string("aBcd")
    #time.sleep(10)

"""Używa się tak:
import keystroke
keystroke.wait(['a','b'])
i to wisi dopóki ktoś nie wciśnie a lub b.

Z innymi klawiszam jest tak, że trzeba ich numery podać, np żeby czakać 
na SPACE można zrobić
keystroke.wait(['32'])

w razie czego funkcja się 'odwiesza' jak się wciśnie ESC."""