scikit-opt 粒子群算法(PSO)动态可视化实战:用 record_mode 记录轨迹并生成 GIF 动画
【免费下载链接】scikit-optGenetic Algorithm, Particle Swarm Optimization, Simulated Annealing, Ant Colony Optimization Algorithm,Immune Algorithm, Artificial Fish Swarm Algorithm, Differential Evolution and TSP(Traveling salesman)项目地址: https://gitcode.com/GitHub_Trending/sci/scikit-opt
本指南以 scikit-opt 开源库中粒子群算法(PSO)的动画展示为核心,讲解如何在带非线性约束的二维寻优问题上运行 PSO,并借助record_mode记录每一代粒子的位置与速度,最后用 Matplotlib 的FuncAnimation将整个收敛过程渲染为动态 GIF。读完本文,你将掌握从「跑算法」到「做动画」的完整链路,并理解 scikit-opt 中记录机制与粒子群迭代的底层实现。
一、整体思路:两个 step 完成 PSO 可视化
粒子群算法的动态展示分为两个阶段,对应仓库中的示例文件 examples/demo_pso_ani.py:
- 做 PSO:构造目标函数与约束,初始化
PSO对象并开启record_mode,运行迭代; - 画动画:读取
pso.record_value中记录的每一代粒子位置X与速度V,用matplotlib.animation.FuncAnimation逐帧播放,导出为pso.gif。
之所以要先开record_mode,是因为 PSO 类在默认情况下只保留最终结果(如gbest_x、gbest_y),并不会保留每一代的全量粒子快照。动画素材正是依赖这个开关在迭代过程中持续"录影"。
二、Step 1:运行带约束的 PSO 并开启记录模式
首先定义目标函数。示例使用的是 Ackley 函数风格的二维测试函数(多峰、存在大量局部最优,适合演示粒子群的收敛行为):
import numpy as np from sko.PSO import PSO def demo_func(x): x1, x2 = x return -20 * np.exp(-0.2 * np.sqrt(0.5 * (x1 ** 2 + x2 ** 2))) - np.exp( 0.5 * (np.cos(2 * np.pi * x1) + np.cos(2 * np.pi * x2))) + 20 + np.e constraint_ueq = ( lambda x: (x[0] - 1) ** 2 + (x[1] - 0) ** 2 - 0.5 ** 2 , ) max_iter = 50 pso = PSO(func=demo_func, n_dim=2, pop=40, max_iter=max_iter, lb=[-2, -2], ub=[2, 2] , constraint_ueq=constraint_ueq) pso.record_mode = True pso.run() print('best_x is ', pso.gbest_x, 'best_y is', pso.gbest_y)关键参数与含义(对应 PSO 构造函数)
| 参数 | 示例取值 | 说明 |
|---|---|---|
func | demo_func | 待优化的目标函数,输入为决策向量x,输出为标量适应值 |
n_dim | 2 | 决策变量维度,即x中分量的个数 |
pop | 40 | 粒子群规模(particle 数量),与 GA 中的种群规模口径一致 |
max_iter | 50 | 最大迭代代数 |
lb/ub | [-2, -2]/[2, 2] | 每个维度的下界/上界;构造时会被广播为n_dim维,并断言dim == len(lb) == len(ub)且ub > lb(见 PSO.py#L99-L101) |
constraint_ueq | 一个元组 | 不等式约束(≤ 0),可传多个,约束函数返回大于 0 视为不可行 |
w | 0.8(默认) | 惯性权重,控制粒子速度的惯性 |
c1/c2 | 0.5(默认) | 认知系数与群体系数,分别控制跟随个体历史最优pbest与全局最优gbest的力度 |
constraint_ueq中传入的约束(x[0] - 1) ** 2 + (x[1] - 0) ** 2 - 0.5 ** 2表示:可行域是圆心在(1, 0)、半径为0.5的圆内(≤ 0)。在 check_constraint 中,只要任一约束函数返回值大于 0,该粒子即被判为不可行;这一判定在更新个体最优update_pbest时生效——即使适应值更优,若违反约束也不会被采纳为新的pbest。
record_mode 与 run() 的底层逻辑
pso.record_mode = True是关键一步。在 PSO 类 中,默认record_mode = False,而record_value初始化为{'X': [], 'V': [], 'Y': []},分别用于存放每一代粒子的位置矩阵(pop × n_dim)、速度矩阵与适应值向量。
调用run()后,每次迭代的主循环为:update_V()→recorder()→update_X()→cal_y()→update_pbest()→update_gbest()(见 run 方法)。recorder()(见 PSO.py#L169-L174)在每次迭代时把当前的X、V、Y追加进record_value,从而形成完整的历史轨迹。迭代结束后:
pso.gbest_x/pso.gbest_y:全局历史最优位置与适应值(示例中的print输出);pso.record_value['X']:长度为max_iter的列表,第i项是第i代的全部粒子位置;pso.record_value['V']:同样的结构,保存每代速度,是动画插值的关键数据。
需要说明的是,位置更新公式遵循标准 PSO:V = w·V + c1·r1·(pbest_x - X) + c2·r2·(gbest_x - X),随后X = clip(X + V, lb, ub)(见 update_V / update_X)。粒子的速度被限制在ub - lb的范围内初始化,保证搜索步长不至于失控。
三、Step 2:读取记录数据并渲染逐帧动画
运行结束后,把record_value交给 Matplotlib 做动画:
import matplotlib.pyplot as plt from matplotlib.animation import FuncAnimation record_value = pso.record_value X_list, V_list = record_value['X'], record_value['V'] fig, ax = plt.subplots(1, 1) ax.set_title('title', loc='center') line = ax.plot([], [], 'b.') X_grid, Y_grid = np.meshgrid(np.linspace(-2.0, 2.0, 40), np.linspace(-2.0, 2.0, 40)) Z_grid = demo_func((X_grid, Y_grid)) ax.contour(X_grid, Y_grid, Z_grid, 30) ax.set_xlim(-2, 2) ax.set_ylim(-2, 2) t = np.linspace(0, 2 * np.pi, 40) ax.plot(0.5 * np.cos(t) + 1, 0.5 * np.sin(t), color='r') plt.ion() p = plt.show() def update_scatter(frame): i, j = frame // 10, frame % 10 ax.set_title('iter = ' + str(i)) X_tmp = X_list[i] + V_list[i] * j / 10.0 plt.setp(line, 'xdata', X_tmp[:, 0], 'ydata', X_tmp[:, 1]) return line ani = FuncAnimation(fig, update_scatter, blit=True, interval=25, frames=max_iter * 10) plt.show() ani.save('pso.gif', writer='pillow')画面构成拆解
- 背景等值线:
np.meshgrid生成[-2, 2] × [-2, 2]的 40×40 网格,代入demo_func得到适应值曲面,ax.contour(..., 30)画出 30 层等高线,直观呈现目标函数的多峰地貌; - 红色约束圆:用参数方程
(0.5·cos(t) + 1, 0.5·sin(t))描绘约束边界,动画中可以看到粒子被限制在圆内搜索; - 蓝色散点:
line = ax.plot([], [], 'b.')先声明一个空的蓝色散点对象,每帧用plt.setp更新其坐标数据。
update_scatter 的插值机制:为什么 frames 是 max_iter × 10
FuncAnimation(fig, update_scatter, blit=True, interval=25, frames=max_iter * 10)共播放50 × 10 = 500帧。每帧的回调函数把帧号frame拆成i, j = frame // 10, frame % 10:
i是迭代代数(0~49),决定取哪一帧记录数据;j是帧内子步(0~9),用X_list[i] + V_list[i] * j / 10.0在相邻两代之间做线性插值——把粒子从第i代位置按速度方向均匀推进 10 步,得到平滑的运动轨迹,避免逐代跳变造成的画面闪烁。
也就是说,record_value里记录的 50 代数据被展开为 500 帧动画,帧间隔interval=25ms,总时长约 12.5 秒,能清晰看到粒子群从随机初始分布逐渐聚拢到约束圆内最优解附近的完整过程。ax.set_title('iter = ' + str(i))同步显示当前帧对应的迭代代数。
保存 GIF
ani.save('pso.gif', writer='pillow')writer='pillow'指定使用 Pillow 作为 GIF 编码器(需已安装pillow库)。保存前请确保当前进程持有plt.show()打开的窗口,或在非交互环境下直接执行ani.save(...)。
四、从动画反推迭代质量:配合 gbest_y_hist 做收敛分析
动画展示的是"空间"维度的收敛过程,若要定量评估优化效果,可结合 PSO 类的gbest_y_hist属性——它在 run 方法 中每轮迭代追加当前gbest_y,构成全局最优适应值的历史曲线。示例中print('best_x is ', pso.gbest_x, 'best_y is', pso.gbest_y)输出的即是最终解:
import matplotlib.pyplot as plt plt.plot(pso.gbest_y_hist) plt.show()将动画与收敛曲线对照使用,可以直观判断:等高线图上粒子是否陷入局部最优、约束圆边界附近粒子是否被罚函数机制拒之门外,以及w、c1、c2参数调整后收敛速度的变化。
五、把动画方案迁移到自己的问题上
本示例是一个完整的可复用模板,迁移到自定义问题时只需替换三处:
- 目标函数:保持
def demo_func(x): x1, x2 = x的签名,n_dim与lb/ub随之调整。注意:若函数只接受单参数向量(func(x)形式),scikit-opt 会通过 func_transformer 自动包装为逐行调用;多入参旧式写法会触发弃用警告。 - 约束元组:
constraint_ueq里可添加多个lambda x: ...不等式,全部 ≤ 0 才视为可行;等式约束constraint_eq在当前 PSO 实现中尚不可用(见 PSO.py#L60-L63 的注释说明)。 - 绘图坐标:把
np.meshgrid(np.linspace(-2.0, 2.0, 40), ...)的网格范围改为与lb/ub一致,约束边界曲线换成你自己的可行域几何形状。
如果不需要插值平滑,也可以让frames=max_iter,在update_scatter中直接取X_list[i],生成逐代跳变的帧序列。
六、小结
本文围绕 scikit-opt 的 PSO 动画示例,完整覆盖了从参数配置(w、c1、c2、pop、lb/ub、constraint_ueq)、record_mode记录机制(PSO.py 的 recorder/record_value)、迭代更新原理(速度-位置公式),到FuncAnimation逐帧渲染与 GIF 导出的全部环节。核心示例代码见 examples/demo_pso_ani.py,基础用法与非线性约束的更多说明见 docs/zh/README.md 第 3 节,相关收敛曲线示例可参考 examples/demo_pso.py。掌握了这套模板,你就可以为自己的优化问题快速生成可复现、可展示的粒子群搜索过程动画。
【免费下载链接】scikit-optGenetic Algorithm, Particle Swarm Optimization, Simulated Annealing, Ant Colony Optimization Algorithm,Immune Algorithm, Artificial Fish Swarm Algorithm, Differential Evolution and TSP(Traveling salesman)项目地址: https://gitcode.com/GitHub_Trending/sci/scikit-opt
创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考