简介：面向信号处理、时序分解方向的课程设计与毕业设计场景，这份资料围绕鲸鱼优化算法WOA自动寻优变分模态分解VMD的模态数K与惩罚因子α展开，给出可直接运行的Matlab实现，并配套边际谱、中心频率图、收敛曲线及超参数迭代过程图，便于观察每次迭代对分解效果的影响。压缩包共1个PDF文件，约600KB，主要承载代码说明、程序设计与结果展示内容，结构紧凑，适合快速查阅整体思路。已有1336人学习，说明该主题在相关学习者中有一定关注度。代码采用参数化编程，参数修改方便，注释较完整，数据以Excel形式组织，替换数据后运行主程序即可复现实验，可用于对比不同分解结果、理解WOA寻优流程，并为论文或项目中的创新点设计提供参考。

1. 拎着一列 Excel 数据来找我调 VMD 的人，十个里有九个栽在 K 上

上周有人拿一段轴承振动信号跑 VMD，K 随手填了 8，分解出来前四个分量几乎一模一样，频谱图上几条中心频率挤在一起，他把结果图放进报告里，被问了一句「这四个模态有什么区别」就答不上来。这就是WOA-VMD 鲸鱼算法优化 VMD 变分模态分解要解决的核心痛点：VMD 的模态数 K 和惩罚因子 alpha 靠手调，调不准就是过分解或者欠分解，而 Matlab 这套实现把 K 和 alpha 交给鲸鱼算法去搜，用包络熵当适应度，迭代几十代给你一组相对靠谱的组合。适合手上有一列时序数据、想做分解‑特征提取链条的电子信息、计算机、数学方向的同学，也适合把它当课程设计、毕业设计的创新点来改。

2. 项目结构与技术选型：main.m 之外那几个文件各管什么

2.1 目录与文件角色拆解

这份资源是一套平铺式 Matlab 工程，没有包管理、没有类，全是函数文件加脚本，好处是丢进当前文件夹就能跑，坏处是函数名冲突时排查起来很烦。按职责拆开看大概是这么几层：

新手最容易忽略的是 fun.m 和熵函数这两层。很多人以为「换了优化算法就完事」，实际上适应度函数才决定优化器往哪走，优化器只是走得更快或更慢。你把包络熵换成排列熵，同样的 K 范围能搜出完全不同的最优解，这个差别比换算法大得多。

VMD.m 的函数签名在源码注释里写得很清楚，我按注释整理一遍，这就是整套代码的参数契约：

% Input and Parameters: % --------------------- % signal - 待分解的时域信号（行向量，长度 N） % alpha - 惩罚因子（带宽约束强度） % tau - time-step of the dual ascent ( pick 0 for noise-slack ) % K - the number of modes to be recovered（模态个数） % DC - true if the first mode is put and kept at DC (0-freq) % init - 0 = all omegas start at 0 % 1 = all omegas start uniformly distributed % 2 = all omegas initialized randomly % tol - tolerance of convergence criterion; typically around 1e-6 % % Output: % ------- % u - 分解变量的集合（K x N，每行一个模态） % u_hat - 变量的频谱（K x N 复数） % omega - 变量中心频率（迭代记录 x K）

这张表里真正影响结果的就三个：alpha、K、init。tau设 0 是给含噪信号留松弛项的常见做法，设大于 0 会强制模态更窄、更容易把噪声也拆成一个模态；DC设 1 会把 0 Hz 那个直流分量单独拎出来，如果你的信号已经去过均值，设 0 和设 1 差别不大；tol一般 1e-6，调到 1e-3 会明显快但中心频率可能没收敛干净。init我一般用 1，均匀分布初始化，比全 0 收敛稳一些，随机初始化虽然有时能跳出局部极小，但代价是同一份数据跑两次结果不一样，写论文的时候会被追问可复现性。

2.2 为什么挑鲸鱼算法去调 K 和 alpha

