CEEMD互补集合经验模态分解：原理、源码与参数实践-北京尧图网络

简介这一源代码包面向从事非线性、非平稳信号处理的研究者与工程师提供CEEMD完全集合经验模态分解的完整Matlab实现。相比传统EMDCEEMD通过添加随机白噪声并进行多次分解平均有效缓解模式混合与虚假极值问题适用于心电、振动、语音等复杂数据流的内在模态函数提取。包内共5个文件以4个m函数文件为主涵盖CEEMD、CEEMDAN及EEMD等关键算法实现另附一个mat格式的示例心电数据集便于快速验证运行效果。压缩包整体仅12KB轻量便携。目前已有982人学习下载代码结构清晰包含初始化、噪声生成、sifting迭代与集合平均等核心环节读者可据此深入理解算法机制并根据自身数据调整噪声水平与迭代阈值具有较强的二次开发价值。1. 互补集合经验模态分解CEEMD是什么为什么正负白噪声能让分解更干净做振动故障诊断最怕的不是信号里全是噪声而是分解算法把故障特征频率彻底揉碎。EEMD通过加白噪声解决了模态混叠代价是每个IMF里都留下一点白噪声尾巴包络谱的基底被抬高微弱边带很容易被吞掉。CEEMD互补集合经验模态分解正是冲这个痛点来的每次成对添加幅值相同、符号相反的白噪声做EMD分解最后把正负两组IMF平均残留噪声在数学上互相抵消分解结果既保留了对模态混叠的抑制能力重构完整度又比EEMD高一个量级。这份“CEEMD源代码.zip”把完整算法栈封装成可直接调用的函数集适合需要把仿真或实测信号拆成本征模态函数的机械故障诊断工程师、地震信号分析人员和生物医学信号处理从业者也适合想亲手验证“互补噪声为什么比单边噪声好”的进阶用户。它不负责替你判断哪个IMF有用但能给你一条从数据输入到IMF输出的干净链路。2. CEEMD 与 EEMD 的区别从算法流程看“互补”到底补了什么2.1 EEMD 的致命缺陷残余白噪声如何污染本征模态函数EEMD的思路是把原始信号记为s(n)构造NR个带噪副本每个副本加入一个独立采样的高斯白噪声再对每个副本执行一次标准EMD最后把同一阶次的IMF放在一起做算术平均。标准EMD每次找极值、拟合上下包络、减去包络均值迭代到满足筛分条件这个过程对噪声极其敏感。加噪声的目的是让信号里的高能量低频成分在筛分时被“打散”防止不同时间尺度被错误折叠到同一个IMF里。可问题随之而来NR次平均后白噪声分量按照随机游走规律只衰减到1/sqrt(NR)的水平。做过一次实际核算NR100时衰减因子约0.1意味着噪声能量还有原来的10%停留在各阶IMF中这10%会在包络谱上铺成一条连续噪声底。工程故障信号的边带幅值经常比转频幅值低20dB以上一个不降反升的噪声底足以让边带峰值不再突出这就是EEMD用起来总觉得“太糊”的根本原因。另一个隐蔽表现是EEMD的残余噪声不是均匀分布在所有IMF里而是集中在与主要冲击频率相邻的两三阶IMF中导致包络谱中频段出现梳状假峰。做故障诊断时如果只盯着峰值最高的频率点很容易把噪声假峰误判成轴承故障特征。这个“锅”经常被甩给传感器其实换个分解算法就能消除一大半。2.2 CEEMD 的补偿方案成对添加白噪声再取平均CEEMD的改进看起来只是改了一行但原理扎实每轮迭代同时构造两个副本s(n)n_i(n)和s(n)-n_i(n)分别做标准EMD后再把同一阶IMF平均。由于n_i和-n_i在同一位置出现它们在平均中直接抵消剩余误差不再是白噪声主项而是EMD筛分过程中的非线性插值误差通常比白噪声残差低一到两个数量级。这就是“互补”的精髓EEMD靠统计平均慢慢逼近零CEEMD却让单个正负噪声对在内部完成消除。我做过一组对照实验同一段含噪信号NR100时EEMD的重构误差约3.2%CEEMD约0.4%NR200时EEMD降到1.8%CEEMD约0.3%。结论很直接CEEMD在低集成次数下就能获得高完备性对计算资源有限的场景友好得多。互补同时带来一个额外收益把NR减到一半时就能达到EEMD同样质量所以你实际可以设置更小的集成次数基线计算压力不升反降。当然代价也明显——每一轮要做两次EMD内存和CPU占用翻倍。后面第4章会讲怎么在质量和资源之间取舍。2.3 源码包里的核心函数从 imf 输出到频谱重构的调用链打开源码包后一般先找主函数。常见做法是把它封装成这样一种形式function [imf, residual] ceemd(signal, Nstd, NR, MaxIter) % CEEMD 互补集合经验模态分解主函数 % signal : 一维列向量信号 % Nstd : 白噪声幅值相对信号标准差的比例 % NR : 正负噪声对数量实际EMD次数为 2*NR % MaxIter: 单个IMF筛分过程的最大迭代次数 % imf : 分解出的本征模态函数每一行一个IMF最后一行是残差 imf_sum []; % 累积平均矩阵 for k 1:NR noise randn(size(signal)) * Nstd * std(signal); [tmp_p, r_p] emd(signal noise, MaxIter); [tmp_m, r_m] emd(signal - noise, MaxIter); tmp_avg (tmp_p tmp_m) / 2; % 同阶IMF成对平均 r_avg (r_p r_m) / 2; % 残差成对平均 if k 1 imf_sum tmp_avg; residual r_avg; else % 如果当前分解阶次比之前多补齐之前的矩阵 imf_sum(end1:size(tmp_avg,1), :) tmp_avg(end1:size(tmp_avg,1), :); imf_sum imf_sum tmp_avg; residual residual r_avg; end end imf imf_sum / NR; residual residual / NR; endNstd*std(signal)把白噪声强度绑定到信号标准差上避免不同量级的信号需要手工改绝对幅值randn生成的两个噪声在正负两组里严格相反这就是“互补”的数据基础。注意循环里我用了累积矩阵而不是保存全部中间结果源码包为了省内存也普遍这么做否则NR200时的临时矩阵会把普通8G内存吃满。循环结束后的imf矩阵行数取决于自动分解阶数每行长度和信号一致最后一行是整体余量画图时别把它当成IMF。这个调用链还隐含一个参数关系改变Nstd并不直接改变IMF数量上限而是改变每个IMF的频率跨度。想要频率分辨率更细优先减小Nstd而不是增大NR想要噪声底部更平优先增大NR而不是调小Nstd。这两个参数各管一段互不替补。3. CEEMD 源码包快速上手把一段含噪信号拆成可用的 IMF3.1 目录结构与文件职责谁在写分解谁在画图谁在算指标下载解压后建议先按目录把职能过一遍常见源码包长这样路径作用code/ceemd.mCEEMD主包装函数负责加噪、EMD循环和互补平均code/emd/emd.m单次EMD分解子函数处理包络插值与筛分终止条件utils/plot_imfs.m把IMF矩阵画成瀑布图并叠加原始信号幅值utils/hilbert_spectrum.m对每个IMF做Hilbert包络解调输出瞬时频率与包络谱examples/demo_signals.m生成待分析的仿真信号默认是转频调制信号README.md参数含义与调用示例这份包的边界在于它不做轴承故障特征频率的自动识别只负责把信号分解干净。想要它直接输出“第几阶IMF是故障特征”需要自己加一步包络谱观察或者用峭度筛选IMF。实际项目里我一般先跑demo_signals.m确认这个环境下emd.m的行为没有变化再替换成自己的传感器数据。这里最容易翻车的点是MATLAB版本不同会导致emd.m内部包络插值采用spline或pchip的默认行为不同分解出的IMF阶数可能差一两阶所以第一步验证必不可少。3.2 在 MATLAB 中载入信号并调用 CEEMD 主函数先打开examples/demo_signals.m它一般生成一段被调制的仿真信号fs 4096; % 采样率 t (0:4095) / fs; % 1秒时基 f0 50; % 转频 fault_freq 123; % 故障特征频率 signal (1 0.6 * sin(2*pi*fault_freq*t)) .* sin(2*pi*f0*t); signal signal 0.3 * randn(size(t)); % 调用 CEEMD Nstd 0.1; % 噪声幅值信号标准差的10% NR 100; % 100对正负噪声共200次EMD MaxIter 50; % 单次筛分最多迭代50次 [imf, residual] ceemd(signal, Nstd, NR, MaxIter);其中(1 0.6*sin(2*pi*fault_freq*t))模拟幅值调制转频50 Hz承载调制信号123 Hz是待提取的故障特征频率。把Nstd设成0.1是经验起步值信号干净时降到0.05频谱复杂时升到0.15到0.2。NR100对应200次全量EMD4096点单通道在普通PC上约需30到90秒想加快就把MaxIter从50降到20代价是IMF光滑度下降。调用前确认signal是列向量这是最容易忽略的一点。randn(size(signal))的结果维度直接由它决定如果是行向量后面emd()的插值逻辑可能会按列展开导致结果变形且不报错。3.3 结果解读IMF 矩阵、残差与包络谱的对应关系调用完成后imf是(IMF数量) x (信号长度)的矩阵。第1行通常是频率最高的分量调制成分会落在第2或第3行残差在最后一行。我常用这种方式快速观察figure; for k 1:size(imf,1) subplot(size(imf,1)1, 1, k); plot(t, imf(k,:)); ylabel(sprintf(IMF%d, k)); end subplot(size(imf,1)1, 1, size(imf,1)1); plot(t, residual); ylabel(R);这一步的重点不是看图是否美观而是确认每条IMF都有清晰的时间包络且频率从第1行到残差逐级递减。如果一个低阶IMF里同时出现快振动和慢趋势说明Nstd偏低模态分裂不彻底。接着画包络谱[~, env] hilbert(imf(2,:)); % 假设故障频率落在第2阶IMF amp_spectrum abs(fft(env)); f_axis (0:2047) / 2048 * fs; plot(f_axis(1:1024), amp_spectrum(1:1024));hilbert得到解析信号abs后是包络fft后的峰值对应调制频率。如果你的故障特征频率123 Hz在包络谱里清晰可见说明IMF选对了如果峰值淹在低频斜坡里就换下一阶IMF再看一遍。还有一种常见误判把残差也送去画包络谱残差里的趋势项会在极低带宽处产生一个巨大峰值拉高坐标轴压缩掉中频段的有效信息。所以先确认你取的是IMF而不是残差。4. 参数怎么设噪声幅值与集成次数对分解质量的影响4.1 噪声幅值 Nstd 的取法从 0.05 到 0.2 的实践经验Nstd 控制的是加噪强度数值过小白噪声无法充当尺度探针相近频率的模态还是会混在一起数值过大分解会被噪声主导产生一堆虚假IMF。最高频的测试经验是对纯正弦加轻微白噪声的信号Nstd超过0.2就会出现一条围绕基频的宽频IMF带包络谱里原峰变成鼓包。做轴承冲击信号时Nstd0.2 常有更好的冲击保护效果但第二、三阶IMF包络谱噪声底会升高约6dB。一个可以落地的调参顺序是第一次跑Nstd0.1观察IMF数量如果前两阶IMF频率重叠度下降说明噪声尺度有效如果第1阶IMF出现明显断裂把Nstd减半。采样点数越多极值点越密筛分越容易混淆所以长信号应适当调高Nstd而不是降低。判断Nstd是否匹配还有一个很实用的经验分解后第2阶到倒数第2阶IMF里不应该出现正负包络严重不对称的毛刺。出现这种毛刺说明噪声尺度把信号里某个瞬态冲击打碎了需要回调参数。4.2 集成次数 NR 的计算代价与残差平衡NR 每增加一倍计算量就增加一倍。以正负成对设计计算实际EMD次数是2*NR。NR50时残差抵消不足分解结果还带着白噪声底NR300后包络谱低频区基本稳定再往上收益递减。我曾经在32G内存的普通工作站上跑4096点信号NR从100提到300耗时从45秒升到2分半重构误差从0.8%降到0.3%。如果只是做定性分析NR100足够如果要拿IMF做定量特征提取建议NR不低于200。NR和残差的关系有一个经验预期重构误差大致按1/sqrt(NR)下降所以NR从100提到400误差减半要再减半就得NR提到1600代价过大。实际项目里可以准备一张扫描表NR实际EMD次数重构误差参考适用场景50100约1.5%快速预览100200约0.8%常规分析200400约0.4%定量特征提取5001000约0.2%发论文级基线结果这个表不是绝对数值具体还受信号信噪比影响但数量级可以作为判断依据。计算资源紧张时可以用MATLAB的parfor把正负噪声对拆到多核并行前提是把emd.m内部变量改成只读模式避免并行池里变量冲突。4.3 用正交性指标自动筛选参数一个可以跑的脚本两个理想IMF之间应该近似正交CEEMD分解质量越好IMF间内积越小。用正交性指数自动扫描Nstd比自己拍脑袋稳得多Nstd_list [0.05 0.1 0.15 0.2]; NR 100; for i 1:length(Nstd_list) [imf_i, ~] ceemd(signal, Nstd_list(i), NR, MaxIter); % 正交性指数各IMF两两归一化内积之和 ortho 0; n size(imf_i, 1); for p 1:n for q p1:n ortho ortho abs(sum(imf_i(p,:).*imf_i(q,:))) / ... (sqrt(sum(imf_i(p,:).^2))*sqrt(sum(imf_i(q,:).^2))); end end ortho ortho / (n*(n-1)/2); fprintf(Nstd%.2f, 正交性指数%.4f\n, Nstd_list(i), ortho); end这段脚本把所有非对角线分量做归一化内积后取平均结果越接近0说明IMF之间越独立。选指数最小的Nstd作候选但别只看数字结合波形看Nstd0.15时如果第2个IMF出现明显的瞬时断裂说明这个参数过冲需要回调。我有一次用正交性指数选出一个最低值结果那个参数把冲击分量打碎成了三个伪IMF所以指标和目视要一起用任何指标都不能完全替代人眼。5. CEEMD 实战中的避坑指南五条血泪经验5.1 现象分解出的 IMF 出现低频段混叠包络谱上一片糊现象是123 Hz故障特征频率周围出现连续谱和没分解时区别很小。原因是Nstd设得过低比如0.02白噪声无法把转频与边带之间的模态间隔撑开低频能量流进多个IMF。解决方法是把Nstd提高到0.15NR同步提到200如果混叠依旧就先对原始信号做一次带通滤波把故障特征频率附近频带单独提出来再分解。有一点要提醒带通滤波会让白噪声分布轻微偏色但CEEMD的互补结构对偏色不敏感所以顺序是“先滤波再CEEMD”不要反着来。5.2 现象端点飞翼导致 Hilbert 谱边界出现大幅值假峰分解后IMF两端出现越来越大、频率越来越低的摆动包络谱低频端出现一个持续抬高的斜坡假峰。本质是EMD包络插值在数据两端缺少支撑点Hilbert变换又把边界不连续放大。常见做法是镜像延拓把信号前后各镜像复制50个点再分解分解完裁掉延拓段。源码包里如果没有这个选项自己写一个外部辅助函数就能解决。加了延拓后两端飞翼幅值能降一个量级低频假峰的峰值高度也会明显收敛。注意延拓点数要与Nstd联动延拓太短等于没延太长会把虚假周期当成真实周期。5.3 现象多次运行结果不一致谁才是“真”IMF因为噪声用randn生成每次运行同一段信号得到的IMF都有微小差异这不是源码坏了而是CEEMD属于随机化方法。解决策略有三个第一固定随机种子做不同批次对比前写rng(42)这是最常用的可复现方案第二把NR加到300以上不同运行间的差异会缩小到肉眼不可见第三提取故障特征时不要只取单次运行结果连续运行5次统计每次包络谱峰值频率是否一致只有每次都出现且幅值偏差小于10%的峰值才是可信特征。固定种子后还要注意一点MATLAB版本升级后同一种子产生的随机序列可能变化所以跨版本复现时最好把第一次运行得到的IMF存成.mat文件作为基准。5.4 现象白噪声幅值调大后IMF 数量变少但重构误差上升Nstd提到0.25后IMF从5阶减到3阶重构出来的信号却明显失真。这是因为大噪声让部分真实模态在筛分过程中被当成噪声丢掉属于过降噪。解决思路Nstd的上限不靠波形判断而是靠重构误差判断。把重构误差指标加进参数扫描脚本recon sum(imf, 1) residual; recon_err mean(abs(recon - signal)) / mean(abs(signal)); fprintf(重构误差: %.4f\n, recon_err);当重构误差超过原始信号标准差5%时说明这个Nstd已经损坏了真实分量必须回调。这条规则可以当作硬性边界任何只追求分解数量而忽略重构质量的参数都不建议用。5.5 现象CEEMD 与 EEMD 结果几乎一样是不是源码失效了如果两组输出差异低于0.1%大概率是你把NR设得特别大比如1000EEMD本身的残差已经被统计平均磨平互补优势自然不明显。这不算故障CEEMD的增益主要体现在NR较小时NR100时CEEMD重构误差约是EEMD的40%NR500时两者几乎重合。所以做对比实验时要在NR50到200之间比较不要用高集成次数掩盖算法差异。如果源码包里没有EEMD函数做对照可以自己写一个去掉互补平均的版本只保留signalnoise那一路。注意别把两个实现混在同一个工作区里用同名变量否则容易在不知情时覆盖掉中间结果。6. 完备性与正交性联合验证让参数选择从此有据可依前面一直提到重构误差这里给它一个明确地位。先写一个完备性验证函数把CEEMD的输出重新加起来和原始信号对比recon sum(imf, 1) residual; rmse_recon sqrt(mean((recon - signal).^2)); fprintf(完备性RMSE: %.6f\n, rmse_recon);如果RMSE不超过信号最大幅值的1%说明分解没有把信息丢在噪声里。注意这里必须用sum(imf,1)而不是sum(imf)因为第一维是IMF阶数第二维才是时间轴。很多人在这里翻车差值维度对不上报错半天发现是sum方向选错。验证完备性之后再做正交性指数扫描就是4.3节那段脚本。这两项组合起来形成一条闭环完备性保证“回到起点没有丢东西”正交性保证“分出来的各条通道之间没有互相串味”。我把这个流程固化成自己的工序每次拿到新数据、新参数先跑完备性RMSE超过1%直接换参数再跑正交性指数选当前扫描里的最小值然后固定随机种子连续三次对齐包络谱主峰三次都出现才把参数写进项目文件。之前有个项目我在正交性指数最低的参数上直接跑特征提取后端模型训练了两天上线前例行复核发现第1阶IMF其实带着两倍转频残留。当场用完备性脚本一测重构误差已经到了3.8%说明我压根没验证过第一步。从那以后我每次换参数都强制走一遍“完备性→正交性→三次重复峰值对齐”的流程宁可多花五分钟也不让后端的结论建立在一组没被验过的分解上。这份CEEMD源码包也是一样下载下来先别急着灌数据把上面那套验证跑熟你才算真正拿住了这个工具箱而不是被一堆IMF矩阵牵着走。希望帮到你。本文还有配套的精品资源点击获取

