脑电功率谱与多尺度熵分析的MATLAB实现：从Welch法到MSE特征提取

简介：面向脑电信号分析和生物医学工程研究者的多尺度熵（MSE）与功率谱分析工具包，由经典C代码重写为MATLAB实现，支持从Brain-atlas III等设备读取脑电数据，完成多尺度熵计算、粗粒化、样本熵估计、频谱分析与批次处理。工具包针对20通道、128Hz采样率的脑电数据设计，输出格式为包含20个通道电压的矩阵，符合常见脑电实验的数据结构。压缩包共40个文件，其中26个M脚本/函数构成主体，覆盖数据加载、滤波、熵计算、批处理和可视化等多个分析环节；另有7个文本说明文档、2个C辅助程序以及电极位置等配套文件，整体仅23KB，目录紧凑、便于二次开发与算法迁移。该资源已有342人学习或浏览，具有较好的参考价值。主脚本可一键串联从原始数据读取到熵值导出的完整流程，批处理模块支持多历元、多通道、多尺度因子组合，同时提供多种谱估计与拓扑绘图工具。研究者可直接运行示例复现结果，也可按需修改参数、替换数据接口，快速嵌入自己的脑电分析管线。

1. 脑电功率谱与mse-analysis：先把这个“质谱”的误解说清楚

先说一个很多人刚看到标题时都会踩的坑：“质谱分析”这四个字，放在脑电功率谱的语境里，基本可以判定是误标。脑电信号是头皮上记录到的电压随时间变化的序列，跟质谱仪里那种按质荷比分离离子的分析完全不是一回事。这里真正要处理的是 mse-analysis 的另一种常见展开：Multiscale Entropy Analysis，也就是多尺度熵分析。它与脑电功率谱放在一起，正好构成一套经典的脑电特征组合——功率谱描述信号在频域上的能量分布，多尺度熵描述信号在不同时间尺度上的复杂度与规律性。前者是经典的线性分析，后者是非线性动力学方法，两者互补，是睡眠分期、疲劳识别、情绪分类、麻醉深度监测这些方向的常见特征来源。这篇文章就把这两块分别拆开，给出能直接用在 MATLAB 里的代码、参数配置和踩坑记录，适合手里已经有一份脑电数据、想马上算出可用特征的从业者。

2. 脑电功率谱怎么算：Welch法在MATLAB里的最小实现

2.1 为什么用Welch不用周期图：方差与平滑性的取舍

脑电功率谱估计最朴素的做法是周期图：把整段数据做傅里叶变换，取模的平方得到功率谱。这方法理论最直接，但有个硬伤——估计出的功率值方差非常大，信号里稍微有点噪声或突发干扰，频谱曲线就会剧烈抖动，很难作为稳定特征使用。对于脑电这种强非平稳信号，直接做整段 FFT 基本只适合看个大概。

Welch 法的思路是先把数据分成有重叠的若干段，每段加窗后分别做周期图，再把所有段的功率谱平均。这样做付出的代价是频率分辨率略有下降，但换来了方差显著减小，曲线平滑得多，特征稳定性也更好。脑电分析里 Welch 法几乎是默认选择，几乎没有实际项目会直接用原始周期图。

脑电的采样率和记录时长决定了 Welch 参数怎么设。这里有个基本的折中逻辑：窗长越长，频率分辨率越高，但脑电在长时间窗内已经不平稳，把 10 秒清醒脑电强行当作平稳信号处理，α节律和慢波混在一起，反而抹平了关键频带差异。我一般遵循的默认是窗长 2 秒，这个值对大多数静息态任务态数据都够用。

频率分辨率由窗长决定，公式是 ( \Delta f = 1 / T_{win} )。2 秒窗对应 0.5 Hz 分辨率，能清楚分开 δ、θ、α、β 这些常用频带。窗长定 1 秒也能做，但 1 Hz 的分辨率会把 α 带（8~13 Hz）的边缘频率算得比较粗糙。后面的代码都以 2 秒窗为基础。

2.2 核心参数pwelch代码：窗长、重叠率与FFT点数

MATLAB 里 Welch 法已经封装在 pwelch 函数里，不需要自己写分段加窗的循环。下面是单通道脑电功率谱估计的最小实现：

