简介：这份PDF文档围绕鲸鱼优化算法WOA优化变分模态分解VMD展开，面向计算机、电子信息工程、数学等专业学生及信号处理学习者，帮助理解如何自动寻优VMD的关键参数K与惩罚因子alpha，从而改善复杂信号的分解效果。文档内容涉及边际谱、频率图、收敛曲线及超参数优化迭代过程图，便于观察每次迭代中参数变化与分解性能的对应关系；同时说明参数化编程、注释清晰、数据以Excel形式便于替换、主程序main可直接运行等Matlab实现特点。资源包共1个PDF文件，约600KB，适合作为课程设计、期末大作业或毕业设计中的算法创新点参考。已有1337人学习，说明该方向具有一定关注度。读者可从中获取WOA-VMD的整体思路、参数设置逻辑、结果展示方式和代码组织参考，为信号分解实验与论文写作提供较完整的案例。

1. 一条振动信号引出的问题：为什么要让 WOA 去挑 VMD 的 K 和 alpha

一段实测振动信号直接丢进 VMD，K 和 alpha 全靠手感试，分解出来的模态不是混叠就是被硬切成碎片——这是我第一次调变分模态分解时最真实的体验。WOA-VMD 这套 Matlab 代码解决的问题很具体：把模态数 K 和惩罚因子 alpha 交给鲸鱼算法自己搜，人只负责给出搜索区间和适应度函数。它适合手头有一份 Excel 采样数据、要在课程设计或课题里做出“参数自适应”这个创新点、又不想从零推变分公式的人。跑通之后你能看到收敛曲线、中心频率图和边际谱三张图，也知道结果不对时该动哪一行。

2. 技术选型：WOA 与 VMD 各自解决什么问题，又怎么拼在一起

2.1 VMD 的关键参数：K 与惩罚因子 alpha

变分模态分解把原始信号拆成 K 个带限的模态分量，每个分量围绕自己的中心频率振荡。它求解的是一个约束变分问题：所有模态加起来要能重构原信号，同时每个模态的带宽之和最小。这个优化问题在频域用交替方向乘子法迭代求解，代码里那一串alpha / tau / K / DC / init / tol就是它的全部输入。

问题在于 VMD 对 K 极其敏感。K 给小了，两个频率接近的分量会被强行塞进同一个模态，出现模态混叠；K 给大了，同一个物理分量会被切成两三块，中心频率挤在一起，包络谱上出现一堆没有物理意义的峰。alpha 控制带宽，alpha 太小模态带宽过宽、噪声跟着进来，alpha 太大带宽被压窄，分解出的波形严重失真。人工试参数通常要跑几十次，看着中心频率图判断“这次是不是切多了”，费时而且不可复现。

所以优化的目标就很清楚了：在给定的 alpha、K 搜索区间里，找一个让分解质量指标最好的组合。这本质是一个二维连续-离散混合优化问题，其中 alpha 连续、K 整数，且目标函数没有解析梯度，只能靠启发式算法去试。

2.2 WOA 的三种行为与二维搜索空间

鲸鱼优化算法模拟座头鲸的捕食行为，把候选解当作鲸鱼个体，最优解当作猎物。每代迭代时，每头鲸鱼按概率 p 在两种动作里选一个：p < 0.5 时走包围或随机搜索，p ≥ 0.5 时走螺旋攻击。包围阶段用|A| < 1判定，系数 A 由收敛因子 a 从 2 线性降到 0 控制，前期 |A| 大、探索范围广，后期 |A| 小、收缩到局部。随机搜索阶段则是从种群中随机挑一个个体当参考点，作用是防止整个种群过早塌到同一个局部最优。螺旋攻击用对数螺旋线让个体绕最优解转圈靠近，模拟鲸鱼围捕时的螺旋上升路径。

放到这个项目里，搜索空间的维度只有 2：第一维是 alpha，第二维是 K。维度低是这套方案能跑得快的前提，因为每次适应度评估都要完整跑一遍 VMD，而 VMD 本身是个迭代求解器，单次耗时并不低。如果维度上升到七八维，种群 5 个、迭代 20 次根本不够，收敛曲线会一路抖到底。

