1、粒子群算法介绍优化问题是工业设计中经常遇到的问题,许多问题最后都可以归结为优化问题. 为了解决各种各样的优化问题,人们提出了许多优化算法,比较著名的有爬山法、遗传算法等.优化问题有两个主要问题:一是要求寻找全局最小点,二是要求有较高的收敛速度. 爬山法精度较高,但是易于陷入局部极小. 遗传算法属于进化算法( Evolutionary Algorithms) 的一种,它通过模仿自然界的选择与遗传的机理来寻找最优解. 遗传算法有三个基本算子:选择、交叉和变异. 但是遗传算法的编程实现比较复杂,首先需要对问题进行编码,找到最优解之后还需要对问题进行解码,另外三个算子的实现也有许多参数,如交叉率和变异
2、率,并且这些参数的选择严重影响解的品质,而目前这些参数的选择大部分是依靠经验.1995 年Eberhart 博士和 kennedy 博士提出了一种新的算法;粒子群优化(Particle Swarm Optimization -PSO) 算法 . 这种算法以其实现容易、精度高、收敛快等优点引起了学术界的重视,并且在解决实际问题中展示了其优越性。粒子群优化(Particle Swarm Optimization - PSO) 算法是近年来发展起来的一种新的进化算法( Evolutionary Algorithm - EA) .PSO 算法属于进化算法的一种,和遗传算法相似,它也是从随机解出发,通过
3、迭代寻找最优解,它也是通过适应度来评价解的品质. 但是它比遗传算法规则更为简单,它没有遗传算法的“交叉”(Crossover) 和“变异”(Mutation) 操作. 它通过追随当前搜索到的最优值来寻找全局最优 。1. 引言粒子群优化算法(PSO)是一种进化计算技术(evolutionary computation),由Eberhart 博士和 kennedy 博士提出。源于对鸟群捕食的行为研究 。PSO 同遗传算法类似,是一种基于迭代的优化算法。系统初始化为一组随机解,通过迭代搜寻最优值。但是它没有遗传算法用的交叉(crossover)以及变异(mutation),而是粒子在解空间追随最优的
4、粒子进行搜索。同遗传算法比较,PSO 的优势在于简单容易实现并且没有许多参数需要调整。目前已广泛应用于函数优化,神经网络训练,模糊系统控制以及其他遗传算法的应用领域2. 背景: 人工生命“人工生命“是来研究具有某些生命基本特征的人工系统。人工生命包括两方面的内容 :1. 研究如何利用计算技术研究生物现象2. 研究如何利用生物技术研究计算问题我们现在关注的是第二部分的内容. 现在已经有很多源于生物现象的计算技巧. 例如, 人工神经网络是简化的大脑模型. 遗传算法是模拟基因进化过程的.现在我们讨论另一种生物系统- 社会系统. 更确切的是, 在由简单个体组成的群落与环境以及个体之间的互动行为. 也可
5、称做“群智能“(swarm intelligence). 这些模拟系统利用局部信息从而可能产生不可预测的群体行为例如 floys 和 boids, 他们都用来模拟鱼群和鸟群的运动规律, 主要用于计算机视觉和计算机辅助设计.在计算智能(computational intelligence)领域有两种基于群智能的算法. 蚁群算法(ant colony optimization)和粒子群算法(particle swarm optimization). 前者是对蚂蚁群落食物采集过程的模拟. 已经成功运用在很多离散优化问题上.粒子群优化算法(PSO) 也是起源对简单社会系统的模拟. 最初设想是模拟鸟群觅
6、食的过程. 但后来发现 PSO 是一种很好的优化工具.3. 算法介绍如前所述,PSO 模拟鸟群的捕食行为。设想这样一个场景:一群鸟在随机搜索食物。在这个区域里只有一块食物。所有的鸟都不知道食物在那里。但是他们知道当前的位置离食物还有多远。那么找到食物的最优策略是什么呢。最简单有效的就是搜寻目前离食物最近的鸟的周围区域。PSO 从这种模型中得到启示并用于解决优化问题。PSO 中,每个优化问题的解都是搜索空间中的一只鸟。我们称之为“粒子”。所有的粒子都有一个由被优化的函数决定的适应值(fitness value),每个粒子还有一个速度决定他们飞翔的方向和距离。然后粒子们就追随当前的最优粒子在解空间
7、中搜索。PSO 初始化为一群随机粒子(随机解)。然后通过迭代找到最优解。在每一次迭代中,粒子通过跟踪两个“极值“来更新自己。第一个就是粒子本身所找到的最优解,这个解叫做个体极值 pBest。另一个极值是整个种群目前找到的最优解,这个极值是全局极值 gBest。另外也可以不用整个种群而只是用其中一部分作为粒子的邻居,那么在所有邻居中的极值就是局部极值。在找到这两个最优值时,粒子根据如下的公式来更新自己的速度和新的位置:v = w * v + c1 * rand() * (pbest - present) + c2 * rand() * (gbest - present) (a)present =
8、 persent + v (b)v 是粒子的速度, w 是惯性权重,persent 是当前粒子的位置. pbest and gbest 如前定义 rand () 是介于(0, 1)之间的随机数. c1, c2 是学习因子. 通常 c1 = c2 = 2.程序的伪代码如下For each particle_Initialize particleENDDo_For each particle_Calculate fitness value_If the fitness value is better than the best fitness value (pBest) in history_se
9、t current value as the new pBest_End_Choose the particle with the best fitness value of all the particles as the gBest_For each particle_Calculate particle velocity according equation (a)_Update particle position according equation (b)_EndWhile maximum iterations or minimum error criteria is not att
10、ained在每一维粒子的速度都会被限制在一个最大速度 Vmax,如果某一维更新后的速度超过用户设定的 Vmax,那么这一维的速度就被限定为 Vmax4. 遗传算法和 PSO 的比较大多数演化计算技术都是用同样的过程 :1. 种群随机初始化2. 对种群内的每一个个体计算适应值(fitness value).适应值与最优解的距离直接有关3. 种群根据适应值进行复制4. 如果终止条件满足的话,就停止,否则转步骤 2从以上步骤,我们可以看到 PSO 和 GA 有很多共同之处。两者都随机初始化种群,而且都使用适应值来评价系统,而且都根据适应值来进行一定的随机搜索。两个系统都不是保证一定找到最优解。但是,
11、PSO 没有遗传操作如交叉(crossover)和变异(mutation). 而是根据自己的速度来决定搜索。粒子还有一个重要的特点,就是有记忆。与遗传算法比较, PSO 的信息共享机制是很不同的. 在遗传算法中,染色体(chromosomes) 互相共享信息,所以整个种群的移动是比较均匀的向最优区域移动. 在PSO 中, 只有 gBest (or lBest) 给出信息给其他的粒子,这是单向的信息流动. 整个搜索更新过程是跟随当前最优解的过程. 与遗传算法比较, 在大多数的情况下,所有的粒子可能更快的收敛于最优解5. 人工神经网络 和 PSO人工神经网络(ANN)是模拟大脑分析过程的简单数学模
12、型,反向转播算法是最流行的神经网络训练算法。进来也有很多研究开始利用演化计算(evolutionary computation)技术来研究人工神经网络的各个方面。演化计算可以用来研究神经网络的三个方面:网络连接权重,网络结构(网络拓扑结构,传递函数),网络学习算法。不过大多数这方面的工作都集中在网络连接权重,和网络拓扑结构上。在 GA 中,网络权重和/或拓扑结构一般编码为染色体(Chromosome),适应函数(fitness function)的选择一般根据研究目的确定。例如在分类问题中,错误分类的比率可以用来作为适应值演化计算的优势在于可以处理一些传统方法不能处理的例子例如不可导的节点传递
13、函数或者没有梯度信息存在。但是缺点在于:在某些问题上性能并不是特别好。2. 网络权重的编码而且遗传算子的选择有时比较麻烦最近已经有一些利用 PSO 来代替反向传播算法来训练神经网络的论文。研究表明 PSO 是一种很有潜力的神经网络算法。PSO 速度比较快而且可以得到比较好的结果。而且还没有遗传算法碰到的问题这里用一个简单的例子说明 PSO 训练神经网络的过程。这个例子使用分类问题的基准函数(Benchmark function)IRIS 数据集。(Iris 是一种鸢尾属植物) 在数据记录中,每组数据包含 Iris 花的四种属性:萼片长度,萼片宽度,花瓣长度,和花瓣宽度,三种不同的花各有 50
14、组数据. 这样总共有 150 组数据或模式。我们用 3 层的神经网络来做分类。现在有四个输入和三个输出。所以神经网络的输入层有 4 个节点,输出层有 3 个节点我们也可以动态调节隐含层节点的数目,不过这里我们假定隐含层有 6 个节点。我们也可以训练神经网络中其他的参数。不过这里我们只是来确定网络权重。粒子就表示神经网络的一组权重,应该是 4*6+6*3=42 个参数。权重的范围设定为-100,100 (这只是一个例子,在实际情况中可能需要试验调整).在完成编码以后,我们需要确定适应函数。对于分类问题,我们把所有的数据送入神经网络,网络的权重有粒子的参数决定。然后记录所有的错误分类的数目作为那个
15、粒子的适应值。现在我们就利用 PSO 来训练神经网络来获得尽可能低的错误分类数目。PSO 本身并没有很多的参数需要调整。所以在实验中只需要调整隐含层的节点数目和权重的范围以取得较好的分类效果。6. PSO 的参数设置从上面的例子我们可以看到应用 PSO 解决优化问题的过程中有两个重要的步骤: 问题解的编码和适应度函数PSO 的一个优势就是采用实数编码, 不需要像遗传算法一样是二进制编码(或者采用针对实数的遗传操作.例如对于问题 f(x) = x12 + x22+x32 求解, 粒子可以直接编码为 (x1, x2, x3), 而适应度函数就是 f(x). 接着我们就可以利用前面的过程去寻优.这个
16、寻优过程是一个叠代过程, 中止条件一般为设置为达到最大循环数或者最小错误PSO 中并没有许多需要调节的参数,下面列出了这些参数以及经验设置粒子数: 一般取 20 40. 其实对于大部分的问题 10 个粒子已经足够可以取得好的结果, 不过对于比较难的问题或者特定类别的问题, 粒子数可以取到 100 或 200粒子的长度: 这是由优化问题决定, 就是问题解的长度粒子的范围: 由优化问题决定,每一维可是设定不同的范围Vmax: 最大速度,决定粒子在一个循环中最大的移动距离,通常设定为粒子的范围宽度,例如上面的例子里,粒子 (x1, x2, x3) x1 属于 -10, 10, 那么 Vmax 的大小
17、就是 20学习因子: c1 和 c2 通常等于 2. 不过在文献中也有其他的取值. 但是一般 c1 等于 c2 并且范围在 0 和 4 之间中止条件: 最大循环数以及最小错误要求. 例如, 在上面的神经网络训练例子中, 最小错误可以设定为 1 个错误分类, 最大循环设定为 2000, 这个中止条件由具体的问题确定.全局 PSO 和局部 PSO: 我们介绍了两种版本的粒子群优化算法: 全局版和局部版. 前者速度快不过有时会陷入局部最优. 后者收敛速度慢一点不过很难陷入局部最优. 在实际应用中, 可以先用全局 PSO 找到大致的结果,再有局部 PSO 进行搜索.另外的一个参数是惯性权重, Shi
18、和 Eberhart 指出(A modified particle swarm optimizer,1998):当 Vmax 很小时(对 schaffer 的 f6 函数,Vmax=3),使用权重 w=0.8 较好.如果没有 Vmax 的信息,使用 0.8 作为权重也是一种很好的选择.另外,对于使用时变的权重,结果不清楚,但是预计结果应比较好.附上一个 C+实现的 C+代码:代码来自 2008 年数学建模东北赛区 B 题,原题描述http:/ “stdafx.h“#include #include #include #include using namespace std;int c1=2;
19、/加速因子int c2=2; /加速因子double w=1; /惯性权重double Wmax=1; /最大惯性权重double Wmin=0.6; /最小惯性权重int Kmax=110; /迭代次数int GdsCnt; /物资总数int const Dim=10; /粒子维数int const PNum=50; /粒子个数int GBIndex=0; /最优粒子索引double a=0.6; /适应度调整因子double b=0.5; /适应度调整因子int XupDim; /粒子位置上界数组int XdownDim=; /粒子位置下界数组int ValueDim; /初始急需度数组i
20、nt VmaxDim; /最大速度数组class PARTICLE; /申明粒子节点void Check(PARTICLE /约束函数void Input(ifstream /输入变量void Initial(); /初始化相关变量double GetFit(PARTICLE /计算适应度void CalculateFit(); /计算适应度void BirdsFly(); /粒子飞翔void Run(ofstream /运行函数/微粒类class PARTICLEpublic:int XDim; /微粒的坐标数组int XBestDim; /微粒的最好位置数组int VDim; /粒子速度数
21、组double Fit; /微粒适合度double FitBest; /微粒最好位置适合度;PARTICLE ParrPNum; /粒子数组int main() /主函数ofstream outf(“out.txt“);ifstream inf(“data.txt“); /关联输入文件infGdsCnt; /输入物资总数Input(inf);Initial();Run(outf,100);system(“pause“);return 0;void Check(PARTICLEint sum=0;for (int i=0;iXup)p.X=Xup;else if (p.XVmax)p.V=Vma
22、x;else if (p.Vcount)p.Xrand()%Dim-;sum=0;for (int i=0;iXup)p.X=Xup;else if (p.XVmax)p.V=Vmax;else if (p.VXup;for (int i=0;iValue;void Initial() /初始化数据GBIndex=0;srand(unsigned)time(NULL);/初始化随机函数发生器for (int i=0;iParrGBIndex.Fit)GBIndex=i;double GetFit(PARTICLEfor (int i=0;i=Parr.FitBest)Parr.FitBest=
23、Parr.Fit;for (int j=0;jParrGBIndex.FitBestvoid Run(ofstreami#include #include using namespace std;double c1=0.7; double c2=0.2; double w=1.0; double Wmax=1.0; double Wmin=0.6; int Kmax=50; int const Dim=7; int const PNum=4; int FB=200; int FC=5; int GBIndex=0; class PARTICLE; void Initial(); void Up
24、date(int *x,int *v); void GetDifferent(int *p1,int *p2,int *v); int GetFit(PARTICLE void CalculateFit(); void BirdsFly(); void Run(int num); double GetRandm(); int A77= 0, 2,1000, 4,1000,1000,1000, 2, 0, 3,1000,1000,1000,1000,1000, 3, 0,1000, 3, 6,1000, 4,1000,1000, 0, 5, 2,1000,1000,1000,1000, 5, 0
25、, 1, 2,1000,1000, 6, 2, 1, 0, 4,1000,1000,1000,1000, 2, 4, 0;int B77= 0, 2,1000, 4,1000,1000,1000, 2, 0, 3,1000,1000,1000,1000,1000, 3, 0,1000, 3, 6,1000, 4,1000,1000, 0, 5, 2,1000,1000,1000,1000, 5, 0, 1, 2,1000,1000, 6, 2, 1, 0, 4,1000,1000,1000,1000, 2, 4, 0;int C77= 1000, 20,1000, 10,1000,1000,1
26、000, 20,1000, 5,1000,1000,1000,1000, 1, 7,1000,1000, 7, 36,1000, 40,1000, 2,1000, 6, 52,1000, 3,1000,1000, 9,1000, 11, 12, 20,1000, 6, 6, 8,1000, 14,1000,1000,1000,1000, 8, 24,1000;int First7=0,1,2,3,4,5,6;int Last7=0,1,2,3,4,5,6; int hp17=0,0,0,0,0,0,0; int hp27=0,0,0,0,0,0,0; /class PARTICLEpublic
27、:int XDim; int XBestDim; int VDim; double Fit; double FitBest; ;PARTICLE PtPNum; int Inp47= 0,6,2,4,1,3,5,0,1,2,4,3,5,6,0,3,6,5,4,2,1,0,3,4,6,5,2,1;double GetRandm()return (double)(rand()/(double)RAND_MAX);/void Initial() GBIndex=0;for(int i=0;iFB)fb=1000;elsefb=0;sum+=(fb+fc);return sum;void Calcul
28、ateFit()for (int i=0;iGetRandm()Update(Last,Pti.V);for(int j=0;jGetRandm()GetDifferent(Pti.XBest,Pti.X,Pti.V);Update(Last,Pti.V);for(j=0;jGetRandm()GetDifferent(PtGBIndex.XBest,Pti.X,Pti.V);Update(Last,Pti.V);for(j=0;jDim;j+) Pti.Vj=0;GetDifferent(First,Last,Pti.V); for(j=0;jDim;j+) Lastj=j;Update(P
29、ti.X,Pti.V);CalculateFit();for (i=0;iPNum;i+)if (Pti.FitPti.FitBest)Pti.FitBest=Pti.Fit;for (int j=0;jDim;j+)Pti.XBestj=Pti.Xj;GBIndex=0;for (i=0;iPNum;i+)if (Pti.FitBest=PtGBIndex.FitBest cout“:“PtGBIndex.FitBestendl;void GetDifferent(int *p1,int *p2,int *p)for(int i=0;iDim;i+)hp1i=p1i;hp2i=p2i;for
30、(i=0;iDim;i+)int sw;if(hp2i!=hp1i)for(int j=i+1;jDim;j+)if(hp2j=hp1i)sw=hp2i;hp2i=hp2j;hp2j=sw;pi=j;for(i=0;iDim;i+)hp1i=0;hp2i=0;void Update(int* x,int* v)int s;for(int i=0; iDim; i+)if(vi!=0)s=xi;xi=xvi;xvi=s;void Run(int num)for(int i=0;inum;i+)BirdsFly();for(int j=0;jDim;j+)coutPtGBIndex.XBestjends;if(PtGBIndex.XBestj=6) return;void main() Initial();Run(Kmax);