CEEMDAN算法实战:从EMD到CEEMDAN的信号分解与故障诊断
简介这份资源聚焦CEEMDAN算法即完全自适应噪声集合经验模态分解是EMD与EEMD的进一步改进方案面向从事非线性、非平稳信号处理的研究生、工程师与科研人员可用于地震波解析、心电与脑电分析、金融时序波动识别及机械故障诊断等场景。压缩包共5个文件以4个m脚本和1个mat数据文件为主脚本分别对应EMD、EEMD与CEEMDAN的分解实现及ICASSP2011示例mat文件提供心电信号样本整体约32KB轻量便于快速上手复现。目前已有5933人学习下载说明该实现具备一定参考价值。读者可借助脚本理解自适应噪声的添加方式、迭代分解与平均流程对比三种算法在模态混合与收敛速度上的差异并直接替换数据验证分解效果为后续特征提取与信号分析提供可运行的基础代码。1. 从模态混叠到 CEEMDAN这份算法资源到底解决了什么如果你处理过轴承振动、心电、地震波或者风速这类非平稳信号大概率被 EMD 的模态混叠折磨过。原始 EMD 在信号里混入一点噪声或间歇成分本征模态函数就会断裂、重叠分解出来的 IMF 物理意义直接崩掉。EEMD 靠加白噪声集合平均来压制这个问题但代价是计算量大、重构误差里残留噪声而且不同次分解的 IMF 数量对不齐。CEEMDAN 就是在 EEMD 基础上继续改进的版本通过自适应加入特定噪声、逐阶计算残差把模态混叠和重构误差同时压下去还让 IMF 的阶数更稳定。这份资源围绕 EMD、EEMD、CEEMDAN 三代算法的递进关系展开适合做信号分解、故障诊断、时序特征提取的工程师也适合正在选分解算法做论文复现的人。它不只是一个函数调用而是把三代算法的差异、参数影响和落地细节讲清楚让你知道什么时候该用哪一代、参数怎么设、结果怎么验证。2. EMD、EEMD、CEEMDAN 的递进逻辑为什么不是越新越好2.1 三代算法的核心差异与选型理由EMD 的出发点是任何复杂信号都能拆成有限个本征模态函数每个 IMF 满足两个条件极值点数量和过零点数量相差不超过一个且上下包络均值为零。它的实现靠筛分过程反复找极值、拟合上下包络、减去均值直到满足 IMF 条件。问题在于当信号里存在间歇性高频成分或噪声时极值点分布会跳变同一个时间尺度上的成分被拆到不同 IMF 里这就是模态混叠。EEMD 的思路很直接既然混叠来自极值点不稳定那就往原始信号里加不同 realization 的白噪声噪声把极值点分布“抹匀”每个尺度上的成分在多次分解中自动对齐到同一阶 IMF最后对集合做平均。白噪声在平均过程中互相抵消理论上残留量随集合数增大而趋近于零。但 EEMD 有两个硬伤一是加噪声后每次分解得到的 IMF 数量可能不同需要对齐处理二是重构信号时残留噪声虽然小但在高精度场景下仍不可忽略。CEEMDAN 的改进点在于不再对原始信号加噪声后直接分解而是每一阶都加入自适应噪声并且用上一阶的残差来驱动下一阶的分解。具体说第一阶 IMF 通过对加噪信号集合平均得到然后从原始信号里减去这个 IMF 得到第一阶残差第二阶 IMF 通过对残差继续加噪分解得到依此类推。这样每一阶的噪声都是独立控制的重构误差可以精确到接近零而且 IMF 阶数在不同集合间更一致。代价是计算量比 EEMD 略高因为每阶都要做集合分解。选型上如果信号干净、长度短、只做粗略趋势提取EMD 够用如果模态混叠明显但能接受轻微残留噪声EEMD 是性价比选择如果要做精确重构、多分量分离或者后续做希尔伯特谱分析CEEMDAN 更稳。常见做法是先用 CEEMDAN 跑一遍看 IMF 数量和频谱分布再决定是否降级到 EEMD 省算力。2.2 用 Python 复现三代算法从 EMD 到 CEEMDAN 的代码骨架下面这段代码用 PyEMD 库演示三代算法的调用差异。PyEMD 是常见的 EMD 系列实现CEEMDAN 在PyEMD.CEEMDAN里。先装依赖pip install EMD-signal numpy matplotlib然后跑一个合成信号包含两个不同频率分量加一点噪声import numpy as np from PyEMD import EMD, EEMD, CEEMDAN import matplotlib.pyplot as plt # 构造合成信号低频趋势 中频振荡 高频噪声 t np.linspace(0, 1, 1000) signal (np.sin(2 * np.pi * 5 * t) # 5 Hz 分量 0.5 * np.sin(2 * np.pi * 20 * t) # 20 Hz 分量 0.1 * np.random.randn(len(t))) # 噪声 # EMD 分解 emd EMD() imfs_emd emd.emd(signal) # EEMD 分解 trials 控制集合数noise_width 控制噪声幅度 eemd EEMD() eemd.trials 100 # 集合次数越大越稳但越慢 eemd.noise_width 0.2 # 噪声标准差相对信号幅度的比例 imfs_eemd eemd.eemd(signal) # CEEMDAN 分解 epsilon 控制每阶噪声幅度trials 控制集合数 ceemdan CEEMDAN() ceemdan.trials 100 ceemdan.epsilon 0.2 imfs_ceemdan ceemdan.ceemdan(signal) # 打印各算法得到的 IMF 数量 print(EMD IMFs:, imfs_emd.shape[0]) print(EEMD IMFs:, imfs_eemd.shape[0]) print(CEEMDAN IMFs:, imfs_ceemdan.shape[0])这段代码里trials是集合次数EEMD 和 CEEMDAN 都靠它做平均设太小残留噪声大设太大跑得慢一般 100 到 200 之间够用。noise_width是 EEMD 的噪声幅度通常取信号标准差的 0.1 到 0.3 倍太大引入伪分量太小压不住混叠。CEEMDAN 的epsilon类似但它是逐阶控制的第一阶可以稍大后面阶可以小一点。跑完看 IMF 数量CEEMDAN 一般比 EEMD 少一阶因为重构更干净。2.3 参数怎么设集合数、噪声幅度与 IMF 阶数的取舍参数设置是这类算法最容易翻车的地方。集合数trials直接决定噪声抵消程度理论上残留噪声标准差按1/sqrt(trials)衰减。设 50 时残留还比较明显设 200 基本干净设 500 以上收益递减但耗时线性增长。我一般先跑 100 看结果如果 IMF 波形毛刺多再往上加。噪声幅度noise_width或epsilon更敏感。设太小极值点分布没被抹匀模态混叠还在设太大噪声本身会分解出伪 IMF污染低频分量。经验值是信号标准差的 0.2 倍左右但具体要看信号信噪比。如果原始信号本身噪声就大可以降到 0.1如果信号很干净但混叠严重可以升到 0.3。IMF 阶数不用手动设算法自适应决定但可以通过max_imf限制最大阶数防止分解出太多无意义的高频分量。常见做法是分解后看每个 IMF 的频谱如果某一阶能量极低且频谱平坦基本是噪声可以丢弃。重构时只取前几阶能量占比高的 IMF能进一步降噪。提示EEMD 和 CEEMDAN 的结果有随机性每次跑 IMF 可能略有差异。做对比实验时固定随机种子或者多跑几次取平均。3. 把 CEEMDAN 用到实际信号上从分解到重构的完整链路3.1 读取真实信号并做预处理实际场景里信号往往带趋势项和异常点直接丢给 CEEMDAN 会分解出奇怪的 IMF。常见做法是先做去趋势和归一化。下面以一段模拟的轴承振动信号为例实际用时替换成你的数据加载逻辑import numpy as np from scipy.signal import detrend from PyEMD import CEEMDAN # 假设 signal 是从 CSV 读入的一维振动信号 # signal np.loadtxt(vibration.csv, delimiter,) # 去线性趋势 signal_detrended detrend(signal) # 归一化到零均值单位方差避免噪声幅度参数依赖原始量纲 signal_norm (signal_detrended - np.mean(signal_detrended)) / np.std(signal_detrended) # CEEMDAN 分解 ceemdan CEEMDAN() ceemdan.trials 150 ceemdan.epsilon 0.2 imfs ceemdan.ceemdan(signal_norm) print(IMF 矩阵形状:, imfs.shape) # (阶数, 样本数)去趋势用detrend去掉线性项归一化让epsilon的参数意义稳定。如果不归一化信号幅度是几千epsilon0.2就完全不起作用。这一步很多人忽略结果分解出来全是噪声。3.2 用频谱和能量占比筛选有效 IMF分解完不是所有 IMF 都有用。高频阶往往是噪声低频阶可能是趋势。筛选方法有两种看频谱峰值或者算能量占比。下面代码计算每个 IMF 的主频和能量占比from scipy.fft import rfft, rfftfreq n len(signal_norm) fs 1000 # 采样率按实际数据改 freqs rfftfreq(n, d1/fs) for i, imf in enumerate(imfs): spectrum np.abs(rfft(imf)) dominant_freq freqs[np.argmax(spectrum)] energy_ratio np.sum(imf**2) / np.sum(signal_norm**2) print(fIMF {i1}: 主频 {dominant_freq:.2f} Hz, 能量占比 {energy_ratio:.4f})主频接近工频或故障特征频率的 IMF 保留能量占比低于 0.01 的基本是噪声。常见做法是保留前 3 到 5 阶能量占比高的 IMF 做重构剩下的丢弃。重构就是把保留的 IMF 直接相加# 假设保留前 4 阶 imfs_kept imfs[:4] signal_reconstructed np.sum(imfs_kept, axis0) # 计算重构误差 reconstruction_error np.mean((signal_norm - signal_reconstructed)**2) print(重构均方误差:, reconstruction_error)CEEMDAN 的重构误差通常能到 1e-6 量级EEMD 在 1e-3 左右EMD 如果模态混叠严重可能到 1e-2。这个指标可以直接用来判断分解质量。3.3 希尔伯特谱与边际谱的快速画法分解完常要做时频分析希尔伯特谱是标配。下面代码对每个 IMF 做希尔伯特变换得到瞬时频率和瞬时幅度再画时频图from scipy.signal import hilbert import matplotlib.pyplot as plt plt.figure(figsize(10, 6)) for i, imf in enumerate(imfs[:5]): # 只画前 5 阶 analytic hilbert(imf) amplitude np.abs(analytic) phase np.unwrap(np.angle(analytic)) inst_freq np.diff(phase) / (2 * np.pi) * fs inst_freq np.append(inst_freq, inst_freq[-1]) plt.scatter(t[1:], inst_freq[1:], camplitude[1:], s1, cmapjet, alpha0.5) plt.colorbar(label幅度) plt.xlabel(时间 (s)) plt.ylabel(瞬时频率 (Hz)) plt.title(CEEMDAN 希尔伯特谱) plt.show()瞬时频率出现负值说明该 IMF 不满足窄带条件可能是分解不充分或噪声残留。常见处理是丢弃负频率点或者回到参数调整那一步加大trials。边际谱就是把希尔伯特谱对时间积分看频率能量分布适合做故障特征频率提取。注意希尔伯特变换要求 IMF 是窄带信号如果 IMF 里混了多个尺度瞬时频率会剧烈跳变。CEEMDAN 的 IMF 通常比 EEMD 更窄带但也不能保证每阶都合格。4. 避坑与排查CEEMDAN 落地时最容易翻车的五个点4.1 现象分解结果每次跑都不一样IMF 数量跳变原因EEMD 和 CEEMDAN 内部有随机噪声生成没有固定种子时每次 realization 不同。解决在调用前设np.random.seed(42)或者用库提供的seed参数。如果 IMF 数量仍然跳变说明trials太小加到 200 以上。4.2 现象低频 IMF 出现明显毛刺重构信号失真原因epsilon或noise_width设太大噪声分解出的伪分量混进了低频阶。解决把噪声幅度降到信号标准差的 0.1 倍同时增大trials补偿噪声抵消效果。如果还不行检查信号是否归一化未归一化时噪声幅度参数会失效。4.3 现象分解跑得极慢几分钟出不来结果原因trials设了几百上千信号长度又长每次集合分解都是 O(N log N) 起步。解决先降采样到 1000 点以内试参数确定后再用原始长度跑。或者用max_imf限制阶数减少后续筛分次数。CEEMDAN 比 EEMD 慢如果只是做趋势提取EEMD 够用。4.4 现象IMF 主频和预期对不上故障特征频率消失原因信号预处理时用了错误的滤波或去趋势把特征频率滤掉了。解决检查detrend是否去掉了有用的低频分量检查采样率fs是否设对频谱轴算错会导致主频读数偏移。常见做法是先用原始信号画频谱确认特征频率存在再跑分解。4.5 现象重构误差很大CEEMDAN 还不如 EEMD原因CEEMDAN 的epsilon逐阶控制如果第一阶设太大后续残差里噪声残留反而更多。解决把epsilon降到 0.1 以下或者改用CEEMDAN的beta参数控制每阶噪声比例。另外检查是否把所有 IMF 都用来重构丢弃噪声阶后误差应该显著下降。5. 进阶技巧用 CEEMDAN 做多分量分离与参数自动寻优5.1 用重构误差和 IMF 正交性做参数自动搜索手动调trials和epsilon很费时间可以写个小循环自动搜。指标用重构均方误差和 IMF 正交性指数后者衡量各阶 IMF 之间是否泄漏def orthogonality_index(imfs): n_imfs imfs.shape[0] total_energy np.sum(imfs**2) cross_energy 0 for i in range(n_imfs): for j in range(i1, n_imfs): cross_energy np.sum(imfs[i] * imfs[j]) return np.abs(cross_energy) / total_energy best_score np.inf best_params {} for trials in [50, 100, 150, 200]: for eps in [0.1, 0.2, 0.3]: ceemdan CEEMDAN() ceemdan.trials trials ceemdan.epsilon eps imfs ceemdan.ceemdan(signal_norm) recon np.sum(imfs, axis0) mse np.mean((signal_norm - recon)**2) oi orthogonality_index(imfs) score mse 0.1 * oi # 加权组合权重按场景调 if score best_score: best_score score best_params {trials: trials, epsilon: eps} best_imfs imfs print(最优参数:, best_params) print(最优得分:, best_score)正交性指数越小越好理想情况接近零。加权系数 0.1 是经验值如果更看重重构精度就调小更看重分量分离就调大。这个搜索跑一遍大概几分钟比盲调靠谱。5.2 把 CEEMDAN 嵌入故障诊断流水线的注意点实际故障诊断里CEEMDAN 通常只是特征提取的第一步后面接包络谱或者机器学习分类器。几个经验第一分解前务必确认采样率足够轴承故障特征频率往往在几千赫兹采样率低于 10 kHz 可能根本采不到第二保留的 IMF 阶数不要固定按能量占比动态选不同工况下有效阶数会变第三重构信号做包络谱时用希尔伯特变换取包络再 FFT比直接对 IMF 做 FFT 更能突出冲击成分。我自己的习惯是每次换数据集先跑一遍 CEEMDAN 看 IMF 数量和主频分布确认分解合理再往下走。有一次偷懒跳过这步直接拿 EEMD 的结果做分类准确率比 CEEMDAN 低了 8 个点回头查才发现模态混叠把故障特征频率拆到了两阶 IMF 里。从那以后我每次换信号都强制走一遍分解质量检查看重构误差和正交性指数不达标就调参数重跑。希望这些能帮到你少走点弯路。本文还有配套的精品资源点击获取