WOA 在这个场景下的优势是参数少。除了种群规模和迭代次数，需要调的只有一个螺旋常数 b 和收敛因子的线性下降策略，不像差分进化那样还有交叉率、变异因子。对课程设计和论文复现来说，少一个参数就少一个被问“这个值为什么这么取”的地方。

2.3 整体数据流：优化器在外层，VMD 在内层的嵌套结构

整条链路是嵌套的：外层是 WOA 的迭代，内层是每一头鲸鱼对应的一组 (alpha, K) 去跑一次 VMD，然后由适应度函数把这次分解结果压成一个标量熵值返回给优化器。优化器只认这个标量，越小越好，它不知道 VMD 内部发生了什么。

这个结构带来一个必须提前想清楚的问题：适应度函数的设计直接决定优化器往哪个方向走。如果用“各模态包络熵的平均值”当适应度，模态数越多、每个模态越纯净、平均熵越小，优化器就会天然偏向大的 K，最后 K 顶到上界。这是这套代码最容易翻车的地方，后面第五章会专门讲怎么加惩罚项。

数据流大致是四段：

1. main.m读 Excel，取出待分解的一维信号，去直流，设定固定参数 tau、DC、init、tol；
2. 把适应度函数句柄@(x) funKAlpha(x, f, tau, DC, init, tol)交给优化器，连同搜索上下界和维度；
3. 优化器返回最优位置Best_pos、最优适应度Best_score、收敛曲线curve和每代最优参数process；
4. 对Best_pos取整，回代到 VMD 再跑一次，拿到最终的u、u_hat、omega，画分解图、频率图、边际谱。

第 4 步的“再跑一次”不能省。优化过程中为了省时间，很多实现会允许 K 和 alpha 带小数参与搜索，最后的整数结果必须用真实的 K 重新分解一次，否则你画的图和适应度评估时用的图对不上。这一步是自己加的还是原代码就有，取决于版本，但只要你的适应度函数里做了round，出图前就必须重新调用一次 VMD。

2.4 为什么不用网格搜索

有人会问，K 从 2 到 10、alpha 从 800 到 2500，网格搜索扫一遍不行吗。可以，但代价摆在那：alpha 即使按 100 步长、K 全扫，也有 18 × 9 = 162 次 VMD 调用；如果想要 alpha 精度到 10，就是 1700 多次。而 WOA 用 5 个个体 20 代，理论上限 100 次评估，实际因为每代只评估新位置，量级差不多。

真正拉开差距的是自适应能力。网格搜索的结果是一张表，你得自己看表选；WOA 直接给一个最优点，而且带收敛曲线，论文里那张“迭代过程图”就是靠process变量画出来的。这也是这套代码被当作创新点用的原因——它不只给结果，还给过程。

3. 目录结构与运行链路：从 main.m 到分解结果图的完整走向

3.1 文件清单与各自职责

拿到源码包先别急着点运行，把文件按职责过一遍，后面改起来才不会改错地方。这类项目的目录通常很扁平，所有 .m 文件和数据放在同一层，常见构成如下。

提示：如果你的包里优化器文件名是别的算法名（比如SSA.m、GWO.m），但主程序注释写的是鲸鱼算法，先看主程序实际调用的是哪个函数。这套代码在不同版本之间换过优化器，函数名和注释没同步是常事。

文件职责清楚之后，运行顺序就明确了：main.m是唯一入口，直接 F5 运行，不需要先跑别的脚本。如果报“未定义函数或变量”，八成是当前工作目录不对——Matlab 只在当前文件夹和搜索路径里找函数，双击 .m 文件打开并不会自动把所在文件夹设为当前目录，这是新手第一个卡点。

3.2 数据入口：Excel 怎么读、列怎么选、去直流为什么必要

数据替换是这套代码最常被用到的功能，也是最容易出错的地方。原始代码通常用xlsread或readmatrix读表，取某一列做信号。这里有两个细节必须处理：转成行向量，以及去直流。