VMD 参数寻优这件事本质是个黑箱优化：目标函数没有解析梯度，因为 K 是整数、分解过程本身还带迭代收敛判断，你对 K 求导毫无意义。这类问题能用的工具就那么几类——网格搜索、粒子群、遗传、模拟退火、鲸鱼，都行。选鲸鱼的实际理由很朴素：

一是参数少。WOA 除了种群数和迭代次数，只有一个从 2 线性降到 0 的控制参数 a，不需要调惯性权重、不需要设交叉变异概率，对写课程设计的人来说省事。二是它把探索和开发塞在一个框架里：|A| >= 1时做随机搜索保证全局性，|A| < 1时向当前最优包围，另外还有一条按概率触发的螺旋路径，用p < 0.5分流。三是代码短，一个主循环三十来行就能写完，学生改起来能看懂每一行在干什么。

不过有一点必须说清楚：WOA 收敛快是真的，早熟也是真的。种群只有 5、迭代只有 20 代的时候，它大概率在第一轮就锁死在某个局部最优附近，然后后面十九代都在原地打转，收敛曲线从第二代开始就是一条水平线。这不是 bug，是参数设置的问题，我后面在避坑章会专门讲。

2.3 适应度函数：包络熵为什么能当裁判

优化器需要一个能比较大小的分数，VMD 参数寻优里最常用的是包络熵。逻辑链是这样的：如果一个模态分量是干净的、有明确物理意义的成分，它的包络应该平滑、能量集中，熵值小；如果这个模态是噪声或者几个成分混在一起，包络起伏杂乱，能量分散，熵值大。所以「让某个分量的包络熵最小」等价于「让分解结果里至少有一个足够干净的成分」。

代码实现上分两步。第一步用希尔伯特变换取包络，第二步把包络归一化成概率分布再算香农熵：

function E = envelope_entropy(u) % u : K x N，VMD 输出的模态矩阵，每行一个分量 [K, N] = size(u); E = zeros(K, 1); for k = 1:K y = abs(hilbert(u(k,:))); % 希尔伯特变换取解析信号，取模得包络 y = y / sum(y); % 归一化，使其满足概率分布条件 y(y == 0) = eps; % 防止 log(0) 出现 -Inf E(k) = -sum(y .* log(y)); % 香农熵 end E = E(:); end

参数上没什么可调的，真正要决策的是「怎么把 K 个熵压成一个数」。三种常见做法：取最小值、取平均值、取所有分量的熵之和。取最小值是这套代码的思路——它鼓励优化器找到一个特别干净的分量，但副作用是它可能纵容其余 K-1 个分量全是垃圾。取平均值更均衡，但容易在 K 偏大时被拉低分数，导致优化器倾向于选小 K。取和则相反，会推着 K 往大走。我在自己项目里一般用「最小包络熵 + 对 K 加一个软惩罚」的组合，兼顾干净度和模态数合理性。

fit 函数的骨架就是把这个熵接到 VMD 上：

function fitness = fun(x, f, tau, DC, init, tol) % x(1) = alpha, x(2) = K，由优化器在 [lb, ub] 内给出 alpha = round(x(1)); % 惩罚因子取整，避免传小数进 VMD K = round(x(2)); % 模态数必须是正整数 K = max(2, min(K, 12)); % 硬性护栏，防优化器越界飞出 alpha = max(100, min(alpha, 5000)); [u, ~, ~] = VMD(f, alpha, tau, K, DC, init, tol); fitness = min(envelope_entropy(u)); % 最小包络熵作为适应度 end

注意round这一行。优化器在连续空间里搜，K 会拿到 6.37 这种值，不取整直接送进 VMD 会在内部循环里炸掉；取整之后又带来一个新问题——6.37 和 5.99 取整都是 6，等于大片搜索空间被映射到同一个点，优化器感觉不到梯度，浪费迭代次数。这是连续优化器解离散参数的通病，后面避坑章也会提。

3. 核心代码拆解：从一列 Excel 到 K 个模态的完整链路

3.1 数据读入与预处理

主程序第一步就是把 Excel 里那列数据拉进来，转成 VMD 能吃的行向量。这一步看着简单，翻车率不低：

