粒子群算法

粒子群算法（particle swarm optimization，PSO）是计算智能领域中的一种生物启发式方法，属于群体智能优化算法的一种，常见的群体智能优化算法主要有如下几类：

除了上述几种常见的群体智能算法以外，还有一些并不是广泛应用的群体智能算法，比如萤火虫算法、布谷鸟算法、蝙蝠算法以及磷虾群算法等等。

而其中的粒子群优化算法（PSO）源于对鸟类捕食行为的研究，鸟类捕食时，找到食物最简单有限的策略就是搜寻当前距离食物最近的鸟的周围。

设想这样一个场景：一群鸟在随机的搜索食物。在这个区域里只有一块食物，所有的鸟都不知道食物在哪。但是它们知道自己当前的位置距离食物还有多远。那么找到食物的最优策略是什么？最简单有效的就是搜寻目前离食物最近的鸟的周围区域。

Step1:确定一个粒子的运动状态是利用位置和速度两个参数描述的，因此初始化的也是这两个参数；

Step2:每次搜寻的结果（函数值）即为粒子适应度，然后记录每个粒子的个体历史最优位置和群体的历史最优位置；

Step3:个体历史最优位置和群体的历史最优位置相当于产生了两个力，结合粒子本身的惯性共同影响粒子的运动状态，由此来更新粒子的位置和速度。

位置和速度的初始化即在位置和速度限制内随机生成一个N x d 的矩阵，而对于速度则不用考虑约束，一般直接在0~1内随机生成一个50x1的数据矩阵。

此处的位置约束也可以理解为位置限制，而速度限制是保证粒子步长不超限制的，一般设置速度限制为[-1,1]。

粒子群的另一个特点就是记录每个个体的历史最优和种群的历史最优，因此而二者对应的最优位置和最优值也需要初始化。其中每个个体的历史最优位置可以先初始化为当前位置,而种群的历史最优位置则可初始化为原点。对于最优值，如果求最大值则初始化为负无穷，相反地初始化为正无穷。

每次搜寻都需要将当前的适应度和最优解同历史的记录值进行对比，如果超过历史最优值，则更新个体和种群的历史最优位置和最优解。

速度和位置更新是粒子群算法的核心，其原理表达式和更新方式：

每次更新完速度和位置都需要考虑速度和位置的限制，需要将其限制在规定范围内，此处仅举出一个常规方法，即将超约束的数据约束到边界（当位置或者速度超出初始化限制时，将其拉回靠近的边界处）。当然，你不用担心他会停住不动，因为每个粒子还有惯性和其他两个参数的影响。

粒子群算法求平方和函数最小值，由于没有特意指定函数自变量量纲，不进行数据归一化。

