From 9b76bc092d34cddd02cc179bb4da4b75f937faa9 Mon Sep 17 00:00:00 2001 From: KALLA GANASEKHAR Date: Tue, 15 Apr 2025 14:36:01 +0530 Subject: [PATCH 1/5] update ica_comparison.py Added visualization plots for different ICA algorithms for clean and noisy data --- examples/preprocessing/ica_comparison.py | 141 ++++++++++++++++++----- 1 file changed, 114 insertions(+), 27 deletions(-) diff --git a/examples/preprocessing/ica_comparison.py b/examples/preprocessing/ica_comparison.py index d4246b80362..cef36fc6d8d 100644 --- a/examples/preprocessing/ica_comparison.py +++ b/examples/preprocessing/ica_comparison.py @@ -1,50 +1,70 @@ """ .. _ex-ica-comp: -=========================================== -Compare the different ICA algorithms in MNE -=========================================== +=========================================================== +Compare the performance of different ICA algorithms in MNE +=========================================================== -Different ICA algorithms are fit to raw MEG data, and the corresponding maps -are displayed. +This example compares various ICA algorithms (FastICA, Picard, Infomax, +Extended Infomax) on the same raw MEG data. For each algorithm: +- The ICA fit time (speed) is shown +- All components (up to 20) are visualized +- The EOG-related component from each method is detected and compared +- Comparison on clean vs noisy data is done + +Note: In typical preprocessing, only one ICA algorithm is used. +This example is for educational purposes. """ -# Authors: Pierre Ablin + +# authors : Ganasekhar Kalla # # License: BSD-3-Clause # Copyright the MNE-Python contributors. # %% -from time import time - import mne -from mne.datasets import sample from mne.preprocessing import ICA +from mne.datasets import sample +import numpy as np +from time import time +from pathlib import Path print(__doc__) # %% + # Read and preprocess the data. Preprocessing consists of: # # - MEG channel selection # - 1-30 Hz band-pass filter -data_path = sample.data_path() -meg_path = data_path / "MEG" / "sample" -raw_fname = meg_path / "sample_audvis_filt-0-40_raw.fif" +# Load sample dataset +data_path = Path(sample.data_path()) +raw_file = data_path / 'MEG' / 'sample' / 'sample_audvis_raw.fif' +raw = mne.io.read_raw_fif(raw_file, preload=True) +raw.pick_types(meg=True, eeg=False, eog=True) +raw.crop(0, 60) -raw = mne.io.read_raw_fif(raw_fname).crop(0, 60).pick("meg").load_data() +# %% -reject = dict(mag=5e-12, grad=4000e-13) -raw.filter(1, 30, fir_design="firwin") +# Copy for clean and noisy +raw_clean = raw.copy() +raw_noisy = raw_clean.copy() +noise = 1e-12 * np.random.randn(*raw_noisy._data.shape) +raw_noisy._data += noise +# Rejection thresholds +reject_clean = dict(mag=5e-12, grad=4000e-13) +reject_noisy = dict(mag=1e-11, grad=8000e-13) # %% -# Define a function that runs ICA on the raw MEG data and plots the components -def run_ica(method, fit_params=None): +# Run ICA +def run_ica(raw_input, method, fit_params=None, reject=None): + print(f"\nRunning ICA with: {method}") ica = ICA( n_components=20, method=method, @@ -53,24 +73,91 @@ def run_ica(method, fit_params=None): random_state=0, ) t0 = time() - ica.fit(raw, reject=reject) + ica.fit(raw_input, reject=reject) fit_time = time() - t0 - title = f"ICA decomposition using {method} (took {fit_time:.1f}s)" + print(f"Fitting ICA took {fit_time:.1f}s.") + + # Updated code with broken long line + title = ( + f"ICA decomposition using {method} on " + f"{'noisy' if raw_input is raw_noisy else 'clean'} data\n" + f"(took {fit_time:.1f}s)" + ) ica.plot_components(title=title) + return ica, fit_time # %% -# FastICA -run_ica("fastica") + + +# Run all ICA methods +def run_all_ica(raw_input, label, reject): + icas = {} + fit_times = {} + eog_components = {} + + + for method, params in [ + ("fastica", None), + ("picard", None), + ("infomax", None), + ("infomax", {"extended": True}), + ]: + name = f"{method}" if not params else f"{method}_extended" + full_label = f"{label}_{name}" + ica, t = run_ica(raw_input, method, params, reject) + icas[full_label] = ica + fit_times[full_label] = t + + eog_inds, _ = ica.find_bads_eog( + raw_input, threshold=3.0, verbose='ERROR' + ) + if eog_inds: + eog_components[full_label] = eog_inds[0] + print(f"{full_label}:Detected EOG comp at index {eog_inds[0]}") + else: + eog_components[full_label] = None + print(f"{full_label}: No EOG component detected") + + return icas, fit_times, eog_components # %% -# Picard -run_ica("picard") + +# Run on both raw versions +icas_clean, times_clean, eog_clean = run_all_ica( + raw_clean, "clean", reject_clean +) +icas_noisy, times_noisy, eog_noisy = run_all_ica( + raw_noisy, "noisy", reject_noisy +) + +# Combine results +icas = {**icas_clean, **icas_noisy} +times = {**times_clean, **times_noisy} +eog_comps = {**eog_clean, **eog_noisy} # %% -# Infomax -run_ica("infomax") + +# Clean EOG components for each algorithm (Column 1) +for method in ["fastica", "picard", "infomax", "infomax_extended"]: + key = f"clean_{method}" + comp = eog_comps.get(key) + if comp is not None: + icas[key].plot_components( + picks=[comp], + title=f"{key} - EOG Component (Clean Data)", + show=True + ) # %% -# Extended Infomax -run_ica("infomax", fit_params=dict(extended=True)) + +# Noisy EOG components for each algorithm (Column 2) +for method in ["fastica", "picard", "infomax", "infomax_extended"]: + key = f"noisy_{method}" + comp = eog_comps.get(key) + if comp is not None: + icas[key].plot_components( + picks=[comp], + title=f"{key} - EOG Component (Noisy Data)", + show=True + ) From dc72b25abe4d2957ad0f0bebe2f6c8bcba0228b2 Mon Sep 17 00:00:00 2001 From: KALLA GANASEKHAR Date: Tue, 15 Apr 2025 15:15:15 +0530 Subject: [PATCH 2/5] create newfeature.rst --- doc/changes/devel/0000.newfeature.rst | 1 + 1 file changed, 1 insertion(+) create mode 100644 doc/changes/devel/0000.newfeature.rst diff --git a/doc/changes/devel/0000.newfeature.rst b/doc/changes/devel/0000.newfeature.rst new file mode 100644 index 00000000000..7b9675e8dd1 --- /dev/null +++ b/doc/changes/devel/0000.newfeature.rst @@ -0,0 +1 @@ +Add new example :ref:`ex-ica-comp` that compares performance of multiple ICA algorithms (FastICA, Picard, Infomax, Extended Infomax) on clean vs noisy MEG data showing differences in runtime, component detection and EOG-related ICs, by `Ganasekhar Kalla`_. From f4b0bd8bd3fb270d540283db8bd8bb0b045f41f3 Mon Sep 17 00:00:00 2001 From: KALLA GANASEKHAR Date: Tue, 15 Apr 2025 20:49:51 +0530 Subject: [PATCH 3/5] update ica_comparison.py added visualization plots for different ICA algorithms and added fit times for both clean and noisy data. --- doc/changes/devel/0000.enhancement.rst | 1 + doc/changes/names.inc | 1 + examples/preprocessing/ica_comparison.py | 3 +-- 3 files changed, 3 insertions(+), 2 deletions(-) create mode 100644 doc/changes/devel/0000.enhancement.rst diff --git a/doc/changes/devel/0000.enhancement.rst b/doc/changes/devel/0000.enhancement.rst new file mode 100644 index 00000000000..ff668050076 --- /dev/null +++ b/doc/changes/devel/0000.enhancement.rst @@ -0,0 +1 @@ +Add example :ref:`ex-ica-comp` comparing ICA algorithms on clean vs noisy MEG data, by :newcontrib:`Ganasekhar Kalla`. diff --git a/doc/changes/names.inc b/doc/changes/names.inc index 0d5ee6a5c73..c4ee9923948 100644 --- a/doc/changes/names.inc +++ b/doc/changes/names.inc @@ -96,6 +96,7 @@ .. _Florin Pop: https://github.com/florin-pop .. _Frederik Weber: https://github.com/Frederik-D-Weber .. _Fu-Te Wong: https://github.com/zuxfoucault +.. _Ganasekhar Kalla: https://github.com/Ganasekhar-gif .. _Gennadiy Belonosov: https://github.com/Genuster .. _Geoff Brookshire: https://github.com/gbrookshire .. _George O'Neill: https://georgeoneill.github.io diff --git a/examples/preprocessing/ica_comparison.py b/examples/preprocessing/ica_comparison.py index cef36fc6d8d..90d212c8231 100644 --- a/examples/preprocessing/ica_comparison.py +++ b/examples/preprocessing/ica_comparison.py @@ -95,8 +95,6 @@ def run_all_ica(raw_input, label, reject): icas = {} fit_times = {} eog_components = {} - - for method, params in [ ("fastica", None), ("picard", None), @@ -123,6 +121,7 @@ def run_all_ica(raw_input, label, reject): # %% + # Run on both raw versions icas_clean, times_clean, eog_clean = run_all_ica( raw_clean, "clean", reject_clean From feaa0b6bd3401b6c8909cf8e8d56ac3917496b09 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 15 Apr 2025 15:22:29 +0000 Subject: [PATCH 4/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- examples/preprocessing/ica_comparison.py | 34 ++++++++++-------------- 1 file changed, 14 insertions(+), 20 deletions(-) diff --git a/examples/preprocessing/ica_comparison.py b/examples/preprocessing/ica_comparison.py index 90d212c8231..3cf083573ef 100644 --- a/examples/preprocessing/ica_comparison.py +++ b/examples/preprocessing/ica_comparison.py @@ -24,12 +24,14 @@ # %% +from pathlib import Path +from time import time + +import numpy as np + import mne -from mne.preprocessing import ICA from mne.datasets import sample -import numpy as np -from time import time -from pathlib import Path +from mne.preprocessing import ICA print(__doc__) @@ -42,7 +44,7 @@ # Load sample dataset data_path = Path(sample.data_path()) -raw_file = data_path / 'MEG' / 'sample' / 'sample_audvis_raw.fif' +raw_file = data_path / "MEG" / "sample" / "sample_audvis_raw.fif" raw = mne.io.read_raw_fif(raw_file, preload=True) raw.pick_types(meg=True, eeg=False, eog=True) raw.crop(0, 60) @@ -87,6 +89,7 @@ def run_ica(raw_input, method, fit_params=None, reject=None): return ica, fit_time + # %% @@ -107,9 +110,7 @@ def run_all_ica(raw_input, label, reject): icas[full_label] = ica fit_times[full_label] = t - eog_inds, _ = ica.find_bads_eog( - raw_input, threshold=3.0, verbose='ERROR' - ) + eog_inds, _ = ica.find_bads_eog(raw_input, threshold=3.0, verbose="ERROR") if eog_inds: eog_components[full_label] = eog_inds[0] print(f"{full_label}:Detected EOG comp at index {eog_inds[0]}") @@ -119,16 +120,13 @@ def run_all_ica(raw_input, label, reject): return icas, fit_times, eog_components + # %% # Run on both raw versions -icas_clean, times_clean, eog_clean = run_all_ica( - raw_clean, "clean", reject_clean -) -icas_noisy, times_noisy, eog_noisy = run_all_ica( - raw_noisy, "noisy", reject_noisy -) +icas_clean, times_clean, eog_clean = run_all_ica(raw_clean, "clean", reject_clean) +icas_noisy, times_noisy, eog_noisy = run_all_ica(raw_noisy, "noisy", reject_noisy) # Combine results icas = {**icas_clean, **icas_noisy} @@ -143,9 +141,7 @@ def run_all_ica(raw_input, label, reject): comp = eog_comps.get(key) if comp is not None: icas[key].plot_components( - picks=[comp], - title=f"{key} - EOG Component (Clean Data)", - show=True + picks=[comp], title=f"{key} - EOG Component (Clean Data)", show=True ) # %% @@ -156,7 +152,5 @@ def run_all_ica(raw_input, label, reject): comp = eog_comps.get(key) if comp is not None: icas[key].plot_components( - picks=[comp], - title=f"{key} - EOG Component (Noisy Data)", - show=True + picks=[comp], title=f"{key} - EOG Component (Noisy Data)", show=True ) From ad7b1a16d6172e9ff8235a138c69134c574b790c Mon Sep 17 00:00:00 2001 From: KALLA GANASEKHAR Date: Tue, 15 Apr 2025 20:53:29 +0530 Subject: [PATCH 5/5] create 13215.enhancement.rst created 13215.enhancement.rst file --- doc/changes/devel/13215.enhancement.rst | 1 + 1 file changed, 1 insertion(+) create mode 100644 doc/changes/devel/13215.enhancement.rst diff --git a/doc/changes/devel/13215.enhancement.rst b/doc/changes/devel/13215.enhancement.rst new file mode 100644 index 00000000000..ff668050076 --- /dev/null +++ b/doc/changes/devel/13215.enhancement.rst @@ -0,0 +1 @@ +Add example :ref:`ex-ica-comp` comparing ICA algorithms on clean vs noisy MEG data, by :newcontrib:`Ganasekhar Kalla`.