VMD 的内部实现大多按行向量处理信号，如果你的数据读进来是列向量，维度对不上会直接报错或者算出一堆无意义的结果。去直流的理由是：信号的均值对应 0 Hz 分量，如果不去掉，第一个模态会被这个直流分量整个占住，中心频率贴着 0，剩下的模态分配就乱了。

%% 数据入口：换成自己的 Excel 只需要改文件名和列号 raw = readmatrix('data.xlsx'); % R2019a 及以上；老版本用 xlsread % 如果第一行是表头，readmatrix 会给出 NaN，这里剔掉 raw = raw(~isnan(raw(:, end)), :); f = raw(:, 2)'; % 取第 2 列，转置成 1×N 行向量 f = f - mean(f); % 去直流，避免 0 Hz 分量吃掉第一个模态 f = f / max(abs(f)); % 归一化，防止幅值量级影响熵值比较 fs = 1000; % 采样率，画频率轴和边际谱时要用 N = numel(f);

参数说明：raw是原始表格矩阵，行是采样点、列是通道；raw(:, 2)'里的列号按你自己的表改，单通道数据就是第 1 列；mean(f)减掉的是全局均值，如果信号有缓慢漂移，可以先用detrend或高通滤波；max(abs(f))归一化不是必须的，但如果你要横向比较不同信号的包络熵，就必须统一幅值尺度，否则熵值不可比。

fs这一行很容易被忽略。它不参与分解，只影响后面画频率轴和计算瞬时频率。如果fs设错，中心频率图的横轴会整体缩放，你会误以为分解出的频率成分不对。拿到数据先确认采样率，这是血泪经验。

3.3 参数区逐行：从 pop 到 lb/ub 怎么改

主程序里最该动手改的是优化器参数区。原始设置一般是种群 5、迭代 20、alpha 区间 [800, 2000]、K 区间 [4, 10]。这组值能跑出图，但从优化角度看偏保守，具体怎么调要看你的目的。

%% 优化器参数区 pop = 10; % 种群数量，原版 5 偏小，2 维问题 10 比较稳 Max_iter = 30; % 最大迭代次数，看收敛曲线平了就可以停 lb = [800, 3]; % 下边界：[alpha 下界, K 下界] ub = [2500, 10]; % 上边界：[alpha 上界, K 上界] dim = 2; % 维度固定为 2，顺序不能反 fobj = @(x) funKAlpha(x, f, tau, DC, init, tol); % 固定随机种子，保证每次跑出来的 K、alpha 一致，写论文必须加这句 rng(2024); [Best_score, Best_pos, curve, process] = WOA(pop, Max_iter, lb, ub, dim, fobj); % 取整后回代，用真实整数参数再分解一次 Best_pos = round(Best_pos); alpha_opt = Best_pos(1); K_opt = Best_pos(2);

参数说明与调整建议：pop和Max_iter是乘法关系，总评估次数大致是pop × Max_iter，每次评估跑一遍 VMD，所以先拿小数据试跑，测出单次 VMD 耗时再决定这两个值。lb和ub的第一位是 alpha、第二位是 K，顺序和dim = 2严格对应，写反了优化器会在 K 的区间里搜 alpha，结果荒谬但不会报错，属于静默失败。rng这行原代码里往往没有，加了之后结果可复现，答辩时被要求“再跑一次”不至于给出两个不同的 K。

round(Best_pos)对两个参数都取整，alpha 取整影响不大，但如果你希望 alpha 有更细的分辨率，可以只对 K 取整，alpha 保留小数传进 VMD。这属于可选改造，第五章会讲为什么 alpha 取整在某些数据上会带来台阶效应。

4. 核心代码拆解：优化循环、适应度函数与 VMD 调用三段

4.1 优化器主循环：WOA 一代是怎么算完的

优化器的骨架是标准的启发式算法结构：初始化种群，评估，找最优，进入迭代，每代更新位置、越界截断、重新评估、更新全局最优。真正决定收敛质量的是位置更新那一小段。

