初窥PSO粒子群算法

碎碎念:终于搞完了毕设写好了论文,现在坐等答辩,可以有大把时间去做自己一直很想做的一些事情了。恰好这几天一个朋友找我帮他写程序,提出想寻找一个数值最优解,这激起了我的兴趣,查了查资料,决定用PSO粒子群算法来做这件事,下面总结下解决这个问题的思路。

PSO算法简介

粒子群算法(Particle Swarm Optimization,PSO)属于现代优化算法,适用于求解没有目标函数解析式或者解析函数很复杂的优化问题,简单地说,就是我们没办法得到诸如 f ( x ) = x 2 + 1 f(x)=x^{2}+1 f(x)=x2+1 这样可以求导的解析式,甚至根本就得不到一个具体的函数表达式,与此同时,我们还需要求出最小(大)值,在这个背景下,就出现了诸如模拟退火算法、遗传算法、粒子群算法等各种各样的优化算法。这些优化算法其实都是在模拟自然界的某些过程,并尝试着去解决问题。

粒子群算法最初是从研究鸟群捕食时得到的灵感,对于鸟群而言,如何在一大片区域中用有限的时间,最快地找到食物丰盛的领地,是一个决定鸟群生存繁衍的大问题。科学家们发现,单只鸟所能搜索的范围十分有限,但整个鸟群却像是被控制一般,很快聚集在食物最丰盛的地方。我们将每只鸟抽象成一个“粒子”(particle),每个粒子都有自己的位置和速度,并记忆自己曾经飞过的最好位置和群体曾经飞过的最好位置,那么在不断迭代的过程中,整体鸟群就很可能向最优目标出飞去。

概况起来说,一个群体中有 m m m个粒子,在一个 D D D维的空间中以一定速度飞行,每个粒子在飞行时,都能感知自己曾飞过的最好位置和整体群体的最好位置,并依据作为依据来变化。对于粒子群的第 i i i个粒子而言,它的信息如下

  • 当前速度: v i = ( v i 1 , v i 2 , . . . . v i D ) v_{i}=(v_{i1},v_{i2},....v_{iD}) vi=(vi1,vi2,....viD)
  • 当前位置: x i = ( x i 1 , x i 2 , . . . . x i D ) x_{i}=(x_{i1},x_{i2},....x_{iD}) xi=(xi1,xi2,....xiD)
  • 粒子 i i i经历的最好位置: p i , b e s t = ( p i 1 , p i 2 , . . . . p i D ) p_{i,best}=(p_{i1},p_{i2},....p_{iD}) pi,best=(pi1,pi2,....piD)
  • 粒子群经历的最好位置: p g , b e s t = ( p g 1 , p g 2 , . . . . p g D ) p_{g,best}=(p_{g1},p_{g2},....p_{gD}) pg,best=(pg1,pg2,....pgD)

PSO算法分为全局PSO和局部PSO,其中全局PSO指的是使用整个群体找到的最好解作为 p g , b e s t p_{g,best} pg,best,而局部PSO令部分成员成为邻域,并取邻域中的最好点作为 p g , b e s t p_{g,best} pg,best。在实际使用中,可以先用全局PSO找到大概位置,然后用局部PSO进行小范围内的搜索。

粒子每一个维度的位置其实就是这个问题中所包含的独立变量。在每一次的迭代过程中,粒子的位置 x i x_{i} xi都会被用来计算适应度,也就是计算目标函数的值,并依据计算的结果来决定下一次迭代过程中粒子的位置和方向。令 l ≤ d ≤ D l\leq d \leq D ldD,那么每一个粒子的第 d d d维的位置和速度按照如下的式子变化

  • v i d = w ⋅ v i d + c 1 ⋅ r a n d ⋅ ( p i d − x i d ) + c 2 ⋅ r a n d ⋅ ( p g d − x i d ) v_{id}=w\cdot v_{id}+c_{1}\cdot rand\cdot (p_{id}-x_{id})+c_{2}\cdot rand\cdot (p_{gd}-x_{id}) vid=wvid+c1rand(pidxid)+c2rand(pgdxid)
  • x i d = x i d + v i d x_{id}=x_{id}+v_{id} xid=xid+vid