clc; clear; close all; %% 1. 读入数据 data = xlsread('data.xlsx'); % 单列时序，第一行建议就是数据不要表头 f = data(:,1)'; % 转行向量：VMD 要求 1 x N f = f - mean(f); % 去均值，避免 0Hz 直流占掉一个模态 %% 2. 归一化（可选，但强烈建议） f = f / max(abs(f));

xlsread在较新的 Matlab 版本里会提示不推荐使用，改用readmatrix更稳，但老版本、老教程里都是xlsread，两者对表头的处理不一样：xlsread遇到表头会返回 NaN 行，readmatrix会直接跳过。你拿到手第一件事就是size(data)看一眼行列数对不对，如果第一行是 NaN，删掉就行。

去均值那行是我自己加的，源码不一定有。原因很实际：VMD 的DC参数如果设 1，第一个模态会固定在 0 Hz 附近，如果信号本身有较大的直流偏置，这个模态会吃掉大量能量，剩下的模态分配到的带宽就少了，看起来像是「分解不出来东西」。先减掉均值，再设DC = 0，是更省心的组合。

归一化也是同理。如果你的信号幅值跨度大（比如传感器原始 ADC 值），包络熵的数值范围会跟着放大，虽然熵本身对尺度不敏感（归一化成概率分布之后），但数值稳定性上归一化更保险。

3.2 优化器主循环：包围、螺旋、随机搜索三条路径

WOA 的核心就是把种群里每个个体往三个方向推：朝最优个体靠、沿螺旋线绕、或者随便找个同伴参照。主循环写出来是这样的：

function [Best_score, Best_pos, curve, process] = WOA(pop, Max_iter, lb, ub, dim, fobj) % pop : 种群规模 % Max_iter : 最大迭代代数 % lb / ub : 1 x dim 的上下界，dim = 2 对应 [alpha, K] % fobj : 适应度句柄 % ---- 初始化：在界内均匀撒点 ---- Pos = repmat(lb, pop, 1) + rand(pop, dim) .* repmat(ub - lb, pop, 1); Best_pos = zeros(1, dim); Best_score = inf; curve = zeros(1, Max_iter); process = zeros(Max_iter, pop, dim); % 记录每代全体，用于画参数演化图 for t = 1:Max_iter % ---- 评估当前代 ---- for i = 1:pop Pos(i,:) = max(Pos(i,:), lb); % 越界裁剪，比反弹更稳 Pos(i,:) = min(Pos(i,:), ub); fit = fobj(Pos(i,:)); if fit < Best_score Best_score = fit; Best_pos = Pos(i,:); end end a = 2 - 2 * (t / Max_iter); % 关键：从 2 线性衰减到 0 % ---- 位置更新 ---- for i = 1:pop A = 2 * a * rand() - a; C = 2 * rand(); b = 1; % 螺旋常数，一般固定 1 l = -1 + 2 * rand(); % [-1, 1] p = rand(); for j = 1:dim if p < 0.5 if abs(A) >= 1 idx = randi(pop); % 随机参照个体 D = abs(C * Pos(idx,j) - Pos(i,j)); Pos(i,j) = Pos(idx,j) - A * D; % 随机搜索 else D = abs(C * Best_pos(j) - Pos(i,j)); Pos(i,j) = Best_pos(j) - A * D; % 包围猎物 end else D = abs(Best_pos(j) - Pos(i,j)); Pos(i,j) = D * exp(b * l) * cos(2*pi*l) + Best_pos(j); % 螺旋 end end end curve(t) = Best_score; process(t,:,:) = Pos; % 供 huatu 画迭代过程 end

逐段说明。a从 2 衰减到 0 是整个算法的节拍器：A = 2*a*rand() - a的取值区间是[-a, a]，前期 a 接近 2，|A|很容易大于 1，算法就在全局乱逛；后期 a 趋近 0，|A|基本落在 (-1, 1) 内，算法收束到包围行为，围着当前最优细搜。这个「先探索后开发」的切换是自动的，你不用管。

