天天看點

粒子群算法PSO粒子群算法PSO

粒子群算法PSO

簡介

粒子群優化算法(PSO)是一種進化計算技術(evolutionary computation),1995 年由Eberhart 博士和kennedy 博士提出,源于對鳥群捕食的行為研究 。該算法最初是受到飛鳥叢集活動的規律性啟發,進而利用群體智能建立的一個簡化模型。粒子群算法在對動物叢集活動行為觀察基礎上,利用群體中的個體對資訊的共享使整個群體的運動在問題求解空間中産生從無序到有序的演化過程,進而獲得最優解。

基本思想

粒子群算法是模拟群體智能所建立起來的一種優化算法,粒子群算法可以用鳥類在一個空間内随機覓食為例,所有的鳥都不知道食物具體在哪裡,但是他們知道大概距離多遠,最簡單有效的方法就是搜尋目前離食物最近的鳥的周圍區域。

對粒子群優化算法操作的一個簡單解釋如下:每一個粒子都代表了目前優化任務的一個可能解決方案。在每次疊代過程中,每個粒子都會朝着自己的最優解的方向加速,也會朝着種群中任何粒子迄今為止發現的全局最佳位置的方向加速。這意味着,如果一個粒子發現了一個更好的解,所有其他粒子都會靠近它,在這個過程中不斷地搜尋最優解。可以總結出粒子群算法地三條簡單規則:(1)飛離最近的個體,以避免碰撞;(2)飛向目标;(3)飛向群體的中心。

假設存在一個次元為S的搜尋空間,由 m ​ m​ m​ 個粒子組成粒子種群,其中第 i ​ i​ i​ 個粒子用一個 S ​ S​ S​ 維的向量表示,具體為 X i = ( x i 1 , x i 2 , … , x i S ) ​ Xi=(x_{i1}, x_{i2},…,x_{iS})​ Xi=(xi1​,xi2​,…,xiS​)​。将 X i ​ X_i​ Xi​​代入目标函數就可以算出其适應度值,适應度值的大小對應着粒子位置即可行解的好壞。另外,第 i ​ i​ i​ 個粒子移動的速度是 S ​ S​ S​ 維向量,記為 V i = ( v i 1 , v i 2 , … , v i S ) ​ Vi=(v_{i1}, v_{i2} ,…, v_{iS})​ Vi=(vi1​,vi2​,…,viS​)​。第 i ​ i​ i​個粒子搜尋到的最優位置為 P i = ( p i 1 , p i 2 , … , p i S ) ​ P_i=(p_{i1},p_{i2},…,p_{iS})​ Pi​=(pi1​,pi2​,…,piS​)​,整個粒子群搜尋到的最優位置為 P g = ( p g 1 , p g 2 , … , p g S ) ​ P_g=(p_{g1}, p_{g2}, …, p_{gS})​ Pg​=(pg1​,pg2​,…,pgS​)​。

粒子位置和速度的更新公式如下:

KaTeX parse error: No such environment: eqnarray at position 8: \begin{̲e̲q̲n̲a̲r̲r̲a̲y̲}̲ V_i&=&V_i+c_1*…

i = 1 , 2 , … , m i=1,2,…,m i=1,2,…,m, m m m是該群體中粒子的總數; V i V_i Vi​ 是粒子的速度; r a n d ( ) \mathbb{rand}() rand()是介于 ( 0 , 1 ) (0,1) (0,1)之間的随機數; X i X_i Xi​ 是粒子的目前位置。 c 1 c_1 c1​和 c 2 c_2 c2​是學習因子, 在每一維,粒子都有一個最大限制速度 V m a x V_{max} Vmax​,如果 某一維的速度超過設定的 V m a x V_{max} Vmax​ ,那麼這一維的速度 就被限定為 V m a x ​ V_{max}​ Vmax​​ 。以上面兩個公式為基礎,形成了後來PSO 的标準形式。

一般在實際應用中,會對速度更新公式進行修正,即加入慣性權重 w w w:

(3) V i = w ∗ V i + c 1 ∗ r a n d ( ) ∗ ( P i − X i ) + c 2 ∗ r a n d ( ) ∗ ( P g − X i ) V_i=w*V_i+c_1*\mathbb{rand}()*(P_i-X_i)+c_2*\mathbb{rand}()*(P_g-X_i)\tag{3} Vi​=w∗Vi​+c1​∗rand()∗(Pi​−Xi​)+c2​∗rand()∗(Pg​−Xi​)(3)