function [Best_score, Best_pos, curve, process] = WOA(pop, Max_iter, lb, ub, dim, fobj) % 在上下界内均匀撒点做初始化 Positions = repmat(lb, pop, 1) + rand(pop, dim) .* repmat(ub - lb, pop, 1); Best_score = inf; Best_pos = zeros(1, dim); curve = zeros(1, Max_iter); % 每代最优适应度 process = zeros(Max_iter, dim); % 每代最优参数，用来画参数变化曲线 for i = 1:pop % 初始种群先评估一轮 fit = fobj(Positions(i, :)); if fit < Best_score Best_score = fit; Best_pos = Positions(i, :); end end for t = 1:Max_iter a = 2 - 2 * t / Max_iter; % 收敛因子，2 线性降到 0 a2 = -1 - t / Max_iter; for i = 1:pop r1 = rand(); r2 = rand(); A = 2 * a * r1 - a; % |A|<1 包围，|A|>=1 随机搜索 C = 2 * r2; b = 1; % 螺旋常数 l = (a2 - 1) * rand() + 1; % [-1,1] p = rand(); if p < 0.5 if abs(A) < 1 D = abs(C * Best_pos - Positions(i, :)); Positions(i, :) = Best_pos - A * D; else idx = randi(pop); D = abs(C * Positions(idx, :) - Positions(i, :)); Positions(i, :) = Positions(idx, :) - A * D; end else D = abs(Best_pos - Positions(i, :)); Positions(i, :) = D .* exp(b * l) .* cos(2 * pi * l) + Best_pos; end % 越界截断，防止 alpha 或 K 跑到区间外 Positions(i, :) = min(max(Positions(i, :), lb), ub); fit = fobj(Positions(i, :)); if fit < Best_score Best_score = fit; Best_pos = Positions(i, :); end end curve(t, :) = Best_score; process(t, :) = Best_pos; end end

逐段说明：初始化用lb + rand .* (ub - lb)，这在二维箱式约束下等价于均匀撒点，比全零初始化好得多，因为全零起点会让整代个体挤在同一位置，前期几乎没有探索。收敛因子a从 2 降到 0，A的绝对值随之前减，前半程 |A| 大概率大于 1，走随机搜索，后半程自动切到包围，这个从探索到开发的过渡是 WOA 能用的关键，如果你把a改成常数，收敛曲线会明显变差。

l用(a2-1)*rand()+1生成，a2从 -2 降到 -1 附近，使得l落在 [-1, 1] 区间，这是螺旋项的标准取法。process记录每代最优参数，是画“超参数迭代过程图”的数据来源，没有这个变量就画不出参数随迭代变化的曲线。

评估次数上有个细节：内层每个个体更新后都调用一次fobj，也就是每代pop次评估。这对 VMD 来说负担不小，因为每次评估内部都是一次完整的变分迭代。如果跑得太慢，可以改成“只评估更新幅度大的个体”或者用parfor并行内层循环，前提是 VMD 本身没有全局状态。

4.2 适应度函数 funKAlpha：包络熵怎么把一个分解结果变成一个数

适应度函数是整个项目的技术核心，它决定了什么叫“分解得好”。这类代码常用的指标是包络熵：对每个模态做希尔伯特变换取包络，归一化成概率分布，算信息熵。物理意义是，如果一个模态是被噪声主导的，它的包络起伏杂乱、能量分布均匀，熵就大；如果它对应一个清晰的冲击或周期成分，包络有明显的起伏结构，能量集中，熵就小。

function fit = funKAlpha(x, f, tau, DC, init, tol) % x(1) -> alpha，x(2) -> K，这里取整是为了让 K 合法 alpha = round(x(1)); K = round(x(2)); % 非法参数直接判死，避免优化器浪费一次 VMD 计算 if K < 2 || K > 15 || alpha <= 0 fit = 1e10; return; end % 对当前参数完整跑一次 VMD [u, ~, ~] = VMD(f, alpha, tau, K, DC, init, tol); % 逐模态计算包络熵，取平均作为适应度 e = zeros(1, K); for k = 1:K env = abs(hilbert(u(k, :))); % 解析信号求模，得到包络 env = env / (sum(env) + eps); % 归一化为概率分布 e(k) = -sum(env .* log(env + eps)); % 香农熵，越小越“干净” end fit = mean(e); end

