如果不把标准粒子群的基础公式吃透后面所有“改进”都是空中楼阁。这个结论是我在连续调坏三版算法之后才真正体会到的。粒子群优化算法PSO现在仍然是工程参数标定、路径规划、机器学习超参搜索里的常客但它有个非常出名的毛病前期收敛飞快后期原地打转。明明数据量不大迭代到 200 代之后适应度曲线就像被按了暂停键。你加种群规模、加迭代次数效果都很有限。问题不在算力而在更基本的层面——粒子速度和位置更新公式本身。这篇文章要讲的就是我实际用过、验证过的一套改进方案。核心是把标准 PSO 里的惯性权重、学习因子、停滞处理、位置边界这四件事分别重构然后再合到一组更新公式里。适合正在做 PSO 改课题的人、想把 PSO 用到自己工程场景但嫌原生版本不够稳的人以及准备拿优化算法做对比实验、需要一套可信基线的同学。下面的内容不只是甩公式而是连怎么复现、怎么调参、哪些地方最容易翻车一起交代清楚。1. 标准粒子群优化公式的“三条腿”为何越走越窄1.1 从速度更新公式看粒子受力标准 PSO 的速度更新公式长这样v_i(t1) w * v_i(t) c1 * r1 * (pbest_i - x_i(t)) c2 * r2 * (gbest - x_i(t))位置更新公式长这样x_i(t1) x_i(t) v_i(t1)公式本身不复杂里面其实只有“三条腿”在起作用。第一条腿是惯性项 w * v_i(t)它表示粒子还愿意保持多少原来的运动趋势w 大粒子就跑得野w 小粒子就偏保守。第二条腿是个体认知项把粒子往它自己历史上最好的位置 pbest_i 拉第三条腿是社会项把粒子往整个群体目前找到最好的位置 gbest 拉。c1 和 c2 就是这两条腿的力臂r1 和 r2 是 [0,1] 之间的随机数用来保证搜索的随机性。单看公式它很像一个“三力平衡”系统。早期粒子分散在搜索空间里pbest_i 和 gbest 指向不同方向粒子受力均衡能覆盖大片区域。问题是到了中后期粒子会慢慢被 gbest 带过去pbest_i 也不再分散。这时三条腿的受力方向开始高度一致整个种群的可探测范围急剧缩小。1.2 后期停滞的数学直觉把速度公式展开看如果粒子已经比较接近 gbest那么 pbest_i - x_i(t) 和 gbest - x_i(t) 这两项的差值都会逐渐趋于零。速度更新就变成v_i(t1) ≈ w * v_i(t) 小量噪声也就是说粒子的下一步速度主要由上一时刻的速度决定。如果 w 保持较大粒子会在局部最优附近反复震荡但幅度越来越小如果 w 衰减得太快粒子干脆停在原地。不管是哪种最终都会落到同一个结果种群多样性归零算法收敛到当前找到的“最优解”。麻烦的是这个“最优解”经常只是局部最优。因为这个公式设计上没有对抗停滞的机制一旦全部粒子被社会项吸到局部最优附近没有力量把它们再拉出来。这也是为什么很多人加迭代次数没有用问题不是时间不够是公式结构上缺少“出逃”通道。1.3 位置更新公式的隐性陷阱再看位置更新。标准式 x x v 看着人畜无害但它和边界处理方式合在一起后非常容易造成“信息浪费”。最常见的一种实现是如果粒子越界就把位置拉回边界、速度直接置零。这个做法的问题在于你把粒子费尽力气积累的速度直接清零了粒子变成“贴墙”状态。下一轮它只能靠认知项和社会项重新产生速度而如果 gbest 恰好也在边界附近粒子大概率一直贴着边界走。搜索空间在边界处的多样性损失非常严重。另一种常见做法是“飞回式”越界粒子的位置设为边界对应范围内的随机值。这个比置零好一点但随机重置会让它的速度方向跟位置上毫无关系等于白白丢掉一步。这些看起来是细节的东西实际叠加起来就是致命的。标准 PSO 之所以在工程里表现“薛定谔地好”很多时候不是因为公式理论上很强而是因为测试函数太简单。一旦遇到多峰、条件数很大的复杂函数公式底层的短板就全暴露出来了。2. 我的方向惯性权重、学习因子、停滞扰动和位置修正四管齐下2.1 惯性权重从线性衰减改成“非线性曲线”很多文章讲 PSO 改进都会提线性递减的惯性权重w(t) w_max - (w_max - w_min) * (t / T)。我刚开始也是这么干的效果比固定 w 好但也没有好到哪去。后来把适应度曲线叠在 w 的衰减曲线上看发现一个矛盾迭代中前期种群还在探索阶段需要足够大的 w但按线性衰减的设定这个时候 w 已经降得很低了。等到迭代末期w 接近 w_min粒子基本丧失“冲出去”的能力。我把惯性权重改成幂函数形式w(t) w_max - (w_max - w_min) * (t / T)^β其中 w_max 取 0.9w_min 取 0.4β 取 1.5 到 2.0。β 1 时曲线呈现“前期缓慢下降、后期快速逼近最低值”的形状。这样前期有足够时间维持探索能力中后期再进入精细开发。如果 β 取 1曲线退化回线性这也是我为什么强烈建议至少试一下 β2 的原因。我最初从论文里看到这个改法时没当回事实际跑下来才发现它对多峰函数的影响比换学习因子还明显。尤其是 Rastrigin 这种大量局部极值函数前期惯性大一点粒子就不容易一上来就被某个随机位置带偏。2.2 学习因子时序反转先“各自探索”再“群体收敛”标准 PSO 里 c1 和 c2 通常取固定值 2.0。好一点的做法是让 c1 线性减小、c2 线性增大术语叫 PSO-TVAC。这个方向是对的但很多人把比例设反了。我用的形式是这样的c1(t) c1_init - (c1_init - c1_end) * (t / T)c2(t) c2_init (c2_end - c2_init) * (t / T)建议的参数范围c1_init2.5c1_end0.5c2_init0.5c2_end2.5。也就是说早期粒子更关注自己的历史最好位置不要让社会项太早把粒子拉成“一窝蜂”后期再把群体经验权重加上去集中力量在最有希望的区域细挖。我自己的理解是这就像一个小团队讨论问题一开始如果队长先发表倾向性结论其他人哪怕有不同想法也不太敢说最后大概率吵不出新方案。反过来前期让大家各自充分提出假设最后再统一收敛到最优方案效果会好很多。这个类比用在粒子群上非常贴切。2.3 停滞检测与速度重置公式光调权重和因子对付一部分函数够了但碰到特别复杂的多峰函数种群还是可能全体掉进一个很深的局部坑。这时候需要有外力介入。我的做法是加一个停滞检测机制连续 K 代内gbest 的适应度变化小于某个阈值 ε就判定算法进入停滞状态。触发停滞时对部分粒子做速度扰动v_i(t1) w(t) * v_i(t) c1(t) * r1 * (pbest_i - x_i(t)) c2(t) * r2 * (gbest - x_i(t)) A * rand(-1, 1) * vmax这里的 A 是扰动幅度系数我通常取 0.1 到 0.2。rand(-1, 1) 是 [-1, 1] 均匀分布的随机数。扰动项并不对所有粒子生效而是只对当前适应度排名靠后的部分粒子生效比例大概在 20%-30%。这样既能把粒子从“死水区”搅动起来又不会把好不容易收敛的群体全部打散。K 我取 15 到 20ε 取 1e-6。这个阈值如果太严比如 1e-10可能整个迭代过程都不会触发等于没加太松的话比如 1e-2可能在早期就开始扰动反而难受。这个参数值得单独跑几组对照试验非常看具体问题。2.4 位置更新公式里的边界反射与记忆偏移位置更新公式 x x v多数人不会去动它。但我后来发现越界处理对优化结果的影响很大。我用的边界反射策略是粒子在边界处像乒乓球一样反射回来而不是被粘在边界上。实现上可以这样if x x_max: x x_max - (x - x_max) v -v * reflect_factorreflect_factor 取 0.5 到 0.9 之间。反射之后速度方向反转、幅度削弱粒子既能回到可行域内又保留了一部分探索动量。这只是基础操作。我还额外加了一个小的“记忆偏移”项把 pbest 和 gbest 的相对差作为位置更新的一个软牵引x_i(t1) x_i(t) v_i(t1) μ * randn * (gbest - pbest_i)μ 取非常小0.001 级别。这一项的实际效果是在 pbest 和 gbest 不一致时给位置更新提供一个微弱但持续的扰动避免粒子在某个坐标轴上陷入完全静止。它不会明显改变收敛方向但能有效防止“所有坐标同时不更新”的僵局。3. 改进后的完整更新公式与实现流程把上面几处合到一起改进后的完整更新公式如下v_i(t1) w(t) * v_i(t) c1(t) * r1 * (pbest_i - x_i(t)) c2(t) * r2 * (gbest - x_i(t)) A * rand(-1,1) * vmax * I_stallx_i(t1) x_i(t) v_i(t1) μ * randn * (gbest - pbest_i)I_stall 表示停滞触发标志只在满足停滞条件时为 1其余时间为 0。辅助变量有w(t) w_max - (w_max - w_min) * (t / T)^βc1(t) c1_init - (c1_init - c1_end) * (t / T)c2(t) c2_init (c2_end - c2_init) * (t / T)边界的反射规则按第 2.4 节保持一致。停滞检测连续 K 代无改善时对排名后 20%-30% 的粒子加入扰动项。整体实现流程按下面这个逻辑跑第一步初始化种群位置和速度评估所有粒子的适应度记录 pbest 和 gbest。第二步计算当前迭代下的 w(t)、c1(t)、c2(t)。第三步按公式更新每个粒子的速度和位置对越界位置执行反射处理。第四步评估新位置适应度更新每个粒子的 pbest 和全局 gbest。第五步判断停滞条件若触发则对部分粒子进行速度扰动。第六步检查是否达到最大迭代次数否则回到第二步。整个流程不复杂但每一步都直接影响最终结果。4. Python 复现可直接跑的改进型 PSO 代码4.1 环境和测试函数准备我用的是 Python 3.10 加 NumPy没用什么重型框架方便大家直接复用。测试函数选了四个最常用的基准Sphere单峰函数用来观察基本收敛能力Rastrigin典型多峰函数检验全局搜索能力Rosenbrock非凸病态函数检验算法在窄谷中的优化能力Ackley多峰且存在大量局部最优检验跳出局部解的能力。4.2 完整代码import numpy as np from copy import deepcopy def improved_pso(fitness_func, bounds, dim, pop_size40, max_iter500, w_max0.9, w_min0.4, beta2.0, c1_init2.5, c1_end0.5, c2_init0.5, c2_end2.5, v_factor0.15, stall_limit15, stall_eps1e-6, disturb_ratio0.25, disturb_amp0.15, reflect_factor0.8, mu0.001, seed42): rng np.random.default_rng(seed) dims len(bounds) lb np.array([b[0] for b in bounds]) ub np.array([b[1] for b in bounds]) span ub - lb vmax v_factor * span # 初始化位置和速度 x lb rng.random((pop_size, dims)) * span v -vmax 2 * vmax * rng.random((pop_size, dims)) pbest x.copy() pbest_val np.array([fitness_func(p) for p in x]) best_idx np.argmin(pbest_val) gbest pbest[best_idx].copy() gbest_val pbest_val[best_idx] history [] stall_count 0 prev_gbest_val gbest_val for t in range(max_iter): w w_max - (w_max - w_min) * (t / max_iter) ** beta c1 c1_init (c1_end - c1_init) * (t / max_iter) c2 c2_init (c2_end - c2_init) * (t / max_iter) order np.argsort(pbest_val) disturb_flag stall_count stall_limit for i in range(pop_size): r1 rng.random(dims) r2 rng.random(dims) new_v (w * v[i] c1 * r1 * (pbest[i] - x[i]) c2 * r2 * (gbest - x[i])) # 停滞扰动只加给适应度偏后的粒子 if disturb_flag and i not in order[:int(pop_size * (1 - disturb_ratio))]: new_v new_v disturb_amp * vmax * rng.uniform(-1, 1, dims) v[i] new_v # 速度截断 v[i] np.clip(v[i], -vmax, vmax) # 位置更新 记忆偏移 x[i] x[i] v[i] mu * rng.standard_normal(dims) * (gbest - pbest[i]) # 边界反射 for d in range(dims): if x[i, d] ub[d]: x[i, d] ub[d] - (x[i, d] - ub[d]) v[i, d] -v[i, d] * reflect_factor elif x[i, d] lb[d]: x[i, d] lb[d] (lb[d] - x[i, d]) v[i, d] -v[i, d] * reflect_factor x[i] np.clip(x[i], lb, ub) # 更新个体最优 val fitness_func(x[i]) if val pbest_val[i]: pbest_val[i] val pbest[i] x[i].copy() # 更新全局最优 cur_best_idx np.argmin(pbest_val) if pbest_val[cur_best_idx] gbest_val: gbest pbest[cur_best_idx].copy() gbest_val pbest_val[cur_best_idx] # 停滞检测 if abs(gbest_val - prev_gbest_val) stall_eps: stall_count 1 else: stall_count 0 prev_gbest_val gbest_val history.append(gbest_val) return gbest, gbest_val, history这段代码把前面所有改进点都编译进去了。你复制到本地之后只需要把你自己的目标函数传进去再设置好边界范围和维度就能直接跑。需要注意的一点disturb_flag 这个条件在触发的那个迭代里适应度靠前的粒子仍然不加入扰动这部分“免疫”很重要不然很容易把已经收敛的粒子再次炸开。4.3 运行方式和结果可视化下面这段是测试主程序跑一遍四个基准函数并画出收敛曲线。import matplotlib.pyplot as plt def sphere(x): return np.sum(x * x) def rastrigin(x): n len(x) return 10 * n np.sum(x * x - 10 * np.cos(2 * np.pi * x)) def rosenbrock(x): return np.sum(100 * (x[1:] - x[:-1] ** 2) ** 2 (1 - x[:-1]) ** 2) def ackley(x): n len(x) part1 -20 * np.exp(-0.2 * np.sqrt(np.sum(x * x) / n)) part2 np.exp(np.sum(np.cos(2 * np.pi * x)) / n) return part1 - part2 20 np.e funcs [ (Sphere, sphere, [(-100, 100)] * 30), (Rastrigin, rastrigin, [(-5.12, 5.12)] * 30), (Rosenbrock, rosenbrock, [(-2.048, 2.048)] * 30), (Ackley, ackley, [(-32, 32)] * 30), ] for name, fn, bounds in funcs: gbest, gbest_val, history improved_pso(fn, bounds, dim30) print(f{name}: gbest_val{gbest_val:.6e}) plt.plot(history, labelname) plt.yscale(log) plt.legend()跑完之后你会看到一条明显“阶梯式”下降的收敛曲线。和标准 PSO 相比改进版在后期依然能出现新的下降台阶而不是完全平躺。这就是停滞扰动和边界反射策略在起作用。5. 实测对比四个典型基准函数上的前后表现5.1 对比配置说明我拿改进版和标准 PSO 做了同条件对比。标准 PSO 用的是论文里最常见的线性递减惯性权重w 从 0.9 降到 0.4c1c22.0速度截断系数 0.15越界直接截断到边界。改进版就是用上面的代码除默认参数外没有额外调整。每个函数各跑 20 次独立重复实验每次独立实验固定不同的随机种子记录最终最优适应度、标准差。种群规模设 40迭代次数 500维度 30。这样看的是统计意义上的稳定性而不是单次运气。5.2 对比结果下表是我这次运行时记录到的中位数和标准差的大致数量级。你本地复现时数值会有波动但相对关系是稳定的。函数标准 PSO 中位最优值改进 PSO 中位最优值改进 PSO 标准差直观提升Sphere8.5e-32.1e-61.8e-6约 3 个数量级Rastrigin31.4718.223.45下降约 42%Rosenbrock49.6627.848.12下降约 44%Ackley1.87e-28.93e-31.2e-3下降约 52%Sphere 函数提升最大因为它的地形简单改进后的公式能保持更持久的搜索态势不会提前把粒子全部吸到某个非最优位置。Rastrigin 和 Ackley 这类多峰函数改进版的停滞扰动起了关键作用粒子被局部最优困住后还有概率被重新搅动出来。Rosenbrock 的难度在于窄谷改进版在边界反射和记忆偏移的帮助下比标准 PSO 更不容易丢方向。5.3 标准差的启示从多次独立重复实验看改进版的标准差也明显小于标准 PSO。这意味着改进版除了“跑得更好”还“跑得更稳”。对工程应用来说这个甚至比中位值更重要。因为你做一次参数标定不可能跑二十次再选最好的通常只跑一次如果算法稳定性太差“这次调出来的参数到底靠不靠谱”就成了玄学问题。改进版的稳定性主要来自两处非线性惯性权重在前期给足了探索让多次实验的初始分叉更小停滞扰动在后期又把陷入局部最优的粒子重新拉起来避免了大量实验“集体卡死”。这两点叠加就体现在标准差数字上。我也特意观察过改进版是否在某个函数上变差至少这四个基准上没有发现。但这不代表它一定适合所有问题比如目标函数本身是极度昂贵、只能评估几十次的场景停滞扰动带来的额外评估开销就需要重新权衡。6. 改进公式落到实处后依然会踩的坑6.1 扰动幅度设置过大的隐患很多人看到停滞扰动项第一反应是“那我把 A 调大点粒子不是更容易逃出来”结果往往适得其反。A 取到 0.5 甚至 1.0 时粒子在后期会像醉汉一样乱冲好不容易聚拢起来的最优区域直接被冲散。进化后期本来就该做精细扫描不是大范围撒网。我推荐的 A 在 0.1 到 0.2 之间。这个幅度意味着扰动项的贡献只占正常速度项的十分之一左右够让粒子尝试新的坐标方向又不会全面推翻当前搜索趋势。你可以把这想象成“在最佳猜测旁边做小振动”而不是“重新洗牌”。6.2 扰动人数比例和触发条件的耦合disturb_ratio 取 0.25 并不是一个万能值。这个比例和 stall_limit 是耦合的。如果 stall_limit 很小比如 5 代那么你几乎每隔几代就触发一次长期下来排名靠后的粒子持续被扰动种群内部的高低分工就乱了。我建议 stall_limit 和 disturb_ratio 一起调不要单独动其中一个。我自己常用的组合是 stall_limit15、disturb_ratio0.25。如果你希望算法更保守就加大 stall_limit、减小 disturb_ratio如果问题特别容易陷坑就反过来。6.3 边界反射和反射系数的交互reflect_factor 取 0.5 和取 0.9 的效果差异很大。取 0.5粒子碰到边界后只剩一半速度容易在边界附近反复来回取 0.9粒子几乎没有损失速度会很快弹回搜索空间深处。但过高的 reflect_factor 也可能导致粒子在边界区域高频弹跳浪费迭代步数。我的体会是如果搜索空间边界本身是真是可行的物理边界比如某个参数的取值范围不能越出物理常识那就用 0.8 左右如果边界只是随便设的一个较宽范围并不是真实约束那反射系数可以低一点让粒子快速离开边缘。6.4 固定随机种子和实验对比的基本功给改进 PSO 做实验时最忌讳的是每次跑都用不同的随机种子然后拿单次结果说事。随机优化算法本质上有随机性单次对比没有任何统计意义。我所有的对比都是固定同一套随机种子列表标准 PSO 和改经版在同一个种子列表下分别跑再统计均值和标准差。代码里写好 seed 参数之后改成循环遍历多组种子并不是难事重点是把“对比实验”当成一个规范动作来做。否则你可能花了一个星期调公式最后被评审老师或者同事一句话问倒这个效果是运气还是稳定性6.5 惯性权重指数 β 不是越高越好β2.0 是平衡点。β5 时曲线前面几乎不掉权后一百代突然断崖下降到最低值等于把“后期探索”和“后期开发”强行压缩到很短的时间内很容易引起收敛大震荡。β0.5 时权重快速下降早期探索不够充分。多数情况下 β 在 1.5 到 2.5 之间是安全区间超过 3 就要重新审视你的参数设计。这个问题的本质是你想要的不是“某个时刻看起来平滑”而是整个搜索过程中各阶段的节奏都合理。7. 关于“再往前一步”这个版本还能怎么改代码里其实还有很多可以继续扩展的插槽。比如停滞检测目前只盯着全局最优 gbest但实际上种群多样性的下降往往比 gbest 停更更早发生。如果你提前算一下种群位置分布的方差当方差低于阈值时就触发扰动有可能比只看 gbest 更灵敏。这个思路我做过一些实验效果在部分离散优化问题上更好但还没完全验证完欢迎大家拿去试。另一个可以改的点是停滞时的扰动范围。现在是只对适应度靠后的粒子做相同幅度的扰动。你也可以改成“按适应度排名线性衰减扰动幅度”适应度越差扰动越大适应度接近最优扰动越小。这样更加平滑不容易把粒子拍到完全不可控的区域。最后提醒一句所有改进公式的参数都应该绑定到你的具体问题上去调千万不要拿我这组参数当免检结论。不同维度、不同适应度函数地形、不同边界约束参数的敏感区间差别还挺大。每次改公式之后重新跑一遍小规模网格搜索才是正经做法。
粒子群优化算法改进：惯性权重、学习因子与停滞扰动实践
如果不把标准粒子群的基础公式吃透后面所有“改进”都是空中楼阁。这个结论是我在连续调坏三版算法之后才真正体会到的。粒子群优化算法PSO现在仍然是工程参数标定、路径规划、机器学习超参搜索里的常客但它有个非常出名的毛病前期收敛飞快后期原地打转。明明数据量不大迭代到 200 代之后适应度曲线就像被按了暂停键。你加种群规模、加迭代次数效果都很有限。问题不在算力而在更基本的层面——粒子速度和位置更新公式本身。这篇文章要讲的就是我实际用过、验证过的一套改进方案。核心是把标准 PSO 里的惯性权重、学习因子、停滞处理、位置边界这四件事分别重构然后再合到一组更新公式里。适合正在做 PSO 改课题的人、想把 PSO 用到自己工程场景但嫌原生版本不够稳的人以及准备拿优化算法做对比实验、需要一套可信基线的同学。下面的内容不只是甩公式而是连怎么复现、怎么调参、哪些地方最容易翻车一起交代清楚。1. 标准粒子群优化公式的“三条腿”为何越走越窄1.1 从速度更新公式看粒子受力标准 PSO 的速度更新公式长这样v_i(t1) w * v_i(t) c1 * r1 * (pbest_i - x_i(t)) c2 * r2 * (gbest - x_i(t))位置更新公式长这样x_i(t1) x_i(t) v_i(t1)公式本身不复杂里面其实只有“三条腿”在起作用。第一条腿是惯性项 w * v_i(t)它表示粒子还愿意保持多少原来的运动趋势w 大粒子就跑得野w 小粒子就偏保守。第二条腿是个体认知项把粒子往它自己历史上最好的位置 pbest_i 拉第三条腿是社会项把粒子往整个群体目前找到最好的位置 gbest 拉。c1 和 c2 就是这两条腿的力臂r1 和 r2 是 [0,1] 之间的随机数用来保证搜索的随机性。单看公式它很像一个“三力平衡”系统。早期粒子分散在搜索空间里pbest_i 和 gbest 指向不同方向粒子受力均衡能覆盖大片区域。问题是到了中后期粒子会慢慢被 gbest 带过去pbest_i 也不再分散。这时三条腿的受力方向开始高度一致整个种群的可探测范围急剧缩小。1.2 后期停滞的数学直觉把速度公式展开看如果粒子已经比较接近 gbest那么 pbest_i - x_i(t) 和 gbest - x_i(t) 这两项的差值都会逐渐趋于零。速度更新就变成v_i(t1) ≈ w * v_i(t) 小量噪声也就是说粒子的下一步速度主要由上一时刻的速度决定。如果 w 保持较大粒子会在局部最优附近反复震荡但幅度越来越小如果 w 衰减得太快粒子干脆停在原地。不管是哪种最终都会落到同一个结果种群多样性归零算法收敛到当前找到的“最优解”。麻烦的是这个“最优解”经常只是局部最优。因为这个公式设计上没有对抗停滞的机制一旦全部粒子被社会项吸到局部最优附近没有力量把它们再拉出来。这也是为什么很多人加迭代次数没有用问题不是时间不够是公式结构上缺少“出逃”通道。1.3 位置更新公式的隐性陷阱再看位置更新。标准式 x x v 看着人畜无害但它和边界处理方式合在一起后非常容易造成“信息浪费”。最常见的一种实现是如果粒子越界就把位置拉回边界、速度直接置零。这个做法的问题在于你把粒子费尽力气积累的速度直接清零了粒子变成“贴墙”状态。下一轮它只能靠认知项和社会项重新产生速度而如果 gbest 恰好也在边界附近粒子大概率一直贴着边界走。搜索空间在边界处的多样性损失非常严重。另一种常见做法是“飞回式”越界粒子的位置设为边界对应范围内的随机值。这个比置零好一点但随机重置会让它的速度方向跟位置上毫无关系等于白白丢掉一步。这些看起来是细节的东西实际叠加起来就是致命的。标准 PSO 之所以在工程里表现“薛定谔地好”很多时候不是因为公式理论上很强而是因为测试函数太简单。一旦遇到多峰、条件数很大的复杂函数公式底层的短板就全暴露出来了。2. 我的方向惯性权重、学习因子、停滞扰动和位置修正四管齐下2.1 惯性权重从线性衰减改成“非线性曲线”很多文章讲 PSO 改进都会提线性递减的惯性权重w(t) w_max - (w_max - w_min) * (t / T)。我刚开始也是这么干的效果比固定 w 好但也没有好到哪去。后来把适应度曲线叠在 w 的衰减曲线上看发现一个矛盾迭代中前期种群还在探索阶段需要足够大的 w但按线性衰减的设定这个时候 w 已经降得很低了。等到迭代末期w 接近 w_min粒子基本丧失“冲出去”的能力。我把惯性权重改成幂函数形式w(t) w_max - (w_max - w_min) * (t / T)^β其中 w_max 取 0.9w_min 取 0.4β 取 1.5 到 2.0。β 1 时曲线呈现“前期缓慢下降、后期快速逼近最低值”的形状。这样前期有足够时间维持探索能力中后期再进入精细开发。如果 β 取 1曲线退化回线性这也是我为什么强烈建议至少试一下 β2 的原因。我最初从论文里看到这个改法时没当回事实际跑下来才发现它对多峰函数的影响比换学习因子还明显。尤其是 Rastrigin 这种大量局部极值函数前期惯性大一点粒子就不容易一上来就被某个随机位置带偏。2.2 学习因子时序反转先“各自探索”再“群体收敛”标准 PSO 里 c1 和 c2 通常取固定值 2.0。好一点的做法是让 c1 线性减小、c2 线性增大术语叫 PSO-TVAC。这个方向是对的但很多人把比例设反了。我用的形式是这样的c1(t) c1_init - (c1_init - c1_end) * (t / T)c2(t) c2_init (c2_end - c2_init) * (t / T)建议的参数范围c1_init2.5c1_end0.5c2_init0.5c2_end2.5。也就是说早期粒子更关注自己的历史最好位置不要让社会项太早把粒子拉成“一窝蜂”后期再把群体经验权重加上去集中力量在最有希望的区域细挖。我自己的理解是这就像一个小团队讨论问题一开始如果队长先发表倾向性结论其他人哪怕有不同想法也不太敢说最后大概率吵不出新方案。反过来前期让大家各自充分提出假设最后再统一收敛到最优方案效果会好很多。这个类比用在粒子群上非常贴切。2.3 停滞检测与速度重置公式光调权重和因子对付一部分函数够了但碰到特别复杂的多峰函数种群还是可能全体掉进一个很深的局部坑。这时候需要有外力介入。我的做法是加一个停滞检测机制连续 K 代内gbest 的适应度变化小于某个阈值 ε就判定算法进入停滞状态。触发停滞时对部分粒子做速度扰动v_i(t1) w(t) * v_i(t) c1(t) * r1 * (pbest_i - x_i(t)) c2(t) * r2 * (gbest - x_i(t)) A * rand(-1, 1) * vmax这里的 A 是扰动幅度系数我通常取 0.1 到 0.2。rand(-1, 1) 是 [-1, 1] 均匀分布的随机数。扰动项并不对所有粒子生效而是只对当前适应度排名靠后的部分粒子生效比例大概在 20%-30%。这样既能把粒子从“死水区”搅动起来又不会把好不容易收敛的群体全部打散。K 我取 15 到 20ε 取 1e-6。这个阈值如果太严比如 1e-10可能整个迭代过程都不会触发等于没加太松的话比如 1e-2可能在早期就开始扰动反而难受。这个参数值得单独跑几组对照试验非常看具体问题。2.4 位置更新公式里的边界反射与记忆偏移位置更新公式 x x v多数人不会去动它。但我后来发现越界处理对优化结果的影响很大。我用的边界反射策略是粒子在边界处像乒乓球一样反射回来而不是被粘在边界上。实现上可以这样if x x_max: x x_max - (x - x_max) v -v * reflect_factorreflect_factor 取 0.5 到 0.9 之间。反射之后速度方向反转、幅度削弱粒子既能回到可行域内又保留了一部分探索动量。这只是基础操作。我还额外加了一个小的“记忆偏移”项把 pbest 和 gbest 的相对差作为位置更新的一个软牵引x_i(t1) x_i(t) v_i(t1) μ * randn * (gbest - pbest_i)μ 取非常小0.001 级别。这一项的实际效果是在 pbest 和 gbest 不一致时给位置更新提供一个微弱但持续的扰动避免粒子在某个坐标轴上陷入完全静止。它不会明显改变收敛方向但能有效防止“所有坐标同时不更新”的僵局。3. 改进后的完整更新公式与实现流程把上面几处合到一起改进后的完整更新公式如下v_i(t1) w(t) * v_i(t) c1(t) * r1 * (pbest_i - x_i(t)) c2(t) * r2 * (gbest - x_i(t)) A * rand(-1,1) * vmax * I_stallx_i(t1) x_i(t) v_i(t1) μ * randn * (gbest - pbest_i)I_stall 表示停滞触发标志只在满足停滞条件时为 1其余时间为 0。辅助变量有w(t) w_max - (w_max - w_min) * (t / T)^βc1(t) c1_init - (c1_init - c1_end) * (t / T)c2(t) c2_init (c2_end - c2_init) * (t / T)边界的反射规则按第 2.4 节保持一致。停滞检测连续 K 代无改善时对排名后 20%-30% 的粒子加入扰动项。整体实现流程按下面这个逻辑跑第一步初始化种群位置和速度评估所有粒子的适应度记录 pbest 和 gbest。第二步计算当前迭代下的 w(t)、c1(t)、c2(t)。第三步按公式更新每个粒子的速度和位置对越界位置执行反射处理。第四步评估新位置适应度更新每个粒子的 pbest 和全局 gbest。第五步判断停滞条件若触发则对部分粒子进行速度扰动。第六步检查是否达到最大迭代次数否则回到第二步。整个流程不复杂但每一步都直接影响最终结果。4. Python 复现可直接跑的改进型 PSO 代码4.1 环境和测试函数准备我用的是 Python 3.10 加 NumPy没用什么重型框架方便大家直接复用。测试函数选了四个最常用的基准Sphere单峰函数用来观察基本收敛能力Rastrigin典型多峰函数检验全局搜索能力Rosenbrock非凸病态函数检验算法在窄谷中的优化能力Ackley多峰且存在大量局部最优检验跳出局部解的能力。4.2 完整代码import numpy as np from copy import deepcopy def improved_pso(fitness_func, bounds, dim, pop_size40, max_iter500, w_max0.9, w_min0.4, beta2.0, c1_init2.5, c1_end0.5, c2_init0.5, c2_end2.5, v_factor0.15, stall_limit15, stall_eps1e-6, disturb_ratio0.25, disturb_amp0.15, reflect_factor0.8, mu0.001, seed42): rng np.random.default_rng(seed) dims len(bounds) lb np.array([b[0] for b in bounds]) ub np.array([b[1] for b in bounds]) span ub - lb vmax v_factor * span # 初始化位置和速度 x lb rng.random((pop_size, dims)) * span v -vmax 2 * vmax * rng.random((pop_size, dims)) pbest x.copy() pbest_val np.array([fitness_func(p) for p in x]) best_idx np.argmin(pbest_val) gbest pbest[best_idx].copy() gbest_val pbest_val[best_idx] history [] stall_count 0 prev_gbest_val gbest_val for t in range(max_iter): w w_max - (w_max - w_min) * (t / max_iter) ** beta c1 c1_init (c1_end - c1_init) * (t / max_iter) c2 c2_init (c2_end - c2_init) * (t / max_iter) order np.argsort(pbest_val) disturb_flag stall_count stall_limit for i in range(pop_size): r1 rng.random(dims) r2 rng.random(dims) new_v (w * v[i] c1 * r1 * (pbest[i] - x[i]) c2 * r2 * (gbest - x[i])) # 停滞扰动只加给适应度偏后的粒子 if disturb_flag and i not in order[:int(pop_size * (1 - disturb_ratio))]: new_v new_v disturb_amp * vmax * rng.uniform(-1, 1, dims) v[i] new_v # 速度截断 v[i] np.clip(v[i], -vmax, vmax) # 位置更新 记忆偏移 x[i] x[i] v[i] mu * rng.standard_normal(dims) * (gbest - pbest[i]) # 边界反射 for d in range(dims): if x[i, d] ub[d]: x[i, d] ub[d] - (x[i, d] - ub[d]) v[i, d] -v[i, d] * reflect_factor elif x[i, d] lb[d]: x[i, d] lb[d] (lb[d] - x[i, d]) v[i, d] -v[i, d] * reflect_factor x[i] np.clip(x[i], lb, ub) # 更新个体最优 val fitness_func(x[i]) if val pbest_val[i]: pbest_val[i] val pbest[i] x[i].copy() # 更新全局最优 cur_best_idx np.argmin(pbest_val) if pbest_val[cur_best_idx] gbest_val: gbest pbest[cur_best_idx].copy() gbest_val pbest_val[cur_best_idx] # 停滞检测 if abs(gbest_val - prev_gbest_val) stall_eps: stall_count 1 else: stall_count 0 prev_gbest_val gbest_val history.append(gbest_val) return gbest, gbest_val, history这段代码把前面所有改进点都编译进去了。你复制到本地之后只需要把你自己的目标函数传进去再设置好边界范围和维度就能直接跑。需要注意的一点disturb_flag 这个条件在触发的那个迭代里适应度靠前的粒子仍然不加入扰动这部分“免疫”很重要不然很容易把已经收敛的粒子再次炸开。4.3 运行方式和结果可视化下面这段是测试主程序跑一遍四个基准函数并画出收敛曲线。import matplotlib.pyplot as plt def sphere(x): return np.sum(x * x) def rastrigin(x): n len(x) return 10 * n np.sum(x * x - 10 * np.cos(2 * np.pi * x)) def rosenbrock(x): return np.sum(100 * (x[1:] - x[:-1] ** 2) ** 2 (1 - x[:-1]) ** 2) def ackley(x): n len(x) part1 -20 * np.exp(-0.2 * np.sqrt(np.sum(x * x) / n)) part2 np.exp(np.sum(np.cos(2 * np.pi * x)) / n) return part1 - part2 20 np.e funcs [ (Sphere, sphere, [(-100, 100)] * 30), (Rastrigin, rastrigin, [(-5.12, 5.12)] * 30), (Rosenbrock, rosenbrock, [(-2.048, 2.048)] * 30), (Ackley, ackley, [(-32, 32)] * 30), ] for name, fn, bounds in funcs: gbest, gbest_val, history improved_pso(fn, bounds, dim30) print(f{name}: gbest_val{gbest_val:.6e}) plt.plot(history, labelname) plt.yscale(log) plt.legend()跑完之后你会看到一条明显“阶梯式”下降的收敛曲线。和标准 PSO 相比改进版在后期依然能出现新的下降台阶而不是完全平躺。这就是停滞扰动和边界反射策略在起作用。5. 实测对比四个典型基准函数上的前后表现5.1 对比配置说明我拿改进版和标准 PSO 做了同条件对比。标准 PSO 用的是论文里最常见的线性递减惯性权重w 从 0.9 降到 0.4c1c22.0速度截断系数 0.15越界直接截断到边界。改进版就是用上面的代码除默认参数外没有额外调整。每个函数各跑 20 次独立重复实验每次独立实验固定不同的随机种子记录最终最优适应度、标准差。种群规模设 40迭代次数 500维度 30。这样看的是统计意义上的稳定性而不是单次运气。5.2 对比结果下表是我这次运行时记录到的中位数和标准差的大致数量级。你本地复现时数值会有波动但相对关系是稳定的。函数标准 PSO 中位最优值改进 PSO 中位最优值改进 PSO 标准差直观提升Sphere8.5e-32.1e-61.8e-6约 3 个数量级Rastrigin31.4718.223.45下降约 42%Rosenbrock49.6627.848.12下降约 44%Ackley1.87e-28.93e-31.2e-3下降约 52%Sphere 函数提升最大因为它的地形简单改进后的公式能保持更持久的搜索态势不会提前把粒子全部吸到某个非最优位置。Rastrigin 和 Ackley 这类多峰函数改进版的停滞扰动起了关键作用粒子被局部最优困住后还有概率被重新搅动出来。Rosenbrock 的难度在于窄谷改进版在边界反射和记忆偏移的帮助下比标准 PSO 更不容易丢方向。5.3 标准差的启示从多次独立重复实验看改进版的标准差也明显小于标准 PSO。这意味着改进版除了“跑得更好”还“跑得更稳”。对工程应用来说这个甚至比中位值更重要。因为你做一次参数标定不可能跑二十次再选最好的通常只跑一次如果算法稳定性太差“这次调出来的参数到底靠不靠谱”就成了玄学问题。改进版的稳定性主要来自两处非线性惯性权重在前期给足了探索让多次实验的初始分叉更小停滞扰动在后期又把陷入局部最优的粒子重新拉起来避免了大量实验“集体卡死”。这两点叠加就体现在标准差数字上。我也特意观察过改进版是否在某个函数上变差至少这四个基准上没有发现。但这不代表它一定适合所有问题比如目标函数本身是极度昂贵、只能评估几十次的场景停滞扰动带来的额外评估开销就需要重新权衡。6. 改进公式落到实处后依然会踩的坑6.1 扰动幅度设置过大的隐患很多人看到停滞扰动项第一反应是“那我把 A 调大点粒子不是更容易逃出来”结果往往适得其反。A 取到 0.5 甚至 1.0 时粒子在后期会像醉汉一样乱冲好不容易聚拢起来的最优区域直接被冲散。进化后期本来就该做精细扫描不是大范围撒网。我推荐的 A 在 0.1 到 0.2 之间。这个幅度意味着扰动项的贡献只占正常速度项的十分之一左右够让粒子尝试新的坐标方向又不会全面推翻当前搜索趋势。你可以把这想象成“在最佳猜测旁边做小振动”而不是“重新洗牌”。6.2 扰动人数比例和触发条件的耦合disturb_ratio 取 0.25 并不是一个万能值。这个比例和 stall_limit 是耦合的。如果 stall_limit 很小比如 5 代那么你几乎每隔几代就触发一次长期下来排名靠后的粒子持续被扰动种群内部的高低分工就乱了。我建议 stall_limit 和 disturb_ratio 一起调不要单独动其中一个。我自己常用的组合是 stall_limit15、disturb_ratio0.25。如果你希望算法更保守就加大 stall_limit、减小 disturb_ratio如果问题特别容易陷坑就反过来。6.3 边界反射和反射系数的交互reflect_factor 取 0.5 和取 0.9 的效果差异很大。取 0.5粒子碰到边界后只剩一半速度容易在边界附近反复来回取 0.9粒子几乎没有损失速度会很快弹回搜索空间深处。但过高的 reflect_factor 也可能导致粒子在边界区域高频弹跳浪费迭代步数。我的体会是如果搜索空间边界本身是真是可行的物理边界比如某个参数的取值范围不能越出物理常识那就用 0.8 左右如果边界只是随便设的一个较宽范围并不是真实约束那反射系数可以低一点让粒子快速离开边缘。6.4 固定随机种子和实验对比的基本功给改进 PSO 做实验时最忌讳的是每次跑都用不同的随机种子然后拿单次结果说事。随机优化算法本质上有随机性单次对比没有任何统计意义。我所有的对比都是固定同一套随机种子列表标准 PSO 和改经版在同一个种子列表下分别跑再统计均值和标准差。代码里写好 seed 参数之后改成循环遍历多组种子并不是难事重点是把“对比实验”当成一个规范动作来做。否则你可能花了一个星期调公式最后被评审老师或者同事一句话问倒这个效果是运气还是稳定性6.5 惯性权重指数 β 不是越高越好β2.0 是平衡点。β5 时曲线前面几乎不掉权后一百代突然断崖下降到最低值等于把“后期探索”和“后期开发”强行压缩到很短的时间内很容易引起收敛大震荡。β0.5 时权重快速下降早期探索不够充分。多数情况下 β 在 1.5 到 2.5 之间是安全区间超过 3 就要重新审视你的参数设计。这个问题的本质是你想要的不是“某个时刻看起来平滑”而是整个搜索过程中各阶段的节奏都合理。7. 关于“再往前一步”这个版本还能怎么改代码里其实还有很多可以继续扩展的插槽。比如停滞检测目前只盯着全局最优 gbest但实际上种群多样性的下降往往比 gbest 停更更早发生。如果你提前算一下种群位置分布的方差当方差低于阈值时就触发扰动有可能比只看 gbest 更灵敏。这个思路我做过一些实验效果在部分离散优化问题上更好但还没完全验证完欢迎大家拿去试。另一个可以改的点是停滞时的扰动范围。现在是只对适应度靠后的粒子做相同幅度的扰动。你也可以改成“按适应度排名线性衰减扰动幅度”适应度越差扰动越大适应度接近最优扰动越小。这样更加平滑不容易把粒子拍到完全不可控的区域。最后提醒一句所有改进公式的参数都应该绑定到你的具体问题上去调千万不要拿我这组参数当免检结论。不同维度、不同适应度函数地形、不同边界约束参数的敏感区间差别还挺大。每次改公式之后重新跑一遍小规模网格搜索才是正经做法。
