前言

粒子群算法 (Particle Swarm Optimization, PSO) 是启发式算法的一种, 其原理是在解空间中设置多个可移动的节点, 每轮迭代中, 每个节点根据自身经历过的最优点和全局经历过的最优点计算自身的下一个移动方向和位置. 这个过程与生物群体觅食的行为十分类似, 因此也可以说 PSO 是一种仿生算法. 本文将介绍基于 Python 的 PySwarms 库, 快速实现一个自定义目标函数的解决方案.

首先给出 PySwarms 库的文档链接, 在读完本文后如果还有不理解的地方可以参考该文档: Welcome to PySwarms’s documentation!.

此外, 本文将默认你已经懂得什么是 PSO 及其原理, 如果你尚不了解 PSO 下列是一些参考文章/链接:

PySwarms 的安装及使用

安装

由于 conda 并没有包括 PySwarms, 因此我们使用 pip 安装 PySwarms. 如果你想要从源代码安装, 请参考 Installation — PySwarms documentation.

1
2
pip install pyswarms
# 注意: pyswarm 是另一个库, 不要安装错了

PySwarms 仅依赖于 Numpy 因此请确保你已经安装了 Numpy.

使用基本功能

PySwamrs 支持 PSO 和 离散二进制粒子群算法 (Discrete Binary Particle Swarm Optimization Algorithm, BPSO), 由于这两种算法在 PySwarms 的调用方式十分类似, 因此本文不对后者的使用进行深入的介绍, 而是以前者为主.

首先在你的代码中导入 PySwarms, 同时在使用 PySwarms 的同时往往会使用到 Numpy, 因此也一并导入.

1
2
import pyswarms as ps
import numpy as np

接下来使用一个简单的函数来作为优化对象, 整个优化过程包括三个步骤: 1. 创建目标函数; 2. 根据参数初始化粒子群; 3. 调用优化函数进行求解.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
# 这是 PySwarms 自带的一个目标函数库, 其中包括了很多函数对象, 更多请参考:
# https://pyswarms.readthedocs.io/en/latest/api/pyswarms.utils.functions.html
from pyswarms.utils.functions import single_obj as fx # 1

# 2. 根据参数初始化粒子群
# 参数由一个字典输入, c1 c2 分别是个体/群体加速常数, w 是惯性权重
# 此外还有两个可选参数 k (考虑邻居的数量) 和 p (要使用的闵可夫斯基 p-范数), 具体请参考文档
options = {'c1': 0.5, 'c2': 0.3, 'w':0.9}
dimensions = 2 # 设置粒子的维度, 也就是有几个 feature
# 此外可以设定粒子所在的位置, 例如我们规定所有粒子的所有维度的值都应当位于 (-5, 5) 之间
max_bound = 5 * np.ones(dimensions)
min_bound = - max_bound
bounds = (min_bound, max_bound)

# n_particles: 粒子群中的粒子数量, dimensions: 每个粒子的维度
# options: 参数选项, bounds: 粒子范围限制
optimizer = ps.single.GlobalBestPSO(n_particles=10,
dimensions=dimensions,
options=options,
bounds=bounds)

# 3. 求解
obj_fun = fx.rastrigin_func # 目标函数
# obj_fun: 要优化的目标函数, print_step: 每次结果打印的间隔, 例如这里是每迭代 100 此打印一次结果
# iters: 迭代次数, verbose: 输出信息的冗长形式
# 第一个返回值是最小 cost, 第二个返回值是最小 cost 对应的粒子位置
cost, pos = optimizer.optimize(obj_fun, print_step=100, iters=1000, verbose=3)

接下来你就会看到类似于下文的输出:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
INFO:pyswarms.single.global_best:Iteration 1/1000, cost: 12.243865048066269
INFO:pyswarms.single.global_best:Iteration 101/1000, cost: 1.1759164022634394
INFO:pyswarms.single.global_best:Iteration 201/1000, cost: 0.9949603350768896
INFO:pyswarms.single.global_best:Iteration 301/1000, cost: 0.9949590581556009
INFO:pyswarms.single.global_best:Iteration 401/1000, cost: 0.9949590570934177
INFO:pyswarms.single.global_best:Iteration 501/1000, cost: 0.9949590570932898
INFO:pyswarms.single.global_best:Iteration 601/1000, cost: 0.9949590570932898
INFO:pyswarms.single.global_best:Iteration 701/1000, cost: 0.9949590570932898
INFO:pyswarms.single.global_best:Iteration 801/1000, cost: 0.9949590570932898
INFO:pyswarms.single.global_best:Iteration 901/1000, cost: 0.9949590570932898
INFO:pyswarms.single.global_best:================================
Optimization finished!
Final cost: 0.9950
Best value: [3.5850411183743393e-09, -0.9949586379966202]

其中 Best value 就是算法求得的最优解, Final cost 是对应的最优值.

可视化

绘制 cost history

有时我们需要将算法过程中的 cost history 绘制成图片, PySwarms 已经提供了这一功能, 使用方法十分简单:

1
2
3
# optimizer: 已经完成了求解的粒子群优化实例
lot_cost_history(cost_history=optimizer.cost_history)
plt.show()
Cost History

Cost History

绘制粒子运动过程

此外, 如果你想把粒子运动的过程绘制成动态图片, 可以参考下的代码块 (二维粒子):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
# 这段代码需要在 Jupiter Notebook 中运行, 否则可能无法输出正确的动态图片
# equivalent to rcParams['animation.html'] = 'html5'
# See http://louistiao.me/posts/notebooks/save-matplotlib-animations-as-gifs/
rc('animation', html='html5')

from pyswarms.utils.plotters.formatters import Mesher
# Initialize mesher with sphere function
m = Mesher(func=fx.sphere_func)

