一群"鸟"帮你解数学题：聊聊粒子群算法

[复制链接]

字体大小: 正常放大

2744557306

1192 主题	4 听众	2946 积分

该用户从未签到

电梯直达

1^#

发表于 2026-8-23 17:38 |只看该作者 |倒序浏览

|招呼Ta 关注Ta

算法笔记 · 粒子群优化 · 附完整可运行代码

先讲一个真实的故事。

1995 年，两位研究者在模拟鸟群觅食：一群鸟在一片区域里随机飞，谁也不知道食物在哪。但每只鸟都记得自己飞过的最好位置，也能看到整个鸟群目前找到的最好位置。于是每只鸟的飞行策略变得很简单——一部分凭自己的经验继续飞，一部分朝群里公认的好方向靠拢。

结果让人惊讶：不需要任何"领导"，不需要知道食物在哪，鸟群竟能自发地聚拢到食物周围。

这两位研究者一位是电子工程师 James Kennedy，一位是社会心理学家 Russell Eberhart。他们把这个模拟简化成几个公式，就成了今天要讲的主角——粒子群优化算法（Particle Swarm Optimization，PSO）。顺带一提，这个算法诞生 25 年后，他们拿到了 IEEE 的诺贝尔级别的荣誉——它至今仍是工程优化里最常用的算法之一。

如果你正在准备数学建模竞赛，这篇文章能帮你在"优化类"赛题里多一件趁手的兵器。

一、把优化问题想成"找食物"

先忘掉算法，想一个具体场景：

你要找一个函数的最小值，比如 $$f(x, y) = x² + y²$$，最小值显然在原点。但现实问题里，函数可能长得像一片丘陵——到处是小山包（局部最优），真正的谷底（全局最优）藏在某个角落。梯度下降这类方法从一点出发，很容易卡在最近的小山包里下不来。

PSO 的思路完全不同：与其派一个人去找，不如派一群鸟去找。

[td]

概念	在算法里是什么
一只鸟 / 粒子（Particle）	一个候选解（一组变量值）
位置（Position）	解本身，比如 (x, y) = (3.2, -1.5)
速度（Velocity）	下一步往哪个方向飞、飞多远
食物的位置	函数的最小值点
离食物多近	适应度：把位置代入函数算出的值
自己的经验	个体极值 pbest：我飞过的最好位置
群体的情报	全局极值 gbest：全群找到的最好位置

二、灵魂公式：只有两行

PSO 全部精髓浓缩在两条更新规则里。每一代（每一轮迭代），每只鸟这样飞：

$$v_{k+1} = w \cdot v_k + c_1 r_1 \cdot (pbest - x_k) + c_2 r_2 \cdot (gbest - x_k)$$

$$x_{k+1} = x_k + v_{k+1}$$

第一条是速度更新，即"下一刻往哪飞"。拆开看，一共三股力量在拉扯这只鸟：

[td]

部分	含义	生活翻译
w·vk	惯性：保留原来的飞行方向	"我正往这边飞呢，别让我急转弯"
c₁r₁·(pbest − x)	认知项：朝自己最好位置飞	"我以前在那儿吃过食，回去看看"
c₂r₂·(gbest − x)	社会项：朝全群最好位置飞	"老王说他找到吃的了，跟过去！"

第二条是位置更新，简单到不像话：新位置 = 旧位置 + 速度。就是"按刚才决定的方向飞过去"。

这里的 r₁、r₂ 是 0 到 1 之间的随机数——这点很关键。如果每次都精确地朝 pbest 和 gbest 飞，所有鸟很快会挤成一条线，搜索就死了。加上随机扰动，鸟群才能保持"散着搜索、又大致向好看的方向去"的状态。

一句话记住 PSO： 每只鸟被三股力拉着——惯性（接着飞）、自己的记忆（回自己找到过的好地方）、群体的情报（去大家公认的好地方）。飞到新位置后，更新记忆，再飞。如此反复，鸟群整体就向最优解聚拢。

三、三个参数怎么调[td]

参数	名字	常用值	调大调小会怎样
w	惯性权重	0.6 ~ 0.8	大→搜得广但收敛慢；小→收敛快但易陷局部最优
c₁	认知因子	1.5 ~ 2	大→各飞各的，多样性好
c₂	社会因子	1.5 ~ 2	大→一窝蜂跟风，收敛快但容易早熟