算法流程

标準PSO算法的流程:

Step1:初始化一群微粒(群體規模為m),包括随機位置和速度;

Step2:評價每個微粒的适應度;

Step3:對每個微粒,将其适應值與其經過的最好位置 P i P_i Pi​作比較,如果較好,則将其作為該粒子目前的最好位置 P i P_i Pi​;

Step4:對每個微粒,将其适應值與全局最好位置 P g P_g Pg​作比較,如果較好,則将其作為目前全局最好位置 P g P_g Pg​;

Step5:根據(2)、(3)式調整微粒速度和位置;

Step6:未達到結束條件則轉Step2。

疊代終止條件根據具體問題一般選為最大疊代次數 k k k或(和)微粒群迄今為止搜尋到的最優位置滿足預定最小适應門檻值 ϵ \epsilon ϵ。

算法實作(MATLAB)

一般的PSO實作1:

function [xm,fv,pbest]=stdPSO(fun,sizepop,c1,c2,w,maxg,D)
% Standard PSO
% Input
%fitness: 待優化的目标函數
%sizepop: 粒子數目
%c1: 學習因子1
%c2: 學習因子2
%w: 慣性權重
%maxg: 最大疊代次數
%D: 自變量的次元
%Output 
%xm: 目标函數取最小值時的自變量值
%fv: 目标函數的最小值
%% 參數初始化
% 初始速度和種群上下邊界值
Vmax=1;
Vmin=-1;
popmax=5;
popmin=-5;
%% 産生初始粒子和速度
for i=1:sizepop
    %随機産生一個種群
    pop(i,:)=5*rands(1,D); %初始種群
    V(i,:)=rands(1,D); %初始化速度
    %計算适應度
    fitness(i)=fun(pop(i,:)); %染色體的适應度
end
%找最好的染色體
[bestfitness bestindex]=min(fitness);
zbest=pop(bestindex,:); %全局最佳
gbest=pop; %個體最佳
fitnessgbest=fitness; %個體最佳适應度值
fitnesszbest=bestfitness; %全局最佳适應度值

%% 疊代尋優
for i=1:maxg   %疊代次數
    for j=1:sizepop
        %速度更新
        V(j,:) = w*V(j,:) + c1*rand*(gbest(j,:) - pop(j,:)) + c2*rand*(zbest- pop(j,:));
        V(j,find(V(j,:)>Vmax))=Vmax;
        V(j,find(V(j,:)<Vmin))=Vmin;
        %種群更新
        pop(j,:)=pop(j,:)+V(j,:);
        pop(j,find(pop(j,:)>popmax))=popmax;
        pop(j,find(pop(j,:)<popmin))=popmin;
        %自适應變異
        if rand>0.8
            k=ceil(2*rand);
            pop(j,k)=rand;
        end
        %适應度值
        fitness(j)=fun(pop(j,:));
        %個體最優更新
        if fitness(j) < fitnessgbest(j)
            gbest(j,:) = pop(j,:);
            fitnessgbest(j) = fitness(j);
        end
        %群體最優更新
        if fitness(j) < fitnesszbest
            zbest = pop(j,:);
            fitnesszbest = fitness(j);
        end
    end
    yy(i)=fitnesszbest;
end
% 結果輸出
xm=zbest; %最佳個體值
fv=fitnesszbest;
pbest=yy;
           

測試函數Rastrigin

f ( x ) = A ∗ n + ∑ i = 1 n [ x i 2 − A ∗ cos ⁡ ( 2 π x i ) ] f(x)=A*n+\sum_{i=1}^n[x_i^2-A*\cos(2\pi x_i)] f(x)=A∗n+i=1∑n​[xi2​−A∗cos(2πxi​)]

function F= Rastrigin(x)
%Rastrigin function
% f(x)=An+\sum_{i=1}^n[x_i^2-Acos(2pix_i)]
% when A=10 and x_i\in[-5.12,5.12]. It has a global minimum at x=0 where f(x)=0

[m,n]=size(x);
A=10;
F=zeros(m,1);
for i=1:m
    F(i)=A*n+sum(x(i,:).^2-A*cos(2*pi*x(i,:)));
end
           

繪出測試圖形:

