Interpolating channels marked as 'bad' in the MNE-BIDS-Pipeline

Support & Discussions
,
1

Dear community,

When preparing my EEG data for analysis with the MNE-BIDS-Pipeline, I mark very noisy or flat channels as ‘bad’, following the recommendations I found on the documentation. However, the pipeline appears to exclude these electrodes completely, whereas previously (when I didn’t mark the electrodes as bad) it would interpolate them.

I want to mark the electrodes as ‘bad’ as I feel that this is good practice for the BIDS data structure. However, I want the electrodes to be interpolated in the pipeline. Does anyone know how I can achieve this?

I include my config.py below.

# -*- coding: utf-8 -*-

##############################################################
# Set these values appropriately before running the pipeline 
subjects = 'all'
task = 'vep'                                               
##############################################################

# Sets the appropriate output directories
# bids_root = fr"C:\helios_bids_{task}"
# deriv_root = fr"C:\helios_bids_{task}\derivatives\mne-bids-pipeline-{task}"
bids_root = fr"Z:\Part B\hbd_{task}"
deriv_root = fr"Z:\Part B\hbd_{task}\derivatives\auto-eog-{task}"

subjects_dir = None

# Use to exclude subjects
exclude_subjects = []


ch_types = ["eeg"]
data_type = "eeg"
eeg_reference = "average"  # EEG reference to use
eeg_template_montage = (
    "biosemi64"  # Apply 64-channel Biosemi 10/20 template montage:
)
eeg_bipolar_channels = {'HEOG': ('HEOG-left', 'HEOG-right'),
                        'VEOG': ('VEOG-lower', 'VEOG-upper')}
drop_channels = ["HEOG-left", "HEOG-right", "VEOG-upper", "VEOG-lower"]
eog_channels = ["HEOG", "VEOG"]
analyze_channels = "ch_types"
plot_psd_for_runs = "all"  # For which runs to add a power spectral density (PSD) plot to the generated report.
random_state = (
    42  # Passed to ICA and decoding algos to ensure reproduicibility
)
# Break detection
find_breaks = (
    True  # Automatically find break periods, and annotate them as BAD_break.
)
min_break_duration = 10.0
t_break_annot_start_after_previous_event = 3.0
t_break_annot_stop_before_next_event = 1.5

# Filtering
l_freq = 0.1  # The low-frequency cut-off in the highpass filtering step.
h_freq = 40.0  # The high-frequency cut-off in the highpass filtering step.
notch_freq = (50)  # Notch filter frequency. More than one frequency can be supplied
epochs_decim = 4  # Decimate epochs to 256 Hz
conditions = [
    "Lum/1",
    "Lum/2",
    "Lum/3",
    "Lum/4",
    "LM/1",
    "LM/2",
    "LM/3",
    "LM/4",
    "S/1",
    "S/2",
    "S/3",
    "S/4"
]

# Set the task specific parameters
if task == 'vep':
    epochs_tmin = -0.2  # The beginning of an epoch, relative to the respective event, in seconds.
    epochs_tmax = 0.8  # The end of an epoch, relative to the respective event, in seconds.
elif task == 'ssvep':
    epochs_tmin = -0.2
    epochs_tmax = 2.5
else:
    raise RuntimeError(f"Task {task} not currently supported")
    
baseline = (-0.1, 0)  # Beginning of epoch until time point zero

# Artifact removal
spatial_filter = "ica"  # Use ica
ica_eog_threshold = 3.  # Default
ica_reject = "autoreject_local"  # Find local (per channel) thresholds and repair epochs before fitting ICA
ica_algorithm = "picard-extended_infomax"
ica_l_freq = 1.0
ica_max_iterations = 1000
ica_n_components = None  # 64 - 1
ica_decim = None
reject = "autoreject_local"  # Before and after ICA recommended

# Sensor level analysis
#contrasts = [("Lum", "S"), ("Lum", "LM"), ("LM", "S")]

if task == 'vep':
    decode = False # No decoding for now
    decoding_time_generalization = False  # ?
    decoding_time_generalization_decim = 1
elif task =='ssvep':
    decode = False
else:
    raise RuntimeError(f"Task {task} not currently supported")

# No source estimation
run_source_estimation = False

# Execution
n_jobs = 4

Thank you!

Joel

OS: Windows 11

MNE-Python: 1.10.1

MNE-BIDS-Pipeline: 1.9.0

2

Hello @lightwind, channels marked as globally bad are not interpolated automatically for individual participant processing. The only time interpolation happens is for calculating the grand average signal, to ensure all participants have the same channels available. This was a deliberate decision when we developed the pipeline as we felt that the rank reduction through interpolation would cause more harm than good. However, I can understand that for example for plotting of a topographical map you’d probably want to have all channels present. Is that what you’re trying to do?

Thanks,

Richard

3

Hello @richard,

Thank you, this makes sense to me.

Before I started marking channels as globally bad in the *_channels.tsv, I seem to recall that bad channels were automatically detected and interpolated by the MNE-BIDS-Pipeline. Is that the case? I never seemed to have the issue of missing channels, in any of the outputs.

When you refer to the ‘grand average signal’, do you mean the *_ave.fif files containing the condition averages for each participant?

At present, I am trying to understand how the pipeline handles channels marked as globally bad. Am I right in thinking that, where channels are missing from raw or epochs objects, this can be resolved with a call to .interpolate_bads() on a case by case basis?

Best wishes,

Joel

4

Hello!

Take the following with a slight grain of salt please, as I’m writing this down from memory:

If a channel is marked as “bad” in the BIDS input data, MNE-BIDS-Pipeline ignores it for all participant-level processing.

If a channel is not marked as “bad” in the BIDS data, and MNE-BIDS-Pipeline is running automated bad channel detection, it might detect the channel to be problematic, will mark it as bad and interpolate it.

By Grand Average Signal I mean averaging across participants. By default, at the end of per-participant processing you might end up with one participant with, say, 60 channels and another one with 62 channels, depending on how many were marked as bad in the BIDS input data. In order to average across participants (grand average), all have to have the same number of channels, so we interpolate. But really just for this one process: building the grand average. The setting is interpolate_bads_grand_average.

May I ask why you want to interpolate on a per-participant level? Is it for viz only? Because for stats etc. it won’t make a difference (or cause trouble in the worst case, as you’re introducing a rank deficiency).

Thank you,

Richard

5

Hi @richard,

Thank you, this makes things a lot more clear.

Yes, the reason we want to interpolate at the participant level is for custom group averaging and visualisation after the pipeline. So I think that for this purpose, when loading the *_ave.fif files, we could just say:

evks = [e.interpolate_bads() if e.info[“bads”] else e
for e in evks]

Does that sound like a reasonable approach?

Cheers,

Joel