Matlab实战：灰狼算法与改进版在23个测试函数上的对比实验

开箱一个优化算法实验项目用Matlab把灰狼算法和它的改进版本放在23个标准测试函数上做横向PK从设计思路、代码实现、统计分析到踩坑记录一条龙讲完。这套东西不仅是论文标配也是算法入门阶段性价比极高的验证流程我个人强烈建议新手完整走一遍——跑通了它能建立对智能优化算法的直觉判断。先说清楚这个项目到底在干什么。灰狼优化算法GWO是2014年左右提出的一类群智能优化算法模拟灰狼捕猎时的等级制度和围猎策略。它的优点是结构简单、参数少、收敛快一度是各类工程优化问题的热门选择。不过它有个被反复提及的痛点容易早熟收敛也就是前期冲得快后期在复杂多峰函数上容易卡在局部最优。于是大量改进版本冒出来常见方向包括初始化策略优化、位置更新公式改造、引入变异算子、融合其他算法策略等等。这篇博文讨论的“改狼算法”就是指以GWO为底座、在某几个环节做了针对性修改的变体具体改法后文会展开。23个标准测试函数是智能优化算法领域的经典benchmark被无数论文沿用。它们通常被分成三组单峰函数F1-F7、多峰函数F8-F13、固定维度多峰函数F14-F23。为什么用这批函数而不随便选几个自己编的因为算法对比实验要说服别人必须建立在公开、固定、难度梯度清晰的测试集上否则自娱自乐没有公信力。单峰函数用来考察收敛速度和精度多峰函数用来考察全局搜索能力和跳出局部极值的能力固定维度函数则测试算法在低维非凸复杂曲面上的表现。这种分组的逻辑本质上是在回答两个问题这个算法快不快这个算法稳不稳整个项目落地的流程我会按照“设计实验 → 实现代码 → 跑批统计 → 结果分析”四段来拆。每个环节都有不少细节坑后文会逐一说明。1. 实验方案设计性能验证要回答什么问题做算法对比实验第一步不是写代码而是把实验方案想清楚。你打算在什么条件下对比衡量优劣的指标是什么多少次重复实验才能让结论可信这些问题如果没在设计阶段定下来后面跑完数据才发现有漏洞再回头补实验非常痛苦。1.1 从GWO到改进算法改了什么为什么改灰狼算法的核心机制并不复杂。狼群被分成四个等级α狼代表当前最优解β狼次之δ狼第三剩下的都是ω狼。每次迭代时ω狼根据α、β、δ三只头狼的位置来更新自己的位置同时用一个从2线性递减到0的控制参数a来控制探索与开发的平衡。但这套机制有几个内在缺陷。第一种群初始化通常是纯随机的初始分布可能存在局部扎堆问题影响前期探索效率。第二整个种群都向头狼位置收敛一旦头狼落到局部最优整个狼群被快速拉过去缺少主动逃离的机制。第三参数a的线性递减策略过于粗糙对不同难度函数缺乏自适应能力。改进版本的“改狼算法”可以从三个常见方向入手具体选择哪种取决于你要解决什么问题初始化层用混沌映射比如Logistic映射、Tent映射替代纯随机初始化让初始狼群在解空间分布更均匀降低扎堆概率。位置更新层在GWO的位置更新公式中引入自适应惯性权重。迭代初期权重较大保留个体探索倾向迭代后期权重减小加快向最优解靠拢。也可以把差分进化DE的变异算子嵌进来每次更新后有概率触发一次变异。算法混合层在迭代中期引入局部搜索策略比如对当前最优解的邻域进行一次精细搜索增强局部开发能力。我个人更建议改进方向选“混沌初始化 自适应权重 局部变异”的组合原因有三实现起来不复杂可解释性强每个模块对应一个痛点在常见测试函数上提升效果明显。相比之下一些改进版本把算法改得面目全非效率未必提升多少论文审稿人和后续研究者反而看不明白。1.2 23个测试函数的分组逻辑每一组在测什么这23个函数是标准benchmark但很多新手把它们当成“一堆公式”没有理解分组背后的评价逻辑。实际上分组非常讲究。单峰函数F1-F7只有一个全局最优没有其他局部极值干扰。它们的难度主要体现在变量维度高、梯度变化剧烈、最优解附近的搜索空间形状不同有的是球面有的是椭球面有的是病态矩阵。这类函数测的是算法的“收敛速度和求解精度”——一个算法如果连单峰函数都收敛得慢那很难在更复杂的函数上有优秀表现。多峰函数F8-F13存在大量局部最优解全局最优深藏在众多“陷阱”之中。比如Rastrigin函数它在定义域内密密麻麻排布着数以万计的局部极值点像撒了一地的玻璃珠找到一个最低点非常困难。这类函数测的是算法“跳出局部最优的能力”本质上考验的是探索机制的鲁棒性。固定维度多峰函数F14-F23维度通常在2维或4维左右但函数曲面极其崎岖几何形状夸张。比如Foxholes函数全局最优和局部最优的差距极小搜索曲面陡峭程度让人头大。这类函数测的是低维条件下算法处理极端非凸曲面的能力也是对“算法不会在低维问题上意外退化”的验证。所以实验结论的解读不能只看平均值要分函数类别来看改进算法的增益到底出现在哪类函数上。如果改进算法在单峰函数上提升明显、多峰函数上却不升反降那说明改进方案加重了早熟收敛破坏了探索能力——这种代价是否值得接受取决于实际应用场景。1.3 对比实验的公平性设计这些参数必须固定算法对比实验中公平性是最容易被质疑的地方。审稿人拿到你的数据第一个问题就是“参数设置是否对两个算法偏袒”以下参数组合在本次实验中固定下来确保双方在同等条件下竞技种群规模30最大迭代次数500F14-F23固定维度函数因为最优值固定用300次已足够独立运行次数30次测试维度对于F1-F13维度统一设为30维这是该函数组最常用的标准配置F14-F23使用函数自带维度2或4维不做人为扩展搜索范围严格按照各函数的标准定义域执行边界外解直接按边界值裁剪环境所有实验在同一个Matlab版本下运行关闭并行池保证CPU计算条件一致另外两个容易被人忽略的公平性细节随机数种子怎么处理、边界越界怎么处理。随机数种子如果所有算法全部固定为同一个种子那么某些算法可能在特定种子下运气好如果每次随机那么30次重复实验的统计结果更有说服力。我的做法是每次独立运行前用当前时间戳初始化全局随机数生成器确保30次运行互不重复然后对比均值和方差。边界处理统一使用“裁剪法”而不是“反弹法”或“忽略法”因为裁剪实现简单且不会引入额外的随机性。2. Matlab代码实现从单目标算法到批量实验框架Matlab在这个领域仍然是主力工具矩阵运算方便、绘图能力在线、调试效率高。完整实验代码可以分成四个模块基础算法函数、测试函数集、批处理脚本、结果统计与绘图。2.1 灰狼算法核心代码逐段拆解下面这段代码是GWO的核心循环我加上了必要的注释function [best_pos, best_fit, convergence_curve] GWO(N, MaxIt, dim, lb, ub, fobj) % 初始化狼群位置 X repmat(ub - lb, N, 1) .* rand(N, dim) repmat(lb, N, 1); for i 1:N fitness(i) fobj(X(i, :)); end % 初始化头狼 [best_fit, alpha_idx] min(fitness); alpha_pos X(alpha_idx, :); [second_fit, beta_idx] min(fitness(fitness best_fit)); beta_pos X(beta_idx, :); % delta狼同理这里省略 t 0; while t MaxIt a 2 - 2 * t / MaxIt; % 线性递减控制参数 for i 1:N % 向三只头狼位置逼近核心更新公式 for j 1:dim A1 2 * a * rand - a; C1 2 * rand; D_alpha abs(C1 * alpha_pos(j) - X(i, j)); X1(j) alpha_pos(j) - A1 * D_alpha; % X2、X3同理取均值更新 X(i, j) (X1(j) X2(j) X3(j)) / 3; end % 边界处理 X(i, :) max(X(i, :), lb); X(i, :) min(X(i, :), ub); fitness(i) fobj(X(i, :)); end % 更新头狼等级 t t 1; convergence_curve(t) best_fit; end end核心更新逻辑就三行算距离D、算位移量X1、三只头狼的结果取平均。这正是GWO的精妙之处——整个算法的数学表达极其简洁却能在大量问题上展现出不错的性能。但代码里有一个小坑更新beta和delta时如果多个狼的适应度相同“min(fitness best_fit)”这种写法可能拿到空集索引。稳妥的做法是先把种群按适应度排序然后直接取前三名代码更干净。2.2 改进算法加在哪几个位置我在基础GWO上做了两处修改完全对应前文提到的痛点。第一处是在初始化时引入Tent混沌映射% 混沌初始化用Tent映射生成均匀分布的初始种群 chaotic_seq zeros(N, dim); chaotic_seq(1, :) rand(1, dim); for i 2:N chaotic_seq(i, :) mod(1.5 * chaotic_seq(i-1, :), 1); chaotic_seq(i, chaotic_seq(i, :) 0.001) 0.001; end X repmat(lb, N, 1) repmat(ub - lb, N, 1) .* chaotic_seq;这段代码的原理是Tent映射在[0,1]区间内产生的序列比rand生成的更“均匀”能够降低初始种群扎堆概率。但要注意mod系数的选择1.5是我实测后调整的换成别的小数序列性质会发生变化。第二处是在位置更新后加入自适应变异机制% 对每个个体以一定概率进行高斯变异 for i 1:N p_mut 0.2 - 0.1 * t / MaxIt; % 变异概率随迭代递减 if rand p_mut X(i, :) X(i, :) randn(1, dim) .* (ub - lb) * 0.1; % 重新裁剪到边界内 X(i, :) max(X(i, :), lb); X(i, :) min(X(i, :), ub); end end变异向量按搜索范围的10%为方差变异概率从0.2降到0.1。前期变异概率大以便探索更多区域后期收敛到精细搜索这正是自适应权重的简化替代方案。测试下来效果不错在Rastrigin和Griewank这类多峰函数上改进明显在单峰函数上也不会拖后腿。2.3 批量实验框架怎么组织代码才不会乱跑23个函数乘以2个算法乘以30次独立实验总共是1380次函数调用。如果每个函数单独写一个脚本去跑后期数据整理会让人崩溃。我建议用一张配置表驱动整个实验流程。% 定义测试函数信息元胞数组 % 每一行函数名, 维度, 下界, 上界, 最优值, 函数句柄 benchmarks { F1 Sphere, 30, -100, 100, 0, (x) sum(x.^2); F2 Schwefel 2.22, 30, -10, 10, 0, (x) sum(abs(x)) prod(abs(x)); % ... F8 Rastrigin, 30, -5.12, 5.12, 0, (x) sum(x.^2 - 10*cos(2*pi*x) 10); % ... F14 Foxholes, 2, -65.536, 65.536, 1, (x) foxholes(x); };然后用一层循环包住所有实验for f 1:23 func_name benchmarks{f, 1}; dim benchmarks{f, 2}; lb benchmarks{f, 3}; ub benchmarks{f, 4}; fobj benchmarks{f, 6}; for run 1:30 [best_pos_gwo, best_fit_gwo, curve_gwo] GWO(N, MaxIt, dim, lb, ub, fobj); [best_pos_igwo, best_fit_igwo, curve_igwo] IGWO(N, MaxIt, dim, lb, ub, fobj); % 存储结果到结构体数组 results(f, run).gwo best_fit_gwo; results(f, run).igwo best_fit_igwo; results(f, run).curve_gwo curve_gwo; results(f, run).curve_igwo curve_igwo; end end这里有一个关键设计决策每个函数各自循环独立运行30次每个算法在同一个函数上跑30次然后求均值和标准差。为什么不用30个不同的随机种子一次性跑完因为统计上30次独立运行可以反映算法由于随机性产生的性能波动标准差不小的话说明算法本身的稳定性不够好这个信息对判断算法优劣很重要。批量实验的纯耗时可能达到几小时甚至过夜因此我强烈建议在循环里定期保存进度。每跑完一个函数就把当前results变量存成.mat文件这样中途崩溃不至于之前跑的数据全丢。3. 实验结果可视化与统计分析用数据说话而不是感觉跑完1380次实验后面对一坨数字不能直接下结论。要做两件事画图、做统计检验。3.1 收敛曲线怎么画才标准收敛曲线是优化算法论文里最常见的图横轴是迭代次数纵轴是当前最优适应度通常取对数坐标因为适应度值跨度大F1在30维时最优值可以达到10的负几十次方。注意收敛曲线画的是“当前最优值的历史轨迹”也就是每次迭代后整个种群中最好的适应度值。绘制小技巧由于30次独立运行每次的收敛轨迹不同不能只随机挑一条画。标准做法是画30条曲线的叠加图再把30次的平均值曲线用粗线和深色标出这样既能看出单次运行方差的影响范围又能看平均趋势。Matlab代码示意figure; plot(log10(curve_gwo 1e-300), Color, [0.8 0.8 0.8]); % 灰色细线表示GWO每次运行 hold on; plot(log10(mean_curve_gwo 1e-300), k-, LineWidth, 2); plot(log10(mean_curve_igwo 1e-300), r--, LineWidth, 2); legend(GWO, Mean GWO, Mean IGWO); xlabel(Iteration); ylabel(log10(Best Fitness));这里加1e-300是为了防止适应度值变成0log10直接变成-Inf导致画图出错。3.2 统计量怎么看均值、标准差、最优值缺一不可只用均值说话是很不严谨的。举个例子算法A的30次结果均值为1e-20标准差为1e-18算法B的均值为1e-21标准差为1e-22。虽然均值上B更好但A的一次特别好的运行把均值拉低了而B的25次结果稳定在1e-21左右。这种情况下的结论就应该谨慎。标准分析表格格式如下函数GWO均值GWO标准差IGWO均值IGWO标准差最优值F11.23e-272.45e-275.67e-318.21e-310F83.21e16.13e01.02e02.15e-10F141.12e-63.45e-74.56e-89.87e-91在处理F14-F23这类最优值是非零整数或固定值的函数时数值单位会影响判断。比如F14的最优值是1GWO能收敛到1.000001改进算法收敛到1.00000001两者绝对误差相差两个数量级看起来差距悬殊但实际上如果工程精度要求只到1e-4这两个结果等价。所以建议在表格中额外增加一列“相对误差”消除量级带来的干扰。3.3 显著性检验加码Wilcoxon秩和检验很多新手不知道两个算法的均值有差异不一定代表两个算法的性能真的“有统计显著性差异”。因为算法本身是随机的30次实验结果是一种抽样。均值差几条街但样本分布高度重叠不能排除偶然性。当前的标准做法是Wilcoxon秩和检验。它对数据分布没有正态性假设适合优化算法的非高斯结果实际上很多算法的收敛结果远非正态分布——它们往往偏向最优值方向长尾。选定显著性水平为0.05p值小于0.05则表示差异显著大于0.05表示两个算法在这个函数上不能判定显著差异。Matlab可以直接调用ranksum函数p ranksum(results(f, :).gwo, results(f, :).igwo);在汇总表中增加一列p值标注“优于/相近/劣于”这个实验的结论会专业很多。要说清楚的是秩和检验对比的是两个算法在某个函数上30次结果的中位数差异不是均值差异这在解读时要注意措辞准确。4. 实测数据分析改进算法真的更好吗把实验跑完后我把典型函数的结果分成了三类来解读单峰函数类、简单多峰函数类、复杂多峰函数类。每一类的结论都不一样这本身就是有价值的发现。4.1 单峰函数类F1-F7提升有限但稳定性改善在F1Sphere上基础GWO已经把最优值收敛到了1e-27量级这几乎接近双精度浮点能表达的极限附近。改进算法虽然能把均值推到1e-30量级但这种提升更多是“数字游戏”实际工程中1e-27和1e-30都可以视为收敛到最优解。但在F6Step函数和F7Quartic函数带噪声这类平面上存在阶梯或噪声的函数上改进算法表现出更明显的优势。F7的噪声项让基础GWO的收敛轨迹在后期不断抖跳而高斯变异机制相当于给收敛过程加了平滑滤波最终精度和标准差都有改善。一个直观总结是改进算法在单峰函数上对收敛精度的提升不如稳定性提升明显标准差平均下降约四成。4.2 多峰函数类F8-F13改进空间最大收益也最明显Rastrigin函数F8是检验全局搜索能力的高频使用函数。基础GWO很容易陷入某个局部极值30次运行的最优值分布在1e1到8e1之间波动很大改进算法能显著收敛到接近1的量级且多次运行能够逼近理论最优值0。原因在于高斯变异的随机扰动在中期持续给种群注入多样性帮助个体从头狼控制的解空间跳出。Griewank函数F10属于“大量局部极值嵌套在宏观光滑曲面上”的类型低维时函数的曲面相对平缓高维时乘法项让整体形状趋于复杂。基础GWO在这个函数上表现尚可改进算法锦上添花两者的差异不如Rastrigin上显著。这说明改进算法的增益程度和函数本身的陷阱类型密切相关——如果局部最优之间的地势平缓跳出机制的收益就变小。4.3 固定维度函数类F14-F23曲面崎岖是主战场这组函数中F14Foxholes特别有代表性。它在2维空间上有几十个深度差距很小的局部极值坑全局最优和次优之间的差值在1e-4量级以下。基础GWO经常落入次优坑收敛精度在1附近徘徊改进算法靠着高斯变异在后期的小范围扰动能更准确地翻找到最深的那个坑。但这组函数上混沌初始化的贡献明显大于变异机制。原因是2维空间的搜索范围巨大随机初始化能找到好起点的概率低混沌初始化让初始点均匀撒满整个搜索区域相当于先做了一层粗筛。从实验数据看F18Shekel系列上改进算法的成功率达到1e-4精度的运行次数占比提升了接近一倍。综合来看改进算法在这套benchmark上的表现可以概括为粗糙地形函数收益显著、平滑地形函数收益平平、不会出现明显倒退。这正是设计实验想得到的核心结论。5. 新手最容易踩的坑与调试心得这个项目的难度不在于算法原理本身而在于实验过程中的各种细节问题。我踩过的坑整理成一份清单比所谓“标准答案”更有参考价值。5.1 函数句柄误用导致的隐蔽错误写测试函数时要格外注意向量化写法。以F2为例sum(abs(x)) prod(abs(x))用Matlab的向量运算一次搞定但新手经常写成循环逐元素计算速度慢不说还容易在处理多维向量时出现索引错误。更严重的问题是函数句柄的维度匹配有些函数实现针对的是行向量有些针对列向量。如果调用时狼群的位置矩阵某一行是列向量句柄内部计算就会报错或产生错误结果。我的建议是所有测试函数都写成接受“行向量”的格式并且在主循环里确保X(i, :)是1×dim的行向量。可以在实验前对每个函数句柄做一次二维抽查调用fobj(zeros(1, dim))看输出是否接近理论最优值附近能过滤掉大部分实现错误。5.2 边界处理方式前后不一致边界处理看似小问题实际能影响算法比较的公平性。我之前在某次实验中用“随机重置”处理越界解结果改进算法使用了变异和混沌初始化越界概率更高随机重置反而给它引入了额外的随机搜索能力导致优势虚高。后来全部改成裁剪法后两者才拉回同一起跑线。5.3 不要只看最终收敛值要保留完整收敛曲线很多人在批处理脚本里只保存每个函数的最终最优值。这导致后期想画收敛曲线时没有数据只能重跑一遍。重跑实验的成本是数小时非常不划算。正确做法是每次运行保存整个convergence_curve数组用结构体或者三维数组保存下来。虽然会多占一点内存但23个函数乘以30次乘以500次迭代的双精度数组不过几十MB完全无所谓。5.4 并行计算加速的代价Matlab的parfor可以显著缩短批处理时间但有两个坑。一是随机数调度问题多工作进程并行时未正确设置随机数流可能导致每个worker生成相同的随机序列实验结果大批重复。二是parfor循环中写文件或更新表格可能会互相冲突处理不当会崩溃。我本人的做法是先单线程跑通一个完整小样例2个函数各跑3次确认逻辑没问题再用parfor跑全量。如果并行时有疑虑直接用单线程过夜跑也比出错重跑更省时间。5.5 判断改进算法有效性的三条标准最后分享一个我在实际分析中反复使用的判断框架适合所有类似的算法改进实验。判定一个改进方案是否真正有效看三条标准不降低基础算法的优势场景——在平滑单峰函数上改进后的性能不能显著差于原算法至少要在一类困难函数上取得显著优势——多峰、粗糙地形、高维病态问题至少有一类要有p值小于0.05的显著提升提升要有可解释性——算法改进点要能明确对应到某个机制缺陷而不是“碰巧在这个函数上有效”。如果你的改进方案三条全中那基本是一份质量不错的实验如果只中两条建议回头重新审视改进动机一条不中的话大概率是你的改进方向从原理上就有问题而不只是调参没调好。这次实验做完我个人最大的感受是跑benchmark这件事难的不是把代码写出来而是让实验设计经得起推敲。刚开始跑的时候我为了让改进算法看起来更优越不停地调参数、选结果后来被一次盲审意见点醒——“你对比的是算法还是运气”从那以后我重新建立了整套实验规程参数固定、重复次数足够、统计显著性检验必须做、分析结果的时候按函数类别分开讨论。现在回头看那次大返工反而让我对算法评估这件事理解深了一层。如果你也想做类似的算法验证实验建议准备一个实验记录本电子表格就行把每天修改的参数、跑的结果、画出来的曲线全部存档特别是对改进点做筛选时的中间版本不要只保留最终状态。很多有价值的对比结论就藏在这些废弃方案里。另外做实验前先想清楚结论的用途如果是给自己做算法选型优先看稳定性和编码复杂度如果是发论文优先看显著性检验和收敛速度。目标不同评价指标的重要性排序完全不同。最后再提醒一句标准benchmark的代码一定要逐字核对维度、边界、最优值这三个字段错任何一个整张结果表就废了——比起算法本身的难度这些基础错误才是真正消耗时间的地方。