% eeg_psd_welch.m % 输入: eeg_data 为1xN的单通道脑电序列，fs为采样率 % 输出: f为频率轴，psd_norm为归一化后的相对功率谱 function [f, psd_norm] = eeg_psd_welch(eeg_data, fs) win_len = 2 * fs; % 窗长2秒 noverlap = win_len / 2; % 50%重叠 nfft = 2^nextpow2(win_len); % FFT点数取2的幂，加速计算 eeg_data = eeg_data - mean(eeg_data); % 去直流分量 [psd, f] = pwelch(eeg_data, hann(win_len), noverlap, nfft, fs); psd_norm = psd / sum(psd); % 转为相对功率，便于跨被试比较 end

逻辑说明：去均值这一步很关键，直流分量不除掉，0 Hz 附近会出现一个巨大的能量峰，连带着污染 δ 带。hann 窗在脑电分析里最常见，主瓣较宽但旁瓣衰减好；换成 hamming 窗主瓣稍窄，旁瓣衰减差一些；blackman 窗旁瓣更低但主瓣更宽，频率选择性变差。pwelch 中 nfft 不小于窗长即可，取 2 的幂是为了配合 FFT 加速，不会改变结果本身的频率分辨率。

重叠率这块，50% 是脑电分析里用得最多的；重叠率提高到 75% 曲线更平滑，但相邻段高度相关，计算量上去了，平滑收益却不大。我这里固定用 50%。如果数据本身很短，比如只有 10 秒，2 秒窗、50% 重叠只能分到 9 段，平均效果有限，可以适当把窗长降到 1 秒换取更多段数。

2.3 从频谱到频带特征：相对功率计算的参数边界

功率谱曲线本身很少直接作为特征，一般要压缩成频带特征。脑电研究里最常用的频带划分如下：

频带拆分用逻辑索引从频谱里取值求和即可：

% band_power.m % 输入归一化功率谱psd_norm和频率轴f，输出各频带相对功率 function [delta, theta, alpha, beta] = band_power(psd_norm, f) delta = sum(psd_norm(f >= 0.5 & f < 4)); theta = sum(psd_norm(f >= 4 & f < 8)); alpha = sum(psd_norm(f >= 8 & f < 13)); beta = sum(psd_norm(f >= 13 & f < 30)); end

这里我用的是相对功率而非绝对功率。绝对功率的单位是 μV²/Hz，受颅骨厚度、电极阻抗、放大器增益影响很大，同一个被试今天和明天测出的绝对功率都可能差一倍以上。相对功率做了归一化，反映的是能量占比，跨被试、跨天比较稳定得多。

真正落地时有个细节常被忽略：psd_norm 是对全频段归一化，还是对 0.5~30 Hz 归一化？两者结果差别很大。全频段归一化会把 30 Hz 以上的肌电噪声也计入分母，导致 α、β 相对功率整体被压低。我一般先做 0.5~30 Hz 的截断，再在这个范围内归一化。如果目标特征包含 γ 带，截断上限相应上调到 45 Hz 左右。

3. mse-analysis的MATLAB实现：多尺度熵算法拆解与完整代码

3.1 样本熵之外为什么要加“多尺度”：脑电复杂度的时间尺度

功率谱是线性分析，它假设信号由不同频率的正弦波叠加而成。但脑电本质上是非线性动力系统的输出，包含大量非线性耦合成分。两个功率谱完全一致的信号，时域复杂度可以完全不同；比如一段纯正弦波和一段带有随机调制、间歇性节律的脑电，频谱上可能都有 α 峰，但一个是高度规律的振荡，另一个在时间结构上复杂得多。多尺度熵就是用来量化这种时间结构复杂度的。

多尺度熵建立在样本熵之上。样本熵衡量的是时间序列中产生新模式的可能性大小：序列越规律，匹配模板向量的数量越多，熵值越低；序列越随机，匹配越少，熵值越高。但样本熵只在单一时间尺度上计算，而脑电节律分布在从毫秒到秒级的多个时间尺度上。多尺度熵的做法是先把原始序列做粗粒化——把每 s 个连续点取平均，得到尺度为 s 的粗粒化序列——再对每个尺度的序列分别计算样本熵。

这就解决了一个具体问题：单尺度样本熵会把白噪声与包含多尺度结构的复杂信号混为一谈，因为两者在细粒度上都显得“不规则”。粗粒化之后，真正存在跨尺度结构的信号会在某些尺度上表现出比较高的熵值或特定的熵变化曲线，而纯随机信号在各个尺度上的熵会持续下降。这条曲线的形态，比单个熵值更有判别力。

