MATLAB经验模态分解（EMD）实战：原理、参数调优与工程应用

如果你在MATLAB里对着一堆非线性非平稳信号发愁，FFT看不出门道，小波又拿不准基函数，那么经验模态分解（EMD）大概率是你需要的东西。这几年我在MATLAB里用EMD处理过不少振动和趋势信号，从最初只会调一句emd(x)，到后来被模态混叠和端点飞翼折磨，也算踩出了一条相对稳定的流程。这篇就把经验模态分解的核心原理、MATLAB内置函数的使用方法、参数怎么调、以及那些文档里不会写的坑，一次性说清楚。

我默认你至少用过MATLAB一段时间，知道什么是脚本、什么是命令行。如果你刚接触信号处理也不用紧张，我会把涉及的概念尽量讲成大白话，代码部分可以直接复制跑通。

1. 为什么要用EMD：先弄懂它解决了什么问题

1.1 傅里叶和小波解决不了什么

传统信号分析第一步，绝大多数人都会想到FFT。频谱确实能告诉你信号里有哪些频率成分，但代价是假设信号是平稳的、线性的。实际工程数据很少这么听话：轴承振动里的冲击成分是瞬态的，脑电信号里节律是时变的，风速序列更是随机性拉满。对这些信号做FFT，结果是一堆宽泛的谱峰，频率和时间的关系被抹得一干二净。

短时傅里叶加了个窗口，但窗口长度固定之后，时间分辨率和频率分辨率就互相打架。小波变换比STFT灵活，不过需要你提前选小波基，比如db4、sym8这类，选不同基函数结果差得不是一点半点。我在实际项目里经常遇到这种情况：小波分解完，同一段数据换个基函数，趋势和细节的边界完全变了。说白了，小波是“你要给信号一个预设的形状，它才按你的形状去拆”。EMD最大的不同，就是它不预设任何基函数，完全靠信号本身的局部极值尺度去自适应地拆解。

1.2 EMD的自适应拆解思路

经验模态分解的核心输出叫本征模态函数（IMF），外加一个残差（residual）。你可以把信号想象成一道菜的最终味道，IMF是逐层剥离出的调料层次，残差是最后剩下的汤底。FFT是在频域上找固定的食材清单，EMD则是从成品倒推每一层加料过程，每一层都是在当前剩余信号里找出来的。

这个“自适应”在实际使用中有个很直观的表现：同一段信号里如果有频率漂移，比如40Hz慢慢变成60Hz，EMD会把这个变化作为一个IMF整体提取出来，你要是用FFT看，只会看到一片宽峰，根本说不清频率是怎么变的。这也是我在故障诊断里优先用EMD而不是直接上频谱的原因之一。

1.3 适用场景与不适用场景

先给结论，EMD适合下面这些场景：

• 振动信号故障诊断，尤其是滚动轴承、齿轮箱这类带冲击特征的信号。
• 生物医学信号，比如脑电、心电、肌电的节律分离。
• 气象水文数据，像风速、水位、径流序列的趋势提取。
• 金融时间序列，低频趋势和高频波动拆开看。

不适合的场景也要说清楚。EMD对低信噪比信号很敏感，纯白噪声直接分解会得到一堆无物理意义的IMF；它也没法做真正的实时处理，一段10万点的数据分解下来可能需要几秒甚至更久。所以我的习惯是，拿到数据先做一次FFT粗看频带，再决定要不要上EMD，而不是无脑分解。

2. EMD算法原理逐步拆解

2.1 什么是本征模态函数

要想理解EMD，必须先理解它想得到什么。Huang当年定义IMF时给了两个条件：

• 整个数据段内，极值点个数和过零点个数相等，或者最多相差一个。
• 任意时刻，由局部极大值包络和局部极小值包络定义的均值必须为零。

第一个条件保证IMF是一个窄带信号，这样用希尔伯特变换算瞬时频率才有意义。第二个条件有点像一个局部对称的要求，确保IMF的波形在时间轴上没有大范围的偏置。你不需要死记这些定义，只要抓住一句话：IMF是那些频率成分在时间轴上比较“单纯”的分量，要么是单一调幅调频波，要么是局部窄带信号。