p < 0.5分流那一段，p < 0.5走包围或随机搜索，p >= 0.5走螺旋。螺旋公式里exp(b*l)*cos(2*pi*l)是一条约等于 1 的对数螺旋线，l在 [-1,1] 抖动，效果是让个体绕着最优点螺旋靠近，而不是直线扑上去，能多扫一些邻域。

process(t,:,:)这个三维数组是给画图用的，源码里huatu(process,'VMD')就是吃它，画出「每一代 alpha 和 K 分别落在哪」的散点图或者轨迹图。这个图放进报告里很加分，因为它直观展示了搜索过程，比自己描述「算法收敛良好」有说服力。

参数上要动的是pop和Max_iter。源码给的pop = 5、Max_iter = 20是能跑通的最小配置，代价是搜索很粗糙。我自己的经验值是pop >= 10、Max_iter >= 30，K 的搜索范围比较宽的时候还要再加。计算量大致是pop × Max_iter次 VMD 分解，每次分解的开销跟信号长度和 K 成正比，信号长度 2000 点、K 10 的时候单次大概零点几秒到一秒，pop=10, iter=30就是 300 次，三五分钟能出结果。

3.3 主程序串起来：优化完再用最优参数分解一次

这是整套代码的主干，也是新手最容易写错的地方——优化器跑完之后，必须用最优参数重新调用一次 VMD，不能指望优化器返回模态：

%% 3. 优化 K 与 alpha pop = 10; % 种群数，源码默认 5 偏小 Max_iter = 30; % 迭代次数，源码默认 20 偏小 lb = [800, 4]; % [alpha下限, K下限] ub = [2000, 10]; % [alpha上限, K上限] dim = 2; % 两个决策变量 tau = 0; DC = 0; init = 1; tol = 1e-6; fobj = @(x) fun(x, f, tau, DC, init, tol); [Best_score, Best_pos, curve, process] = WOA(pop, Max_iter, lb, ub, dim, fobj); %% 4. 取整并再做一次干净分解 alpha_opt = round(Best_pos(1)); K_opt = round(Best_pos(2)); fprintf('最优 alpha = %d, 最优 K = %d, 适应度 = %.6f\n', alpha_opt, K_opt, Best_score); [u, u_hat, omega] = VMD(f, alpha_opt, tau, K_opt, DC, init, tol);

lb和ub的顺序一定要和x(1)、x(2)的语义对齐。这里lb(1)=800、ub(1)=2000管的是 alpha，lb(2)=4、ub(2)=10管的是 K。你要是把两个位置写反，K 会被限制在 800 到 2000 之间，第一次迭代就内存爆掉——这是血泪经验，看到 Matlab 卡死在 VMD 第一行、内存条飙红，先回来检查这个。

alpha 的范围 800～2000 偏保守。VMD 论文里常用的默认值是 2000，很多工程场景会试到 5000 甚至更大以压窄带宽。如果你发现分解出来的模态还是糊在一起，先把 ub(1) 抬到 4000 试试。

3.4 出图：收敛曲线、边际谱、频率图三件套

结果图决定了别人信不信你的分解有效。这套资源里能出的图基本是四类，我按用途排个序：

%% 5. 结果可视化 figure; subplot(2,1,1); plot(curve, 'b-o', 'LineWidth', 1.2); xlabel('迭代次数'); ylabel('适应度（最小包络熵）'); title('WOA 收敛曲线'); grid on; subplot(2,1,2); plot(1:length(f), f, 'k'); hold on; for k = 1:K_opt plot(1:length(f), u(k,:) + 3*(k-1), 'LineWidth', 1); end xlabel('采样点'); ylabel('模态 + 偏移'); title(['VMD 分解结果，K = ' num2str(K_opt) ', alpha = ' num2str(alpha_opt)]); grid on; figure; for k = 1:K_opt subplot(K_opt,1,k); plot(abs(fft(u_hat(k,:)))); % 每个模态的单边幅频 ylabel(['u', num2str(k)]); end xlabel('频率 bin'); title('各模态频谱'); figure; bar(mean(omega,1)); % 中心频率的收敛均值 xlabel('模态序号'); ylabel('中心频率'); title('模态中心频率分布');