3.2 可直接复用的mse_analysis函数：代码与逐段说明

下面这份代码实现了从原始脑电序列到多尺度熵曲线的完整流程，包含粗粒化和样本熵两个环节：

% mse_analysis.m % 输入: eeg_data为1xN的脑电序列，S_max为最大尺度因子， % m为嵌入维度，r为相似容限系数(乘以原始序列标准差得到绝对容限) % 输出: mse_values为1xS_max的多尺度熵值，scales为对应尺度 function [mse_values, scales] = mse_analysis(eeg_data, S_max, m, r) N = length(eeg_data); r_abs = r * std(eeg_data); % 固定容限，以原始序列标准差为基准 mse_values = zeros(1, S_max); for s = 1:S_max n = floor(N / s); % 粗粒化后的序列长度 % 粗粒化: 每s个连续点取平均 coarse = mean(reshape(eeg_data(1:n*s), s, n), 1); mse_values(s) = sample_entropy(coarse, m, r_abs); end scales = 1:S_max; end % sample_entropy.m % 计算一维序列的样本熵 function se = sample_entropy(data, m, r) N = length(data); B = 0; % m维匹配对数 A = 0; % m+1维匹配对数 for i = 1:N-m for j = i+1:N-m % 切比雪夫距离: 两向量对应位置差值的最大值 d_m = max(abs(data(i:i+m-1) - data(j:j+m-1))); if d_m < r B = B + 1; % 在m维匹配基础上检查m+1维是否也匹配 d_m1 = max(abs(data(i:i+m) - data(j:j+m))); if d_m1 < r A = A + 1; end end end end if B == 0 || A == 0 se = NaN; % 匹配数为0时熵值无定义，置为NaN else se = -log(A / B); end end

逻辑说明：粗粒化的 reshape 写法比显式循环快得多。reshape(eeg_data(1:n*s), s, n) 把序列重排成 s 行 n 列的矩阵，每列恰好对应一个粗粒化组，mean(..., 1) 沿行方向取平均，得到 1×n 的粗粒化序列。这一步是 mse-analysis 的核心，尺度 s 越大，粗粒化序列越短，熵值对数据长度的敏感度越高。

样本熵部分，切比雪夫距离是标准选择，它对所有位置的差异一视同仁，取最大值作为两个模板向量的距离。r 是相似容限：距离小于 r 就认为两个模板匹配。A 和 B 的比值反映的是“当两个 m 维向量匹配时，它们的 m+1 扩展也匹配的概率”。这个概率越小，说明序列越难预测，熵值越大。

参数说明：m 通常取 2，这是多尺度熵研究里最主流的设置。取 3 对模板匹配数要求更苛刻，需要更长数据才能得到稳定估计。r 通常取 0.10~0.25，代码里 r_abs 是用原始序列标准差计算的绝对容限，也就是说粗粒化后的所有尺度都共用同一个容限值。这保证了不同尺度之间的熵值可比。如果对每个尺度重新计算标准差再乘 r，各尺度容限会随尺度变化，曲线形态会改变，跨尺度比较的意义就打折扣。

3.3 尺度上限和容限r怎么选：常见参数配置与合理性检查

尺度上限 S_max 的选择首要取决于数据长度。粗粒化长度为 floor(N/s)，样本熵对短序列很敏感，经验上粗粒化后的序列长度最好不低于 200 点。所以 S_max ≈ N/200。对于 10 秒、1000 Hz 采样率的数据，N = 10000，S_max 最多取到 50；但实际项目里 30 是更稳妥的上限，因为尺度 30 时粗粒化序列只有 333 点，已经接近样本熵的估计下限。

r 值的选择需要结合具体数据做敏感性检查。我常用的做法是同时跑 r = 0.10、0.15、0.20 三组，看曲线形态是否一致。如果三组曲线趋势相同、只是整体高低不同，说明熵值估计稳定；如果曲线形状随 r 变化很大，比如在某个尺度上出现交叉，说明数据里可能有伪迹或非平稳成分影响了匹配统计，需要先回预处理阶段排查。

另一个检查方法是看熵曲线末尾的走势。健康脑电的多尺度熵曲线通常在中等尺度区间出现峰值，随后缓慢下降。如果曲线从尺度 2 开始就持续走低且没有回升，常见原因有两个：一是数据本身太短，粗粒化序列在尺度大了之后全是均值，信息被抹平；二是数据包含太多伪迹片段，模板匹配被大振幅噪声主导。

