作者：非妃是公主
专栏：《智能优化算法》
博客地址：https://blog.csdn.net/myf_666
个性签：顺境不惰，逆境不馁，以心制境，万事可成。——曾国藩
专栏推荐
| 专栏名称 | 专栏地址 |
|---|---|
| 软件工程 | 专栏——软件工程 |
| 计算机图形学 | 专栏——计算机图形学 |
| 操作系统 | 专栏——操作系统 |
| 软件测试 | 专栏——软件测试 |
| 机器学习 | 专栏——机器学习 |
| 数据库 | 专栏——数据库 |
| 算法 | 专栏——算法 |
序
粒子群算法也是一种智能优化算法，该算法作为一种群体智能算法，主要借鉴了鸟群的思想（那为什么不叫鸟群算法，非要叫粒子群算法呢？我不得而知，感兴趣的朋友可以在评论区告诉我，我会标明出处后进行补充）。
一、概论
粒子群优化算法对鸟群捕食行为进行研究：一群鸟在区域中随机搜索食物，所有鸟知道自己距离食物的距离。
在粒子群算法中，每一个潜在解对应一只小鸟，可能是不考虑这只小鸟的重量和大小等物理特性，就抽象它为一个点，所以叫做粒子。
粒子群算法首先在给定的空间中随机初始化粒子群，待优化问题的变量决定了解空间的维数。每个粒子有了初始位置和速度后，通过迭代，我们就可以确定它的下一次的位置。
在每一次迭代中，每个粒子通过跟踪两个“极值”来更新自己在解空间中的空间位置和飞行速度：一个极值就是单个粒子本身在迭代过程中找到的最优解粒子，这个粒子叫做个体极值。另一个极值是种群所有粒子在迭代过程中找到的最优解粒子，叫做全局极值。这个方法叫做全局粒子群算法，也是粒子群算法的 v1 版本，但是由于它考虑了全局种群，因此非常容易陷入到局部最优。
对应的改进方法是，只用其中一部分作为邻居粒子，这些邻居粒子中的极值就是局部极值，用它来代替上面的全局极值，这种方法就叫做局部粒子群算法，也是粒子群算法的 v2 版本。
二、粒子群算法原理
1. 算法原理
在一个D维的目标搜索空间中，N个粒子组成一个种群，其中第i个粒子表示为一个D维的向量：
X i = ( x i 1 , x i 2 , . . . , x i D ) , i = 1 , 2 , . . . , N X_i=(x_{i1},x_{i2},...,x_{iD}),\quad i=1,2,...,N Xi=(xi1,xi2,...,xiD),i=1,2,...,N
第i个粒子的“飞行”速度也是一个D维的向量，记为：
V i = ( v i 1 , v i 2 , . . . , v i D ) , i = 1 , 2 , . . . N V_i=(v_{i1},v_{i2},...,v_{iD}),\quad i=1,2,...N Vi=(vi1,vi2,...,viD),i=1,2,...N
迄今为止，第i个粒子的个体最优位置：
p b e s t = ( p i 1 , p i 2 , . . . , p i D ) , i = 1 , 2 , . . . N p_{best}=(p_{i1},p_{i2},...,p_{iD}),\quad i=1,2,...N pbest=(pi1,pi2,...,piD),i=1,2,...N
迄今为止，整个粒子群的最优位置：
g b e s t = ( g 1 , g 2 , . . . , g D ) g_{best}=(g_1,g_2,...,g_D) gbest=(g1,g2,...,gD)
对于每一次迭代，粒子根据如下公式进行更新：
v i j ( t + 1 ) = v i j ( t ) + c 1 r 1 ( t ) [ p i j ( t ) − x i j ( t ) ] + c 2 r 2 ( t ) [ g i j ( t ) − x i j ( t ) ] v_{ij}(t+1)=v_{ij}(t)+c_1r_1(t)[p_{ij}(t)-x_{ij}(t)]+c_2r_2(t)[g_{ij}(t)-x_{ij}(t)] vij(t+1)=vij(t)+c1r1(t)[pij(t)−xij(t)]+c2r2(t)[gij(t)−xij(t)]
x i j ( t + 1 ) = x i j ( t ) + v i j ( t + 1 ) x_{ij}(t+1)=x_{ij}(t)+v_{ij}(t+1) xij(t+1)=xij(t)+vij(t+1)
其中， c 1 c_1 c1和 c 2 c_2 c2为学习因子，就是2个常数，自己设定；
r 1 , r 2 r_1,r_2 r1,r2为 [ 0 , 1 ] [0,1] [0,1]范围内的均匀随机数。
i = 1 , 2 , . . . , N P ; i=1,2,...,NP; i=1,2,...,NP;表示种群数量， j = 1 , 2 , . . . , D j=1,2,...,D j=1,2,...,D表示种群维度。
2. 机理分析
观察速度更新的公式，主要分为三个部分，进行机理分析，如下：
第一部分：惯性、动量部分。保持原有的运动趋势。
第二部分：认知部分。根据以往的经验，忘自己的个体最优值附近进行一定的偏移。
第三部分：社会部分。反映了粒子间的合作与知识共享，向全局最优值附近有一定的偏向。
三、算法流程
1. 算法伪代码
- 初始化种群。包括：NP， x i x_i xi和速度 v i v_i vi。
- 计算每个粒子适应度fit[i]
- 对每个粒子，用它的适应度值fit[i]和个体极值 p b e s t ( i ) p_{best}(i) pbest(i)比较，如果 f i t [ i ] < p b e s t fit[i]<p_{best} fit[i]<pbest，那么就更新 p b e s t p_{best} pbest。
- 将fit[i]和全局极值 g b e s t g_{best} gbest比较。如果 f i t [ i ] < g b e s t fit[i]<g_{best} fit[i]<gbest，则更新 g b e s t g_{best} gbest。
- 迭代更新粒子的速度 v i v_i vi得到 v i + 1 v_{i+1} vi+1。
- 迭代更新 x i x_i xi得到 x i + 1 x_{i+1} xi+1
- 对 x i + 1 x_{i+1} xi+1进行边界处理。
- 判断算法是否满足终止条件。如果满足，结束算法。不满足，返回步骤2.
2. 算法流程图
四、参数设置
1. 粒子种群规模 N N N
一般设置粒子数为20~50。
-
对于大部分问题10个粒子即可取得很好的结果；
-
对于复杂的问题，粒子数量可以取到100或200。
另外，粒子数目越大，搜索空间越大，全局搜索能力更强，运行时间也越长。
2. 惯性权重 ω \omega ω
Ⅰ. 固定权重
在进化过程中保持不变，一般取值为 [ 0.8 , 1.2 ] [0.8,1.2] [0.8,1.2].
Ⅱ. 时变权重
粒子在不同阶段拥有不同的搜索和开发能力。
ω = ω m a x − ( ω m a x − ω m i n ) ⋅ t T m a x \omega=\omega_{max}-\frac{(\omega_{max}-\omega_{min})\cdot t}{T_{max}} ω=ωmax−Tmax(ωmax−ωmin)⋅t
其中， T m a x T_{max} Tmax表示最大进化代数； ω m i n \omega_{min} ωmin表示最小惯性权重； ω m a x \omega_{max} ωmax表示最大惯性权重； t t t表示当前迭代次数。在大多数应用中 ω m a x = 0.9 \omega_{max}=0.9 ωmax=0.9， ω m i n = 0.4 \omega_{min}=0.4 ωmin=0.4
3. 加速常数 c 1 c_1 c1和 c 2 c_2 c2
一般设置 c 1 = c 2 = 1.5 c_1=c_2=1.5 c1=c2=1.5.
4. 粒子最大速度 v m a x v_{max} vmax
设定 v m a x v_{max} vmax和调整惯性权重是等效的，所以 v m a x v_{max} vmax一般用于对种群初始化进行设定，即将 v m a x v_{max} vmax设定为每维变量的变化范围，二不再对最大速度进行细致的选择及调节。
5. 停止准则
最大迭代次数、计算精度或最优解的最大停滞部署 △ t \bigtriangleup t △t.
6. 邻域结构的设定
-
全局最优解收敛块，但回陷入局部最优。
-
局部版本粒子群算法，收敛慢，但不容易陷入局部最优。
-
通常的做法是，全局粒子群算法寻找到最优解的大致范围，然后采用局部粒子群算法在最优点附近进行精细搜索。
7. 边界条件处理
当朝贡国最大位置限制 x m a x x_{max} xmax和最大速度限制 v m a x v_{max} vmax时，在范围内随机产生一个数值代替，或者将其设置为最大值，即边界吸收。
五、仿真实例
1. 问题
求下列函数的最小值，其中x的维数n=10。
f ( x ) = ∑ i = 1 n x i 2 , ( − 20 ≤ x i ≤ 20 ) f(x)=\sum_{i=1}^{n}x_i^2,\quad(-20\leq x_i \leq 20) f(x)=i=1∑nxi2,(−20≤xi≤20)
2. 分析
函数只有一个极小点 x = ( 0 , 0 , . . . 0 ) x=(0,0,...0) x=(0,0,...0)，理论最小值 f ( 0 , 0 , . . . 0 ) = 0 f(0,0,...0)=0 f(0,0,...0)=0。
- 初始化种群。N=100，粒子维数10，最大迭代次数200
参数设置如下表所示：
| 参数名称 | 参数大小 |
|---|---|
| 种群大小 N N N | 100 |
| 维度 D D D | 10 |
| 最大迭代次数 | 200 |
| 学习因子 c 1 = c 2 c_1=c_2 c1=c2 | 1.5 |
| 惯性权重 ω \omega ω | 0.8 |
| X m a x X_{max} Xmax | 20 |
| X m i n X_{min} Xmin | 20 |
| V m a x V_{max} Vmax | 10 |
- 利用上述推导描述的算法框架进行实现。
3. matlab代码实现
核心代码如下，我已经逐行注释，可以自行阅读代码流程。
%%%%%%%%%%%%%%%%%粒子群算法求函数极值%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; %清除所有变量
close all; %清图
clc; %清屏
N=100; %群体粒子个数
D=10; %粒子维数
T=200; %最大迭代次数
c1=1.5; %学习因子1
c2=1.5; %学习因子2
w=0.8; %惯性权重
Xmax=20; %位置最大值
Xmin=-20; %位置最小值
Vmax=10; %速度最大值
Vmin=-10; %速度最小值
%%%%%%%%%%%%%%%%初始化种群个体（限定位置和速度）%%%%%%%%%%%%%%%%
x=rand(N,D) * (Xmax-Xmin)+Xmin; % 初始化粒子位置
v=rand(N,D) * (Vmax-Vmin)+Vmin; % 初始化粒子速度
%%%%%%%%%%%%%%%%%%初始化个体最优位置和最优值%%%%%%%%%%%%%%%%%%%
p=x; % 初始化个体最优位置
pbest=ones(N,1); % 初始化个体最优值
for i=1:N
pbest(i)=func1(x(i,:)); % 计算个体适应度（最优）值
end
%%%%%%%%%%%%%%%%%%%初始化全局最优位置和最优值%%%%%%%%%%%%%%%%%%
g=ones(1,D); % 初始化全局最优位置
gbest=inf; % 将全局最好点设置为无穷大
for i=1:N % 更新全局最优解
if(pbest(i)<gbest) % 如果找到
g=p(i,:); % 全局最优位置
gbest=pbest(i); % 全局最优值
end
end
gb=ones(1,T); % 追踪 T 次全局最优解
%%%%%%%%%%%按照公式依次迭代直到满足精度或者迭代次数%%%%%%%%%%%%%
for i=1:T % 开始迭代
for j=1:N % 遍历所有个体
%%%%%%%%%%%%%%更新个体最优位置和最优值%%%%%%%%%%%%%%%%%
if (func1(x(j,:))<pbest(j))
p(j,:)=x(j,:); % 更新个体最优位置
pbest(j)=func1(x(j,:)); % 更新个体最优值
end
%%%%%%%%%%%%%%%%更新全局最优位置和最优值%%%%%%%%%%%%%%%
if(pbest(j)<gbest)
g=p(j,:); % 更新全局最优位置
gbest=pbest(j); % 更新全局最优值
end
%%%%%%%%%%%%%%%%%跟新位置和速度值%%%%%%%%%%%%%%%%%%%%%
v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))...
+c2*rand*(g-x(j,:)); % 更新速度
x(j,:)=x(j,:)+v(j,:); % 更新位置
%%%%%%%%%%%%%%%%%%%%边界条件处理%%%%%%%%%%%%%%%%%%%%%%
for ii=1:D % 遍历所有维度
if (v(j,ii)>Vmax) || (v(j,ii)< Vmin)
v(j,ii)=rand * (Vmax-Vmin)+Vmin; % 处理速度边界
end
if (x(j,ii)>Xmax) || (x(j,ii)< Xmin)
x(j,ii)=rand * (Xmax-Xmin)+Xmin; % 处理位置边界
end
end
end
%%%%%%%%%%%%%%%%%%%%记录历代全局最优值%%%%%%%%%%%%%%%%%%%%%
gb(i)=gbest;
end
disp('最优个体：')
disp(g); % 最优个体
disp('最优值：')
disp(gb(end)); % 最优值
figure
plot(gb)
xlabel('迭代次数');
ylabel('适应度值');
title('适应度进化曲线')
4. 效果展示
算法都收敛性是较好的，最优值变化曲线如下：
最优个体及最优值如下：
六、算法改进
1. 压缩因子粒子群算法
Clerc 等人提出利用 约束因子 控制系统行为的最终收敛，该方法可以有效搜索不同的区域，并且得到高质量的解。压缩因子的速度更新方法如下：
v i j ( t + 1 ) = λ ⋅ v i j ( t ) + c 1 r 1 ( t ) [ p i j ( t ) − x i j ( t ) ] + c 2 r 2 ( t ) [ p g j ( t ) − x i j ( t ) ] v_{ij}(t+1)=\lambda \cdot v_{ij}(t)+c_1r_1(t)[p_{ij}(t)-x_{ij}(t)]+c_2r_2(t)[p_{gj}(t)-x_{ij}(t)] vij(t+1)=λ⋅vij(t)+c1r1(t)[pij(t)−xij(t)]+c2r2(t)[pgj(t)−xij(t)]
式中， λ \lambda λ为压缩因子：
λ = 2 ∣ 2 − φ − ( φ 2 − 4 φ ) ∣ \lambda=\frac{2}{|2-\varphi-\sqrt{(\varphi^2-4\varphi)}|} λ=∣2−φ−(φ2−4φ)∣2
其中，
φ = c 1 + c 2 \varphi=c_1+c_2 φ=c1+c2
使用具有约束因子的粒子群算法具有更快的收敛速度。
2. 离散粒子群算法
主要针对速度，将连续的速度转化为离散的值，进而进行位置的更新，主要公式如下：
s ( v i , j ) = 1 1 + e − v i , j s(v_{i,j})=\frac{1}{1+e^{-v_{i,j}}} s(vi,j)=1+e−vi,j1
x i j = { 1 , r < s ( v i , j ) 0 , 其它 x_{ij}=\begin{cases} 1,& r< s(v_{i,j})\\ 0,& 其它 \end{cases} xij={1,0,r<s(vi,j)其它
其中， r r r是 U ( 0 , 1 ) U(0,1) U(0,1)分布中产生的随机数。1
the end……
粒子群算法到这里就要结束啦~~到此既是缘分，欢迎您的点赞、评论、收藏！关注我，不迷路，我们下期再见！！
😘😘😘 我是Cherries，一位计算机科班在校大学生，写博客用来记录自己平时的所思所想！
💞💞💞 内容繁杂，又才疏学浅，难免存在错误，欢迎各位大佬的批评指正！
👋👋👋 我们相互交流，共同进步！
注：本文由
非妃是公主发布于https://blog.csdn.net/myf_666，转载请务必标明原文链接：https://blog.csdn.net/myf_666/article/details/129405198文章来源：https://www.toymoban.com/news/detail-491467.html
-
包子阳,智能优化算法及其matlab实例.电子工业出版社. ↩︎文章来源地址https://www.toymoban.com/news/detail-491467.html
到了这里，关于粒子群算法（Particle Swarm Optimization（PSO）附简单案例及详细matlab源码）的文章就介绍完了。如果您还想了解更多内容，请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章，希望大家以后多多支持TOY模板网！