% %% 經典測試函數
function DrawRastrigin()
% % finding the mimimum value.
% clc %清屏
% clear all; %删除workplace 變量
% close all; %關掉顯示圖形視窗
% 繪制Rastrigin 函數圖形
x = [-5 : 0.05 : 5 ];
y = x;
[X,Y] = meshgrid(x,y);
[row,col] = size(X);
for l = 1 :col
for h = 1 :row
z(h,l) = Rastrigin([X(h,l),Y(h,l)]);
end
end
mesh(X,Y,z);
shading interp
end

           
粒子群算法PSO粒子群算法PSO

利用PSO進行求解:

clear
clc
sizepop=200;
c1=1.49445;
c2=1.49445;
w=0.5;
maxg=500;
D=2;
t1=clock;
[xm,fv,pbest]=stdPSO(@Rastrigin,sizepop,c1,c2,w,maxg,D);
t2=clock;
disp(['PSO時間:',num2str(etime(t2,t1)),'s']);
plot(pbest,'x-');
title(['時間為',num2str(etime(t2,t1)),'s','; 最佳值為:',num2str(fv)])
           
粒子群算法PSO粒子群算法PSO

可看出,該函數的尋優是收斂的,最優個體值接近于理論值(0,0);但是當次元增加時,标準PSO就不是這麼有效了。比如次元增加到10維,最終結果為:

粒子群算法PSO粒子群算法PSO

結果離最優值還有差距,這就需要對标準PSO進行改進了。

PSO算法的向量化實作

在進行MATLAB計算時,應減少循環的出現,盡量用向量或矩陣計算代替循環以此來優化程式運作速度。下面,我将用向量化的方式來重寫标準PSO算法程式。

function [xm,fv,pbest]=stdPSO_v1(fun,sizepop,c1,c2,w,maxg,D)
% Standard PSO
% Input
%fitness: 待優化的目标函數
%sizepop: 粒子數目
%c1: 學習因子1
%c2: 學習因子2
%w: 慣性權重
%maxg: 最大疊代次數
%D: 自變量的次元
%Output 
%xm: 目标函數取最小值時的自變量值
%fv: 目标函數的最小值

%% 參數初始化
% 初始速度和種群上下邊界值
Vmax=1;
Vmin=-1;
popmax=5;
popmin=-5;
%% 産生初始粒子和速度
% for i=1:sizepop
%     %随機産生一個種群
%     pop(i,:)=5*rands(1,D); %初始種群
%     V(i,:)=rands(1,D); %初始化速度
%     %計算适應度
%     fitness(i)=fun(pop(i,:)); %染色體的适應度
% end
pop=5*rands(sizepop,D);
V=rands(sizepop,D);
fitness=fun(pop);

%找最好的染色體
[bestfitness,bestindex]=min(fitness);
zbest=pop(bestindex,:); %全局最佳
gbest=pop; %個體最佳
fitnessgbest=fitness; %個體最佳适應度值
fitnesszbest=bestfitness; %全局最佳适應度值
pbest=zeros(1,maxg);

%% 疊代尋優
for i=1:maxg   %疊代次數
    %速度更新
    V=w*V+c1*rand(sizepop,1).*(gbest - pop)+c2*rand(sizepop,1).*(zbest- pop);
    V(V>Vmax)=Vmax;
    V(V<Vmin)=Vmin;
    %種群更新
    pop=pop+V;
    pop(pop>popmax)=popmax;
    pop(pop<popmin)=popmin;
    %自适應變異
    tmp=rand(sizepop,D)>0.8;
    tmpr=rand(sizepop,D);
    pop(tmp)=tmpr(tmp);
    %适應度值
    fitness=fun(pop);
    %個體最優更新
    tmp=fitness<fitnessgbest;
    if any(tmp) == true
        gbest(tmp,:)=pop(tmp,:);
        fitnessgbest(tmp)=fitness(tmp);
    end
    %群體最優更新
    tmp=fitness<fitnesszbest;
    if any(tmp) == true
        tmpz=pop(tmp,:);
        tmpf=fitness(tmp);
        [fitnesszbest,tmpi]=min(tmpf);
        zbest=tmpz(tmpi,:); %全局最佳
    end
    pbest(i)=fitnesszbest;
end
% 結果輸出
xm=zbest; %最佳個體值
fv=fitnesszbest;
           

重新對上面的最優化問題求解,得到:

粒子群算法PSO粒子群算法PSO

時間對比為0.902s比0.084s,向量化的程式運作時間隻有原始程式的1/10。

reference

  1. MATLAB 優化算法案例分析與應用 ↩︎

繼續閱讀