4. 功率谱与多尺度熵联合分析：特征矩阵构建与参数配置

4.1 联合特征的价值：频域能量与非线性复杂度的信息互补

功率谱和多尺度熵看的是信号的不同侧面。功率谱把信号拆成频率成分，描述的是“哪些频率上有能量”；多尺度熵描述的是“时间序列的规律性如何随尺度变化”。这两类特征经常出现互补的情况：一段处于清醒闭眼状态的脑电，α 波能量占比很高，功率谱特征一眼就能识别；但两段 α 波能量几乎相同的脑电，一段是稳定的正弦样节律，另一段伴随明显的调幅波动，功率谱上看不出差别，多尺度熵却能区分。

实际项目里，我会把两类特征拼成一个特征向量喂给分类器。比如做睡眠分期，单用功率谱特征准确率在小样本上容易过拟合，加入多尺度熵之后鲁棒性明显提升。这里的关键是不要只把熵曲线上的一个点当特征——整条曲线在不同尺度上的形态都有判别信息，拿尺度 1~10 的熵值全拼进去通常比只取尺度 3 的熵值效果好。

4.2 批量特征提取脚本：通道×试次的数据组织与输出格式

实际处理脑电时很少只分析一个通道、一段数据。常见的数据组织方式是通道数×采样点×试次数三维矩阵。下面的脚本把功率谱频带特征和每个尺度的熵值组装成一张特征表，每一行对应一个（试次，通道）组合：

% batch_feature_extract.m % 输入: data为 [通道数 x 采样点数 x 试次数] 三维矩阵 % 输出: features为 (通道数*试次数) x (4+S_max) 的特征矩阵 % feature_names为对应的特征名字符串元胞数组 function [features, feature_names] = batch_feature_extract(data, fs, S_max, m, r) [n_ch, n_pts, n_trials] = size(data); n_features = 4 + S_max; % 4个频带 + S_max个尺度熵值 features = zeros(n_ch * n_trials, n_features); idx = 1; for trial = 1:n_trials for ch = 1:n_ch x = squeeze(data(ch, :, trial)); x = x(:)'; % 转成行向量 [f, psd_norm] = eeg_psd_welch(x, fs); [delta, theta, alpha, beta] = band_power(psd_norm, f); [mse_val, ~] = mse_analysis(x, S_max, m, r); features(idx, :) = [delta, theta, alpha, beta, mse_val]; idx = idx + 1; end end feature_names = [{'delta','theta','alpha','beta'}, ... arrayfun(@(s) sprintf('mse_%d', s), 1:S_max, 'UniformOutput', false)]; end

逻辑说明：脚本里 squeeze 取单通道数据后转成行向量，保证与函数接口一致。每个通道每段数据独立算功率谱和熵值，特征矩阵里行顺序是“试次内循环通道”，也就是先取出试次 1 的通道 1、通道 2……，再取试次 2。特征矩阵可以直接喂给 fitcsvm、fitcknn 这类分类器。feature_names 同时返回，方便训练后看特征重要性。

需要提醒的是，mse_analysis 的计算复杂度是 O(N²/S²)，批量跑几十个通道乘几十个试次时，纯 MATLAB 循环会明显变慢。如果脚本跑一个数据集要两个小时以上，优先缩短 S_max 到 20，或者把数据降采样到 500 Hz，而不是急着改算法。

4.3 推荐参数表：不同采样率下的窗长、尺度上限与容限

不同采样率下各参数需要联动调整。下面的表是我在实际项目里验证过的一套起始配置，兼顾了频率分辨率和计算时间：

这套配置的核心逻辑是：功率谱窗长固定为 2 秒，保证频带分辨率一致；S_max 根据采样率和数据时长调整，使最大尺度下的粗粒化序列不低于 200 点；R 值固定为 0.15 倍原始标准差，保证跨尺度可比。

在低采样率设备上，比如便携式脑电的 128 Hz 采样率，2 秒窗只有 256 点，频率分辨率足够但平滑度差，这种情况下窗建议加到 4 秒；S_max 建议小于 10。这里没有万能的固定参数，拿到数据后先跑一小段看特征分布，再根据分布调整，比硬套参数要可靠得多。

5. 避坑记录：脑电功率谱与多尺度熵分析最容易翻车的六件事

