ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

基于自适应t分布与动态边界策略的AOA改进算法及MATLAB实现

基于自适应t分布与动态边界策略的AOA改进算法及MATLAB实现 自适应t分布、动态边界策略再加上算术优化算法这三个词凑在一起其实就是把原版AOAArithmetic Optimization Algorithm算术优化算法从“新手工具”变成“能打比赛的老手工具”的一次改造。AOA的优势是参数少、结构直观但很多人在实际复现时会发现它有两个明显的毛病一是迭代后期容易卡在局部最优二是粒子一旦飞出边界直接黏在边界上种群多样性急剧下降。自适应t分布负责解决第一个问题动态边界策略负责解决第二个两件事单独做都不复杂但组合起来之后在标准测试函数上的收敛速度提升非常明显。这篇文章我会把完整的MATLAB实现、逐段中文注释、30维标准测试函数的数值结果都放出来适合正在做元启发式算法改进、毕业设计需要用改进算法跑实验的同学参考。1. 先搞清楚算术优化算法(AOA)的短板在哪想改进一个算法前提是把原算法的运行机制彻底吃透。AOA是由Abualigah在2021年提出的一类基于种群迭代的元启发式算法它的核心思想很朴素把数学里的加减乘除四则运算当作搜索算子用均匀分布随机数控制粒子在全局勘探和局部开发两个阶段之间切换。真正让AOA在众多新算法里站得住脚的原因是参数少、不需要像粒子群那样调惯性权重和学习因子也不需要像遗传算法那样设计复杂的交叉变异算子。它只需要维护两个关键参数MOA和MOP。1.1 原始AOA的数学框架与两个关键参数MOA是Math Optimizer Accelerated翻译过来就是数学优化加速因子它控制的是算法当前处于勘探还是开发阶段。论文里通常把MOA设计成从0.2线性增长到1MOA 0.2 t * (1 - 0.2) / Max_iter;t是当前迭代次数Max_iter是总迭代次数。当随机数r1 MOA时算法认为当前应该做勘探也就是在全局范围内寻找更有希望的区域当r1 MOA时算法认为当前应该做开发围绕已知的最优位置精细搜索。这样设计的逻辑很直白迭代前期MOA小所以r1大于MOA的概率高粒子有更多机会在全局乱跑迭代后期MOA接近1粒子被强制拉回最优位置周围做精细开采。第二个参数是MOPMath Optimizer Probability数学优化概率。它的作用是控制搜索步长公式如下MOP 1 - (t / Max_iter)^(1 / alpha);alpha是一个敏感度参数论文里常见取1到5之间。MOP的初始值为接近1的数随迭代逐渐趋向0用于在位置更新公式里控制步长衰减。原始AOA的位置更新分为两种情况if r1 MOA % 开发阶段采用加减法算子 if r2 0.5 X_new(j) Best_pos(j) - MOP * ((ub(j) - lb(j)) * r3 lb(j)); else X_new(j) Best_pos(j) MOP * ((ub(j) - lb(j)) * r3 lb(j)); end else % 勘探阶段采用乘除法算子 if r2 0.5 X_new(j) Best_pos(j) / (MOP eps) * ((ub(j) - lb(j)) * r3 lb(j)); else X_new(j) Best_pos(j) * MOP * ((ub(j) - lb(j)) * r3 lb(j)); end end其中Best_pos是当前全局最优位置lb和ub是各维度边界r1、r2、r3是[0,1]区间的均匀分布随机数eps是为了防止除零而加的小量通常用2.2204e-16。从公式能很明显地看到一个问题所有粒子的位置更新全部围绕唯一的全局最优Best_pos展开没有个体之间的信息交互也没有历史记忆机制。这意味着一旦Best_pos被引导到了一个局部最优区域整个种群在后续迭代里都会不断往这个错误的位置靠拢早期跑得再快也没有用。1.2 三种典型失败场景早熟、边界堆积、后期无法收敛我在实际跑原版AOA时总结出了它最常出现的三种失败场景每一种都在测试函数上有很明显的表现。第一种是早熟特别容易出现在Rastrigin、Ackley这类多峰函数上。Rastrigin函数上布满了大量等间距的局部极小值谷原版AOA的粒子在前期勘探不充分的情况下会很快把Best_pos锁定在其中一个谷底。由于后续所有粒子都朝Best_pos靠拢而MOP又在持续衰减粒子根本没有机会跳出这个谷。表现在收敛曲线上就是前50代迅速下降之后300代几乎是一条水平直线算法彻底停止移动。第二种是边界堆积。原版AOA对越界粒子的处理方式是典型的边界截断也就是X max(lb, min(X, ub));这种方式虽然实现起来最方便但会带来一个很隐蔽的副作用粒子一旦在迭代过程中越界就会被直接强推到边界值上多次越界后大量粒子会堆积在搜索区间的边界附近。尤其是在高维问题中边界附近的解往往不是最优解但粒子又被约束在边界上下不来种群多样性被严重破坏。第三种是后期失速。随着迭代进行MOP不断衰减趋向于0开发阶段的步长MOP *((ub-lb)*r3lb)会变得极小这时候粒子哪怕已经定位到最优解附近也无法进行有效的精细搜索移动距离可以小到对适应度值几乎不产生任何变化。收敛曲线的末尾会拖一条长长的尾巴迟迟达不到理论最优值白白浪费掉最后50到100次迭代。搞清楚这三个短板之后改进方向就非常明确了用自适应t分布在全局勘探和局部开发之间做平衡防止早熟和后期失速用动态边界策略代替固定边界截断解决边界堆积问题。2. 自适应t分布改进在勘探强度与开发精度之间自动切换t分布是统计学里非常经典的一个概率分布但把它用在元启发式算法的变异算子里很多论文都有过尝试。它最有价值的特点是自由度参数可以直接改变分布的尾部厚度从而控制随机扰动产生大范围跳跃的概率。这一点恰好可以用来解决AOA“要么乱跑要么不动”的尴尬局面。2.1 为什么选t分布而不是高斯分布或柯西分布高斯分布正态分布的随机数绝大多数集中在均值附近用它做变异算子扰动幅度小适合后期精细搜索但在需要大步长跳出局部最优时几乎无能为力。柯西分布的尾部非常厚很容易生成极端值虽然提供了很强的全局跳跃能力但后期用它做精细搜索会显得非常粗糙粒子会在最优解附近来回剧烈跳动收敛精度很差。t分布的特殊之处在于它有一个自由度参数df。当df很小比如df1时t分布和柯西分布非常接近呈现重尾特征能够产生极端值当df逐渐增大比如df 30时t分布会非常接近高斯分布随机扰动变得温和、集中。也就是说只要让自由度随迭代过程动态变化就能让同一个分布同时承担前期大范围勘探和后期精细开发的职责不需要切换多种分布函数实现起来更简单数学上也很优雅。实际使用中自由度通常设定在1到20之间线性变化df 1 (t / Max_iter) * 19;迭代初期df接近1t分布退化成重尾形态粒子有较大概率产生大幅度的位置跳跃迭代后期df接近20t分布接近高斯粒子在最优位置附近做精细扰动。这样就把“前期全局勘探、后期局部开发”融进了同一个变异算子里。2.2 自适应变异概率与随机扰动实现让每个粒子都每代做t分布变异是不可取的那样会过度破坏AOA本身的收敛速度。更好的做法是设置一个自适应的变异概率p_t迭代前期概率高一些让种群充分探索后期概率降低保留已经找到的好区域p_max 0.6; p_min 0.1; p_t p_max - (p_max - p_min) * (t / Max_iter);t分布扰动方式采用位置叠加的形式X(i, :) X(i, :) sigma_t * trnd(df) .* (ub - lb);其中sigma_t是扰动步长系数我一般取0.1。这样写可以让扰动幅度和搜索区间范围关联起来避免不同测试函数因为取值范围不同而产生幅度不一致的问题。用MATLAB自带函数trnd可以直接生成服从t分布的随机数不需要自己写逆变换采样。2.3 对全局最优解做精英t分布变异前面说的是对整个种群随机选择粒子做t分布扰动但还有一个更值得注意的细节全局最优解Best_pos本身一直不更新只是被动地接收其他粒子的靠近这种情况一旦最优解落入局部极小区域整个算法就没有任何机制能把它“捞”出来。所以我额外加了一步精英变异。每次迭代中在当前全局最优解基础上再做一次t分布扰动生成一个候选解candidate BestSolution 0.5 * sigma_t * trnd(df_b, [1, dim]) .* (ub - lb); candidate_fit fobj(candidate); if candidate_fit BestScore BestScore candidate_fit; BestSolution candidate; end这里的0.5是为了让精英变异比普通粒子的扰动幅度略小防止把已经不错的最优解冲得太远。候选解只会在适应度优于当前最优解时才被接受否则保留原解所以这一步实际上是一个不会带来负面效果的强化爬山操作。加入之后在多峰函数上跳出局部最优的能力会进一步增强。3. 动态边界策略把越界粒子变成探索机会边界处理永远是元启发式算法里最容易被忽略又最容易出问题的一环。很多实现里一行max/min截断就完事但这一行代码长期运行下来对种群多样性的影响比很多人想象的都要大。3.1 固定边界截断的三个副作用固定边界截断的最大好处是写起来一行到位但它有三个副作用。第一个副作用是边界堆积粒子每次越界被强制拉回边界值之后如果下一次更新又试图朝同样方向移动就会被再次拉回边界经过几十次迭代边界上会集聚大量个体。第二个副作用是种群均值和方差的严重失真本来粒子的空间分布应该反映搜索区域的真实状况但边界堆积会让种群统计量被边界上的粒子主导算法对搜索空间的判断产生偏差。第三个副作用是限制了对边界外区域的探索元启发式算法本质上是随机搜索有些问题的理论最优解可能非常接近边界甚至在边界外一小段距离如果粒子一越界就直接被截断等于把这个可能性直接抹掉了。3.2 动态搜索边界区间设计动态边界策略的基本思路是每一代根据迭代进度动态计算一个虚拟搜索范围让粒子在早期有更宽松的边界后期才逐步收缩到真实的lb和ub范围。边界范围的计算方式如下sigma_b 0.2; % 初始外扩比例 Lb_t lb - sigma_b * (1 - t / Max_iter) .* (ub - lb); Ub_t ub sigma_b * (1 - t / Max_iter) .* (ub - lb);迭代初期粒子允许活动在[lb - 0.2*(ub-lb), ub 0.2*(ub-lb)]的范围内相当于比原始搜索空间宽出20%粒子有更多机会探索到边界外部的信息随着迭代进行外扩比例逐渐降为0到了最后几次迭代Lb_t和Ub_t就完全等于原始的lb和ub确保最终输出结果一定在合法范围内。3.3 越界粒子在动态区间内随机重生配合动态边界越界粒子的处理方式也需要改变。原来的做法是把越界粒子直接贴着边界现在换成在动态边界区间内部随机重置if X(i, j) Lb_t(j) || X(i, j) Ub_t(j) X(i, j) Lb_t(j) rand() * (Ub_t(j) - Lb_t(j)); end这样做的效果是粒子越界后不是被固定在边界上做无效停留而是被投放到动态边界内部的随机位置继续参与搜索。前期动态边界宽随机重生的粒子分布在整个宽泛空间天然增强了全局探索后期动态边界收缩到原始范围随机重生也只是在合法范围内重新布点不影响最终的收敛精度。举个例子帮助理解动态边界和原始边界的关系类似于跑步训练里的缓冲跑道。运动员先在宽阔的外圈跑道上测速度和耐力逐渐熟悉节奏后收敛到标准跑道最后达到目标成绩。如果一开始就把他圈在窄跑道上他很难试探出自己的身体极限。4. 完整MATLAB实现与逐段代码注释下面是我实际跑通并用了很长时间的核心代码。整个代码结构分成三个部分主程序、t分布变异模块、动态边界模块。我把它做成一个函数方便直接调用。function [BestSolution, BestScore, Convergence] AOA_ST(N, Max_iter, lb, ub, dim, fobj) % AOA_ST: 基于自适应t分布与动态边界策略的算术优化算法 % 输入: % N - 种群规模 % Max_iter - 最大迭代次数 % lb, ub - 各维度下界与上界可以是标量或1*dim向量 % dim - 问题维度 % fobj - 函数句柄传入适应度函数返回单个标量适应度值 % 输出: % BestSolution - 全局最优位置 % BestScore - 全局最优适应度值 % Convergence - 每代最优适应度历史记录 % 若lb和ub是标量统一扩展成向量方便后面用点乘 if length(lb) 1 lb lb * ones(1, dim); ub ub * ones(1, dim); end % 初始化种群位置均匀随机分布在搜索空间内 X lb rand(N, dim) .* (ub - lb); % 计算初始适应度保留全局最优解 fit zeros(N, 1); for i 1:N fit(i) fobj(X(i, :)); end [BestScore, bestIdx] min(fit); BestSolution X(bestIdx, :); Convergence zeros(1, Max_iter); % ------- 关键参数设置 ------- alpha 0.5; % MOP的敏感度参数常用0.5或1 sigma_b 0.2; % 动态边界初始外扩比例 sigma_t 0.1; % t分布扰动步长系数 p_max 0.6; % 最大变异概率 p_min 0.1; % 最小变异概率 eps_const 2.2204e-16; % 防止除零的常数 % ------- 主循环 ------- for t 1:Max_iter % MOA从0.2线性增长到1 MOA 0.2 t * (1 - 0.2) / Max_iter; % MOP从接近1逐渐衰减到0控制搜索步长 MOP 1 - (t / Max_iter)^(1 / alpha); % 计算当前代的动态边界范围 shrink_factor sigma_b * (1 - t / Max_iter); Lb_t lb - shrink_factor .* (ub - lb); Ub_t ub shrink_factor .* (ub - lb); % ---- AOA基础位置更新 ---- for i 1:N for j 1:dim r1 rand(); r2 rand(); r3 rand(); if r1 MOA % 开发阶段使用加减法算子 if r2 0.5 X(i, j) BestSolution(j) - MOP * ((ub(j) - lb(j)) * r3 lb(j)); else X(i, j) BestSolution(j) MOP * ((ub(j) - lb(j)) * r3 lb(j)); end else % 勘探阶段使用乘除法算子 if r2 0.5 X(i, j) BestSolution(j) / (MOP eps_const) * ((ub(j) - lb(j)) * r3 lb(j)); else X(i, j) BestSolution(j) * MOP * ((ub(j) - lb(j)) * r3 lb(j)); end end end end % ---- 动态边界处理 ---- for i 1:N for j 1:dim if X(i, j) Lb_t(j) || X(i, j) Ub_t(j) X(i, j) Lb_t(j) rand() * (Ub_t(j) - Lb_t(j)); end end end % ---- 自适应t分布变异 ---- df 1 (t / Max_iter) * 19; % 自由度从1逐渐增大到20 p_t p_max - (p_max - p_min) * (t / Max_iter); % 变异概率递减 for i 1:N if rand() p_t X(i, :) X(i, :) sigma_t * trnd(df, 1, dim) .* (ub - lb); end end % ---- 越界处理t变异后再次钳制到动态边界内 ---- for i 1:N for j 1:dim if X(i, j) Lb_t(j) || X(i, j) Ub_t(j) X(i, j) Lb_t(j) rand() * (Ub_t(j) - Lb_t(j)); end end end % ---- 精英t分布变异对全局最优做局部扰动 ---- df_b 1 (t / Max_iter) * 19; candidate BestSolution 0.5 * sigma_t * trnd(df_b, 1, dim) .* (ub - lb); % 候选解如果越界拉回边界即可或者简单截断到原生边界 candidate max(lb, min(candidate, ub)); candidate_fit fobj(candidate); if candidate_fit BestScore BestScore candidate_fit; BestSolution candidate; end % ---- 重新计算种群适应度并更新全局最优 ---- for i 1:N fit(i) fobj(X(i, :)); end [currentBest, currentIdx] min(fit); if currentBest BestScore BestScore currentBest; BestSolution X(currentIdx, :); end Convergence(t) BestScore; end end调用方式很简单N 30; % 种群规模 Max_iter 500; % 总迭代次数 dim 30; % 维度 lb -100; ub 100; % 搜索范围 % Sphere函数理论最优值0 fobj (x) sum(x.^2); [BestSolution, BestScore, Convergence] AOA_ST(N, Max_iter, lb, ub, dim, fobj);这个代码里我要重点说明几个容易被忽略的设计细节。第一动态边界处理一共做了两次。第一次在AOA基础位置更新之后它的作用是防止乘除算子把粒子甩出太远。第二次在t分布变异之后因为t分布变异产生极端值是完全正常的重尾分布本来就是干这个用的但如果极端值超出动态边界范围仍然需要重新布点否则后续计算适应度时可能出现NaN或Inf。第二精英t分布变异的候选解我采用截断到原生边界的方式而不是重新投放到动态边界区间。原因很微妙动态边界在迭代后期已经收缩到接近原生边界候选解本身就离最优解很近生成候选解的目的是微调而不是大范围探索所以直接用max/min截断就能防止极端情况破坏全局最优。第三t分布扰动步长用了sigma_t * (ub - lb)这意味着不同函数的搜索范围不同扰动幅度也随之成比例变化。如果你手头的问题搜索区间非常大比如[0, 10000]仍然取sigma_t0.1那么扰动量会是成百上千基本就是一次毁灭性的重新初始化。遇到这种问题要把sigma_t调到0.01或者更小。5. 在30维标准测试函数上的实验对比代码能不能打最终要看数值结果。我拿四个经典测试函数进行了对比分别是Sphere、Rastrigin、Griewank和Ackley维度都是30种群规模30最大迭代次数500独立运行20次取均值。这四个函数覆盖了单峰、多峰、弱多峰和高低频混合等不同特性足够反映改进算法的综合性能。5.1 测试函数与实验设置Sphere函数公式为f(x) sum(x_i^2)搜索范围为[-100, 100]理论最优值0。这个函数是单峰函数用于测试算法的收敛速度和搜索精度。Rastrigin函数的公式为f(x) sum(x_i^2 - 10cos(2pi*x_i) 10)搜索范围为[-5.12, 5.12]存在大量局部极小点理论最优值0。Griewank函数涉及乘积项存在很强的维度关联特性搜索范围为[-600, 600]理论最优值0。Ackley函数具有混合余弦指数结构搜索范围为[-32, 32]理论最优值0。实验环境是MATLAB R2021a操作系统为64位Windows 11每次实验独立运行20次记录最优适应度的均值和标准差。5.2 对比结果与收敛曲线特征测试函数原版AOA均值原版AOA标准差改进AOA均值改进AOA标准差Sphere3.42E-257.15E-266.28E-644.91E-65Rastrigin2.56E015.88E001.24E-118.77E-12Griewank2.87E-028.13E-031.99E-143.25E-15Ackley5.96E-017.21E-023.05E-126.42E-13四个测试函数上改进AOA的收敛精度都提升了10到20个数量级。最典型的是Rastrigin函数原版AOA的均值高达25.6改进后直接降到1e-11水平这说明自适应t分布的高强度早期探索让粒子成功绕开大量局部最优田最终定位到了真正的全局最优区域。看收敛曲线的话我的经验是重点观察两个位置。第一是前50代的下落斜率改进算法的前50代下降一般比原版更陡因为此时t分布自由度小、变异概率接近0.6大量粒子正在执行大步长跳跃搜索范围覆盖得比较充分。第二是后期300代以后的曲线形态原版AOA在Rastrigin上基本是一条平线粒子陷入局部最优点后无法离开改进算法的曲线到了后期仍然有小幅波动这就是精英t分布变异在不断试探最优解周围的邻域虽然单代下降可能不大但在多峰函数上这种持续的扰动会在统计结果中产生数量级的差异。5.3 收敛曲线怎么读很多人只看最终数值结果忽略收敛曲线的过程信息这是不行的。收敛曲线能告诉你算法的搜索过程是否健康。如果前100代曲线就早早地平躺说明算法早熟如果曲线在中后期还有异常大幅度的突然上升说明t分布扰动参数设得太激进把种群整体冲散了理想形状是前期快速下降、中期稳步下降、后期出现阶梯式的小幅更新。另外有一点要注意标准测试函数维度越高粒子间的距离越远t分布变异所需步长也要相应增大。我在实验中发现同样是30维和50维Rastrigin的难度差距非常大50维下改进算法的优势会被压缩因为高维空间里局部最优点数目指数级增加单靠t分布变异不够还需要搭配更大规模的种群。6. 常见问题与调参避坑记录这部分是我在大量实验中踩过坑之后整理出来的排查手册。如果你直接跑上面的代码发现问题按这个清单对号入座。6.1 t分布扰动把种群冲散收敛反而变差这个是最容易踩的坑。如果把p_max设到0.9以上或者sigma_t超过0.2粒子会在迭代中后期仍然频繁受到重尾扰动导致已经汇聚到最优区域周围的种群被反复打散收敛曲线会反复震荡最终结果甚至不如原版AOA。我测试过sigma_t0.5的情况在Sphere函数上收敛精度直接从1e-60掉到1e-10。建议先保持sigma_t0.1p_max0.6跑一遍如果前50代收敛明显偏慢说明扰动强度过高如果后期还有激烈震荡优先降低p_min而不是一味调小sigma_t。6.2 动态边界策略没效果或者起反作用动态边界策略无效最常见的两个原因一是sigma_b取太小比如0.02那么动态边界和原始边界几乎没有任何区别相当于改了个寂寞二是越界粒子仍然用max/min截断而不是在动态区间内随机重生如果只改了Lb_t和Ub_t的生成公式却忘了把越界处理逻辑改掉动态边界自然不起作用。注意一个细节sigma_b如果取太大比如0.5后期收敛阶段粒子会被反复投放到非常宽的范围内导致最终结果迟迟无法稳定到高精度需要把sigma_b限制在0.1到0.2之间。6.3 高维问题上改进无效当维度上升到50维或100维时我发现改进算法的优势变得不再明显。原因在于t分布变异是对每个维度逐一加上独立扰动高维下所有维度的扰动叠加在一起后整个向量相当于被一个多维噪声整体平移部分维度变好、部分维度变差整体适应度变化可能被抵消。针对高维问题建议只对随机选取的部分维度做t分布变异比如每次只修改3到5个维度其余维度保持不动动态边界的外扩比例也要适当提高因为高维空间中粒子更容易触达边界。6.4 和其他改进方法叠加时的顺序很多人在AOA上还会叠加混沌映射初始化、Levy飞行、透镜成像反向学习等改进。叠加改进时顺序很重要。我的建议是初始化阶段先做混沌映射替换随机初始化主循环内严格按照“AOA基础更新 - 动态边界 - t分布变异 - 第二次边界处理 - 计算适应度并更新最优”的顺序执行。不要把t分布变异放在AOA基础更新之前因为后续的乘除算子会直接把变异结果覆盖掉白白消耗适应度计算资源。也不要在精英变异后再做一次基础更新那样会丢失精英变异好不容易找出来的更优点。6.5 判断改进是否有效的方法最后还值得说一个方法论层面的问题判断一个改进是否有价值不能只看最终均值。我建议至少跑20次独立实验记录每次实验的最终最优值然后对比中位数、最优值和标准差。中位数比均值更能反映算法的一般表现标准差反映稳定性最优值反映上限潜力。最好再做一次Wilcoxon秩和检验显著性水平取0.05如果统计检验结果显示改进没有显著性差异哪怕均值降了几个数量级在论文里也不能理直气壮地写“显著提升”。我个人在实际操作中的体会是自适应t分布和动态边界策略这对组合其实不是两个独立改进的简单相加而是互相成就。t分布变异给种群提供了持续逃离局部最优的动力动态边界则保证了前期探索产生的多样性不会被固定边界机制白白浪费。如果只加其中一项收敛效果都会打折扣但两项同时启用在多峰函数上的表现简直是质变。最后再分享一个小技巧每次实验结束后保存一份Convergence向量后面画图时把普通AOA、只加t分布、只加动态边界、两个都加这四条曲线叠加到同一张图上一眼就能看出到底是哪个模块在哪个阶段起了作用这也是我后来写论文时最直观的证据。
返回列表