EEMD与样本熵:非平稳信号自适应分解与重构实战

发布时间:2026/9/26 17:12:10
EEMD与样本熵:非平稳信号自适应分解与重构实战 拿到一段轴承振动信号或者说一段生理信号很多人的第一反应是直接看波形、看频谱然后滤波、取特征。但如果是非平稳、非线性的数据比如设备启停阶段、变转速工况或者肌电、脑电这类本身就“不乖”的信号傅里叶那套全局频域分析就不太够用了。这时候我习惯的做法是先用EEMD做自适应分解把原始序列拆成一组从高频到低频的固有模态函数然后算每个IMF的样本熵用样本熵大小去衡量每个分量的复杂度和信息含量最后按熵值把IMF归并成低频、中频、高频三组信号重新叠加。这条链路在故障诊断、信号去噪、趋势提取里都非常实用这篇就完整记录一下我的实践过程、参数设置逻辑和踩过的坑给准备做类似工作的同学一个可复现的参考。1. 为什么是EEMD加样本熵这套组合1.1 单一工具不够用EEMD相比EMD、小波和FFT的优势先聊最基本的问题做信号分解EMD、小波、FFT都能干为什么最后我选EEMD当主角如果直接对信号做FFT得到的是一整段时间内的平均频谱完全丢失了瞬时频率信息。信号里如果有一个瞬间出现的冲击或者一个慢慢消失的分量FFT只能告诉你“它在某个频段存在过”但说不清楚它是什么时候出现的。小波变换有窗函数能提供时频定位但前提是你得先选好基函数和分解层数基函数选得不合适分析结果就带上人为偏差。对一个未知特性的信号这一步往往是拍脑袋。EMD最大的优势是数据自适应它根据信号本身的极值包络去逐级抽取固有模态函数不需要预设基函数。但EMD有个老毛病叫模态混叠当一个IMF里同时混着不同时间尺度的成分时分解结果会出现“你中有我、我中有你”的情况直接拿去算样本熵或者做特征提取都会出问题。EEMD的解决思路很朴素——把白噪声作为辅助信号加入原数据利用白噪声频谱均匀分布、在各尺度上都有投影的特性把不同尺度的成分“冲散”再通过多次集总平均抵消噪声影响。简单说EEMD是“以噪声治模态混叠”代价是多算了N次EMD计算量变大。在实测中我更看重的是它让分解结果稳定可复现。EMD对极值点的采样顺序特别敏感信号里只要有一丁点噪声分解出的IMF就可能完全不同EEMD因为做了多次带噪声分解再平均结果对噪声不那么敏感。对后续要做样本熵统计、想比较不同工况下熵值差异的场景来说稳定是底线所以EEMD成了首选。1.2 样本熵才是分组的关键标尺分解完之后我们会得到若干个IMF从高频的细节分量到低频的趋势分量。但具体哪几个归入高频组、哪几个归入中频组、哪几个归入低频组按什么分直接看中心频率是最直觉的做法问题是有些IMF的频谱是重叠的光看主峰频率区分起来容易打架。按能量占比分呢又只反映了幅度大小反映不出信号结构的复杂程度。样本熵解决的就是“这个分量到底有多复杂”的问题。样本熵的物理含义是信号中出现新模式的概率有多大。一个信号如果规律性强、可预测比如标准正弦波那它产生新模式的概率很小样本熵值就低一个信号如果是杂乱无章的随机噪声那几乎每个时刻都在冒出新模式样本熵值就高。机械故障冲击信号里的高频分量往往既有确定性成分又有大量类噪声细节样本熵值会明显高于纯净正弦分量。这就是“按样本熵重构”背后的逻辑高熵的IMF信息结构复杂通常对应信号里的快变细节和噪声低熵的IMF规律性强通常对应平稳趋势和主要能量成分。把IMF按样本熵值升序排列低熵的几层重构出低频趋势信号高熵的几层重构出高频细节信号中间的归为中频主成分信号。这样得到的低频、中频、高频三组重构信号每一组反映的是“结构复杂度相近的成分”而不是单纯“频率相近的成分”对后续分析更有意义。1.3 这条技术链能解决什么实际问题我把这套流程用在过轴承故障数据的降噪和特征分离上。原始加速度信号里往往混着转频的低频振动、轴承元件通过频率的中频成分以及测点附近环境噪声和冲击细节的高频成分。直接用高通或者低通滤波器去分边界频率不好定而且滤波器本身会引入相位畸变。用EEMD加样本熵分组重构之后每组信号都是原信号成分的自然聚类重建出的信号不会再混入滤波器的振铃效应。在去噪这个方向上这套思路的优势尤其突出传统去噪是对所有IMF一刀切比如把前两层高频全扔掉但实际有些有效冲击就落在前两层里按样本熵判断则灵活得多如果前两层里某个IMF样本熵值并不明显高于其它层说明它还有明显结构不该整层删掉。这套方法本质上是把“去噪决策”从拍脑袋变成看指标可解释性很强。2. EEMD分解和样本熵计算的核心参数解析2.1 EEMD的关键参数该怎么取EEMD落地的时候有两个参数决定了分解质量白噪声幅值系数和集总平均次数。这两个参数我在PyEMD里是这样设的from PyEMD import EEMD eemd EEMD(trials200, noise_width0.2) imfs eemd(signal)noise_width是白噪声幅值系数代表添加噪声的标准差占原信号标准差的比值。经验值通常在0.1到0.4之间我一般先取0.2。这个值的物理意义很直观噪声太小时起不到破坏模态混叠的作用太大时会把真实的信号成分淹没导致分解出的IMF里面混入大量白噪声伪分量。有人喜欢在高采样率下用0.1在低采样率或信噪比低的场景下用0.3到0.4我个人的习惯是先看原信号的波形如果噪声本身就比较明显我会适当加大噪声幅值同时把集总次数也提上去。trials是集总平均次数这个参数直接决定残余噪声水平。根据EEMD的统计规律分解结果中残留的噪声大约等于白噪声幅值除以集总次数的平方根也就是残余噪声和trials的平方根成反比。举例来说噪声幅值系数0.2、集成100次时残留噪声大约是0.02倍原信号标准差如果集成400次残余噪声能降到0.01倍。但代价是计算时间几乎线性增加。实测下来200次已经能稳定复现IMF形态再往上加集成次数结果变化非常微小曲线几乎重叠这时就说明分解结果已经收敛。这两个参数其实是联动的噪声幅值设得大就需要更多次集成来压平残余噪声噪声幅值小可以适当减少集成次数。如果有人上来就抱怨EEMD重构误差大八成是noise_width设了0.5却只做了几十次集成。2.2 样本熵的两个关键参数嵌入维数和相似容限样本熵计算里两个参数维护不好整个分组就站不住脚嵌入维数m和相似容限r。我实现的样本熵函数是这样写的def sample_entropy(x, m2, r0.2*np.std(x)): x np.asarray(x, dtypefloat) N len(x) # 构造m维和m1维的相空间向量 def create_templates(data, dim): n N - dim templates np.array([data[i:idim] for i in range(n)]) return templates def count_matches(templates, r): n len(templates) # 向量两两之间用切比雪夫距离 count 0 for i in range(n): # 避免重复计算只和后面的向量比 diff np.max(np.abs(templates[i1:] - templates[i]), axis1) count np.sum(diff r) return count B count_matches(create_templates(x, m), r) A count_matches(create_templates(x, m1), r) if A 0 or B 0: return np.nan # 数据过短或容差过小 return -np.log(A / B)嵌入维数m我基本固定取2这是文献里验证过的默认值也符合工程直觉。太小的m会丢失序列连续性信息太大的m则要求数据得足够长否则模板匹配数量太少熵值估计方差变大。相似容限r则是整个计算里最敏感的参数一般取原IMF标准差的0.1到0.25倍。注意这里的标准差应该用当前IMF自身的标准差因为r的物理含义是“在多大误差范围内就算模式匹配”它应当与该分量的波动幅度成比例。如果对所有IMF用一个固定的绝对r幅度大的IMF匹配很容易熵值偏小幅度小的IMF匹配很困难熵值会出现许多NaN分组就废了。我在处理轴承振动数据时有个习惯不管IMF本身幅度大小先把每个IMF做z-score标准化也就是减去均值、除以标准差再在标准化后的序列上算样本熵。这时r取0.2倍标准化后的标准差也就是0.2对所有IMF都统一公平。标准化不会改变序列的时间结构只改变幅值尺度所以样本熵反映的复杂度信息是不变的。2.3 低中高频IMF分组策略拿到所有IMF的样本熵之后我会把熵值从小到大排一列。排序之后会发现一个很有意思的现象绝大多数情况下IMF的排列顺序和样本熵排序是单调对应的也就是第12个IMF最后一个接近残余项的分量熵值最低第0个IMF熵值最高。但偶尔也会出现中间某个IMF熵值突然偏高的情况这种异常往往意味着它包含了局部瞬态冲击或残留的模态混叠成分是值得重点检查的信号片段。分组时我通常采用相对比例式划分而不是死板地“前X层为低频后Y层为高频”。一种做法是把样本熵值归一化到0到1设定一个低阈值和一个高阈值比如熵值小于0.15的IMF归入低频组熵值大于0.65的归入高频组中间的归入中频组。另一种做法是对熵值做一个简单的KMeans聚类簇中心按低、中、高排序后归类。我更倾向于先看熵值分布直方图再定阈值因为不同信号源分解出的IMF数量差异很大一个固定比例不通用但直方图上如果出现明显的两到三个簇聚类结果通常很稳定。分好组后重构就简单了把组内所有IMF逐点相加。这里的相加不是简单的算术运算因为IMF本身经过了EEMD的均值化处理可以直接线性叠加。但要注意最后EEMD输出的residue残余项默认不在IMF集合里它代表信号的整体趋势基线我习惯根据目的决定是否把它并入低频组。如果关注的是温度漂移这类慢变趋势那残余项必须并进低频组如果关注的是振动成分的分解残余项往往是直流偏置或极低频漂移可以直接丢弃防止重构信号整体被抬平。3. Python实操从EEMD分解到三频段重构完整流程3.1 准备测试信号与运行环境我建议先别拿真实数据直接上先用一段合成信号把整条链路跑通。我构造了一个包含三个确定性分量、一个瞬态冲击和噪声的混合信号这样EEMD分解后各IMF的熵值差异会比较明显也方便验证重构结果。import numpy as np from PyEMD import EEMD import matplotlib.pyplot as plt from scipy import signal as sig fs 1000 t np.arange(0, 5, 1/fs) # 低频分量5Hz slow 1.5 * np.sin(2 * np.pi * 5 * t) # 中频分量50Hz幅度调制包络 medium 0.8 * np.sin(2 * np.pi * 50 * t) * (1 0.3 * np.sin(2 * np.pi * 2 * t)) # 高频分量200Hz间隔出现的衰减振荡冲击 high_impulse np.zeros_like(t) for start in np.arange(0.5, 5.0, 1.0): idx int(start * fs) n np.arange(0, 100) high_impulse[idx:idx100] 0.5 * np.exp(-n/40) * np.sin(2 * np.pi * 200 * n / fs) # 白噪声 noise 0.15 * np.random.randn(len(t)) signal_mix slow medium high_impulse noise信号三条主线分别是5Hz的低频趋势、50Hz调幅中频成分和200Hz的衰减冲击足够让EEMD的分解结果有区分度。实际处理乳酸菌、气温序列这种非振动信号时也可以按同样的公式构造任意频段成分这套测试方法很通用。3.2 EEMD分解与画图检查对混合信号做EEMD分解并画出每个IMF和残余项eemd EEMD(trials200, noise_width0.2) eemd.noise_seed(42) # 固定随机种子保证结果可复现 imfs eemd(signal_mix) num_imfs imfs.shape[0] fig, axes plt.subplots(num_imfs 2, 1, figsize(12, 2 * (num_imfs 2))) fig.tight_layout() axes[0].plot(t, signal_mix, k, linewidth0.8) axes[0].set_title(Original Signal) for i in range(num_imfs): axes[i1].plot(t, imfs[i], linewidth0.8) axes[i1].set_title(fIMF {i1}) axes[-1].plot(t, eemd.residue, linewidth0.8) axes[-1].set_title(Residue) plt.show()每次运行后我都会先扫一眼IMF波形形态。判断分解是否成功的经验标准一是看高频IMF里是否还有明显的低频包络如果有就是模态混叠残留二是看同一时间尺度是否被拆到多个IMF里如果是说明噪声幅值不够或集成次数偏少。如果合成信号还算干净EEMD分解出的IMF数量一般在7到11个之间越往后的IMF时间尺度越大波形越趋近于缓慢正弦型。有个细节必须提eemd.noise_seed(42)这行在日常分析里很容易被忽略但它特别重要。不固定种子的话每次运行EEMD添加的白噪声不一样虽然大量平均后均值会收敛但在少量集总次数下你会发现每次分解出的IMF有微小差异。固定种子之后实验可复现这在写报告、对比参数时非常有利。我建议所有用EEMD做研究的人都加上这行。3.3 统一尺度计算样本熵对每个IMF计算样本熵之前先做一个标准化处理然后调用前面写的样本熵函数entropies [] for i in range(num_imfs): imf_norm (imfs[i] - np.mean(imfs[i])) / np.std(imfs[i]) se sample_entropy(imf_norm, m2, r0.2) entropies.append(se) entropies np.array(entropies) print(IMF sample entropies:, entropies)这里用r固定为0.2是因为IMF已经标准化成单位标准差序列。如果某个IMF的熵值是NaN大概率是它的长度太短或者数据近乎恒定这在最后几个低频IMF里偶尔发生。处理办法是对这个IMF直接跳过熵值计算把它强行归入低频组。因为一个几乎不波动的低频IMF结构复杂度的确很低归入低频组不会产生争议。3.4 根据熵值归组并重构低频中频高频我按熵值大小设定两个阈值低于0.3归入低频组高于0.7归入高频组中间的归入中频组。当然这个阈值要根据实际熵值分布微调所以我建议先打印出熵值列表观察一下有没有明显分层再决定阈值。下面的代码展示了一种动态阈值方式用熵值的最小值和最大值做归一化norm_entropies (entropies - np.min(entropies)) / (np.max(entropies) - np.min(entropies)) low_idx np.where(norm_entropies 0.25)[0] mid_idx np.where((norm_entropies 0.25) (norm_entropies 0.7))[0] high_idx np.where(norm_entropies 0.7)[0] low_freq np.sum(imfs[low_idx], axis0) mid_freq np.sum(imfs[mid_idx], axis0) high_freq np.sum(imfs[high_idx], axis0)np.sum(imfs[low_idx], axis0)这一步就是组内叠加。由于IMF是按行排列的low_idx选出的若干行在轴0方向上一加就得到了低频重构信号。做完重构后我至少会检查三个指标第一重构信号总和是否逼近原始信号也就是low_freq mid_freq high_freq与原始信号做差差的能量应该很小第二三组信号各自做FFT观察主频是否分别落在低频段、中频段和高频段第三每组重构信号的样本熵应该呈现出低到高的单调性。如果这三个检查都通过分组就是合理的。reconstructed low_freq mid_freq high_freq residual signal_mix - reconstructed print(重构误差能量占比:, np.sum(residual**2) / np.sum(signal_mix**2)) freqs, power_low sig.welch(low_freq, fs, nperseg2048) freqs, power_mid sig.welch(mid_freq, fs, nperseg2048) freqs, power_high sig.welch(high_freq, fs, nperseg2048)这个误差能量占比通常能到10的负5次方以下。如果占比明显偏大说明EEMD分解本身不完整或者你丢掉太多IMF没纳入任何分组。我见过有人把所有中间IMF全扔了只留低和高两组重构误差很大还误以为是EEMD的锅。3.5 重构结果怎么看以我这段合成信号为例理想结果是低频组主要保留5Hz的慢波时间波形非常平滑中频组保留50Hz调幅成分能看到包络起伏高频组保留200Hz衰减冲击轮廓清晰噪声被部分抑制和分离。画图时把三组信号纵向叠在一起就能很直观地看到信号被“解剖”成了结构复杂度不同的三层。这个图拿到汇报里也很好用评审基本一眼就能理解你的思路。4. 常见问题与避坑指南4.1 问题排查速查表我把实操中遇到过的典型问题整理成了一张表每次流程跑出奇怪结果时先对照这张表排查能省不少时间。现象可能原因解决方法IMF数量特别少比如只有2个噪声幅值太小模态混叠没有完全破除调大noise_width到0.3左右同时增加trials高频IMF里叠加着低频包络模态混叠残留调大噪声幅值或减少单次EMD的筛分迭代次数样本熵值几乎全是NaNr取值太小或IMF长度太短将r设为0.2倍IMF标准差或先标准化再算所有IMF熵值差异不大数据本身平稳单一或r选取过大缩小r到0.1倍标准差再看分布重构误差大分组时丢掉了部分IMF检查low_idxmid_idxhigh_idx是否覆盖所有IMF同样的数据两次结果不同未固定随机种子调用eemd.noise_seed(42)分解耗时特别长trials过大且数据点过万先降采样再做分解或减少trials到100这里我想单独强调一下降采样的场景当原始数据采样率是10k甚至更高而且有几分钟长时直接跑EEMD一个成百上千次EMD叠加非常慢。如果目标频段是在低频和中频范围建议先低通滤波加降采样再做EEMD计算量能降一两个数量级而且不会丢失目标频段的信息。4.2 边界效应处理EEMD和EMD一样在数据两端有飞翼现象也就是端点处的包络估计容易出现大幅摆动。样本熵对边界附近的异常值并不敏感但重构信号的两端波形会明显失真。如果直接把重构信号拿去计算峰值或能量两端几毫秒的数据可能会污染整个统计结果。我常用的简单处理方式是先对原始信号做镜像延拓也就是把首尾各100个点左右的数据翻转拼接在序列两端延拓后的序列做EEMD再把分解结果的边缘裁掉。PyEMD里可以通过EMD(extrema_detectionparabol)之类的参数调整极值检测但对边界改善有限镜像延拓是最省事也稳定的做法。重构完成后把两端预留的那段裁掉再算样本熵结果会干净得多。4.3 对不同类型信号的调参建议这套方法不是只能用来处理振动信号。处理脑电这类生理信号时由于信号本身就有强噪声样本熵值整体偏高我会把低阈值和高阈值都相应提高比如将归一化熵小于0.2的视为低频、大于0.75的视为高频。处理温度、水位等缓变环境序列时分解出的IMF数量和样本熵的区分度都很好但要注意把EEMD输出的残差并入低频组否则趋势会缺失重构出的低频信号看起来就像一条绕基线波动的曲线而不是真实趋势线。处理低频缓慢信号时噪声幅值系数我建议取大一些比如0.4因为这类信号能量集中在少数低频IMF中需要更强的噪声来分离它们而处理高频振动信号0.15到0.2就足够噪声太大反而会把高频成分打散成多个伪IMF。5. 我的实操心得与后续扩展思路这套EEMD加样本熵加重构的流程做完一轮之后最值钱的产物其实不是那三组重构波形而是中间那张样本熵排序表。它给出了每个IMF“结构复杂度”的量化描述顺着这个排序你能立刻看出信号能量到底集中在一两个高熵层还是被平均打散到了各个尺度。这个信息对判定信号来源很有帮助比如滚动轴承早期故障的振动信号高频冲击IMF的样本熵会显著高于正常工况下的同一层IMF这就是一个天然特征。再一个心得是不要迷信固定的r0.2。标准化后统一取0.2只是一个相对稳妥的基线在具体项目里如果发现某一类样本的熵值区分度不够适当减小r到0.1往往能把差异拉开。但r也不能一直往小调太小会导致匹配数为0熵值变成NaN。我一般会写一个十几行的小脚本去扫描r从0.1到0.25范围内的熵值并观察体现在分布上是否稳定一致如果熵值的相对排序保持稳定说明模型选得没问题。如果后面想进一步扩展可以在每组重构信号上继续计算时域特征比如均方根值、峰值因子或者再做一次排列熵。多尺度熵也是这条路很好的延伸它把样本熵推广到不同时间尺度上非常适合分析生理信号的多尺度复杂度差异。信号重构部分也一样如果分出的高频组里还有明显可辨的周期成分可以再对它做一次EEMD形成迭代式的细分解。这套过程的灵活性很大但对基础流程的掌握才是根本希望这篇记录能让你少走一些弯路。