Finding index of dropped epoch & automatic rejection for bad epoch

Support & Discussions
,
1

Hi,

Could anyone help me how to find indices of epochs that have been dropped ?

Also, does mne have an automatic rejection for bad epochs ?

Thanks in advanced
IG

2

Sorry for the slow reply. epochs.selection gives you the indices of epochs that have not been dropped. epochs.drop_log records the reasons why epochs have been dropped. So to get the indices of dropped epochs, there are a few options, e.g.,;

[n for n, dl in enumerate(epochs.drop_log) if len(dl)]  # result is a list
np.nonzero(map(len, epochs.drop_log))  # result is a numpy array

For your second question, when epochs are created via mne.Epochs(raw, events, ...) there are a few ways that epochs can get dropped:

  1. if you pass something for the reject parameter. See the docstring for examples.
  2. if your raw object had any annotations that start with bad or edge, and an epoch would overlap that annotation. These are signal values used when concatenating raws (edge) or when annotating artifacts (bad). You can prevent this dropping by passing reject_by_annotation=False
  3. if the events that you’re epoching are too close to the start/end of the raw object, such that they would get truncated. There is no way to prevent this kind of dropping (all epochs must be the same length, so we can’t possibly include an epoch that is shorter because it occurs right at the beginning/end of the recording).
2 Likes
3

Thank you so much for this info. I will try your advice

4

Hey @drammock ,

I tried to use this to identify the indices of the rejected epochs.
np.nonzero(map(len, epochs.drop_log)) # result is a numpy array

This line just gives me an array of [0]. How do I fix it ?

5

What do you want to fix? There can be plenty of reasons why the first epoch was dropped. What makes you believe it needs fixing?

Best wishes,
Richard

6

@richard , It is not showing the correct indices as an array. I need the array of rejected epochs.

So when I print -
[n for n, dl in enumerate(epochs.drop_log) if len(dl)]
Output : [20,
21,
22,
23,
24]

Then,
np.nonzero(map(len, epochs.drop_log))
output : (array([0]),)

7

hmm, when I posted the original answer I thought the list comprehension and the np.nonzero(map(len,...)) approach would yield equivalent results. I’m not sure why you’re seeing different results: maybe I was wrong, or maybe there’s something weird about your data or your code. In any event, if the list comprehension is giving you the correct indices, and you just need them to be an array, then wrap the list comprehension in np.array():

np.array([n for n, dl in enumerate(epochs.drop_log) if len(dl)])
8

I tried that and did not work. It just converted the list as an array with the epoch numbers but not indices for those epoch numbers.
So my output is still the same…
Output : [20,
21,
22,
23,
24]

I wonder if there is a way to get the array by using something like [epochs.drop_log].get_data()(…this does not work) , similar to what “autoreject” has epochs[reject_log.bad_epochs].get_data (…this works).

9

Wait, what is the difference here? I suppose I still haven’t fully understood your question, sorry :frowning:

10

OK, starting with the sample data:

import os
import numpy as np
import mne

sample_data_folder = mne.datasets.sample.data_path()
sample_data_raw_file = os.path.join(sample_data_folder, 'MEG', 'sample',
                                    'sample_audvis_raw.fif')
raw = mne.io.read_raw_fif(sample_data_raw_file, verbose=False, preload=False)
events = mne.find_events(raw, stim_channel='STI 014')
event_dict = {'auditory/left': 1, 'auditory/right': 2, 'visual/left': 3}
epochs = mne.Epochs(raw, events, event_id=event_dict, preload=False,
                    proj=False)

then this works:

np.array([n for n, dl in enumerate(epochs.drop_log) if len(dl)])

and so does this:

np.nonzero(list(map(len, epochs.drop_log)))[0]

(the latter one needed an extra list call in between the map and the np.nonzero, and a [0] index at the end because np.nonzero returns a tuple).

Both of those return:

array([  3,   7,  11,  15,  17,  20,  24,  28,  32,  36,  37,  41,  45,
        49,  53,  57,  59,  60,  62,  66,  70,  74,  78,  79,  80,  83,
        87,  91,  95,  99, 100, 104, 108, 112, 116, 119, 120, 121, 125,
       129, 133, 137, 141, 143, 144, 146, 150, 154, 158, 162, 163, 164,
       167, 171, 175, 179, 183, 185, 188, 192, 196, 200, 203, 204, 205,
       209, 213, 217, 221, 225, 226, 227, 230, 234, 238, 242, 246, 247,
       251, 255, 259, 263, 266, 267, 268, 272, 276, 280, 284, 287, 288,
       289, 293, 297, 301, 305, 309, 311, 312, 314, 318, 319])

Hopefully that is what you want? I’m still a bit confused by this statement:

What exactly is the difference in your mind between “epoch number” and “index”? For me, in most (perhaps all?) cases, those terms are interchangable.

1 Like
11

Sorry for not being clear, what I was looking for are the data points for each indexed epoch.

To give more context I am trying identify the phase-amplitude coupling. This is the workflow:

  1. X= Load the continuous resting state eeg raw data and then filter for beta and gamma and perform Hilbert transform (to get power in signal) , so now I have “power” indices. I epoched this data at 1s.

  2. Y= Load the same resting eeg data and epoch it at 1s, then annotate for “bad/noisy” epochs. (so no filtering here)

  3. Extract rejected indices of bad/rejected epochs from Y and then subtract those indices from X. But I was wondering how to extract power indices after rejecting epochs and get_data for remaining “good epochs” in X? Like I figured how to index the rejected epochs but I could not get the data points (eg. epochs.get_data())

  4. calculate phase-amp coupling on “good (epochs)” power indices.

12

You can make use of the selection attribute to do something like

good_epochs_idx = Y.selection
X_good = X[good_epochs_idx]
2 Likes