第 7 章：ERP 分析¶

前面所有内容，最终都是为了这一章。我们先建立"为什么平均许多 epoch 能揭示出隐藏信号"的直观理解（不需要任何 EEG 数据），然后正式实操：把 Epochs 变成 Evoked 对象（也就是 ERP），把它可视化，并比较不同条件。

第一部分：为什么平均有效（一次模拟）¶

核心问题在于：单次试次中的脑反应，通常远小于信号中持续存在的噪声（其他脑活动、肌肉紧张、微小的身体移动）。在任意单一试次上，你往往用肉眼根本看不出这个反应。

但这个反应与事件保持锁时关系（每次的潜伏期都相同），而噪声是随机的（与事件无关）。把许多试次平均在一起，随机噪声会部分相互抵消，而稳定存在的反应则会留存并累积。这个平均结果，就是 ERP。

import numpy as np
import matplotlib.pyplot as plt

# 设置字体，确保图中的中文能正常显示
plt.rcParams["font.sans-serif"] = ["Microsoft YaHei", "SimHei", "DejaVu Sans"]
plt.rcParams["axes.unicode_minus"] = False

rng = np.random.default_rng(42)

sfreq = 250
tmin, tmax = -0.2, 0.8
times = np.arange(tmin, tmax, 1 / sfreq)

# 一个虚构的"真实"脑反应：一个在 300 ms 处达到峰值的凸起（类似 P300）。
# 现实中我们永远无法直接看到这个反应本身 -- 它被埋没在噪声里。
true_response = 3 * np.exp(-0.5 * ((times - 0.3) / 0.08) ** 2)

n_trials = 40
noise_amplitude = 8  # 故意设得比信号大得多 -- 这更符合单试次 EEG 的真实情况

single_trials = np.array(
    [true_response + rng.normal(0, noise_amplitude, size=times.shape) for _ in range(n_trials)]
)

fig, axes = plt.subplots(1, 2, figsize=(12, 4), sharey=True)

axes[0].plot(times, single_trials.T, color="gray", alpha=0.3)
axes[0].plot(times, true_response, color="black", linewidth=2, label="真实（隐藏）反应")
axes[0].set_title(f"{n_trials} 个单一试次（每一条几乎全是噪声）")
axes[0].legend()

axes[1].plot(times, single_trials.mean(axis=0), color="crimson", linewidth=2, label="跨试次平均")
axes[1].plot(times, true_response, color="black", linewidth=1, linestyle="--", label="真实（隐藏）反应")
axes[1].set_title("试次的平均 = ERP")
axes[1].legend()

for ax in axes:
    ax.set_xlabel("相对事件的时间 (s)")
    ax.axvline(0, color="k", linewidth=0.5)
axes[0].set_ylabel("幅度 (a.u.)")
plt.tight_layout()
plt.show()

png

两张图之间没有做任何滤波或清理 —— 唯一的操作就是取平均。对应到 MNE 的对象上：上面的 single_trials 就是 Epochs 所持有的内容（每行一个试次）；single_trials.mean(axis=0) 正是 epochs.average() 所计算的东西，其结果就是一个 Evoked 对象。第 4–6 章里做的所有事情（滤波、ICA、拒绝）存在的唯一目的，都是让单个 epoch 的噪声更小，从而让这个平均结果收敛得更快、更干净 —— 但真正产生 ERP 的，是平均这一步本身。

第二部分：一个真实的 ERP¶

重新构建第 4 章和第 6 章中清理、切分好的数据（本节内容自成一体，所以这个 notebook 可以独立运行）：

import mne
from pathlib import Path

mne.set_log_level("WARNING")

sample_folder = Path(mne.datasets.sample.data_path()) / "MEG" / "sample"
raw = mne.io.read_raw_fif(sample_folder / "sample_audvis_raw.fif", preload=True)
events = mne.find_events(raw, stim_channel="STI 014")
event_id = {
    "auditory/left": 1,
    "auditory/right": 2,
    "visual/left": 3,
    "visual/right": 4,
}

raw.pick(["eeg", "eog"])
raw.filter(l_freq=1.0, h_freq=40.0)