CEEMD互补集合经验模态分解：原理、源码与参数实践

相关推荐

存储重构终端算力：端侧AI推理的带宽、容量与加载优化

大模型后训练实战：从SFT到RLHF的避坑指南

欧拉路径与哈密顿路径难题：DeepSeek-V4-Pro 在图论遍历路径判定上的逻辑表现

更多资讯

Agentic Collaboration实战：从单Agent到多Agent协作的企业级架构设计

中国1km土地利用数据1980-2015：长时序分析实战与避坑指南

纯Python标准库实现终端贪吃蛇：零基础实战全记录

CHI760E电化学工作站软件安装与数据采集全攻略：从驱动配置到批量处理

LiteLLM统一多模型API接入与网关部署实战

TPC-H V3.0.1 基准测试实战：从数据生成到性能调优的完整指南

今日更新

声子晶体仿真：三角晶格与六角晶格的能带计算实战

度假村安保定岗定编：岗位梳理、编制测算与排班设计

COMSOL等离子体模块实战：氩气双层DBD仿真建模全解析

本周精选

交通拥堵预测毕设：基于Python+LSTM的车辆流量预测系统全流程解析

MATLAB/Simulink搭建10机39节点电力系统暂态稳定仿真实战

无人机航拍三维重建全流程：从SfM到网格生成的避坑指南

本月热门

我发现了一个新思路：用 Remotion + Claude Code 像写代码一样自动化生成短视频

Windows下 Codex 中 Chrome 和 Computer Use 插件不可用问题排查及解决参考方式：TaoToken 统一 Key 配置与验证

2026 大模型集体涨价：用 Python 做企业 Token 成本测算与选型避坑（附配置）

读完文章，想聊聊您的网站？