2.2 筛分过程的核心逻辑

EMD的整个流程可以浓缩成几个步骤，我直接结合代码来拆。

假设你有一段信号x，长度为N。第一步是把局部极大值点和局部极小值点找出来，可以用findpeaks，也可以用极值比较的写法。第二步是用三次样条插值把极大值点连成上包络，极小值点连成下包络。第三步是求上下包络的平均值，得到m。第四步，用x减去m，得到候选分量h。

h = x - m，然后判断h是否满足IMF条件。如果满足，h就是第一个IMF；如果不满足，就把h当成新的信号，重复上述过程。这个反复迭代的过程就是“筛分”。

我放一个简化版包络计算代码，方便理解。

% 简化示意，仅供理解原理 [pks_max, loc_max] = findpeaks(x); [pks_min, loc_min] = findpeaks(-x); pks_min = -pks_min; t = 1:length(x); env_upper = spline(loc_max, pks_max, t); env_lower = spline(loc_min, pks_min, t); mean_env = (env_upper + env_lower) / 2; h = x - mean_env;

注意，这段代码只是让你看清楚包络怎么算。真实的内置函数远不止这么简单，它还要处理迭代停止、边界延拓、模态个数限制等一堆问题，但你理解了这段逻辑，就等于理解了EMD的核心。

2.3 停止条件与收敛门限

筛分不能无限迭代下去。迭代次数太少，IMF没“洗”干净；迭代次数太多，信号会被抹成等幅调频波，原本的幅度调制信息就没了。Huang当年用的是标准差门限，也就是相邻两次筛分结果的差异占比，小于0.2到0.3就停止。MATLAB内置的emd函数里也有对应的停止控制，默认的SiftMaxNumIterations是100，实际跑到二三十次多数情况就收敛了。

我在实际调试里的经验是，除非你面对的是特别复杂的信号，否则默认参数直接跑就够了。真正需要调的是后面要讲的MaxNumIMF和插值方式，而不是死磕迭代门限。

2.4 三次样条包络与过冲问题

上下包络用三次样条插值画出来之后，经常会在信号端点附近出现大幅振荡，这就是所谓的包络过冲。原因很简单，样条要满足二阶导数连续，而数据端点的导数信息是无法凭空产生的，插值结果就会“甩尾巴”。

遇到这种情况，有两个直接有效的办法：一是把插值方式从'spline'换成'pchip'，分段三次Hermite插值不会产生那么剧烈的过冲；二是后续讲端点处理技术时，先对信号做延拓，分解完再截掉延拓部分。

3. MATLAB环境准备与工具箱选择

3.1 不同MATLAB版本的函数差异

先说明一个容易踩坑的地方：MATLAB的EMD可不是每个版本都有。从R2020a开始，Signal Processing Toolbox才正式内置了emd函数，同时配套提供hht、ceemdan等函数。如果你用的是R2018b甚至更老的版本，直接执行emd(x)会报“未定义函数或变量”。

我见过不少朋友从网上下载了一个老EMD工具箱，然后跟新版MATLAB内置函数撞了名字，结果一行代码报出一堆莫名其妙的问题。遇到这种情况，第一步先执行：

which emd -all

这一句能列出当前MATLAB能找到的所有emd函数路径。如果有两个以上，说明冲突了。你自己写的工作目录里如果有emd.m，它的优先级比工具箱内置函数还高，这才是很多“结果不对”的根源。

3.2 直接使用内置emd还是第三方工具箱

我自己评估了一下，整理成一张对比表，你选型的时候直接参考：

我的建议很简单：能用内置就用内置，官方维护的东西出问题好查文档。只有当你需要把EMD的算法逻辑改得面目全非，或者研究端点处理、包络优化时，再去看第三方开源实现。

3.3 安装与路径配置常见坑

第三方工具箱解压以后，放到一个全英文路径下，然后在启动脚本或者命令行里执行：