从上面的例子我们可以看到应用PSO解决优化问题的过程中有两个重要的步骤: 问题解的编码和适应度函数 不需要像遗传算法一样是二进制编码(或者采用针对实数的遗传 *** 作.例如对于问题 f(x) = x1^2 + x2^2+x3^2 求解, 粒子可以直接编码为 (x1, x2, x3), 而适应度函数就是f(x). 接着我们就可以利用前面的过程去寻优.这个寻优过程是一个叠代过程, 中止条件一般为设置为达到最大循环数或者最小错误

PSO中并没有许多需要调节的参数,下面列出了这些参数以及经验设置

粒子数: 一般取 20 – 40. 其实对于大部分的问题10个粒子已经足够可以取得好的结果, 不过对于比较难的问题或者特定类别的问题, 粒子数可以取到100 或 200

粒子的长度: 这是由优化问题决定, 就是问题解的长度

粒子的范围: 由优化问题决定,每一维可以设定不同的范围

Vmax: 最大速度,决定粒子在一个循环中最大的移动距离,通常设定为粒子的范围宽度,例如上面的例子里,粒子 (x1, x2, x3) x1 属于 [-10, 10], 那么 Vmax 的大小就是 20

学习因子: c1 和 c2 通常等于 2. 不过在文献中也有其他的取值. 但是一般 c1 等于 c2 并且范围在0和4之间

中止条件: 最大循环数以及最小错误要求. 例如, 在上面的神经网络训练例子中, 最小错误可以设定为1个错误分类, 最大循环设定为2000, 这个中止条件由具体的问题确定.

全局PSO和局部PSO: 我们介绍了两种版本的粒子群优化算法: 全局版和局部版. 前者速度快不过有时会陷入局部最优. 后者收敛速度慢一点不过很难陷入局部最优. 在实际应用中, 可以先用全局PSO找到大致的结果,再用局部PSO进行搜索. 代码来自2008年数学建模东北赛区B题， #includestdaf(x).h#include<math.h>#include<time.h>#include<iostream>#include<fstream>usingnamespacestdintc1=2//加速因子intc2=2//加速因子doublew=1//惯性权重doubleWmax=1//最大惯性权重doubleWmin=0.6//最小惯性权重intKmax=110//迭代次数intGdsCnt//物资总数intconstDim=10//粒子维数intconstPNum=50//粒子个数intGBIndex=0//最优粒子索引doublea=0.6//适应度调整因子doubleb=0.5//适应度调整因子intXup[Dim]//粒子位置上界数组intXdown[Dim]=//粒子位置下界数组intValue[Dim]//初始急需度数组intVmax[Dim]//最大速度数组classPARTICLE//申明粒子节点voidCheck(PARTICLE&,int)//约束函数voidInput(ifstream&)//输入变量voidInitial()//初始化相关变量doubleGetFit(PARTICLE&)//计算适应度voidCalculateFit()//计算适应度voidBirdsFly()//粒子飞翔voidRun(ofstream&,int=2000)//运行函数classPARTICLE//微粒类{public:intX[Dim]//微粒的坐标数组intXBest[Dim]//微粒的最好位置数组intV[Dim]//粒子速度数组doubleFit//微粒适合度doubleFitBest//微粒最好位置适合度}PARTICLEParr[PNum]//粒子数组intmain()//主函数{ofstreamoutf(out.txt)ifstreaminf(data.txt)//关联输入文件inf>>GdsCnt//输入物资总数Input(inf)Initial()Run(outf,100)system(pause)return0}voidCheck(PARTICLE&p,intcount)//参数:p粒子对象,count物资数量{srand((unsigned)time(NULL))intsum=0for(inti=0i<Dimi++){if(p.X>Xup)p.X=Xupelseif(p.X<Xdown)p.X=Xdownif(p.V>Vmax)p.V=Vmaxelseif(p.V<0)p.V=0sum+=p.X}while(sum>count){p.X[rand()%Dim]--sum=0for(inti=0i<Dimi++){if(p.X>Xup)p.X=Xupelseif(p.X<Xdown)p.X=Xdownif(p.V>Vmax)p.V=Vmaxelseif(p.V<0)p.V=0sum+=p.X}}voidInput(ifstream&inf)//以inf为对象输入数据{for(inti=0i<Dimi++)inf>>Xupfor(inti=0i<Dimi++)inf>>Value}voidInitial()//初始化数据{GBIndex=0srand((unsigned)time(NULL))//初始化随机函数发生器for(inti=0i<Dimi++)Vmax=(int)((Xup-Xdown)*0.035)for(inti=0i{for(intj=0j<Dimj++){Parr.X[j]=(int)(rand()/(double)RAND_MAX*(Xup[j]-Xdown[j])-Xdown[j]+0.5)Parr.XBest[j]=Parr.X[j]Parr.V[j]=(int)(rand()/(double)RAND_MAX*(Vmax[j]-Vmax[j]/2))}Parr.Fit=GetFit(Parr)Parr.FitBest=Parr.Fitif(Parr.Fit>Parr[GBIndex].Fit)GBIndex=i}}doubleGetFit(PARTICLE&p)//计算对象适应度{doublesum=0for(inti=0i<Dimi++)for(intj=1j<=p.Xj++)sum+=(1-(j-1)*a/(Xup-b))*Valuereturnsum}voidCalculateFit()//计算数组内各粒子的适应度{for(inti=0i{Parr.Fit=GetFit(Parr)}}voidBirdsFly()//粒子飞行寻找最优解{srand((unsigned)time(NULL))staticintk=10w=Wmax-k*(Wmax-Wmin)/Kmaxk++for(inti=0i{for(intj=0j<Dimj++){Parr.V[j]=(int)(w*Parr.V[j])Parr.V[j]+=(int)(c1*rand()/(double)RAND_MAX*(Parr.XBest[j]-Parr.X[j])Parr.V[j]+=c2*rand()/(double)RAND_MAX*(Parr[GBIndex].XBest[j]-Parr.X[j]))}}Check(Parr,GdsCnt)for(intj=0j<Dimj++){Parr.X[j]+=Parr.V[j]Check(Parr,GdsCnt)}CalculateFit()for(inti=0i{if(Parr.Fit>=Parr.FitBest){Parr.FitBest=Parr.Fitfor(intj=0j<Dimj++)Parr.XBest[j]=Parr.X[j]}}GBIndex=0for(inti=0i{if(Parr.FitBest>Parr[GBIndex].FitBest&&i!=GBIndex)GBIndex=i}}voidRun(ofstream&outf,intnum)//令粒子以规定次数num飞行{for(inti=0i<numi++){BirdsFly()outf<<(i+1)<<ends<for(intj=0j<Dimj++)outf<outf<<endl}cout<<Done!<<endl}

  本系列文章主要针对粒子群算法进行介绍和运用，并给出粒子群算法的经典案例，从而进一步加深对粒子群算法的了解与运用（预计在一周内完成本系列文章）。主要包括四个部分：

  粒子群算法也称粒子群优化算法（Particle Swarm Optimization, PSO），属于群体智能优化算法，是近年来发展起来的一种新的进化算法（Evolutionary Algorithm, EA）。 群体智能优化算法主要模拟了昆虫、兽群、鸟群和鱼群的群集行为，这些群体按照一种合作的方式寻找食物，群体中的每个成员通过学习它自身的经验和其他成员的经验来不断地改变搜索的方向。 群体智能优化算法的突出特点就是利用了种群的群体智慧进行协同搜索，从而在解空间内找到最优解。

  PSO 算法和模拟退火算法相比，也是 从随机解出发，通过迭代寻找最优解 。它是通过适应度来评价解的品质，但比遗传算法规则更为简单，没有遗传算法的“交叉”和“变异”，它通过追随当前搜索到的最大适应度来寻找全局最优。这种算法以其 容易实现、精度高、收敛快 等优点引起了学术界的重视，并在解决实际问题中展示了其优越性。

  在粒子群算法中，每个优化问题的解被看作搜索空间的一只鸟，即“粒子”。算法开始时首先生成初始解，即在可行解空间中随机初始化 粒子组成的种群 ，其中每个粒子所处的位置 ，都表示问题的一个解，并依据目标函数计算搜索新解。在每次迭代时，粒子将跟踪两个“极值”来更新自己， 一个是粒子本身搜索到的最好解 ，另一个是整个种群目前搜索到的最优解 。 此外每个粒子都有一个速度 ，当两个最优解都找到后，每个粒子根据如下迭代式更新：

  其中参数 称为是 PSO 的 惯性权重(inertia weight) ，它的取值介于[0,1]区间；参数 和 称为是 学习因子(learn factor) ；而 和 为介于[0,1]之间的随机概率值。

  实践证明没有绝对最优的参数，针对不同的问题选取合适的参数才能获得更好的收敛速度和鲁棒性，一般情况下 , 取 1.4961 ，而 采用 自适应的取值方法 ，即一开始令 ， 使得 PSO 全局优化能力较强 ；随着迭代的深入，递减至 ， 从而使得PSO具有较强的局部优化能力 。

  参数 之所以被称之为惯性权重，是因为 实际 反映了粒子过去的运动状态对当前行为的影响，就像是我们物理中提到的惯性。 如果 ，从前的运动状态很少能影响当前的行为，粒子的速度会很快的改变；相反， 较大，虽然会有很大的搜索空间，但是粒子很难改变其运动方向，很难向较优位置收敛，由于算法速度的因素，在实际运用中很少这样设置。也就是说， 较高的 设置促进全局搜索，较低的 设置促进快速的局部搜索。

---