5.1 低频段能量虚高：基线漂移混入了功率谱

现象：功率谱在 0.5~2 Hz 区间能量异常高，α 峰几乎看不见，整条曲线呈现“左高右低”的单边下降形状。

原因：脑电记录里电极与头皮之间的慢电位漂移会形成低频大振幅成分，这段漂移在功率谱里表现为集中在极低频的巨大能量。如果信号里还夹杂了缓慢的呼吸伪迹，低频能量会进一步膨胀，把本来微弱的 δ 带信号完全盖住。

解决：在计算功率谱之前先做高通滤波，我一般用 0.5 Hz 或 1 Hz 截止频率的一阶巴特沃斯滤波器。滤波顺序要在去均值之后、pwelch 之前。另一种补救是用 detrend 函数去掉整段数据的线性趋势，但只能对付缓慢坡形漂移，对呼吸伪迹作用有限。做完高通后再看一眼低频段的功率值，正常静息态下 0.5~1 Hz 的绝对功率不应该比 α 带高出几个数量级。

5.2 尺度大了MSE全趋同：粗粒化序列长度不够

现象：多尺度熵曲线在尺度 1~5 时各组差异明显，尺度超过 15 之后所有被试的曲线都挤到一起，熵值稳定在某个区间不再有区分度。

原因：S_max 选得太大，粗粒化后的序列长度太短，样本熵的模板匹配统计不可靠。长度只有 100 点的序列，m=2 时能配对的模板对本来就少，熵值波动剧烈且普遍偏低，所有数据看起来都“随机且不可预测”，可区分性自然消失。

解决：把 S_max 限制在 N/200 以内。打个比方，1000 Hz 采样、20 秒的数据共 20000 点，S_max 取 20 时粗粒化后还有 1000 点，稳定；取 40 时就只剩 500 点，已经偏少。另一个技巧是看熵曲线在最大尺度处的方差：如果最后 3 个尺度的标准差不随尺度增大而变大，说明还能往上取；如果方差急剧膨胀，说明尺度过头了。

5.3 50Hz出现尖峰：工频干扰与陷波滤波的时机

现象：功率谱在 50 Hz（或 60 Hz）处有一个尖锐陡峭的窄峰，高度明显高于两侧。国内市电是 50 Hz，国外部分设备可能是 60 Hz。

原因：脑电放大器虽然做了一定屏蔽，但市电通过电源线、照明设备耦合进记录系统的现象非常普遍，尤其是电极线靠近电源线的时候。这个尖峰如果不去掉，凡是要用到 γ 带或全频段归一化的算法都会被干扰。

解决：用陷波滤波器（notch filter）在 50 Hz 处做窄带衰减。MATLAB 里可以用 iirnotch 设计陷波器，Q 值取 30 左右比较合适，太窄了抑制不干净，太宽了会伤及 45~55 Hz 的真实脑电成分。但要注意，陷波滤波器不是随便加的——如果分析目标仅在 0.5~30 Hz，直接在功率谱计算后截断频率轴，根本不用陷波；只有需要 40 Hz 以上的 γ 带时才必须处理工频。加上陷波之后要检查 50 Hz 附近有没有出现相位畸变的振铃伪迹，这在小样本分类里很容易被误当作真实信号。

5.4 两次测量结果对不上：参考电极与r值标准化不一致

现象：同一被试隔一天测的静息态脑电，α 波相对功率变化不大，但多尺度熵曲线整体偏移，分类器训练集和测试集特征分布明显错位。

原因：两个隐藏变量在捣乱。一是参考电极位置不一致，从乳突参考换成平均参考之后，所有通道的幅值尺度都会变化，r_abs 是用标准差算的，跟着一起变，熵值自然偏移。二是数据标准化方式不一致，有的代码对原始信号做了 z-score，有的没有，这直接改变 r_abs 的基准。

解决：固定整套处理链路。参考方式固定为平均参考或同一参考，记录测前阻抗值；r 值的计算基准在整批数据里统一，不要对单个被试单独算 r。如果要做跨被试可比的多尺度熵，建议先把每个通道的数据做 z-score 标准化，再用全局标准差乘以 r 系数作为容限。这样能显著减弱个体差异对熵值的影响。

5.5 批量分析慢到没法用：样本熵循环的优化方向

现象：mse_analysis 跑单段 30 秒数据只要几秒，但放到 32 通道 × 40 试次的数据集上，跑一个通宵还没结束。