工程上还有个进阶技巧：惯性权重线性递减——从 0.9 慢慢降到 0.4。前期惯性大，鸟群撒开网全局搜索；后期惯性小，大家稳稳地精修局部。这被称为"先广搜、后精修"，是 PSO 最常用的改良。

四、上代码：70 行实现完整 PSO

下面这段代码零依赖，保存成 pso_demo.py 直接运行。我们拿两个经典函数开刀：Sphere（单峰，考收敛速度）和 Rastrigin（满地局部最优，考跳出陷阱的能力）。

# -*- coding: utf-8 -*-- X8 r# H* s% c; `% S2 v\" }
"""粒子群优化算法（PSO）完整演示，零依赖可直接运行"""1 A( D3 L2 f; b
import random- _2 j5 E* D7 V* |8 ~6 Y7 s
import math
( b6 c5 J/ i) M( g1 _\" g3 ^
; [) K3 J! Q\" y* Z5 z
def sphere(x):
: X& n7 m\" f2 U/ D\" }# O) g% i
"""Sphere 函数：f(0,0,...,0)=0，单峰"""
, w8 h' ]9 f# C$ H
return sum(xi ** 2 for xi in x)4 P+ w- l\" I* B1 b7 [0 }
! T! V+ n o0 w2 y1 L5 w1 @1 M
def rastrigin(x):\" s) R# F# N8 c+ a' p. d% |% k
"""Rastrigin 函数：大量局部最优点，全局最优仍在原点"""
% l/ m7 \% c9 Y, J9 v4 T
return sum(xi**2 - 10*math.cos(2*math.pi*xi) + 10 for xi in x)7 m+ J5 Z% A7 v- W/ x
; `9 o3 }( w- B+ @8 ~ F/ @& }
class Particle:
' o) I. S- s\" _$ s* d
"""一只『鸟』：有位置、速度，记得自己飞过的最好地方"""4 D3 B7 ^- f2 S# A4 h+ R; R3 ^
def __init__(self, dim, bounds):/ |* n4 M; n7 `0 w
self.pos = [random.uniform(bounds[i][0], bounds[i][1]) for i in range(dim)]
+ s8 c1 R' X, _/ z7 l5 I
self.vel = [random.uniform(-1, 1) for _ in range(dim)]
/ K# K5 w# z\" s; @; J
self.best_pos = list(self.pos) # 个体极值 pbest4 u$ i3 {) A. q; K {8 a) ^/ |
self.best_val = float("inf")
# R( a- I! k- V: H7 y z8 N
/ j; a' Y) p\" s7 i/ B9 F
def update_velocity(self, gbest, w, c1, c2):9 G4 e! B% _4 Q a
# 速度更新：v = w*v + c1*r1*(pbest-x) + c2*r2*(gbest-x)
- w$ G# W6 F. U( P- c; }! U: V
for i in range(len(self.vel)):
' d% ?6 q; m) G- z, {& [
r1, r2 = random.random(), random.random()+ i+ s5 |/ `4 ?' A O) n
cognitive = c1 * r1 * (self.best_pos[i] - self.pos[i]) # 自己的经验8 n; J% K! Z5 W; T9 ]& S
social = c2 * r2 * (gbest[i] - self.pos[i]) # 群体的情报
' H. m5 D9 x8 E, _8 y\" I( D
self.vel[i] = w * self.vel[i] + cognitive + social
3 C5 ]5 Q5 @4 \
, s p ^$ r$ u& M% _9 i
def update_position(self, bounds):
$ A\" J9 }- b! N$ t
# 位置更新：x = x + v，并夹回边界
1 j/ X( p4 {$ `! R
for i in range(len(self.pos)):/ S% t# F/ D\" n5 o! ^) f/ R
self.pos[i] += self.vel[i]& p3 E' j: l3 i. q* x
lo, hi = bounds[i]
x2 r: n0 Z8 x/ m) p
if self.pos[i] < lo: self.pos[i] = lo, @5 G0 \( A) g: z% M) U' L+ D3 O3 I
elif self.pos[i] > hi: self.pos[i] = hi4 z( |* Q\" x\" B+ w2 ?2 d
$ m. y7 V+ K2 k# D6 t) ~4 r! f ?
def pso(func, dim=2, bounds=None, n_particles=30, max_iter=100,
3 C% F\" p\" Q. y4 A! |\" N) W
w=0.7, c1=1.5, c2=1.5, verbose_every=10):
! E0 m8 [\" Q& M J
"""标准 PSO 主流程"""
8 Z, ~& n) O# i8 w
if bounds is None:& n: P3 o2 a: r5 m( M: J' P
bounds = [(-5.12, 5.12)] * dim8 i5 G4 w( q, z B5 S4 u: ?2 U! f- ]( Q
particles = [Particle(dim, bounds) for _ in range(n_particles)], h: U. w: A& g7 ?
gbest, gbest_val = list(particles[0].pos), float("inf")6 \; V( S) ^- I6 S5 N
history = [], W3 G. o9 s4 c3 [& m1 @4 h9 \
* L\" ]* @: I0 z+ C9 k; W- C
for it in range(max_iter):
9 j8 W\" ~ g! ^2 d
# 第 1 步：评估每只鸟，更新 pbest 和 gbest: Z1 M! G( f4 d8 z; q
for p in particles:$ Z\" D& j; q9 C1 z4 Q
val = func(p.pos)
% z. ?8 |( v2 T
if val < p.best_val:
2 g2 G. ], ]- w( h* ~
p.best_val, p.best_pos = val, list(p.pos)7 d, ^$ m2 u2 }\" O. b; ~
if val < gbest_val:6 B+ L% d) O; v, Q4 ^5 G- _
gbest_val, gbest = val, list(p.pos)
- N+ U3 [- w* @
history.append(gbest_val)\" q3 k: N, a2 B4 a, U8 m# w4 `& B
$ M4 K% `! {# ]* G I
# 第 2 步：更新速度和位置\" O+ C4 {3 f7 f0 G6 e
for p in particles:9 p x; ~. ^$ ]3 p- D
p.update_velocity(gbest, w, c1, c2)
; z0 h$ F$ d& }; a# W! w. F6 s; A1 J
p.update_position(bounds)0 _5 ?+ I- n0 Y1 F( K# `
$ K: H0 N5 X8 f! @
if verbose_every and (it + 1) % verbose_every == 0:
Z\" o1 W# i+ K
print(f" 第 {it+1:3d} 代 | 全局最优值 = {gbest_val:.6f}")
, v# |# @2 A) [7 B. [4 n) u
return gbest, gbest_val, history5 f4 B! r# w& p( H$ j0 E
: C7 h( T# b3 w. B% l
if __name__ == "__main__":
0 U\" s9 _1 O# _' l
# 实验 1：单峰函数，看收敛速度$ O' V. X$ ^, X8 a% l) j B6 j- F
print("[实验 1] Sphere 函数（理论最优值 = 0）")+ J. _7 g) J4 V) m
best_pos, best_val, _ = pso(sphere, dim=2, bounds=[(-10, 10)]*2,
0 |. m' d( Z7 P5 |2 D/ s/ [
n_particles=30, max_iter=100)# h9 B. c( D/ X* U' H\" `7 a7 j3 t
print(f" >>> 最优位置: ({best_pos[0]:.4f}, {best_pos[1]:.4f})")
1 ^( s, k7 [% { P
print(f" >>> 最优值: {best_val:.8f}")& G3 d# {1 H, O' p2 F
$ K; q/ E# z3 e* b1 P2 K0 @
# 实验 2：多峰函数，看跳出局部最优
! t# m% z9 ]1 y3 D
print("[实验 2] Rastrigin 函数（局部最优遍地都是）")
1 x7 e( O' R/ A+ P! B# U
best_pos, best_val, _ = pso(rastrigin, dim=2,
+ C M4 e8 L7 M1 s4 y1 Q$ V
n_particles=50, max_iter=200)' {' T4 q\" W$ i0 U8 M$ P$ W
print(f" >>> 最优值: {best_val:.6f}")
+ ^; |+ s6 S4 E# g, n# u* p
) E) k* B+ R3 G
# 实验 3：惯性权重对比（10 维）
\" B5 K1 E& [/ ]8 C( s. t
print("[实验 3] 惯性权重 w 的影响（10 维 Rastrigin）")
: D& ]( g) W% x' [
for w in (0.2, 0.7, 1.2):
\" _5 M& n# s1 j( j( M1 d& D8 ~
_, val, hist = pso(rastrigin, dim=10, bounds=[(-5.12, 5.12)]*10,4 h# f( P; E! b% S, m3 @
n_particles=50, max_iter=200, w=w, verbose_every=0)9 n: J Y, f8 w5 y$ S$ d& h
print(f" w = {w:<4} | 最终最优值 = {val:10.5f}")

复制代码

三个实验分别说明：
实验 1：单峰函数太简单，50 只鸟迭代 50 代就精确命中原点，误差小于百万分之一。
实验 2：Rastrigin 函数像个布满小坑的月球表面，梯度下降几乎必卡坑里，PSO 靠"一群鸟分散搜索"硬是跳了出来。
实验 3最有意思——升到 10 维后参数的差别立刻显现：w=1.2 时惯性太大，鸟群一直"飘着"落不下来（101.6，基本没搜到）；w=0.2 时收敛太快，早熟地卡在局部最优（24.9）；w=0.7 恰到好处（7.0）。这就是参数敏感性的现场教学。
五、什么时候用 PSO，什么时候别用适合的场景[td]
场景
举例
函数不可导 / 不光滑目标函数是查表、仿真得出的，没有解析式
多峰、怕陷入局部最优神经网络调参、TSP 路径优化、portfolio 投资组合
混合整数 / 离散问题稍加改造即可处理 0-1 变量（粒子位置取整）
需要"能跑就行"的快速方案数模竞赛限时三天，PSO 半小时就能调通
要小心的地方
1. 不是全局最优保证。 PSO 是启发式算法，理论上可能错过全局最优。论文里别写"证明了最优解"，要写"搜索到的近似最优解"。
2. 高维会退化。 维度上百时收敛变慢，可考虑增加粒子数或换差分进化（DE）。
3. 结果有随机性。 每次运行结果略有不同。写论文时要多次运行取平均/最优，并报告标准差——评委很看重这一点。
4. 约束处理要自己写。 常见做法是罚函数：越界就往目标函数上加一个很大的惩罚。
和其他算法的关系[td]
算法
思想
和 PSO 的区别
梯度下降沿坡往下走单点出发、需可导；PSO 群体出发、不需要导数
遗传算法 GA 选择、交叉、变异有"淘汰"机制；PSO 不淘汰，所有粒子活到最后、共享信息
模拟退火 SA 以一定概率接受坏解单点搜索；PSO 是群体协作
差分进化 DE 向量差扰动与 PSO 各有胜负，高维时 DE 常更稳
六、写给数模竞赛的你
回到竞赛场景。PSO 在数模里的出场率非常高，典型的用法：
题型一：优化决策。"如何分配 X 使总成本最小/总收益最大"，目标函数复杂、变量多、不可导——PSO 直接上。
题型二：模型调参。 先用回归/机器学习建了模型，再有一堆超参数要调。可以用 PSO 包在外面做参数寻优，论文里就是个漂亮的"双层优化"结构。
题型三：拟合非线性模型。 最小二乘拟合非线性曲线时，初值选不好就拟合歪。用 PSO 先粗搜一轮找初值，再交给 LM 算法精修——这是很扎实的组合拳写法。
论文里怎么写才加分： 不要只丢一个结果。写清楚三件事——① 粒子数、迭代数、参数 w/c₁/c₂ 怎么定的（最好说明试过几组）；② 收敛曲线（画"迭代代数 vs 最优值"的图，评委一眼看到收敛行为）；③ 多次运行的结果稳定性（跑 10 次报告均值和标准差）。做到这三点，PSO 就从"套模型"升级成了"用模型"。
写在最后
PSO 之所以流传三十年，恰恰因为它简单得"不像话"——两条公式，没有 leaders，没有淘汰，只有一群朴素的个体，各自带着一点点记忆和一点点从众，互相分享情报，就把一片茫茫的搜索空间梳了个遍。
有人说这是对"群体智慧"最好的计算演绎：没有一只鸟知道答案，但鸟群知道。
下次遇到解不出的优化题，不妨也撒一把粒子出去——记得先广搜，后精修。
如果这篇文章帮到了你，点个「在看」让更多备赛的朋友看到
+ J! ^( e, X7 v' q

, l5 K& G, x# k

zan