对上式进行解读,得到的几个结论如下:

  • w w w代表的是惯性权重,代表了粒子上一次迭代时的速度对当前速度的影响,起到了平衡算法全局搜索能力和局部搜索能力的作用,一般令 w w w的取值范围是[0.9, 1.2]。 w w w越大,粒子飞得越快,对全局的搜索能力也越强,而 w w w越小,那么粒子飞得越慢,对局部的搜索能力就越强。
  • c 1 ⋅ r a n d ⋅ ( p i d − x i d ) c_{1}\cdot rand\cdot (p_{id}-x_{id}) c1rand(pidxid) 代表粒子向自我历史最优点进行飞行, c 2 ⋅ r a n d ⋅ ( p g d − x i d ) c_{2}\cdot rand\cdot (p_{gd}-x_{id}) c2rand(pgdxid) 代表粒子向群体最优点进行飞行,其中rand是一个在(0,1)间取值的函数,一般用均分分布的函数即可。 c 1 c_{1} c1 c 2 c_{2} c2分别代表自我历史最优和群体最优点的学习因子,也起到了平衡算法全局搜索能力和局部搜索能力的作用, c 1 c_{1} c1 c 2 c_{2} c2通常均为2,若取其他值,一般令 c 1 = c 2 c_{1}=c_{2} c1=c2且范围在[0,4]之间。这里的学习因子越大,粒子越容易在局部收敛,这与权重 w w w恰好是相反的。
  • 为了防止粒子飞离搜索空间,我们需要对粒子的速度做一个限制,令 v i d v_{id} vid的范围在[ − v m a x -v_{max} vmax, v m a x v_{max} vmax],且通常令 v m a x = k ⋅ x d , m a x v_{max}=k\cdot x_{d,max} vmax=kxd,max 0.1 ≤ k ≤ 1 0.1\leq k\leq 1 0.1k1,每一维均如此设置。

算法的整个流程图如下:

Created with Raphaël 2.3.0 开始 使用随机分布,初始化化粒子的位置和速度 计算所有粒子的适应度 更新个体最优位置和群体最优位置 计算粒子的下一个位置和速度 达到? 结束 继续下一次迭代 yes no

实例分析

这次需要解决的是一个效率优化的问题,需要求出如何分配三个暖通机的输出功率,使得整栋楼的供暖系统的热量损失率最低,下面给出关键部分的伪代码和结果图供大家参考。

#! python3
# time代表总的迭代次数, N为群体数目
for i in time:
	for j in range(N):
		result = fit(x[j])			#fit函数计算第j个粒子的适应度
		if result < pbest[j]:		#依据条件更新第j个粒子的个体历史最优点
			p[j] = x[j]
			pbest[j] = x[j]
		if pbest[j] < gbest:		#依据条件更新群体最优点
			g = p[j]
			gbest = pbest[j]
		while True:
			# p[j]为第j个粒子的个体历史最好点,g为粒子群的最优点
			v_test = w * v[j] + c1 * random.random() * (p[j] - x[j]) + c2 * random.random() * (g - x[j])
			x_test = x[j] + v_test
			if (监测是否溢出边界):
				v[j] = v_test			#若不溢出,则更新位置和速度
				x[j] = x[j] + v[j]
				break				#若不溢出,则跳出循环,更新下一个粒子的位置和速度
			else
				v[j] = random.uniform(0, 0.1)		#若溢出了,则重新分配此粒子的位置和速度
				x[j] = random.uniform(0, 1)

初始位置

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

初始速度

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传
经过50次迭代后

最终位置

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

最终速度

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

收敛曲线

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

可以很明显地看到粒子收敛在一个小区间内。

总结

综上所述,这次上手PSO算法的几个想法如下:

  • PSO算法模拟了自然界鸟群在一个大区域内进行觅食的情况,适用于没有解析解的优化问题。
  • PSO算法对搜索空间进行了全局搜索,且与传统的遗传算法不同的是,PSO且将搜索过程中发现的好的解作为知识进行保存,且分享给所有粒子,这在一定程序上提高了搜索效率。
  • PSO算法主要依靠调节粒子的速度来控制算法,且需要重点注意粒子的越界问题。

附录

《粒子群算法的基本理论及其改进研究》