epochs = mne.Epochs(
    raw,
    events,
    event_id=event_id,
    tmin=-0.2,
    tmax=0.5,
    baseline=(-0.2, 0.0),
    preload=True,
    reject=dict(eeg=150e-6, eog=250e-6),
)
epochs

	General
	MNE object type	Epochs
	Measurement date	2002-12-03 at 19:01:10 UTC
	Participant	Unknown
	Experimenter	MEG
	Acquisition
	Total number of events	269
	Events counts	auditory/left: 63 auditory/right: 69 visual/left: 72 visual/right: 65
	Time range	-0.200 – 0.499 s
	Baseline	-0.200 – 0.000 s
	Sampling frequency	600.61 Hz
	Time points	421
	Metadata	No metadata set
	Channels
	EEG	and
	EOG
	Head & sensor digitization	146 points
	Filters
	Highpass	1.00 Hz
	Lowpass	40.00 Hz

`epochs.average()` → `Evoked`¶

这就是上面那次模拟的"单行版"，只不过这次是在真实 epoch 上运行，而不是虚构的试次。

evoked_aud_left = epochs["auditory/left"].average()
print(evoked_aud_left)
evoked_aud_left.plot(picks="eeg")
plt.show()

<Evoked | 'auditory/left' (average, N=63), -0.1998 – 0.49949 s, baseline -0.2 – 0 s, 60 ch, ~3.1 MiB>

png

解读这张图：什么样的成分才算是"真实"的 ERP 成分？¶

你应该能看到，在点击声播放后大约 100 ms 处有一个明显的负向下凹 —— 经典的听觉 N100 反应。plot_joint 会把关键潜伏期处的头皮地形图（各电极间的电压分布）和波形叠加展示在一起 —— 这是一个很好的合理性检查，能确认某个偏转是空间上一致的脑反应，而不是噪声。真正的成分在头皮地形图上应该看起来平滑而聚焦，而不是散乱/随机的。

evoked_aud_left.plot_joint(picks="eeg")
plt.show()

png

量化一个成分：峰值幅度与潜伏期¶

除了肉眼观察，get_peak() 还能在一个时间窗口内，找出最大偏转所在的具体通道、时间和幅度 —— 这正是你真正会拿去报告或做统计的那种数值。

ch_name, latency, amplitude = evoked_aud_left.get_peak(
    tmin=0.05, tmax=0.15, mode="neg", return_amplitude=True
)
print(f"50-150 ms 之间最大的负向偏转: {amplitude * 1e6:.2f} uV，位于 {ch_name}，{latency * 1000:.0f} ms")

50-150 ms 之间最大的负向偏转: -12.00 uV，位于 EEG 014，90 ms

比较条件¶

大多数 ERP 研究真正关心的是比较不同条件，而不只是看单一条件。这里比较：左耳 vs. 右耳听觉点击。

evoked_aud_right = epochs["auditory/right"].average()

mne.viz.plot_compare_evokeds(
    dict(left=evoked_aud_left, right=evoked_aud_right),
    picks="eeg",
    combine="mean",
)
plt.show()

png

关于统计检验的说明¶

用肉眼比较两条条件曲线只是个起点，不是结论 —— 真正的 ERP 研究会检验条件之间的差异是否在统计上可靠，通常采用基于聚类的置换检验（cluster-based permutation test）（mne.stats.permutation_cluster_test 及相关函数），它能妥善处理"你同时在检验许多个相互关联的时间点/通道，而不是单一一次比较"这一事实。这超出了本教程的范围，但值得记住这个名词，等你准备好从上面这样的比较中得出真正的结论时会用到。

小结¶

epochs.average() 产生一个 Evoked —— 也就是你的 ERP。最重要的核心思想是：平均会抵消随机噪声，而稳定的、锁时的反应则会留存下来。
用 evoked.plot() / evoked.plot_joint() 做可视化；用 evoked.get_peak() 量化一个成分。
用 mne.viz.plot_compare_evokeds() 比较条件；正式的比较需要基于聚类的置换检验统计。

下一章：第 8 章 —— 应用到你自己的数据