addpath(genpath('D:\Toolboxes\emd_toolbox'));

注意不要放在带中文或者空格过多的路径里。有些旧工具箱是.m文件带GUI的，执行前还要先运行它自带的初始化脚本。

还有一个小技巧，建议在项目开头写明：

clear; clc; close all; addpath(genpath('D:\Toolboxes\emd_toolbox')); which emd -all

这样一启动就能确认当前用的是哪个版本，省得后期排查半天发现调用错了函数。

4. EMD算法在MATLAB中的实现与参数调优

4.1 最小可运行示例

网上很多教程一上来就贴几万字代码，反而把人搞晕。我先给你一个最小可运行示例，跑通了再去理解细节。

fs = 1000; t = (0:999)' / fs; x = sin(2*pi*50*t) + 0.5*sin(2*pi*120*t) + 0.2*randn(size(t)); % 内置EMD分解 [imf, residual, info] = emd(x); % 绘图 figure; subplot(size(imf,2)+2,1,1); plot(t, x); title('原始信号'); for k = 1:size(imf,2) subplot(size(imf,2)+2,1,k+1); plot(t, imf(:,k)); title(['IMF ', num2str(k)]); end subplot(size(imf,2)+2,1,size(imf,2)+2); plot(t, residual); title('残差');

运行完你会看到，50Hz的正弦分量被拆到了某个IMF里，120Hz的分量在另一个IMF里，随机噪声基本被压到高阶IMF中。这就是EMD最直观的用处：不用先验信息，把叠加在一起的几种频率分量分开。

有一点要注意，内置emd返回的imf是矩阵形式，每一列是一个IMF，而不是细胞数组。如果你想单独画第二个IMF，用的是imf(:,2)，不是imf{2}。这个坑我在同事的代码里看过好多次，他自己定义了个新变量，把函数返回值赋错了。

4.2 核心参数详解与建议取值

内置emd函数支持不少Name-Value参数，下面几个是我实际用得最多的：

实际经验是，不要一开始就调一堆参数。先用默认参数跑一遍，然后看info结构体里的信息，比如info.NumIMF到底给了几个分量，info.EnergyRatio看各IMF能量占比，再针对性地去改。

4.3 边界效应与端点处理

EMD一个老毛病就是端点飞翼，信号两端分解结果经常明显发散。这是因为包络样条在端点处缺少约束，极值点没法延伸到边界外。

一个实用的处理思路是镜像延拓：先取信号开头和结尾各一段，做镜像翻转，拼接成更长的信号，分解后再裁掉两端。我平时会在正式分析前写一个简单延拓函数，类似这样：

K = 50; % 延拓长度，视信号周期而定 x_ext = [flipud(x(1:K)); x; flipud(x(end-K+1:end))]; % 对x_ext做EMD分解，然后取中间原始长度部分

这个方法不是万灵的，但对付大多数情况下已经够用。如果你处理的信号本身很长，边界影响范围很小，直接忽略也没问题。

5. 实测案例：从振动信号到趋势提取

5.1 用两个正弦叠加噪声验证分解质量

我先用一个仿真信号说明EMD怎么跟FFT配合看结果。假设信号是50Hz和120Hz两个正弦叠加，再加一点噪声，这就是前面那段代码。跑完EMD后，你不仅能看到两个IMF分别对应两个频率，还能从info里看到每个IMF的能量占比。

这时再画FFT：

X = fft(x); f = (0:length(x)-1) * fs / length(x); plot(f, abs(X));

你会发现FFT只能告诉你“存在这两个频率”，但EMD能告诉你这两个频率成分在时间域上是如何叠加和分布的。对于故障诊断，我通常先看FFT确定大致频带，再用EMD把特定频带对应的冲击分量单独抠出来。

5.2 从IMF里提取故障特征

做轴承故障诊断的时候，光把IMF画出来没用，你得提取特征。我常用的是每阶IMF的峭度和能量占比。

energy = sum(imf.^2, 1); energy_ratio = energy / sum(energy); for k = 1:size(imf, 2) kurt_value(k) = kurtosis(imf(:, k)); end