四张图各自的用处：收敛曲线用来证明优化过程有效，但要注意它只反映适应度下降，不反映分解质量；时域分解图是给人看「成分分得开不开」最直接的一张，模态之间加偏移量（那个3*(k-1)）是必须的，不然 K 条曲线叠在一起什么都看不见；频谱图看的是每个模态是不是窄带、有没有重叠，VMD 的理想结果是 K 个互不重叠的窄带，如果两条频谱明显交叠，说明 K 取大了或者 alpha 取小了；中心频率柱状图看的是分布是否均匀，如果多个模态挤在同一个频率附近，那就是过分解的典型信号。

边际谱在这套里通常会单独出一个文件或一段，用希尔伯特变换对各模态做时频分析再对时间积分，得到整体的频率‑能量分布。它的价值在于跟原始信号的傅里叶谱对比——如果边际谱在几个主要频段上出现了清晰的峰，而原始谱上这几段本来就模糊，说明分解确实分离出了东西。

4. 跑通流程：从零到出图的完整步骤

4.1 环境准备与文件放置

Matlab 版本上，R2016b 之后基本都能跑，因为核心只用了hilbert、fft、xlsread这些老函数，没依赖工具箱里太新的东西。唯一要注意的是hilbert属于 Signal Processing Toolbox，如果装的是精简版 Matlab 没这个工具箱，envelope_entropy第一行就报错，这时候得自己用fft手动实现希尔伯特变换，或者换成不依赖解析信号的熵指标（比如排列熵）。

文件放置上，最省事的做法是把所有 .m 文件和 data.xlsx 放在同一个工作目录，然后cd进去或者右键文件夹「添加到路径」。Matlab 对函数文件的搜索是看当前路径和搜索路径的，如果主程序能跑、一到VMD(...)就提示「未定义函数」，八成是路径没加或者文件名大小写不对——Windows 上不敏感，Linux 和 Mac 上vmd.m和VMD.m是两个东西。

Excel 那列数据的要求写成清单更清楚：

• 单列，前面不要有表头行，有就用data = data(2:end,:)删掉
• 长度不低于 200 点，太短的信号 VMD 分解不出有意义的东西
• 采样率不参与计算，但画频率轴的时候需要，所以心里要有数，方便后面把 bin 换算成 Hz
• 数值不要有 NaN 和 Inf，any(isnan(f))提前检查一遍

4.2 参数逐项设置与含义

按下表把 main.m 顶部的参数过一遍，每一项我都写了推荐的起手值和改动方向：

K 的上限怎么估？看原始信号的傅里叶谱有几个明显的峰簇，峰簇数加一到二就是 K 的合理上限。比如一段轴承故障信号，谱上能看到转频基频、二倍频、共振带三个区域，那 K 大概 4～6 就够，上限设 8 已经宽了。设到 10 以上通常是在给优化器挖坑，它很容易挑一个大 K 然后靠其中一个分量拿到很小的包络熵，其余的都被牺牲掉。

跑之前把随机数种子固定一下，这一行能省掉后面写报告时「为什么两次结果不一样」的麻烦：

rng(2024); % 固定随机种子，保证 WOA 初始化可复现

4.3 怎么判断结果是对的

跑出图不等于跑对了。我一般按这个顺序自查：

第一看收敛曲线。健康的曲线应该是前期陡降、中间缓慢下降、最后趋于平缓。如果第一代就是最低点然后一动不动，说明种群太小或者初值正好落在最优附近，加大 pop 重跑；如果曲线到最后一代还在往下掉，说明迭代次数不够，加大 Max_iter。

第二看 K 的取值。如果最优 K 顶到了 ub(2) 的上界，说明你给的搜索空间不够，真正的解可能更大；如果 K 掉到了 lb(2) 的下界，通常意味着信号本身就很简单，或者适应度函数在惩罚大 K。两种情况都要回头调范围而不是直接接受结果。

