简介面向信号处理、时序分解方向的课程设计与毕业设计场景这份资料围绕鲸鱼优化算法WOA自动寻优变分模态分解VMD的模态数K与惩罚因子α展开给出可直接运行的Matlab实现并配套边际谱、中心频率图、收敛曲线及超参数迭代过程图便于观察每次迭代对分解效果的影响。压缩包共1个PDF文件约600KB主要承载代码说明、程序设计与结果展示内容结构紧凑适合快速查阅整体思路。已有1336人学习说明该主题在相关学习者中有一定关注度。代码采用参数化编程参数修改方便注释较完整数据以Excel形式组织替换数据后运行主程序即可复现实验可用于对比不同分解结果、理解WOA寻优流程并为论文或项目中的创新点设计提供参考。1. 拎着一列 Excel 数据来找我调 VMD 的人十个里有九个栽在 K 上上周有人拿一段轴承振动信号跑 VMDK 随手填了 8分解出来前四个分量几乎一模一样频谱图上几条中心频率挤在一起他把结果图放进报告里被问了一句「这四个模态有什么区别」就答不上来。这就是WOA-VMD 鲸鱼算法优化 VMD 变分模态分解要解决的核心痛点VMD 的模态数 K 和惩罚因子 alpha 靠手调调不准就是过分解或者欠分解而 Matlab 这套实现把 K 和 alpha 交给鲸鱼算法去搜用包络熵当适应度迭代几十代给你一组相对靠谱的组合。适合手上有一列时序数据、想做分解‑特征提取链条的电子信息、计算机、数学方向的同学也适合把它当课程设计、毕业设计的创新点来改。2. 项目结构与技术选型main.m 之外那几个文件各管什么2.1 目录与文件角色拆解这份资源是一套平铺式 Matlab 工程没有包管理、没有类全是函数文件加脚本好处是丢进当前文件夹就能跑坏处是函数名冲突时排查起来很烦。按职责拆开看大概是这么几层文件/层次职责你会改动它的频率main.m主程序读数据、设参数、调优化器、画图高天天动VMD.m分解核心求解变分问题返回 u、u_hat、omega低基本不动优化器文件WOA/SSA 之类种群迭代、包围、螺旋、随机搜索中换算法时动fun.m适应度函数把一组 (alpha, K) 映射成一个标量分数高换评价指标就动熵值计算函数包络熵 / 排列熵 / 样本熵中想换指标时动绘图脚本如 huatu.m收敛曲线、参数变化曲线中data.xlsx待分解时序单列高换数据就动新手最容易忽略的是 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 和 alphaVMD 参数寻优这件事本质是个黑箱优化目标函数没有解析梯度因为 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 NVMD 输出的模态矩阵每行一个分量 [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 30K 的搜索范围比较宽的时候还要再加。计算量大致是pop × Max_iter次 VMD 分解每次分解的开销跟信号长度和 K 成正比信号长度 2000 点、K 10 的时候单次大概零点几秒到一秒pop10, iter30就是 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管的是 alphalb(2)4、ub(2)10管的是 K。你要是把两个位置写反K 会被限制在 800 到 2000 之间第一次迭代就内存爆掉——这是血泪经验看到 Matlab 卡死在 VMD 第一行、内存条飙红先回来检查这个。alpha 的范围 8002000 偏保守。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 和 Infany(isnan(f))提前检查一遍4.2 参数逐项设置与含义按下表把 main.m 顶部的参数过一遍每一项我都写了推荐的起手值和改动方向参数源码值我的起手值改动方向pop51015结果波动大就加大Max_iter203050收敛曲线还在降就加大lb(1)/ub(1)alpha800/2000500/4000模态糊在一起就抬上限lb(2)/ub(2)K4/103/8按信号的物理成分数估计tau00含噪严重可试 1e-3DCfalsefalse未去均值时改 trueinit11要可复现就固定 1tol1e-61e-6想跑快改 1e-5K 的上限怎么估看原始信号的傅里叶谱有几个明显的峰簇峰簇数加一到二就是 K 的合理上限。比如一段轴承故障信号谱上能看到转频基频、二倍频、共振带三个区域那 K 大概 46 就够上限设 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)显示的数值对不上比如显示 K6.83实际出来 7 个模态。原因优化器在连续空间搜索Best_pos是浮点数。VMD 内部用 K 做循环上界Matlab 会自动向上取整所以不报错但结果错位。更麻烦的是适应度函数里如果已经round过优化器感知的和你最终用的是两个不同的 K。解决在 main.m 里显式round并且把取整后的值打印出来用这个值重新调 VMD见 3.3 节代码。同时检查fun.m里的取整逻辑和 main 里是否一致两处都用round而不是一处round一处floor否则优化出来的最优点和你实际用的点会差一个。5.3 K 搜索范围设太大导致过分解现象优化器兴高采烈地返回 K10适应度也确实低但画出来一看十個模态里有五六个是几乎一样的窄带噪声中心频率挤在一起。原因最小包络熵这个指标只奖励「有一个分量特别干净」不管其余分量死活。K 越大抽到一个干净分量的概率越高适应度自然越低。优化器没有理由不选大 K。解决三条路选一条。一是压紧 K 的上界按信号实际的物理成分数来定别超过 8二是给适应度加惩罚项fitness min(E) 0.01*K让大 K 付出代价三是把适应度换成所有分量包络熵的均值逼优化器对整体负责。我一般用第一条加第二条的组合简单直接。5.4 同一份数据两次运行结果不同现象今天跑出 K6明天同样的代码跑出 K7alpha 也变了报告里的图对不上。原因WOA 的初始化用rand螺旋和随机搜索也都在用rand、randi没有固定种子每次启动都是新的随机序列。加上init 2随机初始化中心频率的话VMD 内部也不确定。解决main.m 开头加rng(2024)init固定成 1。做完这两件事同一台机器同一个 Matlab 版本上结果完全一致。换机器或换版本可能有微小差异因为rand的实现细节会变但 K 的取值一般稳定。文献里的可复现性要求就是这么满足的。5.5 VMD 内部迭代不收敛却不报错现象程序正常跑完没有警告但中心频率图画出来乱七八糟几个模态的频率在整个时间轴上来回跳重构误差也偏大。原因tol设得太大比如 1e-3或者Max_iterVMD 内部的不是 WOA 的设得太小VMD 在还没收敛的时候就退出了返回的是中间状态。这个函数不会告诉你「我没收敛」它只是安静地返回。解决把tol拉回 1e-6 甚至 1e-7VMD 内部的迭代上限调到 500代价是慢一些。判断是否收敛的土办法是看omega的最后几行是不是稳定了如果最后两行的差还在 1e-3 量级那就是没收敛。另外信号长度太短低于 200 点也会导致收敛困难这时候不是调参能解决的得换数据段。6. 二次开发把裁判换掉比把算法换掉更有用改这套代码的方向大部分人第一反应是「换个优化算法」把 WOA 换成 GWO、PSO、蜣螂、麻雀然后在报告里写「提出了一种改进的 XX-VMD 方法」。这个路子能走通但性价比不高——所有元启发式算法在只有两个决策变量、适应度评估又贵的场景下表现差异其实没那么大你换成 PSO 大概率 K 还是那个 K。真正能拉开差距的是动适应度函数。包络熵是 2016 年前后 VMD 调参的标配但它只看单分量的稀疏性。几个更实用的替代思路适应度设计关注点适合场景最小包络熵单分量最干净基线用来对齐别人结果平均包络熵 K 惩罚整体质量想避免过分解模态间相关系数惩罚分量不冗余模态长得太像时重构误差 稀疏度加权完备性与稀疏性折中写创新点时包络熵 中心频率间距频率不重叠频谱挤在一起时改起来就在 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 对一下如果曲线上有明显的拐点而优化器没落在拐点上那就是适应度函数或者搜索范围的问题不是算法的问题。从那以后我每次改完优化器的任何一行代码都会先固定种子跑三遍看结果稳不稳再去动报告里的图——不然图形对不上数字的时候回头查是代码问题还是画图问题能耗掉一整个晚上。希望帮到你。本文还有配套的精品资源点击获取