冲击类故障在IMF上的表现是峭度明显偏高。你只需要扫一遍kurt_value，哪个IMF峭度异常，就重点看哪个。这个方法我在实际项目里验证过，比直接对原始信号求峭度更可靠，因为原始信号里正常振动成分会把冲击特征稀释掉。

5.3 EMD、EEMD、CEEMDAN怎么选

很多文章会把这三个名词混在一起说，其实它们的关系是递进的。EMD容易产生模态混叠，也就是不同频带的成分出现在同一个IMF里。混叠严重时，用EEMD，思路是多次给原始信号加白噪声，再对分解结果取平均，噪声辅助下的极值分布会更稳定。但EEMD的残留噪声比较明显，于是又有了CEEMDAN，它每次加的噪声是自适应产生的，最终IMF更干净。

MATLAB里直接用ceemdan就可以：

[imf_ceemdan, residual_ceemdan] = ceemdan(x, 'MaxNumIMF', 6);

我的选型经验是：快速预览、看趋势，用emd就够了；要做科研出图或者特征提取，优先ceemdan，计算量大一些但结果干净；EEMD除非你有特别要求，否则日常完全可以用CEEMDAN替代。

5.4 趋势提取案例

有一类需求很常见：从带有波动的数据里把长期趋势提出来。比如设备缓慢劣化变量上叠加了周期性波动。EMD处理这种问题非常顺手，因为残差本身就是趋势。

t = (0:999)' / 100; x = 0.02 * t + sin(2*pi*5*t) + 0.3*randn(size(t)); [imf, residual] = emd(x, 'MaxNumIMF', 3); plot(t, x, 'b'); hold on; plot(t, residual, 'r', 'LineWidth', 2); legend('原始信号', '趋势');

趋势线就是这个红色残差，你不需要自定义滤波器系数，EMD自动就把趋势剥出来了。我第一次用这个功能的时候，感觉比滑动平均省心多了，因为不用选窗口长度。

6. 常见问题与排查技巧实录

6.1 模态混叠怎么处理

模态混叠是EMD被讨论最多的问题。它通常出现在信号里有间歇性高频成分的时候，比如一段平稳振动中突然出现一个冲击，EMD会把冲击的一部分“借”到相邻的IMF里，导致一个IMF包含两种不相关频带。

我的处理顺序是：先看时域波形和频谱，确认混叠的大致位置；然后对信号做带通滤波预处理，把明显不相干的频带先分开；再试一次EMD，如果还不行就直接换CEEMDAN。长信号建议分段分解，不要指望一次把整段数据全处理完。

6.2 端点飞翼效应

前面提过镜像延拓，这里再补充一个排查技巧。如果你发现某个IMF两端明显发散，先用肉眼对比一下原始信号的端点幅值。很多时候是数据本身两端就存在突变，这时候无论怎么延拓都白搭，不如直接裁掉两端各2%的数据再分解。

另外一个经验是，把Interpolation设为'pchip'可以减轻飞翼，但不是所有情况都有效。我在实测里，周期信号用'pchip'效果很好，带冲击的信号还是会飞，这时候只能延拓加裁剪双管齐下。

6.3 分解速度慢、IMF数量超预期

信号太长和噪声太大都会导致筛分次数激增。遇到这种情况，先降采样，但注意降采样前要低通滤波，否则会混叠出假频率。

然后是限制IMF数量：

[imf, residual, info] = emd(x, 'MaxNumIMF', 5, 'SiftMaxNumIterations', 50);

如果分解出的IMF数量还是超过预期，看info，通常会发现问题出在某个IMF上，它很可能是一个没有物理意义的噪声分量。我在项目里有个习惯，所有IMF都要算一遍能量占比，低于0.1%的直接不参与后续分析。

6.4 版本、路径、编码相关报错

这里汇总几个我见过的高频报错：

