.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/zapline/plot_05_adaptive_advanced.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_zapline_plot_05_adaptive_advanced.py: ZapLine-plus: Advanced Settings and Features. ============================================= This example focuses on the extra adaptive features around the core ZapLine step: harmonic handling, chunk metadata, and the hybrid fallback path when a plain component projection is not enough. Authors: Sina Esmaeili (sina.esmaeili@umontreal.ca) .. GENERATED FROM PYTHON SOURCE LINES 11-20 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from scipy import signal from mne_denoise.viz import plot_psd_comparison from mne_denoise.zapline import ZapLine from mne_denoise.zapline.adaptive import detect_harmonics .. GENERATED FROM PYTHON SOURCE LINES 21-24 Simulate Data with Harmonics ---------------------------- Create a signal with fundamental (50 Hz) plus two harmonics (100, 150 Hz). .. GENERATED FROM PYTHON SOURCE LINES 24-63 .. code-block:: Python print("=== ZapLine-plus Advanced Features Demo ===\n") sfreq = 500 # High enough to capture harmonics n_ch = 32 duration = 60 n_times = int(duration * sfreq) times = np.arange(n_times) / sfreq rng = np.random.RandomState(42) # Brain signal (pink noise) brain = rng.randn(n_ch, n_times) * 0.5 # Line noise with harmonics fundamental = 50.0 harmonics = [100.0, 150.0] # 2f, 3f # Create spatial topographies for each frequency topos = [] for i in range(3): topo = rng.randn(n_ch) topo /= np.linalg.norm(topo) topos.append(topo) # Add noise at each frequency noise = np.zeros((n_ch, n_times)) freqs = [fundamental] + harmonics amplitudes = [5.0, 2.0, 1.0] # Decreasing amplitude for harmonics for freq, amp, topo in zip(freqs, amplitudes, topos): ts = np.sin(2 * np.pi * freq * times) * amp noise += np.outer(topo, ts) data = brain + noise print(f"Data: {n_ch} channels × {n_times} samples ({duration}s)") print(f"Frequencies: {freqs} Hz with amplitudes {amplitudes}") .. rst-class:: sphx-glr-script-out .. code-block:: none === ZapLine-plus Advanced Features Demo === Data: 32 channels × 30000 samples (60s) Frequencies: [50.0, 100.0, 150.0] Hz with amplitudes [5.0, 2.0, 1.0] .. GENERATED FROM PYTHON SOURCE LINES 64-67 Step 1: Detect Harmonics ------------------------ The `detect_harmonics` function explicitly finds harmonics present in data. .. GENERATED FROM PYTHON SOURCE LINES 67-73 .. code-block:: Python print("\n--- Step 1: Harmonic Detection ---") detected = detect_harmonics(data, sfreq, fundamental) print(f"Fundamental: {fundamental} Hz") print(f"Detected harmonics: {detected}") .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 1: Harmonic Detection --- Fundamental: 50.0 Hz Detected harmonics: [100.0, 150.0] .. GENERATED FROM PYTHON SOURCE LINES 74-77 Step 2: Apply ZapLine-plus with Harmonic Processing ---------------------------------------------------- Enable `process_harmonics=True` to automatically clean harmonics. .. GENERATED FROM PYTHON SOURCE LINES 77-99 .. code-block:: Python print("\n--- Step 2: ZapLine-plus with Harmonics ---") zapline = ZapLine( sfreq=sfreq, line_freq=None, # Auto-detect adaptive=True, adaptive_params={ "process_harmonics": True, "max_harmonics": 3, "hybrid_fallback": True, "fmin": 45, "fmax": 55, }, ) cleaned = zapline.fit_transform(data) result = zapline.adaptive_results_ data_clean = result["cleaned"] print("Cleaning complete!") print(f"Chunk info entries: {len(result['chunk_info'])}") # print(f"Topography entries: {len(result['removed_topographies'])}") .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 2: ZapLine-plus with Harmonics --- Cleaning complete! Chunk info entries: 3 .. GENERATED FROM PYTHON SOURCE LINES 100-103 Step 3: Examine Per-Chunk Metadata ---------------------------------- The `chunk_info` contains detailed information for each segment processed. .. GENERATED FROM PYTHON SOURCE LINES 103-112 .. code-block:: Python print("\n--- Step 3: Chunk Metadata ---") for i, info in enumerate(result["chunk_info"]): print(f"Chunk {i + 1}:") print(f" Frequency: {info['fine_freq']:.2f} Hz") print(f" Time range: {info['start'] / sfreq:.1f}s - {info['end'] / sfreq:.1f}s") print(f" Components removed: {info['n_removed']}") print(f" Artifact present: {info['artifact_present']}") .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 3: Chunk Metadata --- Chunk 1: Frequency: 50.00 Hz Time range: 0.0s - 60.0s Components removed: 3 Artifact present: True Chunk 2: Frequency: 99.95 Hz Time range: 0.0s - 60.0s Components removed: 0 Artifact present: False Chunk 3: Frequency: 149.95 Hz Time range: 0.0s - 60.0s Components removed: 0 Artifact present: False .. GENERATED FROM PYTHON SOURCE LINES 113-117 Step 4: Visualize Removed Topographies -------------------------------------- Chunk-wise metadata are often more informative than a single global summary when the contamination changes over time. .. GENERATED FROM PYTHON SOURCE LINES 117-121 .. code-block:: Python print("\n--- Step 4: Topography Outputs ---") # (To support this, we would need to store `est.patterns_` from each chunk pass) .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 4: Topography Outputs --- .. GENERATED FROM PYTHON SOURCE LINES 122-125 Step 5: Compare Spectra ----------------------- Visualize noise reduction at fundamental and harmonics. .. GENERATED FROM PYTHON SOURCE LINES 125-133 .. code-block:: Python print("\n--- Step 5: Spectral Comparison ---") # Use our reusable viz function for a quick overview plot_psd_comparison( data, data_clean, sfreq=sfreq, line_freq=fundamental, fmax=180, show=True ) .. image-sg:: /auto_examples/zapline/images/sphx_glr_plot_05_adaptive_advanced_001.png :alt: PSD Comparison :srcset: /auto_examples/zapline/images/sphx_glr_plot_05_adaptive_advanced_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 5: Spectral Comparison ---
.. GENERATED FROM PYTHON SOURCE LINES 134-138 Detailed Harmonic Comparison ---------------------------- Detailed per-frequency comparison makes it easier to check that harmonics drop together rather than leaving isolated residual peaks. .. GENERATED FROM PYTHON SOURCE LINES 138-161 .. code-block:: Python fig, axes = plt.subplots(1, 3, figsize=(14, 4)) for i, (freq, ax) in enumerate(zip(freqs, axes)): f_o, psd_o = signal.welch(data, fs=sfreq, nperseg=int(sfreq * 2), axis=-1) f_c, psd_c = signal.welch(data_clean, fs=sfreq, nperseg=int(sfreq * 2), axis=-1) psd_o_mean = np.mean(10 * np.log10(psd_o + 1e-20), axis=0) psd_c_mean = np.mean(10 * np.log10(psd_c + 1e-20), axis=0) ax.plot(f_o, psd_o_mean, "k-", label="Original", lw=1) ax.plot(f_c, psd_c_mean, "g-", label="Cleaned", lw=1.5) ax.set_xlim(freq - 5, freq + 5) ax.axvline(freq, color="r", ls="--", alpha=0.5) ax.set_xlabel("Frequency (Hz)") ax.set_ylabel("Power (dB)") ax.set_title(f"{int(freq)} Hz {'(Fundamental)' if i == 0 else '(Harmonic)'}") ax.legend(fontsize=8) ax.grid(True, alpha=0.3) plt.suptitle("Noise Reduction at Fundamental and Harmonics") plt.tight_layout() plt.show() .. image-sg:: /auto_examples/zapline/images/sphx_glr_plot_05_adaptive_advanced_002.png :alt: Noise Reduction at Fundamental and Harmonics, 50 Hz (Fundamental), 100 Hz (Harmonic), 150 Hz (Harmonic) :srcset: /auto_examples/zapline/images/sphx_glr_plot_05_adaptive_advanced_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 162-165 Step 6: Hybrid Fallback Demo ---------------------------- When ZapLine alone can't clean fully, hybrid_fallback applies a notch. .. GENERATED FROM PYTHON SOURCE LINES 165-176 .. code-block:: Python print("\n--- Step 6: Quantitative Results ---") # Calculate noise power reduction for freq in freqs: mask = (f_o > freq - 0.5) & (f_o < freq + 0.5) pow_o = np.mean(psd_o[:, mask]) pow_c = np.mean(psd_c[:, mask]) reduction_db = 10 * np.log10(pow_o / max(pow_c, 1e-20)) print(f"{int(freq):>3} Hz: {reduction_db:.1f} dB reduction") print("\n[OK] Advanced features demo complete!") .. rst-class:: sphx-glr-script-out .. code-block:: none --- Step 6: Quantitative Results --- 50 Hz: 27.7 dB reduction 100 Hz: 19.7 dB reduction 150 Hz: 14.0 dB reduction [OK] Advanced features demo complete! .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 2.473 seconds) .. _sphx_glr_download_auto_examples_zapline_plot_05_adaptive_advanced.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_05_adaptive_advanced.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_05_adaptive_advanced.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_05_adaptive_advanced.zip `