Hello everyone,
for my EEG data I want to artifact reject based on an amplitude threshold (-100 to +100 ”V) and a joint probability threshold (SD = 5), but only for the first second after stimulus onset. My trials are 9 seconds long, and of course I can just epoch the data to the first second and run the annotation on it. But after this procedure I need my 9 sec epochs back. I am not sure how to tell python to integrate the corrected first second to the rest 8 seconds.
Thank you!
Franzi
You could:
annotate_amplitude() only works with continuous (ârawâ) data.
right, thanks for the reminder Richard.
I think in that case Iâd recommend using .get_data, numpy, and custom code.
Something like:
time_of_interest = epochs.time_as_index([0, 1]) # from 0 to 1 seconds
data = epochs.get_data() # (n_epochs, n_channels, n_times)
for epoch in data:
# epoch is a numpy array (n_channels, n_times)
epoch_window = epoch[:, time_of_interest[0]:time_of_interest[1]]
# epoch_window is the first second of your data
... # and so on
Probably should use epochs.time_as_index() here to get the correct time indices!
Thank you for your help
I tried to run the following code
time_of_interest = data.time_as_index([0, 1]) # from 0 to 1 seconds
data1 = data.get_data() # (n_epochs, n_channels, n_times)
for epoch in data1:
# epoch is a numpy array (n_channels, n_times)
epoch_window = data1[:, time_of_interest[0]:time_of_interest[1]]
reject_criteria = dict(eeg=100e-6)
flat_criteria = dict(
eeg=-100e-6)
data_annot = mne.Epochs(epoch_window, events, tmin=0, tmax=1, reject_tmax=0,
reject=reject_criteria, flat=flat_criteria,
reject_by_annotation=False, preload=True)
data_annot.plot_drop_log()
It says undefined âeventsâ. But when I do define the events it says "âEpochsFIFâ object has no attribute âannotationsâ
time_of_interest = data.time_as_index([0, 1]) # from 0 to 1 seconds
data1 = data.get_data() # (n_epochs, n_channels, n_times)
for epoch in data1:
# epoch is a numpy array (n_channels, n_times)
epoch_window = data1[:, time_of_interest[0]:time_of_interest[1]]
all_events, all_event_id = mne.events_from_annotations(data)
event_id = slice(all_event_id, 'Stimulus/S ')
events, _ = mne.events_from_annotations(data, event_id=event_id)
reject_criteria = dict(eeg=100e-6)
flat_criteria = dict(
eeg=-100e-6)
data_annot = mne.Epochs(epoch_window, events, tmin=0, tmax=1, reject_tmax=0,
reject=reject_criteria, flat=flat_criteria,
reject_by_annotation=False, preload=True)
data_annot.plot_drop_log()
Any idea how to handle this problem?
If I read your code correctly, you pass data to events_from_annotations(); however, this function expects a Raw and not an Epochs object.
Since you already managed to create epochs, there must be events somewhere??! You should just re-use them.
Best,
Richard
Yes I defined the events before running the ICA
I already tried to include the code here:
#define functions
def slice(sourcedict, string):
newdict = {}
for key in sourcedict.keys():
if key.startswith(string):
newdict[key] = sourcedict[key]
return newdict
def flux_eegEpochs (raw_load):
reject_criteria = dict(eeg=100e-6)
flat_criteria = dict(
eeg=-100e-6)
all_events, all_event_id = mne.events_from_annotations(raw_load)
print(all_event_id)
print(all_events[:5])
#Reload events chooseing only the ones that start with 'Stimulus/S '
event_id = slice(all_event_id, 'Stimulus/S ')
events, _ = mne.events_from_annotations(raw_load, event_id=event_id)
# Epoch data
epochs = mne.Epochs(raw_load,
events,
event_id=event_id,
tmin=-1,
tmax=8,
detrend=1,
reject=reject_criteria,
flat=flat_criteria,
reject_by_annotation=False,
#baseline=(-0.2,0),
proj=True,
preload=True)
len(epochs.events)
return epochs
epochs.plot_drop_log()
Doesnt work unfortunatly
Since your event names all start with Stimulus/, we can abuse MNEâs event name grouping capability: you can trivially select all Stimulus epochs via
epochs_stim = epochs['Stimulus']
No need to slice etc.
Does that help solve your problem?
Unfortunatly not. It does help to get rid of the âeventâ error but then it says "Flat value must be a number >= 0, not â-0.0001â, if I change it to â0â what I actually dont want because I want a -100 flat and +100 peak, it gives me the following error: âNo appropriate channels found for the given picks (array(, dtype=int32))â
epochs = mne.Epochs(raw_load,
events,
event_id=event_id,
tmin=-1,
tmax=8,
detrend=1,
reject=reject_criteria,
flat=flat_criteria,
reject_by_annotation=False,
#baseline=(-0.2,0),
proj=True,
preload=True)
len(epochs.events)
return epochs
epochs.plot_drop_log()
def identify_ica(raw):
""" Use ICA to reject artifacts
"""
# Perform ICA
ica = mne.preprocessing.ICA(
n_components=20, # Number of components to return
max_pca_components=None, # Don't reduce dimensionality too much
random_state=0,
max_iter=800,
verbose='INFO')
ica.fit(raw, reject_by_annotation=True)
# Plot ICA results
ica.plot_components(inst=raw) # Scalp topographies - Click for more info
ica.plot_sources(raw) # Time-courses - click on the ones to exclude
return ica
raw_load=mne.io.read_raw_fif('C:\\Users\Admin\Desktop\Projekt_Franziska\Studie 3\Analyse\EEG\BIDS_mast\sub-56\ses-Expo\eeg/sub-56_ses-Expo_run01_data_interpolated.fif', verbose=False)
mne_data_folder ='C:\\Users\Admin\Desktop\Projekt_Franziska\Studie 3\Analyse\EEG\BIDS_mast\sub-56\ses-Expo\eeg/'
print(raw_load.annotations)
print(len(raw_load.annotations))
print(set(raw_load.annotations.duration))
print(set(raw_load.annotations.description))
print(raw_load.annotations.onset[0])
dataepoch=flux_eegEpochs (raw_load)
print(dataepoch.info)
filt_data = dataepoch.copy()
filt_data.load_data().filter(l_freq=1., h_freq=None)
raw_ica=identify_ica(filt_data);
Hello @FranziMNE, itâs become difficult for me to figure out which line exactly is causing the error youâre reporting. Could you please try to reduce your code snippet to the relevant lines only, and also share the entire traceback when you post error messages? This should make it easier for others to understand where things are going wrong.
Best wishes,
Richard
Got it, no epochs were left when I set the reject_criteria to the stringent threshold I used, so no epochs were left in the ICA, resulting in that error.
So I think everything can be traced back to the reject_criteria.
For the critera: Why cant I set the flat criteria to -100 microVolt? I actually want to have those criteria: amplitude threshold (-100 to +100 ”V) and a joint probability threshold (SD = 5)
The rejection thresholds are specified as so-called âpeak-to-peak amplitudeâ, PTP.
This is basically the signal change from the lowest to the highest point. Therefore, the value can never be negative.
What youâre asking for is, to my knowledge, not supported. But Iâm not sure if this is actually a problem â I believe the concept of PTP is actually better suited for artifact rejection, as itâs independent of signal baseline and DC offset.
Best wishes,
Richard
Alright! 
The backround ist just, that for my current paper, reviewers asked to redo our EEG analysis since the preprocessing pipeline differs between ERP and TFR analysis in the paper. This is due the fact that my co-author ran the ERP analysis with matlab and I did the TFR analysis in python. I now set up the ERP analysis in python and try to build up the pipeline my co-author used. He used the mentioned artifaction correction criteria. This is why I wanted to redo them in python.
Thank you very much for your help again!
Best, Franzi
I think you mis-understand what flat means here. both reject and flat must be positive, and reflect absolute values of signal amplitudes. If you want to reject epochs that have values above 100 microvolts and also below -100 microvolts, then all you need to specify is reject=100 (actually, 100e-6, since MNE-Python stores data in Volts, not microvolts). flat is used to reject epochs where the signal stays below a given absolute value for the entire epoch (as an indication that something went wrong with the connection / the active electrodeâs amplifier isnât working / whatever).
Objection 
Itâs not signal amplitudes, itâs the peak-to-peak amplitudes! That is, the difference between the highest and the lowest value.
But your suggestion does do the trick for the scenario pointed out here, I believe, yes! So my previous statement:
is actually incorrect.
But I just wanted to clarify the nomenclature here a little.
ah, of course you are right. The main point stands (that flat is not appropriate here) but reject will not be sufficient by itself. I wonder if we should allow reject to take a callable?
Actually ⊠no, this is not the case. Say if you specify reject={'eeg': 100e-6} and within a given epoch, the minimum signal amplitude is 40 ”V and the maximum is 120 ”V, the absolute amplitude exceeds 100 ”V but the PTP does not, so the epoch is retained.
We cannot âfixâ this by specifying an additional flat parameter either, because this one also only looks at PTP.
I believe this would be the best approach that grants us greatest flexibility without forcing us to introduce a new parameter. @FranziMNE could then just use a function like:
def custom_reject(single_epoch_data):
if (
single_epoch_data().min() < -100e-6 or
single_epoch_data().max() > +100e-6
):
return False # reject
else:
return True # keep
Hi Thanks for the idea
# Load data (after ICA)
subject = 'sub-52'
session = 'ses-Expo'
bids_root ='C:\\Users\\Admin\\Desktop\\Projekt_Franziska\Studie 3\\Analyse\\EEG\\BIDS_mast'
mne_data_folder=bids_root + '/' + subject + '/' + session + '/' + 'eeg/'
data=mne.read_epochs(mne_data_folder + subject + '_' + session + '_run01_' + 'data_ica.fif', verbose=False)
def custom_reject(data):
if (
data().min() < -100e-6 or
data().max() > +100e-6
):
return False # reject
else:
return True # keep
data.plot_drop_log()
It sadly doesnt reject anything, even when changing the criteria to very unreasonable low values