MNE-Python 第10天学习笔记:结果报告与可视化
一、为什么需要报告和可视化?
1.1 数据分析的"最后一公里"
数据分析的完整流程: 原始数据 → 预处理 → 分段 → 分析 → 📊 报告/图表 ↑ 这是别人看到的! 前面的工作做得再好,如果图和报告不好看: - 审稿人看不懂 → 被拒稿 - 导师不满意 → 返工 - 自己过段时间也忘了分析过程1.2 今天要完成的事
1. 回顾前9天的分析成果 2. 用 MNE Report 生成交互式 HTML 报告 3. 创建发表级别的对比图 4. 保存高质量图片二、环境准备
2.1 导入库
# ========== 环境设置 ========== # 设置 matplotlib 后端 # 今天以保存图片为主,用 Agg 后端即可 import matplotlib matplotlib.use('Agg') # pyplot:matplotlib 的绘图接口 import matplotlib.pyplot as plt # mne:脑电分析核心库 import mne # numpy:科学计算库 import numpy as np # os:操作系统路径处理 import os # warnings:警告控制 import warnings warnings.filterwarnings('ignore') # ========== 中文字体设置 ========== plt.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei'] plt.rcParams['axes.unicode_minus'] = False print("="*60) print("MNE-Python 第10天:结果报告与可视化") print("="*60) print("\n🎉 今天是最后一天!") print("我们将整合前9天的所有知识,生成专业报告!")2.2 加载预处理好的数据
# ---------- 1. 加载数据并快速预处理 ---------- print("\n加载数据...") sample_data_folder = mne.datasets.sample.data_path() raw_fname = os.path.join(sample_data_folder, 'MEG', 'sample', 'sample_audvis_raw.fif') raw = mne.io.read_raw_fif(raw_fname, preload=False) # 提取 EEG 通道 raw_eeg = raw.copy().pick_types(eeg=True, eog=True, stim=True) # 重命名通道 eeg_names = [ch for ch in raw_eeg.ch_names if ch.startswith('EEG')] montage = mne.channels.make_standard_montage('standard_1020') standard_names = montage.ch_names[:len(eeg_names)] raw_eeg.rename_channels(dict(zip(eeg_names, standard_names))) raw_eeg.set_montage(montage) # 加载到内存并滤波 raw_eeg.load_data() raw_eeg.notch_filter(freqs=60, picks='eeg', verbose=False) raw_eeg.filter(l_freq=1, h_freq=40, picks='eeg', verbose=False) raw_eeg.set_eeg_reference('average', verbose=False) print("✅ 数据预处理完成") # 提取事件 events = mne.find_events(raw_eeg, stim_channel='STI 014') event_id = {'听觉/左耳': 1, '听觉/右耳': 2, '视觉/左眼': 3, '视觉/右眼': 4}三、创建 Epochs 和 ERP
# ---------- 2. 创建 Epochs 和计算 ERP ---------- print("\n创建 Epochs...") epochs = mne.Epochs( raw_eeg, events, event_id=event_id, tmin=-0.2, tmax=0.5, baseline=(-0.2, 0), reject=dict(eeg=150e-6), preload=True, verbose=False ) print(f"✅ Epochs: {len(epochs)} 个分段") # 计算各条件的 ERP print("计算 ERP...") evoked_dict = {} for cond_name in event_id.keys(): evoked_dict[cond_name] = epochs[cond_name].average() print(f" {cond_name}: {len(epochs[cond_name])} 个分段") # 合并听觉和视觉条件 evoked_auditory = mne.combine_evoked( [evoked_dict['听觉/左耳'], evoked_dict['听觉/右耳']], weights=[0.5, 0.5]) evoked_auditory.comment = '听觉' evoked_visual = mne.combine_evoked( [evoked_dict['视觉/左眼'], evoked_dict['视觉/右眼']], weights=[0.5, 0.5]) evoked_visual.comment = '视觉' # 差异波 evoked_diff = mne.combine_evoked( [evoked_auditory, evoked_visual], weights=[1, -1]) evoked_diff.comment = '听觉-视觉' print("✅ 差异波计算完成")四、时频分析
# ---------- 3. 时频分析 ---------- print("\n计算时频分析...") # 创建更长的 Epochs 用于时频分析 epochs_tfr = mne.Epochs( raw_eeg, events, event_id=event_id, tmin=-1.0, tmax=2.0, baseline=(-0.5, 0), reject=dict(eeg=150e-6), preload=True, verbose=False ) # 选择 Cz 通道 epochs_tfr.pick(['Cz']) # Morlet 小波变换 frequencies = np.arange(4, 31, 1) n_cycles = frequencies / 2 power = mne.time_frequency.tfr_morlet( epochs_tfr['听觉/左耳'], freqs=frequencies, n_cycles=n_cycles, return_itc=False, average=True, decim=2, verbose=False ) print("✅ 时频分析完成")五、MNE Report:自动生成 HTML 报告
5.1 什么是 MNE Report?
MNE Report = 一个交互式 HTML 网页 里面可以包含: - 原始数据预览 - 功率谱密度图 - 事件标记 - ERP 波形 - 脑地形图 - 时频图 - 源定位结果 好处: - 所有分析结果集中在一个文件 - 可以交互(缩放、切换时间点) - 方便分享给导师/合作者 - 方便自己回顾分析过程5.2 创建报告
# ---------- 4. 创建 MNE Report ---------- print("\n" + "="*60) print("生成 MNE HTML 报告") print("="*60) # Report():创建报告对象 # title:报告标题 # subjects_dir:MRI 数据目录(源定位可视化需要) report = mne.Report( title='脑电数据分析报告 - 听觉/视觉实验', subjects_dir=os.path.join(sample_data_folder, 'subjects'), verbose=False ) print("报告中添加内容...") # ----- 4.1 添加原始数据 ----- print(" + 原始数据预览") report.add_raw( raw=raw_eeg, title='原始数据(滤波后)', psd=True, # 同时显示 PSD butterfly=True # 蝴蝶图模式 ) # ----- 4.2 添加事件 ----- print(" + 事件标记") report.add_events( events=events, event_id=event_id, sfreq=raw_eeg.info['sfreq'], title='实验事件' ) # ----- 4.3 添加 Epochs ----- print(" + 分段数据") report.add_epochs( epochs=epochs, title='数据分段(Epochs)' ) # ----- 4.4 添加 ERP ----- print(" + ERP 波形") for cond_name, evoked in evoked_dict.items(): safe_name = cond_name.replace('/', '_') report.add_evokeds( evokeds=evoked, titles=cond_name, n_time_points=10, # 地形图的时间点数 ) # ----- 4.5 添加 ERP 对比 ----- print(" + ERP 对比") report.add_evokeds( evokeds=[evoked_dict['听觉/左耳'], evoked_dict['听觉/右耳'], evoked_dict['视觉/左眼'], evoked_dict['视觉/右眼']], titles=['听觉/左耳', '听觉/右耳', '视觉/左眼', '视觉/右眼'], n_time_points=8, ) # ----- 4.6 添加时频图 ----- print(" + 时频图") # 将时频图保存为图片,再添加到报告 fig_tfr = power.plot( picks='Cz', baseline=(-0.5, 0), mode='percent', vmin=-50, vmax=50, show=False ) fig_tfr.suptitle('Cz ERD/ERS', fontsize=12) report.add_figure( fig=fig_tfr, title='Cz通道 ERD/ERS', image_format='png' ) plt.close(fig_tfr) # ----- 4.7 添加差异波 ----- print(" + 差异波") report.add_evokeds( evokeds=evoked_diff, titles='听觉 - 视觉', n_time_points=8, ) # ----- 保存报告 ----- report_fname = 'day10_eeg_report.html' report.save(report_fname, overwrite=True) print(f"\n✅ 报告已保存为 {report_fname}") print(f" 用浏览器打开即可查看交互式报告!")add_raw()参数详解:
add_evokeds()参数详解:
六、创建发表级图表
6.1 多子图综合布局
# ---------- 5. 发表级综合图表 ---------- print("\n" + "="*60) print("创建发表级图表") print("="*60) # 创建 3×3 的大图布局 fig = plt.figure(figsize=(18, 14)) fig.suptitle('听觉 vs 视觉 脑电分析综合图', fontsize=18, fontweight='bold') # 选择 Cz 通道(头顶中央,ERP 分析最常用) plot_ch = 'Cz' # ----- 子图1:听觉左耳 ERP ----- ax1 = fig.add_subplot(3, 3, 1) evoked_dict['听觉/左耳'].plot( picks=plot_ch, axes=ax1, show=False, spatial_colors=False, time_unit='ms' ) ax1.set_title('听觉/左耳 ERP', fontsize=12, fontweight='bold') ax1.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax1.grid(True, alpha=0.3) # ----- 子图2:听觉右耳 ERP ----- ax2 = fig.add_subplot(3, 3, 2) evoked_dict['听觉/右耳'].plot( picks=plot_ch, axes=ax2, show=False, spatial_colors=False, time_unit='ms' ) ax2.set_title('听觉/右耳 ERP', fontsize=12, fontweight='bold') ax2.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax2.grid(True, alpha=0.3) # ----- 子图3:听觉地形图(150ms) ----- ax3 = fig.add_subplot(3, 3, 3) evoked_dict['听觉/左耳'].plot_topomap( times=[0.15], ch_type='eeg', axes=ax3, show=False, average=0.05, time_unit='ms' ) ax3.set_title('听觉地形图 (150ms)', fontsize=12, fontweight='bold') # ----- 子图4:视觉左眼 ERP ----- ax4 = fig.add_subplot(3, 3, 4) evoked_dict['视觉/左眼'].plot( picks=plot_ch, axes=ax4, show=False, spatial_colors=False, time_unit='ms' ) ax4.set_title('视觉/左眼 ERP', fontsize=12, fontweight='bold') ax4.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax4.grid(True, alpha=0.3) # ----- 子图5:视觉右眼 ERP ----- ax5 = fig.add_subplot(3, 3, 5) evoked_dict['视觉/右眼'].plot( picks=plot_ch, axes=ax5, show=False, spatial_colors=False, time_unit='ms' ) ax5.set_title('视觉/右眼 ERP', fontsize=12, fontweight='bold') ax5.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax5.grid(True, alpha=0.3) # ----- 子图6:视觉地形图(150ms) ----- ax6 = fig.add_subplot(3, 3, 6) evoked_dict['视觉/左眼'].plot_topomap( times=[0.15], ch_type='eeg', axes=ax6, show=False, average=0.05, time_unit='ms' ) ax6.set_title('视觉地形图 (150ms)', fontsize=12, fontweight='bold') # ----- 子图7:听觉 vs 视觉对比 ----- ax7 = fig.add_subplot(3, 3, 7) mne.viz.plot_compare_evokeds( {'听觉': evoked_auditory, '视觉': evoked_visual}, picks=plot_ch, axes=ax7, show=False, legend='upper right' ) ax7.set_title('听觉 vs 视觉', fontsize=12, fontweight='bold') ax7.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax7.grid(True, alpha=0.3) # ----- 子图8:差异波 ----- ax8 = fig.add_subplot(3, 3, 8) evoked_diff.plot( picks=plot_ch, axes=ax8, show=False, spatial_colors=False, time_unit='ms' ) ax8.set_title('差异波(听觉 - 视觉)', fontsize=12, fontweight='bold') ax8.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax8.axhline(y=0, color='black', linestyle='-', alpha=0.3) ax8.grid(True, alpha=0.3) # ----- 子图9:时频图 ----- ax9 = fig.add_subplot(3, 3, 9) power.plot( picks='Cz', baseline=(-0.5, 0), mode='percent', vmin=-50, vmax=50, axes=ax9, show=False, colorbar=True ) ax9.set_title('Cz ERD/ERS', fontsize=12, fontweight='bold') # 调整布局并保存 plt.tight_layout() fig.savefig('day10_comprehensive_figure.png', dpi=300, bbox_inches='tight', facecolor='white', edgecolor='none') print("✅ 综合图表已保存为 day10_comprehensive_figure.png (300 DPI)") plt.close(fig)6.2 不同条件的地形图时间序列
# ---------- 6. 地形图时间序列 ---------- print("\n创建地形图时间序列...") # 创建 2×4 的地形图对比 # 上行:听觉条件,下行:视觉条件 # 列:50, 100, 150, 200 毫秒 times_topo = [0.05, 0.10, 0.15, 0.20] # 秒 fig, axes = plt.subplots(2, 4, figsize=(20, 10)) fig.suptitle('听觉 vs 视觉 脑地形图时间序列', fontsize=16, fontweight='bold') for col, t in enumerate(times_topo): # 上行:听觉 evoked_dict['听觉/左耳'].plot_topomap( times=[t], ch_type='eeg', axes=axes[0, col], show=False, average=0.05, time_unit='ms', vmin=-3, vmax=3 # 统一颜色范围便于比较 ) axes[0, col].set_title(f'听觉 {t*1000:.0f}ms', fontsize=12, fontweight='bold') # 下行:视觉 evoked_dict['视觉/左眼'].plot_topomap( times=[t], ch_type='eeg', axes=axes[1, col], show=False, average=0.05, time_unit='ms', vmin=-3, vmax=3 ) axes[1, col].set_title(f'视觉 {t*1000:.0f}ms', fontsize=12, fontweight='bold') plt.tight_layout() fig.savefig('day10_topomap_timeseries.png', dpi=300, bbox_inches='tight', facecolor='white', edgecolor='none') print("✅ 地形图时间序列已保存为 day10_topomap_timeseries.png (300 DPI)") plt.close(fig)七、高质量图片保存
7.1 图片格式选择
print("\n" + "="*60) print("高质量图片保存") print("="*60) print(""" 不同图片格式的用途: ┌─────────┬──────────┬──────────────────────┐ │ 格式 │ 文件大小 │ 最适合 │ ├─────────┼──────────┼──────────────────────┤ │ PNG │ 中等 │ 日常使用、PPT、博客 │ │ TIFF │ 较大 │ 论文投稿(位图) │ │ PDF │ 小 │ 论文投稿(矢量图) │ │ SVG │ 小 │ 可编辑矢量图 │ │ EPS │ 小 │ LaTeX 论文 │ └─────────┴──────────┴──────────────────────┘ """) # 保存多种格式 print("保存多种格式的图片...") # 创建示例图 fig, ax = plt.subplots(figsize=(8, 6)) evoked_dict['听觉/左耳'].plot( picks='Cz', axes=ax, show=False, spatial_colors=False, time_unit='ms' ) ax.set_title('听觉/左耳 ERP (Cz)', fontsize=14, fontweight='bold') ax.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax.grid(True, alpha=0.3) # PNG - 日常使用 fig.savefig('day10_erp.png', dpi=150, bbox_inches='tight', facecolor='white', edgecolor='none') print(" ✅ day10_erp.png (150 DPI - 日常使用)") # TIFF - 论文投稿 fig.savefig('day10_erp.tiff', dpi=300, bbox_inches='tight', facecolor='white', edgecolor='none', pil_kwargs={'compression': 'tiff_lzw'}) print(" ✅ day10_erp.tiff (300 DPI - 论文投稿)") # PDF - 矢量图 fig.savefig('day10_erp.pdf', bbox_inches='tight', facecolor='white', edgecolor='none') print(" ✅ day10_erp.pdf (矢量图 - LaTeX 论文)") plt.close(fig)7.2 savefig 参数详解
# savefig() 的参数说明: fig.savefig( 'filename.png', # 文件名(扩展名决定格式) dpi=300, # 分辨率(dots per inch) # 150 = 屏幕显示 # 300 = 论文投稿标准 # 600 = 印刷出版 bbox_inches='tight', # 自动裁剪多余白边 facecolor='white', # 图片背景色 edgecolor='none', # 图片边框颜色 # TIFF 额外参数: pil_kwargs={'compression': 'tiff_lzw'} # LZW 无损压缩 )八、项目总结
# ---------- 7. 生成项目总结 ---------- print("\n" + "="*60) print("🎉 项目总结") print("="*60) print(""" ┌────────────────────────────────────────────┐ │ 10天 MNE-Python 学习之旅 │ ├────────────────────────────────────────────┤ │ │ │ 第1天:环境搭建与数据初探 │ │ 第2天:Montage与通道信息管理 │ │ 第3天:事件与标记处理 │ │ 第4天:数据预处理 - 滤波与重参考 │ │ 第5天:数据预处理 - ICA伪迹处理 │ │ 第6天:分段(Epoching)与基线校正 │ │ 第7天:事件相关电位(ERP)分析 │ │ 第8天:时频分析(ERD/ERS) │ │ 第9天:源定位基础 │ │ 第10天:结果报告与可视化 │ │ │ │ 你已掌握: │ │ ✅ 脑电数据加载和预处理 │ │ ✅ 事件提取和分段 │ │ ✅ ERP 分析和可视化 │ │ ✅ 时频分析和 ERD/ERS │ │ ✅ 源定位基础 │ │ ✅ 生成专业报告和图表 │ │ │ └────────────────────────────────────────────┘ """) print("保存的文件:") print(" 📄 day10_eeg_report.html - 交互式 HTML 报告") print(" 🖼️ day10_comprehensive_figure.png - 综合图表") print(" 🖼️ day10_topomap_timeseries.png - 地形图时间序列") print(" 🖼️ day10_erp.png - ERP 示例图 (PNG)") print(" 🖼️ day10_erp.tiff - ERP 示例图 (TIFF)") print(" 🖼️ day10_erp.pdf - ERP 示例图 (PDF)")九、第10天完整代码
# ========== 环境设置 ========== import matplotlib matplotlib.use('Agg') import matplotlib.pyplot as plt import mne import numpy as np import os import warnings warnings.filterwarnings('ignore') plt.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei'] plt.rcParams['axes.unicode_minus'] = False print("="*60) print("MNE-Python 第10天:结果报告与可视化") print("="*60) print("\n🎉 今天是最后一天!") # ---------- 1. 加载数据并快速预处理 ---------- print("\n加载数据...") sample_data_folder = mne.datasets.sample.data_path() raw_fname = os.path.join(sample_data_folder, 'MEG', 'sample', 'sample_audvis_raw.fif') raw = mne.io.read_raw_fif(raw_fname, preload=False) raw_eeg = raw.copy().pick_types(eeg=True, eog=True, stim=True) eeg_names = [ch for ch in raw_eeg.ch_names if ch.startswith('EEG')] montage = mne.channels.make_standard_montage('standard_1020') standard_names = montage.ch_names[:len(eeg_names)] raw_eeg.rename_channels(dict(zip(eeg_names, standard_names))) raw_eeg.set_montage(montage) raw_eeg.load_data() raw_eeg.notch_filter(freqs=60, picks='eeg', verbose=False) raw_eeg.filter(l_freq=1, h_freq=40, picks='eeg', verbose=False) raw_eeg.set_eeg_reference('average', verbose=False) events = mne.find_events(raw_eeg, stim_channel='STI 014') event_id = {'听觉/左耳': 1, '听觉/右耳': 2, '视觉/左眼': 3, '视觉/右眼': 4} print("✅ 数据预处理完成") # ---------- 2. 创建 Epochs 和计算 ERP ---------- epochs = mne.Epochs(raw_eeg, events, event_id=event_id, tmin=-0.2, tmax=0.5, baseline=(-0.2, 0), reject=dict(eeg=150e-6), preload=True, verbose=False) evoked_dict = {} for cond_name in event_id.keys(): evoked_dict[cond_name] = epochs[cond_name].average() evoked_auditory = mne.combine_evoked( [evoked_dict['听觉/左耳'], evoked_dict['听觉/右耳']], weights=[0.5, 0.5]) evoked_visual = mne.combine_evoked( [evoked_dict['视觉/左眼'], evoked_dict['视觉/右眼']], weights=[0.5, 0.5]) evoked_diff = mne.combine_evoked([evoked_auditory, evoked_visual], weights=[1, -1]) print("✅ ERP 计算完成") # ---------- 3. 创建 MNE Report ---------- print("\n生成 HTML 报告...") report = mne.Report(title='脑电数据分析报告', verbose=False) report.add_raw(raw=raw_eeg, title='原始数据', psd=True) report.add_events(events=events, event_id=event_id, sfreq=raw_eeg.info['sfreq'], title='实验事件') report.add_epochs(epochs=epochs, title='数据分段') for cond_name, evoked in evoked_dict.items(): report.add_evokeds(evokeds=evoked, titles=cond_name, n_time_points=8) report.add_evokeds(evokeds=evoked_diff, titles='听觉-视觉', n_time_points=8) report.save('day10_eeg_report.html', overwrite=True) print("✅ 报告已保存为 day10_eeg_report.html") # ---------- 4. 创建发表级综合图表 ---------- print("\n创建综合图表...") fig = plt.figure(figsize=(18, 14)) fig.suptitle('听觉 vs 视觉 脑电分析综合图', fontsize=18, fontweight='bold') plot_ch = 'Cz' # 9个子图(简化版) for idx, (cond_name, pos) in enumerate(zip( ['听觉/左耳', '听觉/右耳', '视觉/左眼', '视觉/右耳'], [331, 332, 334, 335])): ax = fig.add_subplot(pos) evoked_dict[cond_name].plot(picks=plot_ch, axes=ax, show=False, spatial_colors=False, time_unit='ms') ax.set_title(cond_name, fontsize=12, fontweight='bold') ax.axvline(x=0, color='red', linestyle='--', alpha=0.5) ax.grid(True, alpha=0.3) # 地形图 ax_topo1 = fig.add_subplot(333) evoked_dict['听觉/左耳'].plot_topomap( times=[0.15], ch_type='eeg', axes=ax_topo1, show=False, average=0.05) ax_topo1.set_title('听觉地形图', fontsize=12) ax_topo2 = fig.add_subplot(336) evoked_dict['视觉/左眼'].plot_topomap( times=[0.15], ch_type='eeg', axes=ax_topo2, show=False, average=0.05) ax_topo2.set_title('视觉地形图', fontsize=12) # 对比 ax_comp = fig.add_subplot(337) mne.viz.plot_compare_evokeds( {'听觉': evoked_auditory, '视觉': evoked_visual}, picks=plot_ch, axes=ax_comp, show=False, legend='upper right') ax_comp.set_title('听觉 vs 视觉', fontsize=12) ax_comp.axvline(x=0, color='red', linestyle='--', alpha=0.5) # 差异波 ax_diff = fig.add_subplot(338) evoked_diff.plot(picks=plot_ch, axes=ax_diff, show=False, spatial_colors=False, time_unit='ms') ax_diff.set_title('差异波', fontsize=12) ax_diff.axhline(y=0, color='black', linestyle='-', alpha=0.3) plt.tight_layout() fig.savefig('day10_comprehensive_figure.png', dpi=300, bbox_inches='tight', facecolor='white') print("✅ 综合图表已保存") plt.close(fig) # ---------- 5. 项目总结 ---------- print("\n" + "="*60) print("🎉 10天 MNE-Python 学习之旅完成!") print("="*60) print("\n你已掌握:") print(" ✅ 脑电数据加载和预处理") print(" ✅ 事件提取和分段") print(" ✅ ERP 分析和可视化") print(" ✅ 时频分析和 ERD/ERS") print(" ✅ 源定位基础") print(" ✅ 生成专业报告和图表") print("\n保存的文件:") print(" 📄 day10_eeg_report.html") print(" 🖼️ day10_comprehensive_figure.png") print("\n🚀 恭喜完成!你已经是一名合格的脑电数据分析者了!")十、今日总结
📝 核心概念
🛠️ 掌握的技能
🔑 图片保存速查
# 日常使用 fig.savefig('figure.png', dpi=150, bbox_inches='tight') # 论文投稿 fig.savefig('figure.tiff', dpi=300, bbox_inches='tight', pil_kwargs={'compression': 'tiff_lzw'}) # LaTeX 论文 fig.savefig('figure.pdf', bbox_inches='tight')