原因：样本熵的双重循环是 O(N²) 复杂度，嵌套在外层尺度循环里。N=30000 时，内层循环要比较约 4.5 亿对模板，纯 MATLAB 解释执行根本跑不动。

解决：基本思路是向量化内层距离计算。把 m 维模板向量排列成矩阵后，用 pdist2 一次性计算所有模板对的切比雪夫距离，再统计小于 r 的数量。内存允许时这种方法能提速数十倍。另一个是降低数据量——把 1000 Hz 数据降采样到 500 Hz，粗粒化的计算量直接降到四分之一。如果数据量实在太大，就把尺度循环改成 parfor 并行，但要注意随机数流和输出顺序的处理。最实用的一招是先跑一小段验证参数，确认熵曲线形态没问题再全量计算，别一上来就对着整个数据集硬跑。

5.6 特征全被眼电伪迹主导：预处理没做就提特征

现象：特征矩阵里前几个主成分载荷几乎全落在某个通道上，分类器学到的不是脑电模式，而是眨眼波形。功率谱在 δ 带或 θ 带能量异常高。

原因：眨眼和眼动产生的眼电幅度通常比脑电大一个数量级，混入数据后功率谱低频段被污染，样本熵也会被大幅不规则的瞬态脉冲主导。没有做伪迹去除就直接提特征，特征质量完全取决于伪迹多寡。

解决：提取特征前先做伪迹剔除。幅度超过 ±100 μV 的片段直接剔除是最快的办法；眨眼密集的数据还需要用 ICA 去掉眼电分量。这里有个容易被忽略的细节：伪迹剔除后要重新检查各通道方差是否稳定，如果某个通道剔除的片段占比超过 30%，这个通道的特征可信度就很低了，不建议再进入特征矩阵。

6. 一个验证技巧：用替代数据检验MSE结果有没有统计意义

6.1 替代数据检验解决什么问题

多尺度熵有一个让人心虚的地方：它把信号的“不可预测程度”解释为复杂度，但有些熵值差异可能只是线性相关结构造成的，并非真正的非线性动力学特征。替代数据检验就是来兜这个底的。它的思路是：把原始数据做傅里叶变换，随机打乱相位再逆变换回去，得到一组与原始数据幅值谱完全相同、但相位结构被破坏的替代序列。如果原始序列的多尺度熵显著不同于替代数据，说明该序列确实存在超越线性结构的非线性成分；如果两者差不多，那么你测到的高熵值很可能只是线性滤波或者幅值谱特征带来的幻觉。

6.2 相位随机化的MATLAB实现与判读

相位随机化需要保证逆变换后还是实信号，对实数序列而言，FFT 结果具有共轭对称性，随机化相位时必须保持这种对称性，否则 ifft 出来的是复数——很多第一次写这个函数的人都栽在这里。

% phase_randomized_surrogate.m % 相位随机化替代数据: 保持幅值谱，打乱相位结构 function surr = phase_randomized_surrogate(data) N = length(data); F = fft(data); half = floor(N/2); phases = 2 * pi * rand(1, half-1); % 随机相位 F_surr = F; F_surr(2:half) = F(2:half) .* exp(1i * phases); % 修改正频率相位 F_surr(N:-1:N-half+2) = conj(F_surr(2:half)); % 保持共轭对称 surr = real(ifft(F_surr)); end

判读方法：生成 30~50 组替代数据，分别计算多尺度熵，得到每个尺度上熵值的零分布。如果原始序列的熵值落在替代数据分布的 95% 区间之外，说明该尺度上的复杂度不是线性结构能解释的。我自己的习惯是保留每个尺度的判定结果，单独报告“尺度 3、5、7 显著”，比笼统说“该被试多尺度熵显著”更严谨。

这套检验对论文和课题汇报都很有价值。但要注意，相位随机化对均匀采样且无缺失段的脑电数据才有效，如果数据经过重采样引入了人工平滑，检验的灵敏度会下降。替代数据检验本质上是一个“证伪”工具，它不能证明一段信号一定蕴含非线性结构，只能帮你排除“其实只是线性信号”的平庸解释。我一般把它放在特征提取流程的最后一步，每个被试只跑一次，确认关键尺度上显著之后，再放心地把特征交给下一步建模。希望这个思路能帮到你。

本文还有配套的精品资源，点击获取