• 报错“未定义函数或变量 'emd'”：版本低于R2020a，或者没装Signal Processing Toolbox。
• 报错“多个函数具有相同名称”：自己写的或有第三方emd.m跟官方冲突，执行which emd -all确认。
• 中文注释乱码：MATLAB 2023里脚本编码是GBK，如果你用VS Code或记事本存成UTF-8，打开就会乱。解决方法是把脚本用外部工具统一转成GBK，或者在MATLAB的预设项-常规-文件编码里改成UTF-8，再重新打开文件。
• 激活异常、License Manager Error这类属于环境问题，优先检查许可文件路径和环境变量，跟算法无关。

6.5 分解结果不稳定

EMD对微小扰动很敏感，哪怕两次运行中间加了条randn，结果都可能不一样。这不算bug，是算法本身特性。解决思路是：在代码开头固定随机种子，如果是CEEMDAN，能提高可复现性；或者同一信号分解十次取平均IMF，代价是慢，但结果更稳。

7. 进阶扩展：从EMD走向HHT谱与工程落地

7.1 画希尔伯特谱

EMD只是第一步，分解出的IMF配合希尔伯特变换，就得到了希尔伯特-黄变换HHT，能画出时间-频率-能量三维谱图。

hht(imf, fs);

运行之后会弹出一张谱图，横轴是时间，纵轴是频率，颜色深浅代表能量强弱。这个图特别适合观察频率随时间的变化，比如轴承故障早期，某个频带的能量会周期性增强，时频谱上很容易看出来。

7.2 边际谱与频谱对比

对希尔伯特谱在时间维度上求和，得到的就是边际谱。它的特点是不像FFT那样需要假设信号平稳，分辨率在局部频带上更有优势。

[hs, f] = hht(imf, fs); marginal_spectrum = sum(hs, 2); plot(f, marginal_spectrum);

我一般会把FFT谱和边际谱画在一起对比看。FFT比较直观，边际谱能突出局部窄带成分，两个结合着用，信号的频率结构基本就摸清了。

7.3 结合机器学习做自动诊断

EMD后接分类器，是故障诊断里很常见的工作流。大致流程是：采集数据，分段，每段做EMD或CEEMDAN，然后对每阶IMF提取能量占比、峭度、样本熵这些特征，拼成一个特征向量，最后丢给fitcsvm或者fitcensemble训练分类器。

我这里给一个最朴素的特征提取片段：

features = []; for k = 1:size(imf, 2) features = [features, sum(imf(:,k).^2), kurtosis(imf(:,k))]; end

实际项目里特征可以加很多，但加多了容易过拟合。我的经验是先做特征重要性排序，保留靠前的10个以内，否则小数据集上分类器会飘。

7.4 封装成自己的工具函数

工程落地时，别每次都复制粘贴十几行代码。我会把EMD流程封装成一个函数，输入原始信号和采样率，输出IMF矩阵、残差和一张总览图。这样后来接手项目的人只要调一个接口就行。

function [imf, residual, info] = my_emd_pipeline(x, fs) [imf, residual, info] = emd(x, 'MaxNumIMF', 5, 'Display', 0); t = (0:length(x)-1) / fs; figure; subplot(size(imf,2)+2,1,1); plot(t, x); title('原始信号'); for k = 1:size(imf, 2) subplot(size(imf,2)+2,1,k+1); plot(t, imf(:,k)); title(['IMF ', num2str(k)]); end subplot(size(imf,2)+2,1,size(imf,2)+2); plot(t, residual); title('残差'); end

这样封装之后，不管后面换多少数据，调用方式都不变，同事拿去用也说省事。

我个人在实际调试中最大的体会是：别把EMD当成黑盒。它的数学形式不如小波漂亮，但只要理解筛分和包络这两件事，很多问题自己就能定位。我目前习惯的工作流是：先FFT看频带，再EMD快速分解看层数，混叠严重就直接切CEEMDAN，最后用hht看时频联合分布。另外建议大家在项目代码里加一个判断，检查info里的能量占比是否达标，防止后端程序拿到一个奇怪的分解结果。如果你也在MATLAB里折腾EMD，卡在哪个环节了，可以先按“检查版本→检查路径→重置参数”的顺序排查，八成能解决。