1 粒子群算法
1 粒子群算法
1.1 概述
粒子群算法(Particle Swarm Optimization,PSO)由Kennedy和Eberhart于1995年提出。该算法的思想来源于对鸟类捕食行为的研究,鸟之间通过集体的协作使得群体能够找到最多的食物,PSO便是通过模拟鸟群飞行觅食的行为,来寻找最优解的算法,这是一种基于群体智能(Swarm Intelligence)的优化方法。
在粒子群算法中,我们将鸟群抽象成粒子群,用一个粒子来代表一只鸟。PSO目标就是:通过一群粒子在解空间中进行搜索,找到使得适应度函数(fitness function)取得最大值(或最小值)的解。
1.2 算法介绍
在正式介绍算法前,先约定好一些符号。假设粒子在\(d\)维空间中进行搜索,种群大小为N(一共有N个粒子),适应度函数为\(f(\cdot)\)。
- \(x_i^{(k)}\in R^d\)表示第\(i\)个粒子在迭代\(k\)次后的位置
- \(v_i^{(k)}\in R^d\)表示第\(i\)个粒子的在迭代\(k\)次后的速度
- \(pbest_i\in R^d\)表示第\(i\)个粒子所经过的历史最优位置
- \(gbest\in R^d\)表示整个种群所找到的历史最优位置
那么第\(i\)个粒子位置更新公式为:
其中
- \(w\)是惯性系数,\(w\)越大说明粒子更倾向于保持原来的运动状态
- \(c_1\)是自我认知系数,\(c_1\)越大说明粒子更倾向于相信自己的经验(认知)
- \(c_2\)是全局认知系数,\(c_2\)越大说明粒子更倾向于相信整个群体的经验(认知)
- \(r_1,r_2\)是随机数,通常是\([0,1]\)上均匀产生的随机数
从更新公式可以看出,粒子下一时刻的运动速度收到三个因素的影响:当前时刻的速度、当前自己找到的最优解的位置以及当前全局最优解的位置。可以用一个图表示粒子位置的更新过程:
基本算法描述如下:
- 首先,人为设置粒子的种群大小\(N\),适应度函数\(f(\cdot)\),惯性系数\(\omega\),以及\(c_1,c_2\);
- 随机初始化每个粒子的位置\(x_i^{(0)}\),初始速度\(v_i^{(0)}\);然后记录当前每个粒子的历史最优解和全局最优解。
- 不断循环迭代,每循环一次就更新种群中每个粒子的状态,直至达到算法的停止条件。最后,输出所找到的最优解。
算法的停止条件:
- 最简单的就是直接设置一个最大迭代上限,超出最大迭代次数后直接退出。
- 设置一个计数器,并设置一个阈值\(C\)。如果,PSO在搜索过程中,连续\(C\)次循环所找到的全局最优解都没有发生变化,就停止。
- 设置一个阈值\(T\),如果PSO所找到的全局最优的解满足\(f(gbest)\ge T\),那么停止算法。
1.3 算法实现与测试
这里我们PSO算法来搜索ackley函数的最小值,ackley函数是一个具有非常多个局部极小值的函数,具体表达式如下:
通常取\(a=20,b=0.2,c=2\pi\)。其中\(d\)表示空间的维数,即\(x\in\mathbb{R}^d\)。该函数具有全局最小值\(f(x^*)=0,x^*=0\)。函数图像如下(三维空间):
实现代码:
import numpy as np # v(k+1) = w*v(k) + c1*r1*(gbest-x) + c2*r2*(pbest-x) # x(k+1) = x(k) + v(k+1) class Particle: """ st_x:位置约束,粒子每个维度上的坐标范围必须处于[st_x[0], st_x[1]]之间 st_v:速度约束,粒子每个维度上的速度范围必须处于[st_v[0], st_v[1]]之间 position:粒子的当前位置 velocity:粒子的当前速度 pbest:粒子自身历史记录的最佳位置 pvalue:粒子自身历史记录的最佳值 """ st_x = None st_v = None def __init__(self, x, v, pbest, pvalue): self.position = x self.velocity = v self.pbest = pbest self.pvalue = pvalue def update_velocity(self, v): v[v < self.st_v[0]] = self.st_v[0] v[v > self.st_v[1]] = self.st_v[1] self.velocity = v def update_position(self, func): # 更新粒子自身的位置,以及判断是否要更新pbest self.position += self.velocity self.position[self.position < self.st_x[0]] = self.st_x[0] self.position[self.position > self.st_x[1]] = self.st_x[1] if self.pvalue > func(self.position): self.pvalue = func(self.position) self.pbest = self.position.copy() class PSO: """ gbest:粒子群历史记录的最佳位置 gvalue:粒子群历史记录的最佳值 """ gbest = None gvalue = np.inf def __init__(self, n_dims, n_particles, st_x, st_v, w, c1, c2, num_iter, func): # 初始化空间维度 self.n_dims = n_dims # 初始化粒子群数目 self.n_particles = n_particles # 目标函数 self.func = func # 粒子惯性权重 self.w = w # 全局部分学习率 self.c1 = c1 # 自我认知部分学习率 self.c2 = c2 # 迭代次数 self.num_iter = num_iter # 存放粒子群的列表(容器) self.particles = [] # 初始化粒子的位置和速度约束 Particle.st_x = st_x Particle.st_v = st_v # 初始化粒子群 for _ in range(n_particles): # 初始化粒子的随机位置在 st_x[0]~st_x[1] x = (st_x[1] - st_x[0]) * np.random.rand(n_dims) + st_x[0] # 计算当前评估值 pvalue = func(x) # 初始化一个粒子 self.particles.append( Particle( x=x, v=(st_v[1] - st_v[0]) * np.random.rand(n_dims) + st_v[0], pbest=x.copy(), pvalue=pvalue ) ) if self.gvalue > pvalue: self.gvalue = pvalue self.gbest = x.copy() def solve(self): # 开始迭代 for index in range(1, self.num_iter + 1): for particle in self.particles: v = self.w * particle.velocity + self.c1 * np.random.rand() * (self.gbest - particle.position) + self.c2 * np.random.rand() * (particle.pbest - particle.position) particle.update_velocity(v) particle.update_position(self.func) for particle in self.particles: if particle.pvalue < self.gvalue: self.gvalue = particle.pvalue self.gbest = particle.pbest.copy() return self.gbest, self.gvalue def ackley(x): return - 20 * np.exp(-0.2 * np.sqrt((x * x).mean())) - np.exp(np.cos(2 * np.pi * x).mean()) + 20 + np.exp(1) if __name__ == "__main__": # 测试不同迭代次数搜索出来的最优解情况 for i in range(11): pso = PSO( n_dims=10, n_particles=100, st_x=(-20, 20), st_v=(-1, 1), w=0.8, c1=2, c2=2, num_iter=50 * i, func=ackley ) gbest, gvalue = pso.solve() print("number of interations:%d\tgvalue:%f" % (50 * i, gvalue))
测试结果:
可以看到随着迭代次数的增加,所搜索到的最优解基本上都是越来越小的。大多数情况下都是这样的,不排除个别情况下,迭代次数多,反而找到更差的解。还需要注意的是当搜索空间维数太高,或者搜索区域很大的时候,粒子数太少往往很难找到全局最优解。但是如果增加粒子数,增加迭代次数提高算法搜索能力,又需要很大的时间开销。