優化演算法——粒子群演算法(PSO)
一、粒子群演算法的概述
粒子群演算法(PSO)屬於群智慧演算法的一種,是通過模擬鳥群捕食行為設計的。假設區域裡就只有一塊食物(即通常優化問題中所講的最優解),鳥群的任務是找到這個食物源。鳥群在整個搜尋的過程中,通過相互傳遞各自的資訊,讓其他的鳥知道自己的位置,通過這樣的協作,來判斷自己找到的是不是最優解,同時也將最優解的資訊傳遞給整個鳥群,最終,整個鳥群都能聚集在食物源周圍,即我們所說的找到了最優解,即問題收斂。
二、粒子群演算法的流程
粒子群演算法通過設計一種無質量的粒子來模擬鳥群中的鳥,粒子僅具有兩個屬性:速度和位置,速度代表移動的快慢,位置代表移動的方向。每個粒子在搜尋空間中單獨的搜尋最優解,並將其記為當前個體極值
(PSO流程)
下面我們具體解釋下流程圖裡面的每一個步驟:
1、初始化
首先,我們需要設定最大的速度區間,防止超出最大的區間。位置資訊即為整個搜尋空間,我們在速度區間和搜尋空間上隨機初始化速度和位置。設定群體規模
2、個體極值與全域性最優解
個體極值為每個粒子找到的歷史上最優的位置資訊,並從這些個體歷史最優解中找到一個全域性最優解,並與歷史最優解比較,選出最佳的作為當前的歷史最優解。
3、更新速度和位置的公式
更新公式為:
其中,稱為慣性因子,和稱為加速常數,一般取。表示區間上的隨機數。表示第個變數的個體極值的第維。表示全域性最優解的第維。
4、終止條件
有兩種終止條件可以選擇,一是最大代數:;二是相鄰兩代之間的偏差在一個指定的範圍內即停止。我們在實驗中選擇第一種。
三、實驗
我們選擇的測試函式是:Griewank。其基本形式如下:
影象為:
(Griewank函式影象)
在實驗中我們選擇的維數是20;MATLAB程式程式碼如下:
主程式:
-
c1=2;%學習因子
-
c2=2;%學習因子
-
Dimension=20;
-
Size=30;
-
Tmax=500;
-
Velocity_max=1200;%粒子最大速度
-
F_n=2;%測試函式名
-
Fun_Ub=600;%函式上下界
-
Fun_Lb=-600;
-
Position=zeros(Dimension,Size);%粒子位置
-
Velocity=zeros(Dimension,Size);%粒子速度
-
Vmax(1:Dimension)=Velocity_max;%粒子速度上下界
-
Vmin(1:Dimension)=-Velocity_max;
-
Xmax(1:Dimension)=Fun_Ub;%粒子位置上下界,即函式自變數的上下界
-
Xmin(1:Dimension)=Fun_Lb;
-
[Position,Velocity]=Initial_position_velocity(Dimension,Size,Xmax,Xmin,Vmax,Vmin);
-
Pbest_position=Position;%粒子的歷史最優位置,初始值為粒子的起始位置,儲存每個粒子的歷史最優位置
-
Gbest_position=zeros(Dimension,1);%全域性最優的那個粒子所在位置,初始值認為是第1個粒子
-
for j=1:Size
-
Pos=Position(:,j);%取第j列,即第j個粒子的位置
-
fz(j)=Fitness_Function(Pos,F_n,Dimension);%計算第j個粒子的適應值
-
end
-
[Gbest_Fitness,I]=min(fz);%求出所有適應值中最小的那個適應值,並獲得該粒子的位置
-
Gbest_position=Position(:,I);%取最小適應值的那個粒子的位置,即I列
-
for itrtn=1:Tmax
-
time(itrtn)=itrtn;
-
Weight=1;
-
r1=rand(1);
-
r2=rand(1);
-
for i=1:Size
-
Velocity(:,i)=Weight*Velocity(:,i)+c1*r1*(Pbest_position(:,i)-Position(:,i))+c2*r2*(Gbest_position-Position(:,i));
-
end
-
%限制速度邊界
-
for i=1:Size
-
for row=1:Dimension
-
if Velocity(row,i)>Vmax(row)
-
Veloctity(row,i)=Vmax(row);
-
elseif Velocity(row,i)<Vmin(row)
-
Veloctity(row,i)=Vmin(row);
-
else
-
end
-
end
-
end
-
Position=Position+Velocity;
-
%限制位置邊界
-
for i=1:Size
-
for row=1:Dimension
-
if Position(row,i)>Xmax(row)
-
Position(row,i)=Xmax(row);
-
elseif Position(row,i)<Xmin(row)
-
Position(row,i)=Xmin(row);
-
else
-
end
-
end
-
end
-
for j=1:Size
-
P_position=Position(:,j)';%取一個粒子的位置
-
fitness_p(j)=Fitness_Function(P_position,F_n,Dimension);
-
if fitness_p(j)< fz(j) %粒子的適應值比運動之前的適應值要好,更新原來的適應值
-
Pbest_position(:,j)=Position(:,j);
-
fz(j)=fitness_p(j);
-
end
-
if fitness_p(j)<Gbest_Fitness
-
Gbest_Fitness=fitness_p(j);
-
end
-
end
-
[Gbest_Fitness_new,I]=min(fz);%更新後的所有粒子的適應值,取最小的那個,並求出其編號
-
Best_fitness(itrtn)=Gbest_Fitness_new; %記錄每一代的最好適應值
-
Gbest_position=Pbest_position(:,I);%最好適應值對應的個體所在位置
-
end
-
plot(time,Best_fitness);
-
xlabel('迭代的次數');ylabel('適應度值P_g');
初始化:
-
function [Position,Velocity] = Initial_position_velocity(Dimension,Size,Xmax,Xmin,Vmax,Vmin)
-
for i=1:Dimension
-
Position(i,:)=Xmin(i)+(Xmax(i)-Xmin(i))*rand(1,Size); % 產生合理範圍內的隨機位置,rand(1,Size)用於產生一行Size個隨機數
-
Velocity(i,:)=Vmin(i)+(Vmax(i)-Vmin(i))*rand(1,Size);
-
end
-
end
適應值計算:
-
function Fitness=Fitness_Function(Pos,F_n,Dimension)
-
switch F_n
-
case 1
-
Func_Sphere=Pos(:)'*Pos(:);
-
Fitness=Func_Sphere;
-
case 2
-
res1=Pos(:)'*Pos(:)/4000;
-
res2=1;
-
for row=1:Dimension
-
res2=res2*cos(Pos(row)/sqrt(row));
-
end
-
Func_Griewank=res1-res2+1;
-
Fitness=Func_Griewank;
-
end
最終的收斂曲線:
(收斂曲線)