### 粒子群优化(PSO)算法的流程图解析 粒子群优化(Particle Swarm Optimization, PSO)是一种基于种群智能的优化方法,其核心思想来源于鸟类群体觅食行为的研究。以下是关于PSO算法的标准流程及其对应的图解说明。 #### 1. 始化阶段 在始化过程中,需要定义种群规模、搜索空间维度以及参数设置(如惯性权重 \(w\)、加速常数 \(c_1\) 和 \(c_2\))。每个粒子的位置和速度均被随机始化于指定范围内[^2]。 ```plaintext for i=1 to N: Initialize particle[i].position randomly within the search space. Initialize particle[i].velocity with a small random value. ``` #### 2. 计算适应度函数 对于每一个粒子,计算其当前位置所对应的目标函数值作为适应度评估依据。此过程用于衡量当前解的质量水平[^3]。 #### 3. 更新个体极值和个人最佳位置(pbest) 如果某个粒子的新位置产生了更优的结果,则更新该粒子的历史最优记录 `pbest`[^4]。 #### 4. 更新全局最佳位置(gbest) 在整个种群中比较各个粒子的最佳历史表现,选取其中效果最好的一个设为整体的 `gbest` 值。 #### 5. 调整速度与位置向量 利用下面两个公式调整粒子的速度和位置: \[ v_{id}^{t+1}=wv_{id}^t+c_1r_1\left(x_{id}^\text{best}-x_{id}^t\right)+c_2r_2\left(y_d-x_{id}^t\right)\] \[ x_{id}^{t+1}=x_{id}^t+v_{id}^{t+1}\] 这里 \( r_1 \),\( r_2 \) 是独立均匀分布于区间 [0,1] 的随机变量;其他符号含义见前述解释。 #### 6. 判断终止条件 当满足预设的最大迭代次数或者达到足够的精确程度时停止循环运算,否则返回第2步继续执行直至完成全部周期。 --- ### 图形化表示 虽然无法直接展示图像文件,但可以概述如下结构来帮助理解完整的PSO工作流: 1. **Start**: 开始程序运行。 2. **Initialize Parameters & Population**: 设置始参数并创建粒子集合。 3. **Evaluate Fitness Function**: 对每只鸟儿/颗粒进行目标测试打分。 4. **Update Personal Best Position (PBest)**: 如果发现更好地点就记住它。 5. **Determine Global Best Position (GBest)**: 找出谁飞得最高最远。 6. **Adjust Velocity and Move Particles Accordingly**: 根据规则改变方向前进一小段距离。 7. **Check Stopping Criteria**: 查看是否结束?否 -> 回到步骤3; 是 -> 输出最终答案。 8. **End**: 完成整个处理序列。 --- ### 示例代码片段 以下是一个简单的Python实现版本供参考学习之用: ```python import numpy as np def pso(objective_func, bounds, num_particles, max_iter): n_dim = len(bounds) particles_pos = [] particles_vel = [] for _ in range(num_particles): pos = [np.random.uniform(bound[0], bound[1]) for bound in bounds] vel = [np.random.uniform(-abs(bound[1]-bound[0])/10, abs(bound[1]-bound[0])/10) for bound in bounds] particles_pos.append(pos) particles_vel.append(vel) personal_best_positions = list(map(list,particles_pos)) global_best_position = min(personal_best_positions,key=lambda x:objective_func(x)) c1,c2,w = 2,2,.9 history=[] for t in range(max_iter): for idx,(pos,vel,pbest) in enumerate(zip(particles_pos,particles_vel,personal_best_positions)): fitness_candidate = objective_func(pos) if(fitness_candidate<fitness_pbest:=objective_func(pbest)): personal_best_positions[idx]=list(pos) gbest_fitnesstmp=objective_func(global_best_position) candidate_fitness=objective_func(pos) if(candidate_fitness<gbest_fitnesstmp):global_best_position=list(pos) new_velocity=[wi+(c1*np.random.rand()*(pb-pi))+(c2*np.random.rand()*(gbi-pi)) for pi, pb,gbi, wi in zip(pos,pbest,global_best_position,vel)] new_position=[pi+vi for pi, vi in zip(pos,new_velocity)] particles_pos[idx]=new_position particles_vel[idx]=new_velocity history.append([*map(lambda z:list(z),zip(*[(objfctn(ij),ij)for ij in particles_pos]))]) return global_best_position,min(history[-1][0]),history ``` ---
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值