粒子群优化算法 PSO

1、基本内容

1.1 概念

在动物的群体行为中，科学家们很早就发现了自然界的鸟群、兽群、鱼群等在其迁徙、捕食过程中，往往表现出高度的组织性和规律性。这些现象受到了高度的重视和广泛的关注，吸引着大批生物学家、动物学家、计算机科学家、行为学家和社会心理学家等的深入研究。

粒子群优化算法(Particle Swarm optimization,PSO)又翻译为粒子群算法、微粒群算法或微粒群优化算法，是通过模拟鸟群觅食行为而发展起来的一种基于群体协作的随机搜索算法。 通常认为粒子群优化算法是群集智能 (Swarm intelligence, SI) 的一种，它可以被纳入多主体优化系统(Multiagent Optimization System, MAOS)。

粒子群优化算法是由Eberhart博士和Kennedy博士于1995年提出的一种全局搜索算法。

1.2 基本原理

事实上，捕食的鸟群都是通过各自的探索与群体的合作最终发现食物所在的位置的。

• 一群分散的鸟在随机地飞行觅食，它们不知道食物所在的具体位置。
• 但各个小鸟会在飞行的过程中不断地记录和更新它曾经到达的离食物最近位置，同时通过信息交流的方式比较大家所找到的最好位置，得到一个当前整个群体已经找到的最佳位置。
• 各个小鸟在飞行过程中结合自身的经验和整个群体的经验，调整自己的飞行速度和所在位置，不断地寻找更加接近食物的位置，最终使得群体聚集到食物位置。

因此，在粒子群优化算法中，

• 每只鸟抽象为一个无质量、无体积的 “ 粒子 ”
• 每个粒子有一个适应度函数，以模拟每只鸟与食物的丰盛程度
• 每个粒子有一个速度决定它的飞行方向和距离，初始值可以随机确定
• 每一次单位时间的飞行后，所有粒子分享信息，下一步将飞向自身最佳位置和全局最优位置的加权中心

1.3 基本步骤

初始化一群随机粒子，通过迭代找到最优。每次迭代中，粒子通过跟踪 “ 个体极值(pbest) ” 和 “ 全局极值(gbest) ” 来更新自己的速度和位置。

假设在D维搜索空间中，有m个粒子，各粒子速度和位置的更新如下：

（1）其中第i个粒子的位置为矢量：

（2）其飞行速度也是一个矢量， 记为：

（3）第i个粒子搜索到的最优位置为：

（4）整个粒子群搜索到的最优位置为

（5）第i个粒子的位置和速度更新为

2、代码分析

2.1 Rastrigrin()函数

数学表达式：
代码：

function y = Rastrigin(x)
% Rastrigin函数
% 输入x,给出相应的y值,在x = ( 0 , 0 ,…, 0 )处有全局极小点0.
[row,col] = size(x);
if  row > 1 
    error( ' 输入的参数错误 ' );
end
y =sum(x.^2-10*cos(2*pi*x)+10);
%y =-y;

2.2 Schaffer()函数

数学表达式：

代码：

function y=Schaffer(x)
[row,col]=size(x);
if row>1
    error('输入的参数错误');
end
y1=x(1,1);
y2=x(1,2);
temp=y1^2+y2^2;
y=0.5+(sin(sqrt(temp))^2-0.5)/(1+0.001*temp)^2;

2.3 Griewank()函数

数学表达式：

代码：

function y=Griewank(x)
%Griewan函数
%输入x,给出相应的y值,在x=(0,0,…,0)处有全局极小点0.
[row,col]=size(x);
if row>1
    error('输入的参数错误');
end
y1=1/4000*sum(x.^2);
y2=1;
for h=1:col
    y2=y2*cos(x(h)/sqrt(h));
end
y=y1-y2+1;
%y=-y;

2.4 fun()函数

作用：根据label的值，选择不同的算法计算粒子适应度值。
代码：

function y = fun(x,label)
%函数用于计算粒子适应度值
%x           input           输入粒子 
%y           output          粒子适应度值 
if label==1
    y=Rastrigin(x);
elseif label==2
    y=Schaffer(x);
else
    y= Griewank(x);
end

2.5 fun1()函数

作用：通过函数进行粒子适应度值计算。
代码：

function y = fun1(x)
%函数用于计算粒子适应度值
%x           input           输入粒子 
%y           output          粒子适应度值 
y=-20*exp(-0.2*sqrt((x(1)^2+x(2)^2)/2))-exp((cos(2*pi*x(1))+cos(2*pi*x(2)))/2)+20+exp(1);
%y=x(1)^2-10*cos(2*pi*x(1))+10+x(2)^2-10*cos(2*pi*x(2))+10;

2.6 Drawfunc()函数

作用：根据label的值绘制不同的三维曲面网格图
代码：

function Drawfunc(label)
%该函数用于绘制图像
x=-5:0.05:5;%1行201列的向量

if label==1
    y = x;
    [X,Y] = meshgrid(x,y);
 %meshgrid函数生成的X，Y是大小相等的矩阵，x，y是两个网格矢量，x，y都是行向量。
 %X：通过将x复制length(y)行（严格意义上是length(y)-1行）得到
 %Y：首先对y进行转置得到y'，将y'复制（length(x)-1）次得到
    [row,col] = size(X);
    for  l = 1 :col
         for  h = 1 :row
            z(h,l) = Rastrigrin([X(h,l),Y(h,l)]);
        end
    end
    surf(X,Y,z);%surf()用于绘制比较光滑的三维曲面网格图
    shading interp%在flat的基础上进行色彩的插值处理，使色彩平滑过渡
    xlabel('x1-axis'),ylabel('x2-axis'),zlabel('f-axis'); 
    title('Rastrigrin mesh'); 
end

if label==2
    y = x;
    [X,Y] = meshgrid(x,y);
    [row,col] = size(X);
    for  l = 1 :col
         for  h = 1 :row
            z(h,l) = Schaffer([X(h,l),Y(h,l)]);
        end
    end
    surf(X,Y,z);