参数说明：x是优化器传进来的候选解，长度等于dim；f是主程序里预处理过的行向量信号；tau取 0 表示不加噪声松弛，带噪数据可以考虑给一个小的tau；DC为 1 时第一个模态被固定在 0 Hz，前面去了直流的话保持 0 即可；init控制中心频率初值分布，0 是全零、1 是均匀分布、2 是随机，实测init = 1收敛最稳；tol一般取 1e-6 到 1e-7，太小会拖慢每次评估。

hilbert需要 Signal Processing Toolbox，如果环境里没有这个工具箱，可以自己用 FFT 实现解析信号：补零到偶数长度、只保留正频率、乘 2 再逆变换。eps加在分母和 log 里是为了防止某个包络点为 0 导致log(0)得到 -Inf，一旦出现 -Inf 或 NaN，优化器的比较逻辑就崩了，表现出来的现象是整条收敛曲线都是 NaN。

注意：用mean(e)作适应度存在单调偏移风险，K 越大平均熵往往越小，优化器会一路把 K 顶到上界。这不是代码写错了，是指标设计的问题，解决办法见 6.2。

4.3 最优参数回代与三张图的画法

拿到最优参数之后，要重新跑一次 VMD，拿到完整的u、u_hat、omega用来出图。这一步的代码很短，但参数必须用取整后的值，不能直接把优化器返回的浮点解丢进去。

%% 用最优参数重新分解 [u, u_hat, omega] = VMD(f, alpha_opt, tau, K_opt, DC, init, tol); %% 图 1：收敛曲线 + 两个参数随迭代的变化 figure; subplot(3,1,1); plot(curve, 'b-o', 'LineWidth', 1.2); xlabel('迭代次数'); ylabel('适应度'); title('WOA 收敛曲线'); grid on; subplot(3,1,2); plot(process(:,1), 'r-s', 'LineWidth', 1.2); xlabel('迭代次数'); ylabel('alpha'); title('惩罚因子 alpha 变化'); grid on; subplot(3,1,3); plot(process(:,2), 'k-^', 'LineWidth', 1.2); xlabel('迭代次数'); ylabel('K'); title('模态数 K 变化'); grid on;

收敛曲线看两件事：一是是否单调不增（优化器每代都保留历史最优，理论上必须单调，如果出现上升说明你的实现里最优解没被正确保存）；二是尾部是否已经压平，没有压平就说明迭代次数不够，加大Max_iter再看。参数变化曲线看的是 K 在什么位置稳定下来，如果 K 一直在 9、10 之间跳，那就是过分解的信号，回头去改适应度。

频率图直接用 VMD 返回的omega，它是每个模态的中心频率轨迹，通常只画最后一行，也就是收敛后的中心频率，横轴是模态序号。如果两个模态的中心频率几乎重合，说明 K 给多了或者 alpha 太小，模态没分开。

边际谱稍微麻烦一点，标准做法是对每个模态做希尔伯特变换得到瞬时幅值和瞬时频率，然后把幅值按频率分箱累加。简化实现如下。

%% 图 2：边际谱（简化版，按频点累加瞬时幅值） MS = zeros(1, floor(N/2)); for k = 1:K_opt h = hilbert(u(k, :)); A = abs(h); % 瞬时幅值 fi = [0 diff(unwrap(angle(h)))] * fs / (2*pi); % 瞬时频率 for n = 1:N idx = round(fi(n) * N / fs) + 1; % 频率映射到频点下标 if idx >= 1 && idx <= numel(MS) MS(idx) = MS(idx) + A(n); end end end fre = (0:numel(MS)-1) * fs / N; figure; plot(fre, MS, 'LineWidth', 1.2); xlabel('频率 (Hz)'); ylabel('幅值'); title('边际谱'); grid on;

