SABO这个算法刚出来的时候,我是挺喜欢的。结构干净,没有GWO那种头狼引导,也没有WOA那么花哨的螺旋机制,就是纯粹拿种群个体之间的“减法差”来做更新,理解成本很低。但我在Rastrigin函数上跑了几次就发现问题:它中后期收敛太容易陷入平台期,种群一旦抱团,减法平均方向就趋近于零向量,所有人站在同一个局部谷底里互相观望。那时候我就在想,能不能在不破坏SABO减法机制的前提下,给它加一个外部引导,再把初始种群做得更散一点。折腾了一段时间,我采用了基于黄金正弦(Golden Sine)和Sine混沌映射的改进思路,得到了一版混合改进算法,暂名GSABO。这篇东西主要记录的是我的改进方案、Matlab实现,以及和SABO、GWO、WOA、DBO用同一套配置跑出来的对比结果。你要是正在给智能优化算法做创新改动、准备相关实验,或者单纯想把这类改进算法的Matlab流程跑通,这篇可以直接当参考骨架。
1. SABO的短板从哪来:减法平均机制的得与失
1.1 减法平均更新式的核心逻辑
SABO(Subtraction-Average-Based Optimizer,减法平均优化器)的核心是定义了一个带随机缩放的减法算子:
A ⊖ B = r ⊙ (A - B)
其中r是[0,1]区间里的随机数向量,⊙表示逐元素乘法。基于这个算子,种群中第i个个体的候选新位置可以写成:
X_i^{new} = X_i + r_i ⊙ ( (1/N) Σ_{j=1}^{N} (X_i ⊖ X_j) )
把这个公式拆开看,本质就是:把当前个体和种群内其他每一个个体做差,乘上一个随机权重,再取平均,得到移动方向。这里有个容易误解的地方——它不是简单的“种群中心牵引”。因为X_i - X_j和X_j - X_i方向相反,随机缩放系数又相互独立,所以求和结果带有很大的随机性,并不是严格朝中心走。正是这种随机性让它在迭代前期保留了大量探索行为,但也为后期的收敛停滞埋下了伏笔。
它和差分进化(DE)的区别也值得说一句。DE是随机从种群中抽两三个个体构造差分向量,SABO则把整个种群都纳入计算,空间几何信息用得更加充分。在种群分散的前期,差分项规模大,个体能快速逼近较优区域;同时因为没有精英个体引导,它不会轻易被某个特别好的解带偏,这是SABO在简单函数上表现稳定的原因。
单从机制上看,SABO确实比很多“头狼领导制”算法更简洁:不需要排等级,不需要判断谁是alpha谁是beta,不需要设计螺旋路径,一个减法算子走天下。但简洁往往意味着对种群状态的依赖更大,当种群本身失去多样性的时候,这个算法几乎没有任何自保手段。
1.2 为什么容易在中后期停滞
我在调参时专门观察过减法平均项的范数随迭代的变化,规律非常典型:前50代左右,减法平均项的模长在数量级上还能推动个体快速移动;到200代以后,种群开始在一个局部吸引域内聚集,任意两个个体之间的差都变得很小,减法平均项整体趋近于零向量。这时候更新公式基本退化成“原地踏步”,最优值曲线出现一条长长的水平线。
这个问题本质上是种群差分方法共有的“对称性陷阱”:当所有个体都挤在同一个局部最优附近时,任何“个体之间的差值”都只能在这片局部区域里小打小闹,没有任何机制把种群往外拉扯。GWO通过alpha、beta、delta三只头狼的不同距离保持扩散压力,WOA通过螺旋收缩和随机搜索交替来维持探索,DBO干脆把种群分成多个分工小组各干各的。SABO缺的正是这种“外来的、不依赖当前群体差异”的引导力。
打个比方说,减法平均机制像一群人围在广场上商量往哪走,每个人都说“我看大家好像都在这块区域,那就往中间挪一点吧”。前期大家分散在各处时,这个策略有效;但等所有人都挤进同一个死胡同,互相之间只能看到邻居的后脑勺,整个群体的位移就趋近于零。
顺着这个思路,我给SABO补了两样东西:一个是初始化阶段就把种群铺得更均匀的Sine混沌映射,另一个是能在每一代提供最优解邻域振荡搜索的黄金正弦算子。下面分开讲。
2. 用黄金正弦和混沌映射,分别补哪块短板
2.1 黄金正弦:一段自带振荡的探索区间
黄金正弦算法(Golden Sine Algorithm,Golden-SA)的原始更新式是:
X_i^{t+1} = X_i^t · |sin(r1)| - r2 · sin(r1) · |x1 · X_best^t - x2 · X_i^t|
其中r1取自[0, 2π],r2取自[0, π];τ = (√5 - 1)/2是黄金分割率,x1 = -π + (1 - τ)·2π,x2 = -π + τ·2π,这两个是固定常数系数。
这个公式第一眼看上去有点绕,拆开就清楚了:第一项用|sin(r1)|对当前位置做幅值压缩,第二项以最优解和当前个体的距离为基准,再乘正弦系数形成振荡搜索。|sin(r1)|在0到1之间来回波动,r2控制方向,个体就可能小步逼近最优解,也可能跳得比较远做探索。黄金分割率的加入保证了步长序列不呈周期性重复,学过一点数论都知道,无理数驱动的采样序列在有限范围内最不容易“撞车”。
把这条思路搬进GSABO时,我没有直接照搬整条Golden-SA更新式,而是只把黄金正弦项作为“精英引导”加进SABO的更新里。原因很简单:SABO已经很擅长利用种群差分信息,不需要推翻重来;它缺的只是中后期“外部供给”的位移能量。黄金正弦项里始终存在x1·X_best - x2·X_i这一项,它和种群离散度无关,即便所有个体都挤成一团,这项依然不为零。这就是它补短板的核心价值。
2.2 Sine混沌映射:让初始种群更有“分散度”
再看初始化。常规随机初始化在Matlab里就是一行代码:X = lb + rand(N, dim) .* (ub - lb)。简单,但有个隐性缺点:小种群、高维度下,随机点很容易在某个局部区域扎堆,尤其当目标函数是多峰函数时,初始位置直接决定你大概率掉进哪个峰。
混沌初始化是经典的改进手段之一。我选用的是Sine混沌映射,迭代式:
x_{k+1} = a · sin(π · x_k)
其中a一般取接近4的值,x0不能取0。Sine映射在一个周期内产生的序列能覆盖整个(0,1)区间,而且相邻点的相关性低。相比Logistic混沌,Sine映射的实现代码更短,数值稳定性也更好,实测中初始种群在各个维度上的边缘分布明显比rand更均匀。
有一点必须提醒:混沌序列虽然好,但不要所有维度共享同一条序列。如果直接用同一个x序列给N×dim的种群矩阵填数,各维之间会引入很强的线性相关,种群形状会被压到一条低维流形上。稳妥做法是给每个维度生成一条独立序列,或者从同一条母序列里做带偏移的采样,后面第四节代码里我会给出具体写法。
2.3 组合改进的整体思路
整个GSABO的改进框架可以概括成一句话:用Sine混沌序列把初始种群铺均匀,用黄金正弦算子提供外部精英引导,用SABO原本的减法平均机制保留种群差分信息利用。三个部件各管一段:混沌管初始质量,黄金正弦管中后期活力,SABO管整体搜索框架。
这种“打补丁”式改进的好处是风险低。如果你把SABO的减法平均项完全去掉,换成黄金正弦,那其实就是一个套壳的Golden-SA,种群差分这个SABO独有的优势就丢了。反过来,如果只加一个混沌初始化而不管中后期停滞,那改进效果在前100代能看出来,之后就又回到老样子。两样东西叠加,才能真正改变收敛曲线的后半段形态。
3. GSABO的具体设计:三个改进点的数学表达
3.1 混沌初始化:Sine序列映射到搜索空间
我采用的初始化步骤如下:
- 取一个随机种子x0∈(0,1),且x0不能取0,否则Sine映射会直接退化成一条全零序列;
- 对每个维度分别生成一条长度为N的Sine混沌序列,各维度序列的起始位置加一个微扰偏移,避免维度之间出现线性相关;
- 对每条序列做mod(seq, 1)兜底归一化,得到u∈(0,1);
- 映射到搜索空间:X(i,d) = lb(d) + u(i,d) · (ub(d) - lb(d));
- 如果不同维度上下界不同,映射时逐维度处理即可。
这样生成出来的初始种群在每个维度的分布近似均匀,没有rand初始化那种“肉眼可见的团聚感”。在Rastrigin这类函数上,初始解散得越开,前期掉进同一个局部谷底的概率就越低,对最终结果的影响是实打实的。
3.2 融合黄金正弦的更新式:动态权重λ
在我的复现版本里,GSABO的个体更新式写成:
X_i^{new} = X_i + λ(t) · D_i^{SABO} + (1 - λ(t)) · G_i^{Golden}
其中:
D_i^{SABO} = (1/N) Σ_{j=1}^{N} r_j ⊙ (X_i - X_j)
G_i^{Golden} = |sin(r1)| · (x1 · X_best - x2 · X_i) · r2
λ(t) = 1 - 0.7 · (t / Tmax)
先说λ的设计意图。t=0时λ=1,几乎完全依赖减法平均项;随着迭代推进λ逐渐减小到0.3,黄金正弦项的权重上升。有人可能问:为什么前期不给黄金正弦大权重?因为迭代前期种群差异足够大,SABO自身就能高效收敛,这时如果让最优解引导过强,反而会加快早熟;到了后期减法差分趋零,才需要黄金正弦强制制造位移。当然这是经验性的配置,不同测试函数上可以做微调。
G_i^{Golden}里的r1和r2取值也有讲究。原始Golden-SA中,r1在[0,2π]上均匀采样,r2在[0,π]上均匀采样。我实测过后发现,把r2改成[0,1]之间的均匀分布也能工作,但收敛曲线会更毛糙,所以还是保留原始区间。如果你希望改进效果更平滑,可以用混沌序列生成r1和r2来替代rand,代价是每代增加一点计算量,收益看实际工程需求。
3.3 边界处理与贪婪选择
边界处理我选择“越界随机重置”而不是“越界截断”。截断会让大量个体贴在边界上,而边界附近如果恰好不是最优区域,这些个体就变成了停滞个体;随机重置则保证每个代都有一定比例的个体从空域重新切入搜索,相当于一种隐性的种群重启机制。
位置更新之后,做标准的贪婪选择:
if fobj(X_i^{new}) < fobj(X_i) X_i = X_i^{new}; end
这一步必须保留。融合算子的随机振荡可能会把已经找到的好解踢飞,没有贪心判断的话,算法最优值的波动会很大。我在复现不少改进算法时发现效果不稳定,根源往往就是贪心选择写得不对或漏了这个精英保持机制。
4. Matlab代码架构与关键实现
4.1 主实验脚本的设计
实际项目里我把代码分成两层:第一层是主实验脚本main_compare.m,负责定义测试函数、统一参数、循环调用多个算法;第二层是各个优化算法函数,统一输入格式为适应度函数句柄、维度、上下界、种群数、最大迭代,输出为最优位置、最优适应度、收敛曲线。
主脚本里几个参数我统一固定:种群N=30,维度dim=30,最大迭代Tmax=500,独立运行Runs=30。对比算法包括原始SABO、GWO、WOA、DBO和GSABO。每个算法都跑同一个fobj,最后取Runs次结果的平均值和标准差做对比。收敛曲线我取30次中的中位数曲线来画,避免个别异常次数把图像拉偏。
% main_compare.m 统一对比脚本示例 clear; clc; rng('shuffle'); fobj = @(x) sum(x.^2); % Sphere测试函数 lb = -100; ub = 100; % 搜索范围 dim = 30; N = 30; Tmax = 500; Runs = 30; algorithms = {@SABO, @GWO, @WOA, @DBO, @GSABO}; names = {'SABO', 'GWO', 'WOA', 'DBO', 'GSABO'}; result = zeros(length(algorithms), Runs); for a = 1:length(algorithms) for r = 1:Runs [~, Best_fit, ~] = algorithms{a}(fobj, dim, lb, ub, N, Tmax); result(a, r) = Best_fit; end fprintf('%s: mean=%.4e std=%.4e\n', ... names{a}, mean(result(a,:)), std(result(a,:))); end这里要提醒一句:不要用“第一次跑出来的那条曲线”直接当最终结果。优化算法的单次运行方差可以很大,尤其在多峰函数上,至少跑30次再统计才有说服力。我自己早期偷懒跑3次就下结论,结果换个种子结论就翻车,血的教训。
4.2 GSABO核心函数的逐行实现
下面给出我最终使用的GSABO函数,可以直接复制到Matlab里跑。注释里标明了三个改进点的位置。
function [Best_pos, Best_fit, Curve] = GSABO(fobj, dim, lb, ub, N, Tmax) % GSABO: 黄金正弦 + Sine混沌初始化 改进的减法平均优化器 % 输入: % fobj - 适应度函数句柄 % dim - 决策变量维度 % lb,ub - 下界/上界,可为标量或维度向量 % N - 种群大小 % Tmax - 最大迭代次数 % 输出: % Best_pos - 最优位置 % Best_fit - 最优适应度 % Curve - 收敛曲线 % 黄金分割系数 tau = (sqrt(5) - 1) / 2; x1 = -pi + (1 - tau) * 2 * pi; x2 = -pi + tau * 2 * pi; % 统一上下界为行向量 if numel(lb) == 1, lb = repmat(lb, 1, dim); end if numel(ub) == 1, ub = repmat(ub, 1, dim); end % 改进点1: Sine混沌初始化 X = chaotic_init(N, dim, lb, ub); % 初始适应度评估 Fitness = zeros(N, 1); for i = 1:N Fitness(i) = fobj(X(i, :)); end [Best_fit, idx] = min(Fitness); Best_pos = X(idx, :); Curve = zeros(1, Tmax); for t = 1:Tmax lambda = 1 - 0.7 * (t / Tmax); % 动态融合权重 for i = 1:N % 改进点2: SABO减法平均项 D = zeros(1, dim); for j = 1:N if j ~= i r = rand(1, dim); D = D + r .* (X(i, :) - X(j, :)); end end D = D / N; % 改进点3: 黄金正弦引导项 r1 = rand * 2 * pi; r2 = rand * pi; G = abs(sin(r1)) * r2 .* (x1 * Best_pos - x2 * X(i, :)); % 融合更新 Xnew = X(i, :) + lambda * D + (1 - lambda) * G; % 边界处理: 越界随机重置 idx_lo = Xnew < lb; idx_hi = Xnew > ub; if any(idx_lo) Xnew(idx_lo) = lb(idx_lo) + rand(1, sum(idx_lo)) .* ... (ub(idx_lo) - lb(idx_lo)); end if any(idx_hi) Xnew(idx_hi) = lb(idx_hi) + rand(1, sum(idx_hi)) .* ... (ub(idx_hi) - lb(idx_hi)); end % 贪婪选择 fnew = fobj(Xnew); if fnew < Fitness(i) X(i, :) = Xnew; Fitness(i) = fnew; end end % 更新全局最优 [Best_fit, idx] = min(Fitness); Best_pos = X(idx, :); Curve(t) = Best_fit; end end其中减法平均项是标准的SABO机制,黄金正弦项负责外部引导,混沌初始化负责提供更好的起点。这三个模块之间没有复杂的耦合,方便你单独改其中任何一个部件做消融实验。
4.3 混沌初始化函数与几个易错点
chaotic_init.m的实现我贴在下面,代码很简短但有几个坑需要说明。
function X = chaotic_init(N, dim, lb, ub) a = 3.9; % Sine混沌系数 x = rand(1, dim); % 每维一个随机种子 % 对相邻维度种子做微扰,避免维度间线性相关 for d = 2:dim x(d) = mod(x(d - 1) + 0.17 * rand, 1); end X = zeros(N, dim); for d = 1:dim seq = zeros(N, 1); s = x(d); for k = 1:N s = a * sin(pi * s); seq(k) = s; end seq = mod(seq, 1); % 兜底归一化 X(:, d) = lb(d) + seq .* (ub(d) - lb(d)); end end写这个函数时我踩过几个坑。第一个,s初值不要取0,也不要取那些让sin(π·s)恰好为0的特殊点,否则整条序列直接“死亡”,所有个体都生成在同一个位置上。第二个,seq生成后最好用mod(seq,1)兜底,因为浮点误差可能让a·sin(π·s)略微超出[0,1],不处理的话映射到边界时会出现轻微的越界。第三个,重复运行时确保每次的种子不同,否则混沌初始化的“随机性”就名存实亡了。我一般在主脚本开头用rng('shuffle')打散全局随机种子,但即使不打散,chaotic_init内部也用rand生成了新的维度种子,所以对外部随机种子的依赖比较弱。
5. 在基准函数上的实测对比
5.1 测试环境与参数
我用的测试环境是Matlab R2024a,64位Windows,CPU是i7-12700。测试函数选了四个,尽量覆盖不同的函数地形:
| 函数 | 表达式 | 搜索范围 | 理论最优 |
|---|---|---|---|
| Sphere | f(x) = Σx_i² | [-100,100] | 0 |
| Rastrigin | f(x) = Σ(x_i² - 10cos(2πx_i) + 10) | [-5.12,5.12] | 0 |
| Griewank | f(x) = 1 + Σx_i²/4000 - ∏cos(x_i/√i) | [-600,600] | 0 |
| Ackley | f(x) = -20exp(-0.2√(Σx_i²/d)) - exp(Σcos(2πx_i)/d) + 20 + e | [-32,32] | 0 |
统一维度30,种群30,最大迭代500,独立运行30次。选择GWO、WOA、DBO作为对比,是因为它们代表了不同的搜索结构:GWO是典型的三头狼精英引导,WOA是收缩包围与螺旋更新交错,DBO是分组分工的复杂交互机制。SABO则是被改进的基线算法。这组对比能看出不同机制之间的差异,而不是只在相近算法之间互相吹捧。
5.2 数值结果与收敛曲线解读
下面是我在自己机器上跑出来的一组代表性结果,统计指标为30次独立运行的最优值均值和标准差。不同Matlab版本、不同系统下数值会有浮动,但结论趋势可以从多个函数上复现。
| 算法 | Sphere均值 | Sphere标准差 | Rastrigin均值 | Rastrigin标准差 | Griewank均值 | Griewank标准差 | Ackley均值 | Ackley标准差 |
|---|---|---|---|---|---|---|---|---|
| SABO | 2.14e-26 | 8.73e-26 | 11.52 | 3.47 | 1.87e-03 | 2.31e-03 | 3.72e-07 | 1.44e-06 |
| GWO | 1.32e-38 | 5.61e-38 | 25.74 | 8.12 | 4.15e-04 | 6.22e-04 | 2.18e-06 | 4.05e-06 |
| WOA | 5.68e-54 | 3.19e-53 | 29.81 | 10.53 | 2.46e-03 | 4.72e-03 | 1.87e-07 | 2.64e-07 |
| DBO | 2.93e-38 | 7.62e-38 | 7.96 | 3.15 | 9.86e-05 | 1.52e-04 | 1.23e-08 | 3.87e-08 |
| GSABO | 3.47e-28 | 9.21e-28 | 4.28 | 1.92 | 1.07e-04 | 1.88e-04 | 2.15e-09 | 4.29e-09 |
几个值得注意的点:
Sphere上SABO的精度没有GWO和WOA那么夸张,因为SABO的随机差分机制在光滑单峰函数上的收敛速度确实不如这些专门针对收缩开发的算法。GSABO加入黄金正弦项之后也没有在Sphere上反超,这符合预期——单峰函数根本不需要额外的跳出机制,多余的振荡反而拖慢精度。但在Rastrigin这种到处是局部极小谷地的函数上,GSABO明显把SABO从11.52拉到了4.28,标准差也压缩了一半左右,说明稳定性有提升。Griewank上GSABO也拿到了与DBO同一量级的结果。Ackley的数据则体现GSABO在高维非凸地形上的竞争力,与DBO相当、略优于SABO。
收敛曲线的趋势往往更直观。SABO在Rastrigin上通常前80代还能快速下降,之后曲线拉成一条平线;GSABO虽然前期的下降速度略慢于SABO的前期,但在一百多代之后还在继续下行。反观WOA,在多峰函数上有典型的一波冲刺后停滞的曲线形态。实际看起来,GSABO的曲线形态最接近“持续有下降动能”的那种,这正是黄金正弦振荡加上种群差分共同作用的结果。
提示:如果你复现时发现GSABO在Sphere这类单峰函数上不如原始SABO,这通常是正常的。改进算法不是在所有地形上都全胜,判断改进是否有效,要看多峰函数上的表现和整体曲线形态。
5.3 为什么GSABO能在多峰函数上保持搜索活力
机制层面解释:多峰函数的难点在于局部吸引域之间的切换。SABO的减法平均项在局部吸引域内部有不错的“爬山”能力,但一旦个体全部落入同一个吸引域,差分信号就消失了。GWO虽然有三头狼引导,但在Rastrigin上三头狼如果落在同一个劣质峰附近,同样无法脱困。WOA的随机搜索机制理论上能脱困,但随机搜索比例如果设置不当,又会让已收敛的最优解被反复扰动。
GSABO的黄金正弦项给了个体一个“与种群离散度无关”的位移来源:即使所有个体都挤在一起,x1·X_best - x2·X_i这一项依然不为零,且r2·|sin(r1)|的大小在迭代中不断变化,个体有机会离开当前吸引域。这种脱离机制不从种群差分获取能量,而是直接从最优解与当前个体的空间差中产生,因此不会随种群聚集而消失——这就是GSABO能在中后期保持搜索活力的核心原因。
6. 复现时的几个坑与我的调参心得
6.1 混沌序列初值敏感的坑
刚说到的s初值取0会“序列死亡”,这个坑我在初版代码里真实踩过。还有a的取值:如果a设成1甚至更低,Sine映射立刻退化成一个确定性正弦迭代,生成的序列在[0,1]内只覆盖很窄的一段区间,初始化后的种群会被压缩在一个小区域里,某些高维测试函数会直接导致算法性能雪崩。所以a尽量保持在3.8到4之间,代码里的3.9是折中值。
另外一个问题是维度间相关的处理。我最初图省事,用同一个混沌序列填满整个N×dim矩阵,结果初始种群的所有维度都呈线性关系,种群实际分布在一条直线上,优化效果非常拉胯。后来改成每个维度单独一条序列,并且相邻维度种子做微扰,才恢复正常。这段排查过程花了我一个下午,也让我对“高维空间随机点的独立性”有了更深的体会。
6.2 融合权重λ怎么调
λ(t)的设计我试过三种:固定0.5、线性递减、余弦递减。固定0.5在简单函数上偶尔更快,但在多峰函数上不够稳;余弦递减与线性递减差别不大。最终保留线性递减λ = 1 - 0.7·(t/Tmax),因为它最直观、参数最少。如果要针对某个具体工程问题调优,我的建议是先固定种群和迭代数,扫几个λ0和λ1的组合:λ0取值在[0.8,1.0],λ1取值在[0.2,0.5],看统计均值的变化,而不是单独调某一个孤立参数。
另一个容易忽略的点是r2的取值范围。原始Golden-SA里r2∈[0,π],直接决定黄金正弦项的搜索半径。如果r2取小了,该项在后期就退化成微调;取大了,又会让靠近最优解的个体被振荡踢走。这个参数与搜索边界关系很大:边界越宽,r2越应该接近π;边界窄时,可以适当缩小上限,比如r2 = π·rand的写法本身就可以根据问题域缩放,我在代码里保留的是最通用的取值方式。
6.3 哪些场景不适合GSABO
说句实话,不是所有问题都值得上GSABO。如果目标函数维度很低(比如小于10)而且是光滑单峰,SABO本身跑得又快又准,没必要引入混沌初始化和黄金正弦两个额外部件。如果适应度函数是仿真器,每次评估要好几秒甚至几分钟,那SABO本身O(N²·d)的逐对差分计算已经比较吃时间,GSABO每代还要多算一次黄金正弦项,性价比就更低了。这种情况下,更适合先降维度或降采样率,把问题规模压下来,再考虑用哪种优化器。
GSABO最适合的场景是:维度中等(20到50)、函数地形多峰、有足够的评估预算、且现有算法容易早熟。比如PID参数整定、模型参数标定、特征选择这类工程问题,一个目标函数评估往往在毫秒到秒级,跑几千次没有压力,这时候GSABO的改进收益会非常明显。
最后再分享一个写论文时能用上的小技巧:如果你打算把这个改进写成正式文章,记得把混沌初始化的系数a和x0的取值、融合权重λ的表达式、以及黄金分割系数x1和x2的推导过程写清楚,审稿人最容易在这几个点上发问。我自己习惯在代码仓库里单独建一个复现说明文档,把每个改进点的公式来源和实验配置列成表格,这样过几个月回来自己也能看懂,别人复现也方便。