第三看频谱重叠度。把各模态的频谱叠在同一张图上，用眼睛判断有没有明显交叠。理想状态是像梳子一样排开，每条带宽有限。交叠严重就调 alpha 上限。

第四看重构误差。把 K 个模态加起来跟原信号比：

err = norm(f - sum(u, 1)) / norm(f); fprintf('重构相对误差 = %.4f%%\n', err * 100);

VMD 本身是完备分解，理论上模态相加能还原原信号，这个误差正常应该在 1% 以内。如果超过 5%，说明分解没收敛，去查 tol 和 Max_iter。

5. 避坑与排查：五个把分解结果带偏的细节

5.1 优化器函数名和算法对不上

现象：标题写的是鲸鱼算法 WOA，主程序里调的却是SSA(...)，注释还写着「麻雀算法参数设置」，跑出来的结果你说不清到底是哪个算法优化的。

原因：这类代码往往是在一个通用模板上换算法改出来的，改了标题和注释，函数名和内部实现忘了同步。SSA 和 WOA 的调用签名恰好一样（pop, Max_iter, lb, ub, dim, fobj），所以能跑通，但内部逻辑完全不同——SSA 有发现者、加入者、警戒者三种角色，WOA 是包围‑螺旋‑随机搜索。

解决：打开优化器文件，看主循环里有没有a = 2 - 2*(t/Max_iter)和exp(b*l)*cos(2*pi*l)这两个特征式，有就是 WOA，没有就是别的。要对得上标题，就照 3.2 节那份实现自己写一个WOA.m替换进去，签名保持一致，主程序一行都不用改。这件事必须做，不然答辩时被问「你的鲸鱼算法三种行为分别在哪几行」直接卡住。

5.2 优化器返回的最优参数没取整就送进 VMD

现象：报错「K 必须为正整数」，或者更隐蔽——程序跑完了，但分解出来的模态数跟Best_pos(2)显示的数值对不上，比如显示 K=6.83，实际出来 7 个模态。

原因：优化器在连续空间搜索，Best_pos是浮点数。VMD 内部用 K 做循环上界，Matlab 会自动向上取整，所以不报错但结果错位。更麻烦的是适应度函数里如果已经round过，优化器感知的和你最终用的是两个不同的 K。

解决：在 main.m 里显式round，并且把取整后的值打印出来，用这个值重新调 VMD（见 3.3 节代码）。同时检查fun.m里的取整逻辑和 main 里是否一致，两处都用round而不是一处round一处floor，否则优化出来的最优点和你实际用的点会差一个。

5.3 K 搜索范围设太大导致过分解

现象：优化器兴高采烈地返回 K=10，适应度也确实低，但画出来一看，十個模态里有五六个是几乎一样的窄带噪声，中心频率挤在一起。

原因：最小包络熵这个指标只奖励「有一个分量特别干净」，不管其余分量死活。K 越大，抽到一个干净分量的概率越高，适应度自然越低。优化器没有理由不选大 K。

解决：三条路选一条。一是压紧 K 的上界，按信号实际的物理成分数来定，别超过 8；二是给适应度加惩罚项，fitness = min(E) + 0.01*K，让大 K 付出代价；三是把适应度换成所有分量包络熵的均值，逼优化器对整体负责。我一般用第一条加第二条的组合，简单直接。

5.4 同一份数据两次运行结果不同

现象：今天跑出 K=6，明天同样的代码跑出 K=7，alpha 也变了，报告里的图对不上。

原因：WOA 的初始化用rand，螺旋和随机搜索也都在用rand、randi，没有固定种子，每次启动都是新的随机序列。加上init = 2（随机初始化中心频率）的话，VMD 内部也不确定。

解决：main.m 开头加rng(2024)，init固定成 1。做完这两件事，同一台机器同一个 Matlab 版本上结果完全一致。换机器或换版本可能有微小差异，因为rand的实现细节会变，但 K 的取值一般稳定。文献里的可复现性要求就是这么满足的。