这段是内层双重循环，N上万时明显慢。常见做法是先对fi做histcounts加权统计，或者把内层换成向量化的accumarray：MS = accumarray(idx', A', [numel(MS), 1])'一行顶掉整个 for 循环。参数上fs必须和信号真实采样率一致，unwrap不能省，否则相位卷绕会让瞬时频率出现锯齿。

5. 避坑与排查：跑不通、收敛曲线抖动、K 顶到上界的 5 个真实坑

5.1 现象：一运行就报未定义函数 VMD

原因基本只有两种：一是当前工作目录不在源码文件夹，Matlab 的搜索路径没包含它；二是包里的函数文件名和调用名不一致，比如文件叫vmd.m而代码里写VMD(...)。第三种少见但存在的情况是函数被放在了子文件夹里而没addpath。

解决：先在命令行敲pwd确认当前目录，不对就cd过去；然后which VMD看能不能定位到文件。如果返回空，用addpath(genpath(pwd))把当前文件夹及子文件夹全加进路径。文件名大小写在 Windows 上不敏感、在 Linux 上敏感，跨平台跑之前统一改成和调用一致的名字，这是最省事的做法。

5.2 现象：每次运行的 K 和 alpha 都不一样，图也对不上

原因是没有固定随机种子，而 WOA 的初始化、rand()、randi()全都是随机的。种群只有 5 个个体、迭代 20 代时，初值的影响被放大得很明显，两次运行给出 K=6 和 K=8 完全正常。

解决：在主程序开头加rng(2024)，放在任何随机调用之前。如果要做统计，就不能只跑一次，正确做法是固定不同的种子跑 10 次，记录每次的最优 K 和 alpha，取众数或者取适应度中位数对应的那一组。把“单次最优”当成结论写进报告是常见误用，评审问一句“重复性如何”就答不上来。

5.3 现象：K 每次都收敛到搜索上界，中心频率图挤成一团

这是这套代码最高频的技术问题。原因是包络熵对 K 单调：模态越多，每个模态的频率成分越单一，包络越规整，平均熵越小。优化器只认数值，自然一路加 K。

解决有两条路。一是改适应度，加过分解惩罚，比如fit = mean(e) + lambda * K，lambda取 0.01 到 0.1 量级，具体值通过观察惩罚前后 K 的落点来定；二是换判据，用“最小包络熵 + 模态间相关系数约束”或者“重构误差 + 中心频率间隔”的组合指标。还有一条更简单的：直接收紧ub的第二位，把 K 上限设成 8，用先验知识兜住。这不是优雅方案，但在明确知道信号有几种成分时最有效。

5.4 现象：优化器前期就收敛，整条曲线一条直线

原因是种群初始化撞在了同一个位置。rand(pop, dim)在 pop 很小时容易生成一堆相近的点，加上 alpha 的量级是几百到几千，位置更新里A * D的步长非常大，一步就飞出边界被截断回lb或ub，整代个体全部贴在边界上，看起来就像“秒收敛”。

解决：给 alpha 做尺度归一化，把搜索区间映射到 [0,1]，在适应度函数里再还原成真实值，这样两个维度的步长量级一致，A * D不会在 alpha 方向过大。具体写法是lb = [0, 3]; ub = [1, 10];，然后在funKAlpha里alpha = round(800 + x(1) * (2500 - 800))。这一个改动通常能让收敛曲线从直线变成正常的阶梯下降。

5.5 现象：读进来的 Excel 数据分解结果全是噪声

原因通常是列选错了或者表头没处理。readmatrix遇到表头会给出 NaN 行，VMD 拿到含 NaN 的信号不会报错，而是算出一堆 NaN 传下去，最后画图空白；另一种是选到了时间列而不是幅值列，时间列是单调递增的，VMD 会把它当成一个超低频趋势分量。

解决：读完之后先size(raw)和raw(1:3, :)看清结构，再决定列号；用any(isnan(f))检查信号里有没有 NaN，有就整行剔除或者用fillmissing插值；最后plot(f(1:min(2000,N)))肉眼确认一下波形像不像振动信号。这一步花十秒，能省掉半小时的“为什么图是空的”。

