Laboratorium EEG/EEGLAB: Różnice pomiędzy wersjami
m (→SKRYPTY) |
|||
Linia 69: | Linia 69: | ||
==== <tt>filtruj_i_tnij.m</tt>==== | ==== <tt>filtruj_i_tnij.m</tt>==== | ||
Musimy teraz napisać funkcję <tt>filtruj_i_tnij</tt>. Musi ona kolejno wczytać nasze zbiory danych, stworzyć montaż kanałów, przefiltrować je, kanał po kanale, wyciąć interesujące nas epoki, wykonać dla nich korektę poziomu "0", i zapisać tak przetworzone dane do nowych "*.set"-ów. Szkilet tej funkcji wygląda następująco, proszę uzupełnić brakujące fragmenty: | Musimy teraz napisać funkcję <tt>filtruj_i_tnij</tt>. Musi ona kolejno wczytać nasze zbiory danych, stworzyć montaż kanałów, przefiltrować je, kanał po kanale, wyciąć interesujące nas epoki, wykonać dla nich korektę poziomu "0", i zapisać tak przetworzone dane do nowych "*.set"-ów. Szkilet tej funkcji wygląda następująco, proszę uzupełnić brakujące fragmenty: | ||
+ | |||
<source lang = matlab> | <source lang = matlab> | ||
− | + | ||
+ | function filtruj_i_tnij(EEG, ALLEEG, CURRENTSET, directoryIn, directoryOut, addName, subjects, LowFreqCut, HighFreqCut, startT, stopT, eegChan, baseline) | ||
+ | |||
+ | % zapamiętujemy stan wejściowy zmiennych globalnych eeglaba | ||
+ | EEG0 = EEG; | ||
+ | ALLEEG0 = ALLEEG; | ||
+ | CURRENTSET0 = CURRENTSET; | ||
+ | events = {'Pseudo','ANeg','ANeu','APos','ONeg','ONeu','OPos','RNeg','RNeu','RPos'}; % typ eventow do wyciecia | ||
+ | |||
+ | for sub = subjects | ||
+ | |||
+ | % tworzymy nazwę setu dla osoby z numerem 'sub' | ||
+ | if sub < 10 | ||
+ | nameOfSubj = ['B0' num2str(sub)]; | ||
+ | else | ||
+ | nameOfSubj = ['B' num2str(sub)]; | ||
+ | end | ||
+ | |||
+ | % wczytujemy set z danymi | ||
+ | EEG = pop_loadset('filename',[nameOfSubj '.set'],'filepath',directoryIn); | ||
+ | |||
+ | % uaktualniamy zestaw danych globalnych dla eeglaba: | ||
+ | [ALLEEG, EEG, CURRENTSET] = eeg_store( ALLEEG, EEG, 0 ); | ||
+ | EEG = eeg_checkset( EEG); | ||
+ | |||
+ | % zostawiamy tylko potrzebne kanaly EEG , te które są wymienione w zmiennej eegChan: | ||
+ | EEG.data = EEG.data(....); % wycianmy dane | ||
+ | EEG.nbchan = .... ; % uaktualniamy liczbę kanałów | ||
+ | EEG.chanlocs = EEG.chanlocs(1,...); % obcianamy lstę położeń kanałów | ||
+ | |||
+ | % referencja - odjecie sredniej ze wszystkich kanalow | ||
+ | EEG = pop_reref(.... ); | ||
+ | EEG = eeg_checkset(EEG); | ||
+ | |||
+ | % filtracja danych - zastosujemy filtr Butterwartha 2-go rzędu | ||
+ | [b1,a1] = butter(2, .... , 'high'); | ||
+ | [b2,a2] = butter(2, ... , 'low'); | ||
+ | [b3,a3] = butter(2, 2*[49.5 50.5]/EEG.srate, 'stop'); % filtr pasmowozaporowy "notch" na sieć | ||
+ | for channel = 1:size(EEG.data,1) | ||
+ | EEG.data(channel,:) = ... % górnoprzepustowy | ||
+ | EEG.data(channel,:) = ... % dolno | ||
+ | EEG.data(channel,:) = ... % notch | ||
+ | end | ||
+ | EEG = eeg_checkset( EEG ); | ||
+ | |||
+ | % ekstrakcja epok: | ||
+ | % zauważmy, że podajemy w tej funkcji listę eventów wokół których mają być robione wycinki i okresy wycinków | ||
+ | % przy okazji nadajemy setowi nazwę z dodanym stringiem addName: | ||
+ | EEG = pop_epoch( EEG, events, [startT, stopT], 'newname', [nameOfSubj addName], 'epochinfo', 'yes'); | ||
+ | |||
+ | [ALLEEG EEG CURRENTSET] = pop_newset(ALLEEG, EEG, 1,'gui','off'); | ||
+ | EEG = eeg_checkset( EEG ); | ||
+ | |||
+ | % usuniecie baseline | ||
+ | EEG = pop_rmbase(...); | ||
+ | EEG = eeg_checkset( EEG ); | ||
+ | |||
+ | % zapisanie setu | ||
+ | EEG = pop_saveset( EEG, 'filename',[nameOfSubj addName],'filepath',directoryOut); | ||
+ | disp(['filtrowanie i cięcie ----------- dane ' nameOfSubj ' : GOTOWE ----------------']) | ||
+ | |||
+ | end | ||
+ | |||
+ | % przywracamy stan wejściowy zmiennych globalnych eeglaba | ||
+ | EEG0 = EEG; | ||
+ | EEG = EEG0; | ||
+ | ALLEEG = ALLEEG0; | ||
+ | CURRENTSET = CURRENTSET0; | ||
+ | end | ||
</source> | </source> | ||
+ | ==== Wracamy na chwilę do <tt>run_me.m</tt> ==== | ||
<source lang = matlab> | <source lang = matlab> | ||
Linia 84: | Linia 154: | ||
%Tools -> Reject data epochs -> Reject data (all methods) | %Tools -> Reject data epochs -> Reject data (all methods) | ||
%link do tutoriala EEGlab: https://sccn.ucsd.edu/wiki/Chapter_01:_Rejecting_Artifacts#Rejecting_artifacts_in_epoched_data | %link do tutoriala EEGlab: https://sccn.ucsd.edu/wiki/Chapter_01:_Rejecting_Artifacts#Rejecting_artifacts_in_epoched_data | ||
− | + | </source> | |
+ | ====Automatyczne usowanie artefaktow==== | ||
+ | <source lang = matlab> | ||
%% Automatyczne usowanie artefaktow | %% Automatyczne usowanie artefaktow | ||
Wersja z 11:50, 14 mar 2018
Laboratorium_EEG/EEGLAB
Zajęcia bazują na podręczniku praktycznym użytkowania pakietu EEGLAB dostępnym na stronach:
- rozdziały 1-8 z części I: http://sccn.ucsd.edu/wiki/EEGLAB#The_EEGLAB_Tutorial_Outline
- rozdziały z części: http://sccn.ucsd.edu/wiki/Chapter_02:_Writing_EEGLAB_Scripts
Spis treści
Analiza potencjałów wywołanych
Wprowadzenie
W tej części zajęć naszym celem jest przećwiczenie pracy z danymi z eksperymentu w którym badamy potencjały wywołane na grupie osób. Dane pochodzą z eksperymentu opisanego w artykule: Effects of Valence and Origin of Emotions in Word Processing Evidenced by Event Related Potential Correlates in a Lexical Decision Task
Skupimy się na badaniu różnic pomiędzy przetwarzaniem słów i pseudo-słów.
Przebieg eksperymentu
Uczestnicy zostali poinformowani o celu eksperymentu i charakterze pomiaru EEG. Zachęciliśmy ich, aby utrzymywali wygodną postawę i kontrolowali mruganie. Protokół zapewniał 3-s przerwy na normalne mruganie co 10 prób, jak również dwie dłuższe przerwy, których czas trwania jest kontrolowany przez uczestnika, dla odpoczynku i poprawy postawy. Długie przerwy występowały co 270 prób.
Zadaniem było odczytanie bodźców, które pojawiły się na środku ekranu, i zaklasyfikowanie ich jako słów lub pseudo-słów przez naciśnięcie oznaczonych klawiszy na klawiaturze. Rejestrowano typ i opóźnienie odpowiedzi. Pojedynczy blok eksperymentalny zawierał 135 słów i 135 pseudo-słów; ten blok powtórzono trzy razy. Słowa i pseudo słowa były wyświetlane w losowej kolejności we wszystkich blokach. Próby przebiegały w następujący sposób: (1) punkt fiksacji wyświetlany przez 500 ms; (2) bodziec wyświetlany, dopóki uczestnik nie odpowie; (3) pusty ekran wyświetlany przez losowo zmienny odstęp czasu między 1000 a 1100 ms.
Przygotujmy katalogi do analizy. Przydadzą się następujące:
- DANE
- SKRYPTY
- RYSUNKI
Dane dostępne są tutaj.
Proszę je ściągnąć i wypakować do katalogu DANE. Powinno pojawić się 66 plików, po 2 na każdą osobę. Pliki *.set to pliki matlabowe zawierające struktury EEG, pliki *fdt zawierają właściwe sygnały.
SKRYPTY
Etapy analizy, które musimy wykonać to:
- FILTROWANIE i CIECIE sygnałów
- DOBÓR PROGÓW DLA ARTEFAKTOW na podstawie oglądania sygnału i prób z różnymi progami do autoamtu (abnormal values / abnorlmal trends)
- Automatyczne usuwanie artefaktów
- Przygotowanie danych do analiz wg. przynależności do czynnika grupującego
- Wybór obszarów zainteresowania (grup elektrod do analizy)
- Analiza wariancji
- Zilustrowanie wyników analizy
Aby nasze analizy były łatwo powtarzalne przygotujemy główny skrypt np: run_me.m, do którego będziemy dodawać stopniowo skrypty i funkcje kolejnych etapów analizy. Dla czytelności kodu i możliwości uruchamiania osobno poszczególnych etapów posłużymy się technologią komórek w Matlabie (komórkę rozpoczynają dwa procenty '%%');
Zacznijmy zatem pisać run_me.m. Ustawiamy zmienne ze ścieżkami, wybieramy, które dane wchodzą do analizy (część została odrzucona ze względu na nadmiar artefaktów), częstości odcięcia filtrów dolno- i górno-przepustowego, czasy wokół momentu wyświetlenia bodźca, informacja, które kanały zawierają EEG, zakres czasu do korekty poziomu "0":
%% USTAWIENIA:
directoryIn = '../DANE/'; % lokalizacja folderu z surowymi danymi
directoryOut = '../DANE/'; % lokalizacja do zapisu przygotowanych danych
addName = '_epoch'; % dodatek do nazwy pliku przy zapisywaniu przygotowanych danych
subjects = [1:6,8:11 ,15:24,26:37]; % lista badanych do wczytania (z poprawnymi danymi)
LowFreqCut = 0.1; % granica odciecia filtru gornoprzepustowego
HighFreqCut = 30; % granica odciecia filtru dolnoprzepustowego
startT = -0.3; % start do wyciecia epok [s]
stopT = 1.0; % stop do wyciecia epok [s]
eegChan = 1:19; % numery kanalow EEG
baseline = [-200, 0]; % baseline [od,do] [ms] do odjecia
%% FILTROWANIE i CIECIE
eeglab
filtruj_i_tnij(EEG , ALLEEG , CURRENTSET, directoryIn, directoryOut, addName, subjects, LowFreqCut, HighFreqCut, startT, stopT, eegChan, baseline)
filtruj_i_tnij.m
Musimy teraz napisać funkcję filtruj_i_tnij. Musi ona kolejno wczytać nasze zbiory danych, stworzyć montaż kanałów, przefiltrować je, kanał po kanale, wyciąć interesujące nas epoki, wykonać dla nich korektę poziomu "0", i zapisać tak przetworzone dane do nowych "*.set"-ów. Szkilet tej funkcji wygląda następująco, proszę uzupełnić brakujące fragmenty:
function filtruj_i_tnij(EEG, ALLEEG, CURRENTSET, directoryIn, directoryOut, addName, subjects, LowFreqCut, HighFreqCut, startT, stopT, eegChan, baseline)
% zapamiętujemy stan wejściowy zmiennych globalnych eeglaba
EEG0 = EEG;
ALLEEG0 = ALLEEG;
CURRENTSET0 = CURRENTSET;
events = {'Pseudo','ANeg','ANeu','APos','ONeg','ONeu','OPos','RNeg','RNeu','RPos'}; % typ eventow do wyciecia
for sub = subjects
% tworzymy nazwę setu dla osoby z numerem 'sub'
if sub < 10
nameOfSubj = ['B0' num2str(sub)];
else
nameOfSubj = ['B' num2str(sub)];
end
% wczytujemy set z danymi
EEG = pop_loadset('filename',[nameOfSubj '.set'],'filepath',directoryIn);
% uaktualniamy zestaw danych globalnych dla eeglaba:
[ALLEEG, EEG, CURRENTSET] = eeg_store( ALLEEG, EEG, 0 );
EEG = eeg_checkset( EEG);
% zostawiamy tylko potrzebne kanaly EEG , te które są wymienione w zmiennej eegChan:
EEG.data = EEG.data(....); % wycianmy dane
EEG.nbchan = .... ; % uaktualniamy liczbę kanałów
EEG.chanlocs = EEG.chanlocs(1,...); % obcianamy lstę położeń kanałów
% referencja - odjecie sredniej ze wszystkich kanalow
EEG = pop_reref(.... );
EEG = eeg_checkset(EEG);
% filtracja danych - zastosujemy filtr Butterwartha 2-go rzędu
[b1,a1] = butter(2, .... , 'high');
[b2,a2] = butter(2, ... , 'low');
[b3,a3] = butter(2, 2*[49.5 50.5]/EEG.srate, 'stop'); % filtr pasmowozaporowy "notch" na sieć
for channel = 1:size(EEG.data,1)
EEG.data(channel,:) = ... % górnoprzepustowy
EEG.data(channel,:) = ... % dolno
EEG.data(channel,:) = ... % notch
end
EEG = eeg_checkset( EEG );
% ekstrakcja epok:
% zauważmy, że podajemy w tej funkcji listę eventów wokół których mają być robione wycinki i okresy wycinków
% przy okazji nadajemy setowi nazwę z dodanym stringiem addName:
EEG = pop_epoch( EEG, events, [startT, stopT], 'newname', [nameOfSubj addName], 'epochinfo', 'yes');
[ALLEEG EEG CURRENTSET] = pop_newset(ALLEEG, EEG, 1,'gui','off');
EEG = eeg_checkset( EEG );
% usuniecie baseline
EEG = pop_rmbase(...);
EEG = eeg_checkset( EEG );
% zapisanie setu
EEG = pop_saveset( EEG, 'filename',[nameOfSubj addName],'filepath',directoryOut);
disp(['filtrowanie i cięcie ----------- dane ' nameOfSubj ' : GOTOWE ----------------'])
end
% przywracamy stan wejściowy zmiennych globalnych eeglaba
EEG0 = EEG;
EEG = EEG0;
ALLEEG = ALLEEG0;
CURRENTSET = CURRENTSET0;
end
Wracamy na chwilę do run_me.m
%% DOBOR PROGOW DLA ARTEFAKTOW na podstawie ogladanie sygnalu i prob z roznymi progami do autoamtu (abnormal values / abnorlmal trends)
EEG = pop_loadset('filename','B01_epoch.set','filepath',directoryOut);
[ALLEEG, EEG, CURRENTSET] = eeg_store(ALLEEG, EEG, 0);
EEG = eeg_checkset( EEG);
eeglab redraw
%Tools -> Reject data epochs -> Reject data (all methods)
%link do tutoriala EEGlab: https://sccn.ucsd.edu/wiki/Chapter_01:_Rejecting_Artifacts#Rejecting_artifacts_in_epoched_data
Automatyczne usowanie artefaktow
%% Automatyczne usowanie artefaktow
abnormalValue = 60; % wartosc [uV] powyzej/ponizej(na minusie) ktorej oznaczane sa artefakty
abnormalTrend = [50, 0.3]; % [max slope [uV/epoch], R squared limit (0-1)]
markArtefacts(EEG , ALLEEG , CURRENTSET, directoryOut, addName, subjects, abnormalValue, abnormalTrend)
%% przygotowanie danych
data = prepareData(EEG , ALLEEG , CURRENTSET, directoryOut, addName, subjects);
% data - macierz z danymi o rozmiarze [2(1.Words|2.Pseudo) x subjects x kanaly x probki]
%% GFP
timeMarks = [];
plotGFP(data,startT,stopT,EEG.srate,timeMarks,EEG.chanlocs)
timeMarks = [70,120,225,290,375,670];
plotGFP(data,startT,stopT,EEG.srate,timeMarks,EEG.chanlocs)
%% ANOVA
ROI.channels = [[18,13];[15,10];[19,17];[9,4];[11,6]];
ROI.labels = {'LF','CF','RF','LP','RP'};
dataANOVA = prepareDataForANOVA(data,timeMarks,ROI,startT,stopT,EEG.srate);
% dataANOVA - [5 okien czasowych x obserwacje x |1.srednia amplituda 2. subject 3.slowo(0)/pseudo(1) 4.ROI|]
for window = 1:(length(timeMarks)-1)
dANOVA = squeeze(dataANOVA(window,:,:));
stats = anova(dANOVA(:,1),dANOVA(:,2),dANOVA(:,3),dANOVA(:,4),{'type','ROI'});
stats2 = PostHocStats(stats,dANOVA,ROI);
disp(['------------- Time window: ' num2str(timeMarks(window)) '-' num2str(timeMarks(window+1)) '-------------'])
disp(stats)
disp(stats2)
end
%% rysowanie ERP w ROI
plotERP_ROI(data,startT,stopT,EEG.srate,timeMarks,ROI)
Laboratorium_EEG/EEGLAB