# Make animation
animation = plot_contour(pos_history=optimizer.pos_history,
mesher=m,
mark=(0,0))
# Enables us to view it in a Jupyter notebook
HTML(animation.to_html5_video())
2D 粒子运动过程

2D 粒子运动过程

上图展示了二维粒子的运动, 下图展示了三维例子的运动.

3D 粒子运动过程

3D 粒子运动过程

一般来说, 我们只能绘制不高于三维空间中的粒子, PySwarms 的文档中已经包括了绘制二维粒子运动和三维粒子运动的代码, 详情请参考 Visualization — PySwarms documentation.

自定义目标函数

普通目标函数

上文中使用的是库自身提供的优化目标函数, 但往往我们需要自己选择优化目标而库中并不提供, 因此下面介绍如何使用自定义的函数来进行粒子群优化.

首先我们介绍一下 obj_fun 的输入和输出 (参考: Feature Subset Selection — PySwarms documentation#Writing the custom-objective function):

  • obj_fun 的输入应当是一个 numpy.ndarray 格式的 m*n 多维列表, 每一行代表一个粒子, 共 m 个粒子, 每一列表示粒子在一个维度上的值, n 等于创建 optimizer 时的 dimensions 参数.
  • obj_fun 输出应当是一个长度为 m 的 array-like 类型, 值为各个粒子的 cost.

由于 Python 的鸭子类型特点, 我们只需要构造一个能够实现正确输入输出的目标函数即可. 下面我们构造一个二次函数优化作为例子:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
import pyswarms as ps
import numpy as np

def obj_fun(x):
# 注意: 这里要选择每个粒子的第一个 feature, 否则 power 函数接收的将会是一个 array
# 选择第一个 feature 的原因是我们的例子中每个粒子只有一个 feature, 请灵活运用
# 如果不会写成 numpy 形式的话用循环或列表生成器来解决也是可以的, 只是可能会影响效率
# 例如: return [np.power(i[0], 2) for i in x]
return np.power(x[:, 0], 2)

max_bound = 5 * np.ones(1)
min_bound = - max_bound

optimizer = ps.single.GlobalBestPSO(n_particles=10, dimensions=1,
options={'c1': 0.5, 'c2': 0.3, 'w':0.9},
bounds=(min_bound, max_bound))

cost, pos = optimizer.optimize(obj_fun, print_step=10, iters=100, verbose=3)

有时在初始化粒子群的时候可能会出现下面这句话, 但是似乎并不影响结果.

1
INFO:pyswarms.backend.topology.base:Running on `dynamic` topology, neighbors are updated regularly.Set `static=True` for fixed neighbors.

运行结束后, 你将得到类似下文的输出:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
INFO:pyswarms.single.global_best:Arguments Passed to Objective Function: {}
INFO:pyswarms.single.global_best:Iteration 1/100, cost: 0.0034462728980131944
INFO:pyswarms.single.global_best:Iteration 11/100, cost: 0.00012490737667841898
INFO:pyswarms.single.global_best:Iteration 21/100, cost: 5.316504393857646e-05
INFO:pyswarms.single.global_best:Iteration 31/100, cost: 2.3769348137312317e-06
INFO:pyswarms.single.global_best:Iteration 41/100, cost: 8.29034246626045e-07
INFO:pyswarms.single.global_best:Iteration 51/100, cost: 8.29034246626045e-07
INFO:pyswarms.single.global_best:Iteration 61/100, cost: 4.565048551676232e-08
INFO:pyswarms.single.global_best:Iteration 71/100, cost: 4.565048551676232e-08
INFO:pyswarms.single.global_best:Iteration 81/100, cost: 3.964141359304464e-08
INFO:pyswarms.single.global_best:Iteration 91/100, cost: 5.7765847812579206e-09
INFO:pyswarms.single.global_best:================================
Optimization finished!
Final cost: 0.0000
Best value: [7.600384714774588e-05]

最后的结果表明, 当粒子位于 [7.600384714774588e-05] 时, 目标函数取得最优值 0.0000.

带参数的目标函数

有时你可能会遇到这样的问题: 目标函数除了粒子的位置, 还有一些额外但是固定的参数 (例如上文中二次函数的对称轴, 截距等等), 此时当然可以将这些参数写死, 或者使用函数式编程中偏函数的思路来解决, 但 PySwarms 本身就提供了传入参数的功能, 在调用 optimizer.optimize() 时可以直接传入 obj_fun:

1
2
3
4
5
6
# 省略了初始化粒子群的部分
def rosenbrock_with_args(x, a, b, c=0):
f = (a - x[:, 0]) ** 2 + b * (x[:, 1] - x[:, 0] ** 2) ** 2 + c
return f
cost, pos = optimizer.optimize(rosenbrock_with_args, 1000, print_step=100, verbose=3,
a=1, b=100, c=0) # 注意这里, 直接传入了 rosenbrock_with_args() 的参数

除了上面的传参方式之外, 还可以传入一个参数字典:

1
2
kwargs={"a": 1.0, "b": 100.0, 'c':0}  # 在 rosenbrock_with_args() 中可以直接调用 a, b, c
cost, pos = optimizer.optimize(rosenbrock_with_args, 1000, print_step=100, verbose=3, **kwargs)

结语

自此, 基于 PySwarms 的自定义目标函数粒子群优化就实现完成了, 如果有更多问题, 欢迎联系我或者查看 PySwarms 的文档. 如果你想要了解 PySwarms 的更多细节, 可以到 PySwarms 的 Github 页面 GitHub - ljvmiranda921/pyswarms: A research toolkit for particle swarm optimization in Python 阅读源码.