6. 二次开发：换信号、换适应度、换优化器的三条改造路线

6.1 换优化器：只要对齐四个返回值

这套代码的优化器和 VMD 是解耦的，只要新算法满足“输入 pop、Max_iter、lb、ub、dim、fobj，输出最优适应度、最优位置、收敛曲线、参数变化矩阵”这个契约，主程序里换一行调用就能跑。这也是它适合做对比实验的原因——同一份数据、同一个适应度，换三种优化器就能出一张对比表。

% 只改这一行，其余代码不动 [Best_score, Best_pos, curve, process] = GWO(pop, Max_iter, lb, ub, dim, fobj); % 或者 [Best_score, Best_pos, curve, process] = SSA(pop, Max_iter, lb, ub, dim, fobj);

需要留意的是返回值顺序。不同来源的优化器函数返回的参数个数和顺序常常不一致，有的是[Best_score, Best_pos, curve]三个，有的把curve放在第一位。直接替换会出现“收敛曲线画出来是二维矩阵”这种诡异现象，原因是把Best_pos当成了curve。稳妥做法是换完之后先whos看一眼各变量的尺寸，curve应该是 1×Max_iter，process应该是 Max_iter×dim。

6.2 换适应度：给 K 加惩罚，避免过分解

针对 5.3 的问题，最值得做的改造是把适应度从“纯熵”改成“熵 + 惩罚”。下面这个版本还顺手加了中心频率间隔项，两个模态靠得太近就扣分。

function fit = funKAlpha_v2(x, f, tau, DC, init, tol, lambda) alpha = round(x(1)); K = round(x(2)); if K < 2 || alpha <= 0, fit = 1e10; return; end [u, ~, omega] = VMD(f, alpha, tau, K, DC, init, tol); e = zeros(1, K); for k = 1:K env = abs(hilbert(u(k, :))); env = env / (sum(env) + eps); e(k) = -sum(env .* log(env + eps)); end % 惩罚 1：模态数越多扣得越多，抑制过分解 pen_K = lambda * K; % 惩罚 2：相邻中心频率间隔过小（<5Hz）时扣分 fc = sort(omega(end, :)); gap = diff(fc); pen_f = 0.5 * sum(max(0, 5 - gap)); fit = mean(e) + pen_K + pen_f; end

参数说明：lambda是惩罚强度，从 0.05 开始试，观察最优 K 是否从 10 降到 4 到 6 的合理区间；omega(end, :)取的是最后一次迭代的中心频率，gap小于阈值说明模态没分开，max(0, 5 - gap)只惩罚过近的部分，不奖励过远。这套改法的代价是引入两个新超参数，答辩时可能被追问，所以最好配一组lambda的敏感性实验，把 K 随lambda变化的曲线也画出来。

6.3 验证方法：三组对照实验

改完之后怎么证明改得对，比改本身更重要。我一般固定做三组对照：

第一组是固定参数 VMD，K 用手工试出来的值、alpha 用默认 2000，作为基线。第二组是原始 WOA-VMD，不做惩罚。第三组是加了惩罚的版本。三组用同一份数据、同一个rng种子，比较四个指标：最优 K、平均包络熵、重构误差 RRMSE，以及单次运行耗时。

RRMSE 会随着 K 增大单调下降，所以它不能单独当判据，必须和 K 一起看——K 加了一大截而 RRMSE 只降了 1%，那就是过分解，不划算。如果三组结果的 K 差异不明显但中心频率间隔差很多，说明你这批数据本身成分就少，优化空间有限，换信号比换算法更有意义。

从那以后我每次改适应度函数，都会先把“固定种子的两次相同运行”和“固定参数基线”这两步强制走一遍，确认改动本身带来的差异大于随机波动，再去看图表。不然很容易把随机噪声当成改进效果写进结论里。这套代码的骨架不大，改起来不费劲，难的是判断改完到底有没有变好，希望帮到你。

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