MNE 溯源实战：最小范数估计
在无穷多个能解释数据的解里，挑"总能量最小"的那个——这就是 MNE。 这一课从直觉走到完整代码，顺便认识它最出名的毛病：偏爱浅部源。
🎯 本节学习目标
- 用自己的话讲清最小范数估计：在无穷多解中挑"总能量最小"的那个
- 理解深度偏置的来源，掌握 lambda 正则化在"拟合数据"与"抑制噪声"之间的折中
- 独立完成 FieldTrip MNE 实战：头模型 → 前置场 →
ft_sourceanalysis - 读懂 source 结构体的关键字段：
pos/trial·avg/mom/pow（激活指数） - 会用
ft_sourceinterpolate+ft_sourceplot把源激活画到 MRI 三视图上 - 理解组级分析为何必须做 MNI 空间归一化，了解源空间配对置换检验的思路
一、最小范数估计：在无穷多解里挑"最省能量"的
MNE（Minimum Norm Estimate，最小范数估计）是分布式溯源的"入门标准"。 它对逆问题的回答朴素得可爱：在所有能解释头皮数据的解里，选总能量（平方和）最小的那个。
借用 7.1 的小例子：电极测到 [3, 3]，候选解 [3,0,0] 与 [1,1,1] 都能解释数据。
它们的总能量（平方和）分别是 9 和 3——MNE 选后者：宁要三个温和的源，不要一个凶猛的源。
这就是"最小范数"四个字的全部含义。
范数是衡量一个向量"有多大"的尺子；L2 范数 = 各元素平方和再开方，对应"总能量"。MNE 用它惩罚"太猛"的解——解越猛罚得越重，最优解自然倾向于温和、分散。
MNE 的解有漂亮的闭式表达（了解即可，不必推导）：
J = Kᵀ (K Kᵀ + λ²C)⁻¹ Y。逐个看角色：
Y 是头皮数据（你要解释的对象）；
K 是把前置场排成的大矩阵（"传播字典"，7.1 讲过）；
C 是噪声协方差（描述噪声长什么样）；
λ 是折中旋钮（下一节讲）。
括号里的 + λ²C 就是正则化——它让求逆这件事从"病态的悬崖"变成"有护栏的缓坡"。
因为解对数据是纯线性的、每个时间点独立计算，MNE 天然适合 ERP 的每个时刻。
二、深度偏置与 lambda 正则化：MNE 的两个内置难题
难题一：深度偏置（欺负"深沉"的源）。
想象麦克风前的两个人：离麦近的人小声说就够响，离麦远的人必须喊。 大脑里也一样：深部源（如扣带回、岛叶）要在头皮上产生同样大小的电位，需要大得多的电流强度。 而 MNE 偏爱"省能量"的解——于是深部活动被系统性低估，浅部（贴近头皮的沟回）活动被高估。 看结果图时皮层表面总是亮、深部总是暗，一部分原因就在这里，不全是生物学。
难题二：噪声放大。逆问题"解不随数据连续变化"，直接对 K 求逆会把噪声放大得离谱。解法是正则化：
在"紧贴数据"与"保持解温和"之间加一个折中系数 λ。λ 越小，解越贴合数据（噪声放大越严重）；λ 越大，解越平滑保守（细节被抹掉）。FieldTrip 里既可给数值，也可用 cfg.mne.snr 按信噪比自动换算——SNR 越高，正则化越少。
lambda 通常在几个数量级内试探；主结果务必做稳健性检查（换一档 λ，结论还在不在），并在方法部分报告取值。深度偏置可用深度加权缓解，但会引入新假设——两难之下，把参数报清楚比追求完美更实际。
三、实战（一）：数据、头模型与前置场
假设你已有一份按模块 2 预处理、模块 3 分好段的数据 data（FieldTrip 格式的 epochs）。三步搭好正问题：
% ---- 第 0 步：锁时平均 + 噪声协方差（MNE 的正则化要用它）----
cfg = [];
cfg.trials = 'all';
cfg.channel = 'EEG';
cfg.covariance = 'yes'; % 估计噪声协方差
cfg.covariancewindow = [-inf 0]; % 用刺激前基线估计"纯噪声"
timelock = ft_timelockanalysis(cfg, data); % timelock.cov 即协方差矩阵
% ---- 第 1 步：头模型（教学用同心球；正式研究换个体/模板 BEM）----
cfg = [];
cfg.method = 'concentricspheres';
headmodel = ft_prepare_headmodel(cfg, timelock.elec);
% ---- 第 2 步：网格 + 前置场（与数据无关，算一次存一次）----
cfg = [];
cfg.resolution = 1; % 1cm 网格：数千个候选偶极子
cfg.headmodel = headmodel;
cfg.elec = timelock.elec;
cfg.channel = timelock.label;
leadfield = ft_prepare_leadfield(cfg); % 每个网格点的"传播字典"
save('leadfield_sphere.mat', 'leadfield', 'headmodel');四、实战（二）：MNE 求解与结果初检
% ---- 第 3 步：最小范数求解 ----
cfg = [];
cfg.method = 'mne';
cfg.sourcemodel = leadfield; % 旧教程写作 cfg.grid，二者等价
cfg.headmodel = headmodel;
cfg.elec = timelock.elec;
cfg.mne.prewhiten = 'yes'; % 用噪声协方差预白化：数值上更稳定
cfg.mne.scalesourcecov = 'yes'; % 源协方差归一，让 lambda 量纲可解释
cfg.mne.lambda = 3; % 正则化强度（官方教程值；见第二节）
source = ft_sourceanalysis(cfg, timelock); % 噪声协方差从 timelock.cov 自动读取
% ---- 第 4 步：结果初检——什么时候、哪里最亮 ----
pow_mean = mean(source.avg.pow, 1); % 每个时间点的全脑平均激活
[~, t_peak] = max(pow_mean); % 最强的时间点
[~, v_peak] = max(source.avg.pow(:, t_peak)); % 该时刻最强的体素编号
fprintf('峰值时刻 %.0f ms，坐标 %s\n', ...
timelock.time(t_peak)*1000, mat2str(source.pos(v_peak, :)));五、读懂输出：source 结构体解剖
跑完之后，source 是一个大结构体。逐字段认识它：
| 字段 | 尺寸/类型 | 含义 |
|---|---|---|
source.pos | 体素数 × 3 | 每个网格点的三维坐标（单位随头模型，常为 mm） |
source.dim | 1 × 3 | 规则网格的维度（用于把向量重排成三维图） |
source.avg / source.trial | 结构 | 平均结果；cfg.keeptrials='yes' 时还有分试次的 trial |
source.avg.mom | 体素数 × 3 × 时间 | 每个体素在三个朝向上的源时间过程（单位 A·m） |
source.avg.pow | 体素数 × 时间 | 激活指数：三个朝向平方和合成一个数（老教程里记作 ai），画图统计都用它 |
source.avg.filter | （keepfilter 时） | 线性逆算子本身，供复用（7.3 的虚拟电极要靠它） |
MNE 结果里皮层表面亮、深部暗是方法天性。报告"某深部结构无激活"是高风险结论，需要方法学证据（如深度加权、与 sLORETA 标准化结果交叉验证）支撑，而不是只看一张 MNE 图。
六、可视化：ft_sourceinterpolate 与 ft_sourceplot
裸网格上的数字不直观，标准做法是插值到 MRI 再切三视图：
% ---- 第 5 步：插值到 MRI + 正交三视图 ----
mri = ft_read_mri('subject01.mri'); % 个体 MRI；没有就用 FieldTrip 模板
cfg = [];
cfg.parameter = 'avg.pow';
interp = ft_sourceinterpolate(cfg, source, mri); % 'avg.' 前缀被剥掉 → interp.pow
cfg = [];
cfg.method = 'ortho'; % 矢状/冠状/横断三视图
cfg.funparameter = 'pow';
cfg.location = source.pos(v_peak, :); % 三条线的交点放在峰值体素上
ft_sourceplot(cfg, interp);
% 想看激活随时间"放电影"：
% cfg = []; cfg.funparameter = 'pow'; ft_sourcemovie(cfg, interp);七、从单被试到组级：MNI 归一化与统计思路
为什么必须归一化：每个人的头形状、脑沟回都不一样，被试 A 的"额上回"和被试 B 的"额上回"
在各自的坐标里根本不在同一个位置。先把所有被试的源图配到同一个标准空间（MNI），才能跨被试平均与统计。
实操上：要么把个体 MRI 用 ft_volumenormalise 归一化到 MNI 后再做插值；
要么直接使用 MNI 空间的模板网格（FieldTrip 自带 standard_sourcemodel3d6mm 这类模板源模型），一步到位省去配准。
单被试内比较两个条件，官方教程的做法是"手工作坊版"：
% ---- 两条件的差值图（官方 MNE 教程的做法）----
cfg = [];
cfg.projectmom = 'yes';
sdA = ft_sourcedescriptives(cfg, sourceA); % sourceA / sourceB：两条件的 MNE 结果
sdB = ft_sourcedescriptives(cfg, sourceB);
sdDIFF = sdA;
sdDIFF.avg.pow = sdA.avg.pow - sdB.avg.pow; % 差值图：哪里条件 A 比条件 B 强
% ---- 组级预告：N 个被试的源图（同一 MNI 网格）进入置换检验 ----
cfg = [];
cfg.method = 'montecarlo'; % 非参数置换（模块 3 的老朋友）
cfg.statistic = 'ft_statfun_depsamplesT'; % 被试内配对 t
cfg.parameter = 'pow'; % 统计对象：每体素的激活
cfg.correctm = 'cluster'; % cluster 校正控制多重比较
cfg.numrandomization = 1000;
% cfg.design = ... % 设计矩阵（写法同模块 3 统计课）
% stat = ft_sourcestatistics(cfg, condA{:}, condB{:}); % condA/condB 为被试列表
% 7.4 会从"统计决策"的角度再回到这个问题一套可复用的正向模型（headmodel + leadfield，存盘）＋ 一张插值到 MRI 的 MNE 激活图。前者换方法不换（7.3 直接复用），后者是你第一条"从头皮走进大脑"的结果。
📌 本节小结
- MNE = 奥卡姆剃刀：在所有能解释数据的解里挑总能量最小的那个，天然偏爱"温和而分散"
- 深度偏置是 MNE 的内置缺陷：深部源被系统性低估——图上皮层表面总亮，不全是生物学
- 正则化 lambda 在噪声放大与细节丢失之间折中：正式分析要报告取值并检验稳健性
- source 结构体：mom 是三朝向时间过程，pow（老教程称 ai）是合成后的激活指数，插值后用于画图与统计
- 组级源分析的前提：所有被试的源图先归一化到同一个 MNI 空间网格
✏️ 课后练习
- 用一份 ERP 数据跑通本课完整流水线（同心球 + 1cm 网格 + MNE），保存
leadfield与source，报告峰值体素坐标与峰值时刻。 - 把 cfg.mne.lambda 换成明显更小与更大的两个值分别重跑，对比激活图的变化，写 3 句观察结论。
- 用
ft_sourceinterpolate把结果插到模板 MRI，分别以峰值体素与两个解剖标志点为中心各出一张正交三视图。