cmpt_ic_art_ftrs() - Computes features of independent components (ICs) that
are useful for identifying ICs that correspond to EEG
artifacts and should be removed the data. These features
are:
**TOPOGRAPHY FEATURES**
1. Absolute cosine similarity to prototypical blink topography
("blink_topo"): Takes a value between 0 and 1 with 1 being
identical to the blink topography prototype.
2. Absolute cosine similarity to prototypical heart artifact
topography ("heart_topo"): Takes a value between 0 and 1
with 1 being identical to the heart artifact topography
prototype.
3. Absolute cosine similarity to prototypical horizontal eye
movement topography ("he_topo"): Takes a value between 0
and 1, with 1 being identical to the horizontal eye movement
topography prototype.
4. Extreme topography channel weight ("xtrm_chan"): The
proportion of a topography's length accounted for by a
single channel. This feature takes a value between 1/n and
1, where n is the number of channels. A value of 1 means
that the IC's topography is zero at all electrodes save one.
This feature tends to be large for ICs that correspond to a
bad channel (e.g., an electrode that changed impedence during
the experiment) or EMG ICs.
**ACTIVATION SPECTRUM FEATURES**
5. EMG-like spectra ("emg_spectra"): Mean spectral power
between 25-55 Hz minus mean spectral power between 0-17 Hz.
This feature is large for ICs composed primarily of EMG.
Spectral power is first transformed to t-scores to
normalize this feature across participants.
6. EOG-like spectra ("eog_spectra"): Mean spectral power
between 0-6 Hz minus mean spectral power between 8-18 Hz.
This feature is large for ICs composed primarily of EOG.
Spectral power is first transformed to t-scores to
normalize this feature across participants.
**ACTIVATION MAGNITUDE FEATURE**
7. Mean standard deviation at each time point ("mn_std"):
This feature is large for ICs characterized by highly
variable activations (typically due to blinks, eye
movements, or extreme EMG). Only computable for epoched
data.
Usage:
>> EEG=cmpt_ic_art_ftrs(EEG_or_setfname,verblevel)
Required Input:
EEG_or_setfname - An EEGLAB EEG struct variable or the filename of an
EEGLAB EEG struct variable (i.e., a *.set file).
Optional Inputs:
'verblevel' - An integer specifiying the amount of information you want
functions to provide about what they are doing during runtime.
Options are:
0 - quiet, only show errors, warnings, and EEGLAB reports
1 - stuff anyone should probably know
2 - stuff you should know the first time you start working
with a data set {default value if not already globally
specified}
3 - stuff that might help you debug (show all reports)
'reset_thresh' - [0 | 1] If "1", the IC feature thresholds are
reset to "NaN" (or added and set to "NaN" if those
fields have not been added yet to the EEG variable).
Output:
EEG - An EEGLAB EEG struct variable. Same as input EEG but with changes
to EEG.icaact, EEG.icawinv, EEG.icasphere, EEG.icaweights, EEG.icspectra
EEG.icfeatures and EEG.reject.thresh* fields. The latermost field
is simply initialized to "NaN" rejection thresholds that have no
effect. The EEG.ica* and EEG.icspectra features are changed only
if the ICs have not been normalized to have topographies of unit
length.
Global Variable:
VERBLEVEL - level of verbosity (i.e., tells functions how much
how much to report about what they're doing during
runtime) set by the optional function argument 'verblevel'
Notes:
-Use the function reject_ics_by_ftr.m to objectively label ICs as
artifacts based on the features computed here. Use set_ic_ftr_thresh.m to
interactively set feature thresholds that these labels will be based on.
-This function needs to load artifact topography prototypes from the file
topo_templates.mat. This file should be in the default Kutaslab MATLAB
path.
-Becuase these features are sensitive to how ICs are scaled, ICs are
normalized to have topographies of unit length before the features are
computed. Doing this means that the IC sum squared activations equals
that of their scalp projections.
Author:
David Groppe
Kutaslab, 12/20090001 function EEG=cmpt_ic_art_ftrs(EEG_or_setfname,verblevel,reset_thresh) 0002 % cmpt_ic_art_ftrs() - Computes features of independent components (ICs) that 0003 % are useful for identifying ICs that correspond to EEG 0004 % artifacts and should be removed the data. These features 0005 % are: 0006 % **TOPOGRAPHY FEATURES** 0007 % 1. Absolute cosine similarity to prototypical blink topography 0008 % ("blink_topo"): Takes a value between 0 and 1 with 1 being 0009 % identical to the blink topography prototype. 0010 % 2. Absolute cosine similarity to prototypical heart artifact 0011 % topography ("heart_topo"): Takes a value between 0 and 1 0012 % with 1 being identical to the heart artifact topography 0013 % prototype. 0014 % 3. Absolute cosine similarity to prototypical horizontal eye 0015 % movement topography ("he_topo"): Takes a value between 0 0016 % and 1, with 1 being identical to the horizontal eye movement 0017 % topography prototype. 0018 % 4. Extreme topography channel weight ("xtrm_chan"): The 0019 % proportion of a topography's length accounted for by a 0020 % single channel. This feature takes a value between 1/n and 0021 % 1, where n is the number of channels. A value of 1 means 0022 % that the IC's topography is zero at all electrodes save one. 0023 % This feature tends to be large for ICs that correspond to a 0024 % bad channel (e.g., an electrode that changed impedence during 0025 % the experiment) or EMG ICs. 0026 % 0027 % **ACTIVATION SPECTRUM FEATURES** 0028 % 5. EMG-like spectra ("emg_spectra"): Mean spectral power 0029 % between 25-55 Hz minus mean spectral power between 0-17 Hz. 0030 % This feature is large for ICs composed primarily of EMG. 0031 % Spectral power is first transformed to t-scores to 0032 % normalize this feature across participants. 0033 % 6. EOG-like spectra ("eog_spectra"): Mean spectral power 0034 % between 0-6 Hz minus mean spectral power between 8-18 Hz. 0035 % This feature is large for ICs composed primarily of EOG. 0036 % Spectral power is first transformed to t-scores to 0037 % normalize this feature across participants. 0038 % 0039 % **ACTIVATION MAGNITUDE FEATURE** 0040 % 7. Mean standard deviation at each time point ("mn_std"): 0041 % This feature is large for ICs characterized by highly 0042 % variable activations (typically due to blinks, eye 0043 % movements, or extreme EMG). Only computable for epoched 0044 % data. 0045 % 0046 % 0047 % Usage: 0048 % >> EEG=cmpt_ic_art_ftrs(EEG_or_setfname,verblevel) 0049 % 0050 % Required Input: 0051 % EEG_or_setfname - An EEGLAB EEG struct variable or the filename of an 0052 % EEGLAB EEG struct variable (i.e., a *.set file). 0053 % 0054 % Optional Inputs: 0055 % 'verblevel' - An integer specifiying the amount of information you want 0056 % functions to provide about what they are doing during runtime. 0057 % Options are: 0058 % 0 - quiet, only show errors, warnings, and EEGLAB reports 0059 % 1 - stuff anyone should probably know 0060 % 2 - stuff you should know the first time you start working 0061 % with a data set {default value if not already globally 0062 % specified} 0063 % 3 - stuff that might help you debug (show all reports) 0064 % 'reset_thresh' - [0 | 1] If "1", the IC feature thresholds are 0065 % reset to "NaN" (or added and set to "NaN" if those 0066 % fields have not been added yet to the EEG variable). 0067 % 0068 % 0069 % Output: 0070 % EEG - An EEGLAB EEG struct variable. Same as input EEG but with changes 0071 % to EEG.icaact, EEG.icawinv, EEG.icasphere, EEG.icaweights, EEG.icspectra 0072 % EEG.icfeatures and EEG.reject.thresh* fields. The latermost field 0073 % is simply initialized to "NaN" rejection thresholds that have no 0074 % effect. The EEG.ica* and EEG.icspectra features are changed only 0075 % if the ICs have not been normalized to have topographies of unit 0076 % length. 0077 % 0078 % 0079 % Global Variable: 0080 % VERBLEVEL - level of verbosity (i.e., tells functions how much 0081 % how much to report about what they're doing during 0082 % runtime) set by the optional function argument 'verblevel' 0083 % 0084 % Notes: 0085 % -Use the function reject_ics_by_ftr.m to objectively label ICs as 0086 % artifacts based on the features computed here. Use set_ic_ftr_thresh.m to 0087 % interactively set feature thresholds that these labels will be based on. 0088 % 0089 % -This function needs to load artifact topography prototypes from the file 0090 % topo_templates.mat. This file should be in the default Kutaslab MATLAB 0091 % path. 0092 % 0093 % -Becuase these features are sensitive to how ICs are scaled, ICs are 0094 % normalized to have topographies of unit length before the features are 0095 % computed. Doing this means that the IC sum squared activations equals 0096 % that of their scalp projections. 0097 % 0098 % Author: 0099 % David Groppe 0100 % Kutaslab, 12/2009 0101 0102 0103 0104 global VERBLEVEL 0105 0106 if ~nargin, 0107 error('cmpt_ic_art_ftrs.m requires at least one input.'); 0108 end 0109 0110 if (nargin>1), 0111 VERBLEVEL=verblevel; 0112 elseif isempty(VERBLEVEL), 0113 VERBLEVEL=2; 0114 end 0115 0116 if (nargin<3), 0117 reset_thresh=0; 0118 elseif ischar(reset_thresh) 0119 if strcmpi(reset_thresh,'on') || strcmpi(reset_thresh,'yes'), 0120 reset_thresh=1; 0121 end 0122 end 0123 0124 if ischar(EEG_or_setfname) 0125 EEG=pop_loadset(fname); 0126 elseif isstruct(EEG_or_setfname), 0127 EEG=EEG_or_setfname; %note this does not create a second copy of EEG in memory since EEG_or_setfname was never modified 0128 %(in other words it's not memory inefficient) 0129 clear EEG_or_setfname; 0130 else 0131 error('Input argument "EEG_or_setfname" needs to be a set filename of an EEG struct variable.'); 0132 end 0133 0134 if isempty(EEG.icawinv), 0135 error('There are no ICA components saved with this dataset. You need to run ICA first.'); 0136 end 0137 n_ic=size(EEG.icawinv,2); 0138 0139 EEG=norm_ics(EEG,'topo length'); 0140 0141 n_chan=size(EEG.data,1); 0142 n_freq=length(EEG.icfreqs); 0143 spectra_tscores=EEG.icspectra-repmat(mean(EEG.icspectra),n_ic,1); 0144 spectra_tscores=spectra_tscores./repmat(std(spectra_tscores),n_ic,1); 0145 0146 % Non-muscle frequency range 0-17 Hz 0147 nnz_freqs(1)=min(find(EEG.icfreqs>=0)); 0148 nnz_freqs(2)=max(find(EEG.icfreqs<=17)); 0149 0150 % Muscle frequency range 25-55 Hz 0151 nz_freqs(1)=min(find(EEG.icfreqs>=25)); 0152 nz_freqs(2)=max(find(EEG.icfreqs<=55)); 0153 0154 % EOG frequency range 0-6 Hz 0155 eog_freqs(1)=min(find(EEG.icfreqs>=0)); 0156 eog_freqs(2)=max(find(EEG.icfreqs<=6)); 0157 0158 % non-EOG frequency range 8-18 Hz 0159 neog_freqs(1)=min(find(EEG.icfreqs>=8)); 0160 neog_freqs(2)=max(find(EEG.icfreqs<=18)); 0161 0162 % Change less conventional EOG channel names to single convention 0163 EEGlocs=EEG.chanlocs; 0164 for c=1:n_chan, 0165 if strcmpi(EEGlocs(c).labels,'lle') 0166 EEGlocs(c).labels='LLEy'; 0167 elseif strcmpi(EEGlocs(c).labels,'rle') 0168 EEGlocs(c).labels='RLEy'; 0169 elseif strcmpi(EEGlocs(c).labels,'lhe') || strcmpi(EEGlocs(c).labels,'lhz') 0170 EEGlocs(c).labels='LHEy'; 0171 elseif strcmpi(EEGlocs(c).labels,'rhe') || strcmpi(EEGlocs(c).labels,'rhz') 0172 EEGlocs(c).labels='RHEy'; 0173 end 0174 end 0175 0176 % Align channels with template 0177 load /usr/local/matlab-toolboxes/matlabmk/topo_templates; 0178 if n_chan<=31, 0179 % Use 31 channel templates 0180 HE_topo_tmplt=HE_topo_tmplt31; 0181 blink_topo_tmplt=blink_topo_tmplt31; 0182 heart_topo_tmplt=heart_topo_tmplt31; 0183 tmplt_use_chans=zeros(1,31); 0184 EEG_crspnd_chans=zeros(1,31); 0185 crspnd_ct=0; 0186 for a=1:31, 0187 fnd=0; 0188 for c=1:n_chan, 0189 if strcmpi(EEGlocs(c).labels,locs31(a).labels), 0190 crspnd_ct=crspnd_ct+1; 0191 fnd=1; 0192 tmplt_use_chans(crspnd_ct)=a; 0193 EEG_crspnd_chans(crspnd_ct)=c; 0194 break; 0195 end 0196 end 0197 if ~fnd 0198 VerbReport(sprintf('Channel %s was not found in this EEG variable. It will be ignored.', ... 0199 locs31(a).labels),1,VERBLEVEL); 0200 end 0201 end 0202 tmplt_use_chans=tmplt_use_chans(1:crspnd_ct); 0203 EEG_crspnd_chans=EEG_crspnd_chans(1:crspnd_ct); 0204 else 0205 % Use 64 channel templates 0206 HE_topo_tmplt=HE_topo_tmplt64; 0207 blink_topo_tmplt=blink_topo_tmplt64; 0208 heart_topo_tmplt=heart_topo_tmplt64; 0209 tmplt_use_chans=zeros(1,64); 0210 EEG_crspnd_chans=zeros(1,64); 0211 crspnd_ct=0; 0212 for a=1:64, 0213 fnd=0; 0214 for c=1:n_chan, 0215 if strcmpi(EEGlocs(c).labels,locs64(a).labels), 0216 crspnd_ct=crspnd_ct+1; 0217 fnd=1; 0218 tmplt_use_chans(crspnd_ct)=a; 0219 EEG_crspnd_chans(crspnd_ct)=c; 0220 break; 0221 end 0222 end 0223 if ~fnd 0224 VerbReport(sprintf('Channel %s was not found in this EEG variable. It will be ignored.', ... 0225 locs64(a).labels),1,VERBLEVEL); 0226 end 0227 end 0228 tmplt_use_chans=tmplt_use_chans(1:crspnd_ct); 0229 EEG_crspnd_chans=EEG_crspnd_chans(1:crspnd_ct); 0230 end 0231 0232 for a=1:n_ic, 0233 EEG.icfeatures(a).xtrm_chan=max(EEG.icawinv(:,a).^2); 0234 EEG.icfeatures(a).emg_spectra=mean(spectra_tscores(a,nz_freqs(1):nz_freqs(2))) ... 0235 -mean(spectra_tscores(a,nnz_freqs(1):nnz_freqs(2))); 0236 EEG.icfeatures(a).eog_spectra=mean(spectra_tscores(a,eog_freqs(1):eog_freqs(2))) ... 0237 -mean(spectra_tscores(a,neog_freqs(1):neog_freqs(2))); 0238 0239 % re-normalize topography templates to unit length in case some 0240 % channels are missing in current EEG variable 0241 tpo=EEG.icawinv(EEG_crspnd_chans,a)/norm(EEG.icawinv(EEG_crspnd_chans,a)); 0242 EEG.icfeatures(a).blink_topo=abs(tpo'*blink_topo_tmplt(tmplt_use_chans)/norm(blink_topo_tmplt(tmplt_use_chans))); 0243 EEG.icfeatures(a).heart_topo=abs(tpo'*heart_topo_tmplt(tmplt_use_chans)); 0244 EEG.icfeatures(a).he_topo=abs(tpo'*HE_topo_tmplt(tmplt_use_chans)); 0245 0246 %ERP mean standard deviation 0247 s=size(EEG.icaact); 0248 if length(s)==3, 0249 ic_act=reshape(EEG.icaact(a,:,:),s(2),s(3))'; 0250 sd=std(ic_act); 0251 EEG.icfeatures(a).mn_std=mean(sd); %average across time points 0252 else 0253 %can't compute erp_nsr for non-epoched data 0254 watchit('Cannot compute IC feature "mean time point standard deviation" for continuous data (i.e., data need to be epoched).'); 0255 end 0256 end 0257 0258 0259 0260 if reset_thresh, 0261 %add fields for possible rejection thresholds 0262 EEG.reject.thresh_ic_xtrm_chan=NaN; 0263 EEG.reject.thresh_ic_mn_std=NaN; 0264 EEG.reject.thresh_ic_emg_spectra=NaN; 0265 EEG.reject.thresh_ic_eog_spectra=NaN; 0266 EEG.reject.thresh_blink_topo=NaN; 0267 EEG.reject.thresh_heart_topo=NaN; 0268 EEG.reject.thresh_he_topo=NaN; 0269 end 0270 0271 % If thresholds haven't been set, add em 0272 thresh{1}='thresh_ic_xtrm_chan'; 0273 thresh{2}='thresh_ic_mn_std'; 0274 thresh{3}='thresh_ic_emg_spectra'; 0275 thresh{4}='thresh_ic_eog_spectra'; 0276 thresh{5}='thresh_blink_topo'; 0277 thresh{6}='thresh_heart_topo'; 0278 thresh{7}='thresh_he_topo'; 0279 0280 fldnms=fieldnames(EEG.reject); 0281 n_ftrs=length(thresh); 0282 found_ftrs=zeros(1,n_ftrs); 0283 for dg=1:length(fldnms), 0284 for dg2=1:n_ftrs, 0285 if strcmpi(thresh{dg2},fldnms{dg}) 0286 found_ftrs(dg2)=1; 0287 break; 0288 end 0289 end 0290 end 0291 for dg=1:n_ftrs, 0292 if ~found_ftrs(dg), 0293 cmd=['EEG.reject.' thresh{dg} '=NaN;']; 0294 eval(cmd); 0295 end 0296 end