-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathPracticalMEEG_Session_4_ERP_Source_Analysis_TOEDIT.m
More file actions
217 lines (179 loc) · 11.9 KB
/
Copy pathPracticalMEEG_Session_4_ERP_Source_Analysis_TOEDIT.m
File metadata and controls
217 lines (179 loc) · 11.9 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
% Practical MEEG 2022
% Wakeman & Henson Data analysis: Ereprocess Data Session #1
% Authors: Ramon Martinez-Cancino, Brain Products, 2022
% Arnaud Delorme, SCCN, 2022
% Johanna Wagner, Zander Labs, 2022
%
% Copyright (C) 2022 Johanna Wagner
%
% This program is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
% Clearing all is recommended to avoid variable not being erased between calls
clear;
clear globals;
% Path to data below. Using relative paths so no need to update.
RootFolder = fileparts(pwd); % Getting root folder
path2data = fullfile(RootFolder,'Data', 'sub-01'); % Path to data
filename = 'wh_S01_run_01_preprocessing_data_session_1_out.set';
% Start EEGLAB
[ALLEEG, EEG, CURRENTSET] = eeglab;
% Loading data
EEG = pop_loadset('filename', filename,'filepath',path2data)
[ALLEEG EEG CURRENTSET] = pop_newset(ALLEEG, EEG, 1);
%% ERP analysis
%% Extract event-locked trials using events listed in 'eventlist'
ALLEEG(2) = pop_epoch( ALLEEG(1), {'Famous'}, [-1 2], 'newname', 'Famous Epoched', 'epochinfo', 'yes');
ALLEEG(3) = pop_epoch( ALLEEG(1), {'Unfamiliar'}, [-1 2], 'newname', 'Unfamiliar Epoched', 'epochinfo', 'yes');
ALLEEG(4) = pop_epoch( ALLEEG(1), {'Scrambled'}, [-1 2], 'newname', 'Scrambled Epoched', 'epochinfo', 'yes');
%% Perform baseline correction
ALLEEG(2) = pop_rmbase(ALLEEG(2), [-1000 0]);
ALLEEG(3) = pop_rmbase(ALLEEG(3), [-1000 0]);
ALLEEG(4) = pop_rmbase(ALLEEG(4), [-1000 0]);
%% Clean data by rejecting epochs.
[ALLEEG(2), rejindx] = pop_eegthresh(ALLEEG(2), 1, 1:ALLEEG(2).nbchan, -400, 400, ALLEEG(2).xmin, ALLEEG(2).xmax, 0, 1);
[ALLEEG(3), rejindx] = pop_eegthresh(ALLEEG(3), 1, 1:ALLEEG(3).nbchan, -400, 400, ALLEEG(3).xmin, ALLEEG(3).xmax, 0, 1);
[ALLEEG(4), rejindx] = pop_eegthresh(ALLEEG(4), 1, 1:ALLEEG(4).nbchan, -400, 400, ALLEEG(4).xmin, ALLEEG(4).xmax, 0, 1);
%% Save dataset
EEG_famous = pop_saveset( ALLEEG(2),'filename', 'wh_S01_run_01_ERP_Analysis_Session_2_famous_out.set','filepath',path2data);
EEG_unfamiliar = pop_saveset( ALLEEG(3),'filename', 'wh_S01_run_01_ERP_Analysis_Session_2_unfamiliar_out.set','filepath',path2data);
EEG_scrambled = pop_saveset( ALLEEG(4),'filename', 'wh_S01_run_01_ERP_Analysis_Session_2_scrambled_out.set','filepath',path2data);
%% ----------------------
%% BELOW IS PLOTTING ONLY
%% ----------------------
%% plot ERP scalp distribution
figure; pop_timtopo(ALLEEG(2), [-100 600], [NaN], 'ERP data and scalp maps of Famous Epoched');
figure; pop_timtopo(ALLEEG(3), [-100 600], [NaN], 'ERP data and scalp maps of Unfamiliar Epoched');
figure; pop_timtopo(ALLEEG(4), [-100 600], [NaN], 'ERP data and scalp maps of Scrambled Epoched');
%% plot Component 5 contributions to ERP scalp envelope with channel projection
pop_envtopo(ALLEEG(2), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'compnums',[5],'title', 'Famous','electrodes','off', 'plotproj', 1);
%% plot Component contributions to ERP scalp envelope
pop_envtopo(ALLEEG(2), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'title', 'Largest ERP components of Famous Epoched','electrodes','off');
pop_envtopo(ALLEEG(3), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'title', 'Largest ERP components of Unfamiliar Epoched','electrodes','off');
pop_envtopo(ALLEEG(4), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'title', 'Largest ERP components of Scrambled Epoched','electrodes','off');
%% plot Component contributions to ERP scalp envelope - removing artefact ICs
pop_envtopo(ALLEEG(2), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Largest ERP components of Famous Epoched','electrodes','off');
pop_envtopo(ALLEEG(3), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Largest ERP components of Unfamiliar Epoched','electrodes','off');
pop_envtopo(ALLEEG(4), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Largest ERP components of Scrambled Epoched','electrodes','off');
%% plot Component contributions to ERP scalp envelope - removing artefact ICs - plot differences between conditions
pop_envtopo(([ALLEEG(2) ALLEEG(4)]), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Famous - Scrambled','electrodes','off');
pop_envtopo(([ALLEEG(2) ALLEEG(3)]), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Famous - Unfamiliar','electrodes','off');
pop_envtopo(([ALLEEG(3) ALLEEG(4)]), [-100 600] ,'limcontrib',[0 400],'compsplot',[7],'subcomps',[1 6],'title', 'Unfamiliar - Scrambled','electrodes','off');
%% Identify Brain ICs using IC Label classification results
if isfield(ALLEEG(1).etc, 'ic_classification')
[M,I] = max(ALLEEG(1).etc.ic_classification.ICLabel.classifications,[],2); % Use max prob for classification
Brain_comps = find(I == find(strcmp(ALLEEG(1).etc.ic_classification.ICLabel.classes, 'Brain')));
%% Subtract artefactual components from the EEG
ALLEEG(2) = pop_subcomp( ALLEEG(2), Brain_comps, 0, 1);
ALLEEG(3) = pop_subcomp( ALLEEG(3), Brain_comps, 0, 1);
ALLEEG(4) = pop_subcomp( ALLEEG(4), Brain_comps, 0, 1);
end
%% Rename datasets
ALLEEG(2) = pop_editset(ALLEEG(2), 'setname', 'Famous', 'run', []);
ALLEEG(3) = pop_editset(ALLEEG(3), 'setname', 'Unfamiliar', 'run', []);
ALLEEG(4) = pop_editset(ALLEEG(4), 'setname', 'Scrambled', 'run', []);
%% plot ERP scalp distribution
figure; pop_timtopo(ALLEEG(2), [-100 600], [NaN], 'Famous');
figure; pop_timtopo(ALLEEG(3), [-100 600], [NaN], 'Unfamiliar');
figure; pop_timtopo(ALLEEG(4), [-100 600], [NaN], 'Scrambled');
%% plot ERP scalp distribution at each ERP peak
figure; pop_timtopo(ALLEEG(2), [-100 600], [120 170 250], 'Famous');
figure; pop_timtopo(ALLEEG(3), [-100 600], [120 170 250], 'Unfamiliar');
figure; pop_timtopo(ALLEEG(4), [-100 600], [120 170 250], 'Scrambled');
%% plot 3 largest contributing ICs to ERP
pop_envtopo(ALLEEG(2), [-100 600] ,'limcontrib',[0 400],'compsplot',[3],'title', 'Famous','electrodes','off');
pop_envtopo(ALLEEG(3), [-100 600] ,'limcontrib',[0 400],'compsplot',[3],'title', 'Unfamiliar','electrodes','off');
pop_envtopo(ALLEEG(4), [-100 600] ,'limcontrib',[0 400],'compsplot',[3],'title', 'Scrambled','electrodes','off');
%% Visualize channel ERPs in 2D
pop_topoplot(ALLEEG(2), 1, [25:25:300] ,'Famous',[3 4] ,0,'electrodes','on');
pop_topoplot(ALLEEG(3), 1, [25:25:300] ,'Unfamiliar',[3 4] ,0,'electrodes','on');
pop_topoplot(ALLEEG(4), 1, [25:25:300] ,'Scrambled',[3 4] ,0,'electrodes','on');
%% Plot channel ERPs in topographic array
figure; pop_plottopo(ALLEEG(2), [1:EEG.nbchan] , 'Famous', 0, 'ydir',1);
figure; pop_plottopo(ALLEEG(3), [1:EEG.nbchan] , 'Unfamiliar', 0, 'ydir',1);
figure; pop_plottopo(ALLEEG(4), [1:EEG.nbchan] , 'Scrambled', 0, 'ydir',1);
%% plot average ERPs for each condition with standard deviation
% find channel index of eeg065
Chanind = find(strcmp({ALLEEG(2).chanlocs.labels},'EEG065'));
if isempty(Chanind)
Chanind = 1;
end
% create timevector for plotting
[val, indL] = min(abs(ALLEEG(2).times+200)); %get timepoints for -200 and 800 Latencies
[val, indU] = min(abs(ALLEEG(2).times-800));
timevec = ALLEEG(2).times(indL:indU); % create timevector
% create datavectors for plotting each condition
av_datavecF = mean(ALLEEG(2).data(Chanind,indL:indU,:),3); % average
std_datavecF = std(ALLEEG(2).data(Chanind,indL:indU,:),1,3); % standard deviation
figure;
X2 = [[timevec],fliplr([timevec])]; %#create continuous x value array for plotting
Y2 = [av_datavecF-std_datavecF,fliplr(av_datavecF+std_datavecF)]; %#create y values for out and then back
fill(X2,Y2,[153/255 204/255 255/255]);
hold on
plot(timevec,av_datavecF, 'b', 'LineWidth',2)
xline(0, 'LineWidth',2)
yline(0, 'LineWidth',2)
xlabel('Latency ms')
ylabel('mu Volt')
title([ 'famous faces channel ' EEG.chanlocs(Chanind).labels ]);
set(gca, 'FontSize', 15)
av_datavecU = mean(ALLEEG(3).data(Chanind,indL:indU,:),3); % average
std_datavecU = std(ALLEEG(3).data(Chanind,indL:indU,:),1,3); % standard deviation
figure;
X2 = [[timevec],fliplr([timevec])]; %#create continuous x value array for plotting
Y2 = [av_datavecU-std_datavecU,fliplr(av_datavecU+std_datavecU)]; %#create y values for out and then back
fill(X2,Y2,[153/255 204/255 255/255]);
hold on
plot(timevec,av_datavecU, 'b', 'LineWidth',2)
xline(0, 'LineWidth',2)
yline(0, 'LineWidth',2)
xlabel('Latency ms')
ylabel('mu Volt')
title([ 'unfamiliar faces channel ' EEG.chanlocs(Chanind).labels ]);
set(gca, 'FontSize', 15)
av_datavecS = mean(ALLEEG(4).data(Chanind,indL:indU,:),3); % average
std_datavecS = std(ALLEEG(4).data(Chanind,indL:indU,:),1,3); % standard deviation
figure;
X2 = [[timevec],fliplr([timevec])]; %#create continuous x value array for plotting
Y2 = [av_datavecS-std_datavecS,fliplr(av_datavecS+std_datavecS)]; %#create y values for out and then back
fill(X2,Y2,[153/255 204/255 255/255]);
hold on
plot(timevec,av_datavecS, 'b', 'LineWidth',2)
xline(0, 'LineWidth',2)
yline(0, 'LineWidth',2)
xlabel('Latency ms')
ylabel('mu Volt')
title([ 'scrambled faces channel ' EEG.chanlocs(Chanind).labels ]);
set(gca, 'FontSize', 15)
%% plot superimposed ERPs
figure;plot(timevec,av_datavecF, 'LineWidth',2, 'color', 'r'); hold on
plot(timevec,av_datavecU, 'LineWidth',2, 'color', 'b')
plot(timevec,av_datavecS, 'LineWidth',2, 'color', 'g')
fillcurves(timevec,av_datavecF-std_datavecF,av_datavecF+std_datavecF, 'r', 0.2);
fillcurves(timevec,av_datavecU-std_datavecU,av_datavecU+std_datavecU, 'b', 0.2);
fillcurves(timevec,av_datavecS-std_datavecS,av_datavecS+std_datavecS, 'g', 0.2);
xline(0, 'LineWidth',2)
yline(0, 'LineWidth',2)
xlabel('Latency ms')
ylabel('mu Volt')
legend('famous', 'unfamiliar', 'scrambled')
set(gca, 'FontSize', 15)
title([ 'face types for channel ' EEG.chanlocs(Chanind).labels ]);
%% ERPimage
figure; pop_erpimage(ALLEEG(2),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{},[],'' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );
figure; pop_erpimage(ALLEEG(3),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{},[],'' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );
figure; pop_erpimage(ALLEEG(4),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{},[],'' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );
% sort by event latency (when there are enough of left_nonsym and right_sym events
% figure; pop_erpimage(ALLEEG(2),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{ 'left_nonsym' 'right_sym'},[],'latency' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );
% figure; pop_erpimage(ALLEEG(3),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{ 'left_nonsym' 'right_sym'},[],'latency' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );
% figure; pop_erpimage(ALLEEG(4),1, [Chanind],[[]],EEG.chanlocs(Chanind).labels,3,1,{ 'left_nonsym' 'right_sym'},[],'latency' ,'yerplabel','\muV','erp','on','limits',[-100 1200 NaN NaN NaN NaN NaN NaN] ,'cbar','on','topo', { [Chanind] EEG.chanlocs EEG.chaninfo } );