5.5 VMD 内部迭代不收敛却不报错

现象：程序正常跑完，没有警告，但中心频率图画出来乱七八糟，几个模态的频率在整个时间轴上来回跳，重构误差也偏大。

原因：tol设得太大（比如 1e-3）或者Max_iter（VMD 内部的，不是 WOA 的）设得太小，VMD 在还没收敛的时候就退出了，返回的是中间状态。这个函数不会告诉你「我没收敛」，它只是安静地返回。

解决：把tol拉回 1e-6 甚至 1e-7，VMD 内部的迭代上限调到 500，代价是慢一些。判断是否收敛的土办法是看omega的最后几行是不是稳定了，如果最后两行的差还在 1e-3 量级，那就是没收敛。另外信号长度太短（低于 200 点）也会导致收敛困难，这时候不是调参能解决的，得换数据段。

6. 二次开发：把裁判换掉，比把算法换掉更有用

改这套代码的方向，大部分人第一反应是「换个优化算法」，把 WOA 换成 GWO、PSO、蜣螂、麻雀，然后在报告里写「提出了一种改进的 XX-VMD 方法」。这个路子能走通，但性价比不高——所有元启发式算法在只有两个决策变量、适应度评估又贵的场景下，表现差异其实没那么大，你换成 PSO 大概率 K 还是那个 K。

真正能拉开差距的是动适应度函数。包络熵是 2016 年前后 VMD 调参的标配，但它只看单分量的稀疏性。几个更实用的替代思路：

改起来就在 fun.m 里动手，比如加上中心频率间距惩罚，防止模态频率挤在一起：

function fitness = fun_v2(x, f, tau, DC, init, tol) alpha = round(x(1)); K = round(x(2)); K = max(2, min(K, 10)); alpha = max(200, min(alpha, 5000)); [u, ~, omega] = VMD(f, alpha, tau, K, DC, init, tol); E = envelope_entropy(u); f_est = mean(omega, 1); % 各模态收敛后的中心频率 if K > 1 gap = min(diff(sort(f_est))); % 相邻中心频率的最小间距 else gap = 0; end % 三项加权：稀疏性 + 模态数惩罚 + 频率分离度（间距越大越好，取负） fitness = mean(E) + 0.02 * K - 0.5 * gap; end

写完先别急着换算法，用同一份数据、同一个rng种子，把 v1 和 v2 各跑三遍，比较 K 的稳定性、重构误差、以及画出来的频谱重叠度。这三个指标都改善了，才算改对了；只有一个变好另外两个变差，那是把问题挪了个位置。

验证的时候还有一个坑我踩过：改完适应度忘了同步改lb、ub的量纲。比如新指标把频率间距也加进去了，适应度的数量级从 0.1 变成 10，那么原来收敛曲线的纵轴范围就完全没法看了，画图脚本里的ylim得跟着调。这种事不影响结果正确性，但图难看，报告里一眼就露怯。

再说一个具体技巧：如果你只想快速验证某个 K 是不是比优化器给的更好，不用重跑整条优化链，直接写个小脚本扫一遍固定 alpha 下的 K，十分钟就能画出一条 K-包络熵曲线，比盯着收敛曲线猜要实在得多。

Ks = 3:10; scores = zeros(size(Ks)); for i = 1:numel(Ks) u = VMD(f, 2000, 0, Ks(i), false, 1, 1e-6); scores(i) = min(envelope_entropy(u)); end plot(Ks, scores, '-o'); xlabel('K'); ylabel('最小包络熵');

把这条曲线和优化器返回的 K 对一下，如果曲线上有明显的拐点而优化器没落在拐点上，那就是适应度函数或者搜索范围的问题，不是算法的问题。

从那以后我每次改完优化器的任何一行代码，都会先固定种子跑三遍看结果稳不稳，再去动报告里的图——不然图形对不上数字的时候，回头查是代码问题还是画图问题，能耗掉一整个晚上。希望帮到你。

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