From 045c65e4cf6ec7544869f3508264517cc3344809 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Tue, 7 Nov 2023 20:49:04 +0100 Subject: [PATCH 01/17] Add support for MNE-ICALabel --- mne_bids_pipeline/_config.py | 6 ++ .../steps/preprocessing/_06a_run_ica.py | 56 ++++++++++++++++--- 2 files changed, 54 insertions(+), 8 deletions(-) mode change 100644 => 100755 mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py diff --git a/mne_bids_pipeline/_config.py b/mne_bids_pipeline/_config.py index 32e0c9735..d5507b756 100644 --- a/mne_bids_pipeline/_config.py +++ b/mne_bids_pipeline/_config.py @@ -1367,6 +1367,12 @@ false-alarm rate increases dramatically. """ +ica_use_icalabel = False +""" +Whether to use MNE-ICALabel to automatically label ICA components. Only available for +EEG data. +""" + # Rejection based on peak-to-peak amplitude # ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py old mode 100644 new mode 100755 index efd8bec84..b827774dd --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -17,6 +17,7 @@ import pandas as pd import numpy as np import autoreject +from mne_icalabel import label_components import mne from mne.report import Report @@ -184,7 +185,7 @@ def make_eog_epochs( return eog_epochs -def detect_bad_components( +def detect_bad_components_mne( *, cfg, which: Literal["eog", "ecg"], @@ -195,7 +196,7 @@ def detect_bad_components( session: Optional[str], ) -> Tuple[List[int], np.ndarray]: artifact = which.upper() - msg = f"Performing automated {artifact} artifact detection …" + msg = f"Performing automated {artifact} artifact detection (MNE) …" logger.info(**gen_log_kwargs(message=msg)) if which == "eog": @@ -395,7 +396,18 @@ def run_ica( # Set an EEG reference if "eeg" in cfg.ch_types: - projection = True if cfg.eeg_reference == "average" else False + if cfg.ica_use_icalabel: + assert cfg.eeg_reference == "average" + projection = False # Avg. ref. needs to be applied for MNE-ICALabel + elif cfg.eeg_reference == "average": + projection = True + else: + projection = False + + if not projection: + msg = "Applying average reference to EEG epochs used for ICA fitting." + logger.info(**gen_log_kwargs(message=msg)) + epochs.set_eeg_reference(cfg.eeg_reference, projection=projection) if cfg.ica_reject == "autoreject_local": @@ -448,7 +460,7 @@ def run_ica( # ECG and EOG component detection if epochs_ecg: - ecg_ics, ecg_scores = detect_bad_components( + ecg_ics, ecg_scores = detect_bad_components_mne( cfg=cfg, which="ecg", epochs=epochs_ecg, @@ -461,7 +473,7 @@ def run_ica( ecg_ics = ecg_scores = [] if epochs_eog: - eog_ics, eog_scores = detect_bad_components( + eog_ics, eog_scores = detect_bad_components_mne( cfg=cfg, which="eog", epochs=epochs_eog, @@ -473,11 +485,30 @@ def run_ica( else: eog_ics = eog_scores = [] + # Run MNE-ICALabel if requested. + if cfg.ica_use_icalabel: + msg = "Performing automated artifact detection (MNE-ICALabel) …" + logger.info(**gen_log_kwargs(message=msg)) + + label_results = label_components(inst=epochs, ica=ica, method="iclabel") + icalabel_ics = [] + icalabel_labels = [] + for idx, label in enumerate(label_results["labels"]): + if label not in ["brain", "other"]: + icalabel_ics.append(idx) + icalabel_labels.append(label) + + msg = f"Detected {len(icalabel_ics)} artifact-related ICs in {len(epochs)} epochs." + logger.info(**gen_log_kwargs(message=msg)) + # Save ICA to disk. # We also store the automatically identified ECG- and EOG-related ICs. msg = "Saving ICA solution and detected artifacts to disk." logger.info(**gen_log_kwargs(message=msg)) - ica.exclude = sorted(set(ecg_ics + eog_ics)) + if cfg.ica_use_icalabel: + ica.exclude = sorted(set(ecg_ics + eog_ics + icalabel_ics)) + else: + ica.exclude = sorted(set(ecg_ics + eog_ics)) ica.save(out_files["ica"], overwrite=True) _update_for_splits(out_files, "ica") @@ -495,12 +526,20 @@ def run_ica( for component in ecg_ics: row_idx = tsv_data["component"] == component tsv_data.loc[row_idx, "status"] = "bad" - tsv_data.loc[row_idx, "status_description"] = "Auto-detected ECG artifact" + tsv_data.loc[row_idx, "status_description"] = "Auto-detected ECG artifact (MNE)" for component in eog_ics: row_idx = tsv_data["component"] == component tsv_data.loc[row_idx, "status"] = "bad" - tsv_data.loc[row_idx, "status_description"] = "Auto-detected EOG artifact" + tsv_data.loc[row_idx, "status_description"] = "Auto-detected EOG artifact (MNE)" + + if cfg.ica_use_icalabel: + for component, label in zip(icalabel_ics, icalabel_labels): + row_idx = tsv_data["component"] == component + tsv_data.loc[row_idx, "status"] = "bad" + tsv_data.loc[ + row_idx, "status_description" + ] = f"Auto-detected {label} (MNE-ICALabel)" tsv_data.to_csv(out_files["components"], sep="\t", index=False) @@ -588,6 +627,7 @@ def get_config( ica_reject=config.ica_reject, ica_eog_threshold=config.ica_eog_threshold, ica_ctps_ecg_threshold=config.ica_ctps_ecg_threshold, + ica_use_icalabel=config.ica_use_icalabel, autoreject_n_interpolate=config.autoreject_n_interpolate, random_state=config.random_state, ch_types=config.ch_types, From 307b9c08aa82b702bc7beb171b4c285e337d3912 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Tue, 7 Nov 2023 20:59:55 +0100 Subject: [PATCH 02/17] Refactor --- .../steps/preprocessing/_06a_run_ica.py | 95 ++++++++++--------- 1 file changed, 52 insertions(+), 43 deletions(-) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index b827774dd..56af5fbc9 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -458,33 +458,6 @@ def run_ica( if cfg.task is not None: title += f", task-{cfg.task}" - # ECG and EOG component detection - if epochs_ecg: - ecg_ics, ecg_scores = detect_bad_components_mne( - cfg=cfg, - which="ecg", - epochs=epochs_ecg, - ica=ica, - ch_names=None, # we currently don't allow for custom channels - subject=subject, - session=session, - ) - else: - ecg_ics = ecg_scores = [] - - if epochs_eog: - eog_ics, eog_scores = detect_bad_components_mne( - cfg=cfg, - which="eog", - epochs=epochs_eog, - ica=ica, - ch_names=cfg.eog_channels, - subject=subject, - session=session, - ) - else: - eog_ics = eog_scores = [] - # Run MNE-ICALabel if requested. if cfg.ica_use_icalabel: msg = "Performing automated artifact detection (MNE-ICALabel) …" @@ -500,15 +473,41 @@ def run_ica( msg = f"Detected {len(icalabel_ics)} artifact-related ICs in {len(epochs)} epochs." logger.info(**gen_log_kwargs(message=msg)) + ica.exclude = sorted(icalabel_ics) + else: + # Run MNE's built-in ECG and EOG component detection + if epochs_ecg: + ecg_ics, ecg_scores = detect_bad_components_mne( + cfg=cfg, + which="ecg", + epochs=epochs_ecg, + ica=ica, + ch_names=None, # we currently don't allow for custom channels + subject=subject, + session=session, + ) + else: + ecg_ics = ecg_scores = [] + + if epochs_eog: + eog_ics, eog_scores = detect_bad_components_mne( + cfg=cfg, + which="eog", + epochs=epochs_eog, + ica=ica, + ch_names=cfg.eog_channels, + subject=subject, + session=session, + ) + else: + eog_ics = eog_scores = [] + + ica.exclude = sorted(set(ecg_ics + eog_ics)) # Save ICA to disk. # We also store the automatically identified ECG- and EOG-related ICs. msg = "Saving ICA solution and detected artifacts to disk." logger.info(**gen_log_kwargs(message=msg)) - if cfg.ica_use_icalabel: - ica.exclude = sorted(set(ecg_ics + eog_ics + icalabel_ics)) - else: - ica.exclude = sorted(set(ecg_ics + eog_ics)) ica.save(out_files["ica"], overwrite=True) _update_for_splits(out_files, "ica") @@ -523,16 +522,6 @@ def run_ica( ) ) - for component in ecg_ics: - row_idx = tsv_data["component"] == component - tsv_data.loc[row_idx, "status"] = "bad" - tsv_data.loc[row_idx, "status_description"] = "Auto-detected ECG artifact (MNE)" - - for component in eog_ics: - row_idx = tsv_data["component"] == component - tsv_data.loc[row_idx, "status"] = "bad" - tsv_data.loc[row_idx, "status_description"] = "Auto-detected EOG artifact (MNE)" - if cfg.ica_use_icalabel: for component, label in zip(icalabel_ics, icalabel_labels): row_idx = tsv_data["component"] == component @@ -540,6 +529,20 @@ def run_ica( tsv_data.loc[ row_idx, "status_description" ] = f"Auto-detected {label} (MNE-ICALabel)" + else: + for component in ecg_ics: + row_idx = tsv_data["component"] == component + tsv_data.loc[row_idx, "status"] = "bad" + tsv_data.loc[ + row_idx, "status_description" + ] = "Auto-detected ECG artifact (MNE)" + + for component in eog_ics: + row_idx = tsv_data["component"] == component + tsv_data.loc[row_idx, "status"] = "bad" + tsv_data.loc[ + row_idx, "status_description" + ] = "Auto-detected EOG artifact (MNE)" tsv_data.to_csv(out_files["components"], sep="\t", index=False) @@ -549,10 +552,16 @@ def run_ica( logger.info(**gen_log_kwargs(message=msg)) report = Report(info_fname=epochs, title=title, verbose=False) + ecg_evoked = None if epochs_ecg is None else epochs_ecg.average() eog_evoked = None if epochs_eog is None else epochs_eog.average() - ecg_scores = None if len(ecg_scores) == 0 else ecg_scores - eog_scores = None if len(eog_scores) == 0 else eog_scores + + if cfg.ica_use_icalabel: + # We didn't run MNE's scoring + ecg_scores = eog_scores = None + else: + ecg_scores = None if len(ecg_scores) == 0 else ecg_scores + eog_scores = None if len(eog_scores) == 0 else eog_scores with _agg_backend(): if cfg.ica_reject == "autoreject_local": From fbfa83e902163261030f62701120f7a16ae8f7cd Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Tue, 7 Nov 2023 21:02:22 +0100 Subject: [PATCH 03/17] Add dependency --- pyproject.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/pyproject.toml b/pyproject.toml index 5d217e36d..9f9fb55aa 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -44,6 +44,7 @@ dependencies = [ "autoreject", "mne[hdf5] >=1.2", "mne-bids[full]", + "mna-icalabel", "filelock", "setuptools >=65", ] From a04dac10a975bc941615539c74b6c52248cbeeb2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Tue, 7 Nov 2023 21:08:47 +0100 Subject: [PATCH 04/17] Sanity check --- mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 56af5fbc9..66c906380 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -79,6 +79,10 @@ def fit_ica( algorithm = "infomax" fit_params = dict(extended=True) + if cfg.ica_use_icalabel: + # The ICALabel network was trained on extended-Infomax ICA decompositions + assert algorithm in ["picard-extended_infomax", "extended_infomax"] + ica = ICA( method=algorithm, random_state=cfg.random_state, From 86ba18d6e671d8b7c074af5e3b1eeb61bc04ec54 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Tue, 7 Nov 2023 22:30:20 +0100 Subject: [PATCH 05/17] Update pyproject.toml --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 9f9fb55aa..9993a034e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -44,7 +44,7 @@ dependencies = [ "autoreject", "mne[hdf5] >=1.2", "mne-bids[full]", - "mna-icalabel", + "mne-icalabel", "filelock", "setuptools >=65", ] From 2e629c00bfcc2bed2a3c4f5bcf79c7eff0a67808 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 11:26:17 +0100 Subject: [PATCH 06/17] Add validation and fix docs --- docs/source/settings/preprocessing/ssp_ica.md | 1 + mne_bids_pipeline/_config.py | 8 ++++++++ mne_bids_pipeline/_config_import.py | 14 ++++++++++++++ .../steps/preprocessing/_06a_run_ica.py | 12 ++++++++---- 4 files changed, 31 insertions(+), 4 deletions(-) diff --git a/docs/source/settings/preprocessing/ssp_ica.md b/docs/source/settings/preprocessing/ssp_ica.md index b132ef4bf..23f510796 100644 --- a/docs/source/settings/preprocessing/ssp_ica.md +++ b/docs/source/settings/preprocessing/ssp_ica.md @@ -24,6 +24,7 @@ tags: - ssp_ecg_channel - ica_reject - ica_algorithm + - ica_use_icalabel - ica_l_freq - ica_max_iterations - ica_n_components diff --git a/mne_bids_pipeline/_config.py b/mne_bids_pipeline/_config.py index d5507b756..7582b7d54 100644 --- a/mne_bids_pipeline/_config.py +++ b/mne_bids_pipeline/_config.py @@ -1371,6 +1371,14 @@ """ Whether to use MNE-ICALabel to automatically label ICA components. Only available for EEG data. + +!!! info + Using MNE-ICALabel mandates that you also set: + ```python + eeg_reference = "average" + ica_l_freq = 1 + h_freq = 100 + ``` """ # Rejection based on peak-to-peak amplitude diff --git a/mne_bids_pipeline/_config_import.py b/mne_bids_pipeline/_config_import.py index 14a55df2e..81a492d41 100644 --- a/mne_bids_pipeline/_config_import.py +++ b/mne_bids_pipeline/_config_import.py @@ -345,6 +345,20 @@ def _check_config(config: SimpleNamespace, config_path: Optional[PathLike]) -> N f"but got shape {destination.shape}" ) + # MNE-ICALabel + if config.ica_use_icalabel: + if config.ica_l_freq != 1.0 or config.h_freq != 100.0: + raise ValueError( + f"When using MNE-ICALabel, you must set ica_l_freq=1 and h_freq=100, " + f"but got: ica_l_freq={config.ica_l_freq} and h_freq={config.h_freq}" + ) + + if config.eeg_reference != "average": + raise ValueError( + f'When using MNE-ICALabel, you must set eeg_reference="average", but ' + f"got: eeg_reference={config.eeg_reference}" + ) + def _default_factory(key, val): # convert a default to a default factory if needed, having an explicit diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 66c906380..3c388ab77 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -79,10 +79,6 @@ def fit_ica( algorithm = "infomax" fit_params = dict(extended=True) - if cfg.ica_use_icalabel: - # The ICALabel network was trained on extended-Infomax ICA decompositions - assert algorithm in ["picard-extended_infomax", "extended_infomax"] - ica = ICA( method=algorithm, random_state=cfg.random_state, @@ -276,6 +272,14 @@ def run_ica( in_files: dict, ) -> dict: """Run ICA.""" + if cfg.ica_use_icalabel: + # The ICALabel network was trained on extended-Infomax ICA decompositions fit + # on data flltered between 1 and 100 Hz. + assert cfg.algorithm in ["picard-extended_infomax", "extended_infomax"] + assert cfg.ica_l_freq == 1.0 + assert cfg.h_freq == 100.0 + assert cfg.eeg_reference == "average" + raw_fnames = [in_files.pop(f"raw_run-{run}") for run in cfg.runs] bids_basename = raw_fnames[0].copy().update(processing=None, split=None, run=None) out_files = dict() From ef07ed4d2aa5b2493602d88705b0b452331c17a9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 11:30:27 +0100 Subject: [PATCH 07/17] Update changelog --- docs/source/v1.5.md.inc | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/docs/source/v1.5.md.inc b/docs/source/v1.5.md.inc index 5522271b1..40a86d747 100644 --- a/docs/source/v1.5.md.inc +++ b/docs/source/v1.5.md.inc @@ -3,6 +3,10 @@ This release contains a number of very important bug fixes that address problems related to decoding, time-frequency analysis, and inverse modeling. All users are encouraged to update. +We also improved logging during parallel processing, added support for finding and repairing bad epochs via +[`autoreject`](https://autoreject.github.io), and included support for automatic labeling of ICA artifacts +via [MNE-ICALabel][https://mne.tools/mne-icalabel]. + ### :new: New features & enhancements - Added `deriv_root` argument to CLI (#773 by @vferat) @@ -22,6 +26,7 @@ All users are encouraged to update. - Added support for "local" [`autoreject`](https://autoreject.github.io) to find (and repair) bad channels on a per-epoch basis before submitting them to ICA fitting. This can be enabled by setting [`ica_reject`][mne_bids_pipeline._config.ica_reject] to `"autoreject_local"`. (#810 by @hoechenberger) +- Added support for automated labeling of ICA components via [MNE-ICALabel][https://mne.tools/mne-icalabel] (#812 by @hoechenberger) - Website documentation tables can now be sorted (e.g., to find examples that use a specific feature) (#808 by @larsoner) [//]: # (### :warning: Behavior changes) From 85dcc8a552ef5d8cac10f9f66e45c430b0d58ec4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 11:39:55 +0100 Subject: [PATCH 08/17] Fix check and add test case --- mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py | 2 +- mne_bids_pipeline/tests/configs/config_ERP_CORE.py | 7 ++++++- 2 files changed, 7 insertions(+), 2 deletions(-) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 3c388ab77..98390d3e6 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -275,7 +275,7 @@ def run_ica( if cfg.ica_use_icalabel: # The ICALabel network was trained on extended-Infomax ICA decompositions fit # on data flltered between 1 and 100 Hz. - assert cfg.algorithm in ["picard-extended_infomax", "extended_infomax"] + assert cfg.ica_algorithm in ["picard-extended_infomax", "extended_infomax"] assert cfg.ica_l_freq == 1.0 assert cfg.h_freq == 100.0 assert cfg.eeg_reference == "average" diff --git a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py index 47fcb5846..ad625fae4 100644 --- a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py +++ b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py @@ -71,12 +71,17 @@ t_break_annot_start_after_previous_event = 3.0 t_break_annot_stop_before_next_event = 1.5 +# Settings for autoreject and ICA if task == "N400": # test autoreject local without ICA spatial_filter = None reject = "autoreject_local" autoreject_n_interpolate = [2, 4] -elif task == "N170": # test autoreject local before ICA +elif task == "N170": # test autoreject local before ICA, and MNE-ICALabel spatial_filter = "ica" + ica_algorithm = "picard-extended_infomax" + ica_use_icalabel = True + ica_l_freq = 1 + h_freq = 100 ica_reject = "autoreject_local" reject = "autoreject_global" autoreject_n_interpolate = [2, 4] From e86f179bb9093beb9801a9540045a1444fa0b871 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 12:21:30 +0100 Subject: [PATCH 09/17] Add cfg.h_freq --- mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py | 1 + 1 file changed, 1 insertion(+) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 98390d3e6..89492da04 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -649,6 +649,7 @@ def get_config( random_state=config.random_state, ch_types=config.ch_types, l_freq=config.l_freq, + h_freq=config.h_freq, epochs_decim=config.epochs_decim, raw_resample_sfreq=config.raw_resample_sfreq, event_repeated=config.event_repeated, From cd405c26f874915fd9763ade98b4dbbb929a75e7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 13:11:27 +0100 Subject: [PATCH 10/17] Add onnxruntime dependency --- pyproject.toml | 1 + 1 file changed, 1 insertion(+) diff --git a/pyproject.toml b/pyproject.toml index 9993a034e..344f0dc8b 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -45,6 +45,7 @@ dependencies = [ "mne[hdf5] >=1.2", "mne-bids[full]", "mne-icalabel", + "onnxruntime", # for mne-icalabel "filelock", "setuptools >=65", ] From 4a85e5130446b8f2ccb5b4d3a89d16896e07fa1c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 13:30:34 +0100 Subject: [PATCH 11/17] Better logging --- mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 89492da04..0622cee4f 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -225,7 +225,7 @@ def detect_bad_components_mne( logger.warning(**gen_log_kwargs(message=warn)) else: msg = ( - f"Detected {len(inds)} {artifact}-related ICs in " + f"Detected {len(inds)} {artifact}-related independent component(s) in " f"{len(epochs)} {artifact} epochs." ) logger.info(**gen_log_kwargs(message=msg)) @@ -479,7 +479,10 @@ def run_ica( icalabel_ics.append(idx) icalabel_labels.append(label) - msg = f"Detected {len(icalabel_ics)} artifact-related ICs in {len(epochs)} epochs." + msg = ( + f"Detected {len(icalabel_ics)} artifact-related independent component(s) " + f"in {len(epochs)} epochs." + ) logger.info(**gen_log_kwargs(message=msg)) ica.exclude = sorted(icalabel_ics) else: From 2d63230164ccccf7de9f0f46cbeb4f1cbdc90328 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 13:49:16 +0100 Subject: [PATCH 12/17] Adjust test config --- mne_bids_pipeline/tests/configs/config_ERP_CORE.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py index ad625fae4..f9c94699a 100644 --- a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py +++ b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py @@ -84,7 +84,7 @@ h_freq = 100 ica_reject = "autoreject_local" reject = "autoreject_global" - autoreject_n_interpolate = [2, 4] + autoreject_n_interpolate = [5, 10] # Only for testing! else: spatial_filter = "ica" ica_reject = dict(eeg=350e-6, eog=500e-6) @@ -254,6 +254,7 @@ baseline = (None, 0) conditions = ["stimulus/face/normal", "stimulus/car/normal"] contrasts = [("stimulus/face/normal", "stimulus/car/normal")] + cluster_forming_t_threshold = 1.5 # Only for testing! elif task == "P3": rename_events = { "response/201": "response/correct", From acc2acdc6f6fcce7c30cc1c69adc747ac28eba4c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 13:49:20 +0100 Subject: [PATCH 13/17] Better doc --- mne_bids_pipeline/_config.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/mne_bids_pipeline/_config.py b/mne_bids_pipeline/_config.py index 7582b7d54..e72f5e222 100644 --- a/mne_bids_pipeline/_config.py +++ b/mne_bids_pipeline/_config.py @@ -1457,7 +1457,8 @@ !!! info This setting only takes effect if [`reject`][mne_bids_pipeline._config.reject] has - been set to `"autoreject_local"`. + been set to `"autoreject_local"`. It is not applied when using + `"autoreject_global"`. !!! info Channels marked as globally bad in the BIDS dataset (in `*_channels.tsv)`) will not From 3032dd18a33d8a9a82295e5a3adc19c96aed2ac0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 14:17:07 +0100 Subject: [PATCH 14/17] Change test config --- mne_bids_pipeline/tests/configs/config_ERP_CORE.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py index f9c94699a..d645adda5 100644 --- a/mne_bids_pipeline/tests/configs/config_ERP_CORE.py +++ b/mne_bids_pipeline/tests/configs/config_ERP_CORE.py @@ -84,7 +84,7 @@ h_freq = 100 ica_reject = "autoreject_local" reject = "autoreject_global" - autoreject_n_interpolate = [5, 10] # Only for testing! + autoreject_n_interpolate = [12] # Only for testing! else: spatial_filter = "ica" ica_reject = dict(eeg=350e-6, eog=500e-6) @@ -254,7 +254,7 @@ baseline = (None, 0) conditions = ["stimulus/face/normal", "stimulus/car/normal"] contrasts = [("stimulus/face/normal", "stimulus/car/normal")] - cluster_forming_t_threshold = 1.5 # Only for testing! + cluster_forming_t_threshold = 1.25 # Only for testing! elif task == "P3": rename_events = { "response/201": "response/correct", From 03a3da9b70857ceb947e43671426cee3d108dafa Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 14:51:19 +0100 Subject: [PATCH 15/17] Try running both MNE and ICALabel artifact detection --- .../steps/preprocessing/_06a_run_ica.py | 67 ++++++++++--------- 1 file changed, 35 insertions(+), 32 deletions(-) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index 0622cee4f..b748a016d 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -136,7 +136,7 @@ def make_ecg_epochs( del raw # Free memory if len(ecg_epochs) == 0: - msg = "No ECG events could be found. Not running ECG artifact " "detection." + msg = "No ECG events could be found. Not running ECG artifact detection." logger.info(**gen_log_kwargs(message=msg)) ecg_epochs = None else: @@ -174,7 +174,7 @@ def make_eog_epochs( eog_epochs = create_eog_epochs(raw, ch_name=ch_names, baseline=(None, -0.2)) if len(eog_epochs) == 0: - msg = "No EOG events could be found. Not running EOG artifact " "detection." + msg = "No EOG events could be found. Not running EOG artifact detection." logger.warning(**gen_log_kwargs(message=msg)) eog_epochs = None else: @@ -466,14 +466,42 @@ def run_ica( if cfg.task is not None: title += f", task-{cfg.task}" + # Run MNE's built-in ECG and EOG component detection + if epochs_ecg: + ecg_ics, ecg_scores = detect_bad_components_mne( + cfg=cfg, + which="ecg", + epochs=epochs_ecg, + ica=ica, + ch_names=None, # we currently don't allow for custom channels + subject=subject, + session=session, + ) + else: + ecg_ics = ecg_scores = [] + + if epochs_eog: + eog_ics, eog_scores = detect_bad_components_mne( + cfg=cfg, + which="eog", + epochs=epochs_eog, + ica=ica, + ch_names=cfg.eog_channels, + subject=subject, + session=session, + ) + else: + eog_ics = eog_scores = [] + # Run MNE-ICALabel if requested. if cfg.ica_use_icalabel: + icalabel_ics = [] + icalabel_labels = [] + msg = "Performing automated artifact detection (MNE-ICALabel) …" logger.info(**gen_log_kwargs(message=msg)) label_results = label_components(inst=epochs, ica=ica, method="iclabel") - icalabel_ics = [] - icalabel_labels = [] for idx, label in enumerate(label_results["labels"]): if label not in ["brain", "other"]: icalabel_ics.append(idx) @@ -484,36 +512,10 @@ def run_ica( f"in {len(epochs)} epochs." ) logger.info(**gen_log_kwargs(message=msg)) - ica.exclude = sorted(icalabel_ics) else: - # Run MNE's built-in ECG and EOG component detection - if epochs_ecg: - ecg_ics, ecg_scores = detect_bad_components_mne( - cfg=cfg, - which="ecg", - epochs=epochs_ecg, - ica=ica, - ch_names=None, # we currently don't allow for custom channels - subject=subject, - session=session, - ) - else: - ecg_ics = ecg_scores = [] - - if epochs_eog: - eog_ics, eog_scores = detect_bad_components_mne( - cfg=cfg, - which="eog", - epochs=epochs_eog, - ica=ica, - ch_names=cfg.eog_channels, - subject=subject, - session=session, - ) - else: - eog_ics = eog_scores = [] + icalabel_ics = [] - ica.exclude = sorted(set(ecg_ics + eog_ics)) + ica.exclude = sorted(set(ecg_ics + eog_ics + icalabel_ics)) # Save ICA to disk. # We also store the automatically identified ECG- and EOG-related ICs. @@ -534,6 +536,7 @@ def run_ica( ) if cfg.ica_use_icalabel: + assert len(icalabel_ics) == (icalabel_labels) for component, label in zip(icalabel_ics, icalabel_labels): row_idx = tsv_data["component"] == component tsv_data.loc[row_idx, "status"] = "bad" From 421ac0f86c274c6bc1b7e422dc2808e9ae3fd567 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Richard=20H=C3=B6chenberger?= Date: Wed, 8 Nov 2023 15:16:11 +0100 Subject: [PATCH 16/17] Typo --- mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py index b748a016d..dd6355c05 100755 --- a/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py +++ b/mne_bids_pipeline/steps/preprocessing/_06a_run_ica.py @@ -536,7 +536,7 @@ def run_ica( ) if cfg.ica_use_icalabel: - assert len(icalabel_ics) == (icalabel_labels) + assert len(icalabel_ics) == len(icalabel_labels) for component, label in zip(icalabel_ics, icalabel_labels): row_idx = tsv_data["component"] == component tsv_data.loc[row_idx, "status"] = "bad" From e5d426ba7b483a7758d60dca22182767e276c6b1 Mon Sep 17 00:00:00 2001 From: Eric Larson Date: Wed, 8 Nov 2023 09:34:35 -0500 Subject: [PATCH 17/17] FIX: Clearer test --- mne_bids_pipeline/_config.py | 8 ++++---- mne_bids_pipeline/tests/test_documented.py | 18 +++++++++++------- 2 files changed, 15 insertions(+), 11 deletions(-) diff --git a/mne_bids_pipeline/_config.py b/mne_bids_pipeline/_config.py index e72f5e222..9d1782585 100644 --- a/mne_bids_pipeline/_config.py +++ b/mne_bids_pipeline/_config.py @@ -1237,7 +1237,7 @@ """ Peak-to-peak amplitude limits to exclude epochs from ICA fitting. This allows you to remove strong transient artifacts from the epochs used for fitting ICA, which could -negatively affect ICA performance. +negatively affect ICA performance. The parameter values are the same as for [`reject`][mne_bids_pipeline._config.reject], but `"autoreject_global"` is not supported. @@ -1262,7 +1262,7 @@ to **not** specify rejection thresholds for EOG and ECG channels here – otherwise, ICA won't be able to "see" these artifacts. -???+ info +???+ info This setting is applied only to the epochs that are used for **fitting** ICA. The goal is to make it easier for ICA to produce a good decomposition. After fitting, ICA is applied to the epochs to be analyzed, usually with one or more components @@ -1367,7 +1367,7 @@ false-alarm rate increases dramatically. """ -ica_use_icalabel = False +ica_use_icalabel: bool = False """ Whether to use MNE-ICALabel to automatically label ICA components. Only available for EEG data. @@ -1398,7 +1398,7 @@ If `None` (default), do not apply artifact rejection. -If a dictionary, manually specify rejection thresholds (see examples). +If a dictionary, manually specify rejection thresholds (see examples). The thresholds provided here must be at least as stringent as those in [`ica_reject`][mne_bids_pipeline._config.ica_reject] if using ICA. In case of `'autoreject_global'`, thresholds for any channel that do not meet this diff --git a/mne_bids_pipeline/tests/test_documented.py b/mne_bids_pipeline/tests/test_documented.py index 097fc1032..7bb622dc4 100644 --- a/mne_bids_pipeline/tests/test_documented.py +++ b/mne_bids_pipeline/tests/test_documented.py @@ -18,15 +18,19 @@ def test_options_documented(): with open(root_path / "_config.py", "r") as fid: contents = fid.read() contents = ast.parse(contents) - in_config = [ + unannotated = [ + item.targets[0].id for item in contents.body if isinstance(item, ast.Assign) + ] + assert unannotated == [] + _config_py = [ item.target.id for item in contents.body if isinstance(item, ast.AnnAssign) ] - assert len(set(in_config)) == len(in_config) - in_config = set(in_config) + assert len(set(_config_py)) == len(_config_py) + _config_py = set(_config_py) # ensure we clean our namespace correctly config = _get_default_config() - config_names = set(d for d in dir(config) if not d.startswith("_")) - assert in_config == config_names + _get_default_config_names = set(d for d in dir(config) if not d.startswith("_")) + assert _config_py == _get_default_config_names settings_path = root_path.parent / "docs" / "source" / "settings" assert settings_path.is_dir() in_doc = set() @@ -51,8 +55,8 @@ def test_options_documented(): assert val not in in_doc, "Duplicate documentation" in_doc.add(val) what = "docs/source/settings doc" - assert in_doc.difference(in_config) == set(), f"Extra values in {what}" - assert in_config.difference(in_doc) == set(), f"Values missing from {what}" + assert in_doc.difference(_config_py) == set(), f"Extra values in {what}" + assert _config_py.difference(in_doc) == set(), f"Values missing from {what}" def test_datasets_in_doc():