粒子群算法簡介

一、粒子群算法的歷史

??????? 粒子群算法源于復雜適應系統（Complex Adaptive System,CAS）。CAS理論于1994年正式提出，CAS中的成員稱為主體。比如研究鳥群系統，每個鳥在這個系統中就稱為主體。主體有適應性，它能夠與環境及其他的主體進行交流，并且根據交流的過程“學習”或“積累經驗”改變自身結構與行為。整個系統的演變或進化包括：新層次的產生（小鳥的出生）；分化和多樣性的出現（鳥群中的鳥分成許多小的群）；新的主題的出現（鳥尋找食物過程中，不斷發現新的食物）。

??? 所以CAS系統中的主體具有4個基本特點（這些特點是粒子群算法發展變化的依據）：

• 首先，主體是主動的、活動的。
• 主體與環境及其他主體是相互影響、相互作用的，這種影響是系統發展變化的主要動力。
• 環境的影響是宏觀的，主體之間的影響是微觀的，宏觀與微觀要有機結合。
• 最后，整個系統可能還要受一些隨機因素的影響。

粒子群算法就是對一個CAS系統－－－鳥群社會系統的研究得出的。

???????粒子群算法（ Particle Swarm Optimization, PSO）最早是由Eberhart和Kennedy于1995年提出，它的基本概念源于對鳥群覓食行為的研究。設想這樣一個場景:一群鳥在隨機搜尋食物，在這個區域里只有一塊食物，所有的鳥都不知道食物在哪里，但是它們知道當前的位置離食物還有多遠。那么找到食物的最優策略是什么呢?最簡單有效的就是搜尋目前離食物最近的鳥的周圍區域。???

????? PSO算法就從這種生物種群行為特性中得到啟發并用于求解優化問題。在PSO中，每個優化問題的潛在解都可以想象成d維搜索空間上的一個點，我們稱之為“粒子”（Particle），所有的粒子都有一個被目標函數決定的適應值(Fitness Value )，每個粒子還有一個速度決定他們飛翔的方向和距離，然后粒子們就追隨當前的最優粒子在解空間中搜索。Reynolds對鳥群飛行的研究發現。鳥僅僅是追蹤它有限數量的鄰居但最終的整體結果是整個鳥群好像在一個中心的控制之下.即復雜的全局行為是由簡單規則的相互作用引起的。

二、粒子群算法的具體表述

????? 上面羅嗦了半天，那些都是科研工作者寫論文的語氣，不過，PSO的歷史就像上面說的那樣。下面通俗的解釋PSO算法。

?????? PSO算法就是模擬一群鳥尋找食物的過程，每個鳥就是PSO中的粒子，也就是我們需要求解問題的可能解，這些鳥在尋找食物的過程中，不停改變自己在空中飛行的位置與速度。大家也可以觀察一下，鳥群在尋找食物的過程中，開始鳥群比較分散，逐漸這些鳥就會聚成一群，這個群忽高忽低、忽左忽右，直到最后找到食物。這個過程我們轉化為一個數學問題。尋找函數? y=1-cos(3*x)*exp(-x)的在[0,4]最大值。該函數的圖形如下：

當x=0.9350-0.9450，達到最大值y=1.3706。為了得到該函數的最大值，我們在[0，4]之間隨機的灑一些點，為了演示，我們放置兩個點，并且計算這兩個點的函數值，同時給這兩個點設置在[0，4]之間的一個速度。下面這些點就會按照一定的公式更改自己的位置，到達新位置后，再計算這兩個點的值，然后再按照一定的公式更新自己的位置。直到最后在y=1.3706這個點停止自己的更新。這個過程與粒子群算法作為對照如下：

• 這兩個點就是粒子群算法中的粒子。
• 該函數的最大值就是鳥群中的食物?
• 計算兩個點函數值就是粒子群算法中的適應值，計算用的函數就是粒子群算法中的適應度函數。
• 更新自己位置的一定公式就是粒子群算法中的位置速度更新公式。

下面演示一下這個算法運行一次的大概過程：

第一次初始化

第一次更新位置

第二次更新位置

第21次更新

最后的結果（30次迭代）

最后所有的點都集中在最大值的地方。

標準的粒子群算法

????????? 在上一節的敘述中，唯一沒有給大家介紹的就是函數的這些隨機的點（粒子）是如何運動的，只是說按照一定的公式更新。這個公式就是粒子群算法中的位置速度更新公式。下面就介紹這個公式是什么。在上一節中我們求取函數y=1-cos(3*x)*exp(-x)的在[0,4]最大值。并在[0,4]之間放置了兩個隨機的點，這些點的坐標假設為x1=1.5； x2=2.5；這里的點是一個標量，但是我們經常遇到的問題可能是更一般的情況－－x為一個矢量的情況，比如二維的情況?z=2*x1+3*x2²的情況。這個時候我們的每個粒子為二維，記粒子P1＝(x11,x12),P2=(x21,x22),P3=(x31,x32)，......Pn=(xn1,xn2)。這里n為粒子群群體的規模，也就是這個群中粒子的個數，每個粒子的維數為2。更一般的是粒子的維數為q，這樣在這個種群中有n個粒子，每個粒子為q 維。

?????? 由n個粒子組成的群體對Q維（就是每個粒子的維數）空間進行搜索。每個粒子表示為：x_i＝（x_i1,x_i2,x_i3,...,x_iQ），每個粒子對應的速度可以表示為v_i=(v_i1,v_i2,v_i3,....,v_iQ)，每個粒子在搜索時要考慮兩個因素：

1。自己搜索到的歷史最優值 p_{i ，}p_i=(p_i1,p_i2,....,p_iQ)，i=1,2,3,....,n。

2。全部粒子搜索到的最優值p_g，p_g=(p_g1,p_g2,....,p_gQ)，注意這里的p_g只有一個。

下面給出粒子群算法的位置速度更新公式：

??????? 這里有幾個重要的參數需要大家記憶，因為在以后的講解中將會經常用到：

它們是：

是保持原來速度的系數，所以叫做慣性權重。

是粒子跟蹤自己歷史最優值的權重系數，它表示粒子自身的認識，所以叫“認知”。通常設置為2。

是粒子跟蹤群體最優值的權重系數，它表示粒子對整個群體知識的認識，所以叫做“社會知識”，經常叫做“社會”。通常設置為2。

是[0,1]區間內均勻分布的隨機數。

是對位置更新的時候，在速度前面加的一個系數，這個系數我們叫做約束因子。通常設置為1。

??? 這樣一個標準的粒子群算法就結束了。

下面對整個基本的粒子群的過程給一個簡單的圖形表示：

判斷終止條件可是設置適應值到達一定的數值或者循環一定的次數。

???????? 注意：這里的粒子是同時跟蹤自己的歷史最優值與全局（群體）最優值來改變自己的位置預速度的，所以又叫做全局版本的標準粒子群優化算法。

標準的粒子群算法(局部版本)

在全局版的標準粒子群算法中，每個粒子的速度的更新是根據兩個因素來變化的，這兩個因素是：1. 粒子自己歷史最優值p_i。2. ?粒子群體的全局最優值p_g。如果改變粒子速度更新公式，讓每個粒子的速度的更新根據以下兩個因素更新，A. 粒子自己歷史最優值p_i。B. 粒子鄰域內粒子的最優值pn_k。其余保持跟全局版的標準粒子群算法一樣，這個算法就變為局部版的粒子群算法。

????? 一般一個粒子i 的鄰域隨著迭代次數的增加而逐漸增加，開始第一次迭代，它的鄰域為0，隨著迭代次數鄰域線性變大，最后鄰域擴展到整個粒子群，這時就變成全局版本的粒子群算法了。經過實踐證明：全局版本的粒子群算法收斂速度快，但是容易陷入局部最優。局部版本的粒子群算法收斂速度慢，但是很難陷入局部最優。現在的粒子群算法大都在收斂速度與擺脫局部最優這兩個方面下功夫。其實這兩個方面是矛盾的。看如何更好的折中了。

???? 根據取鄰域的方式的不同，局部版本的粒子群算法有很多不同的實現方法。

第一種方法：按照粒子的編號取粒子的鄰域，取法有四種：1，環形取法 2，隨機環形取法 3，輪形取法 4，隨機輪形取法。

???? 1? 環形2 隨機環形

???? 3 輪形 4隨機輪形

因為后面有以環形取法實現的算法，對環形取法在這里做一點點說明：以粒子1為例，當鄰域是0的時候，鄰域是它本身，當鄰域是1時，鄰域為2，8；當鄰域是2時，鄰域是2，3，7，8；......，以此類推，一直到鄰域為4，這個時候，鄰域擴展到整個例子群體。據文獻介紹（國外的文獻），采用輪形拓撲結構，PSO的效果很好。

第二種方法：按照粒子的歐式距離取粒子的鄰域

????? ?在第一種方法中，按照粒子的編號來得到粒子的鄰域，但是這些粒子其實可能在實際位置上并不相鄰，于是Suganthan提出基于空間距離的劃分方案，在迭代中計算每一個粒子與群中其他粒子的距離。記錄任何2個粒子間的的最大距離為dm。對每一粒子按照||x_a-x_b||/dm計算一個比值。其中||x_a-x_b||是當前粒子a到b的距離。而選擇閾值frac根據迭代次數而變化。當另一粒子b滿足||x_a-x_b||/dm<frac時，認為b成為當前粒子的鄰域。

????? 這種辦法經過實驗，取得較好的應用效果，但是由于要計算所有粒子之間的距離，計算量大，且需要很大的存儲空間，所以，該方法一般不經常使用。

粒子群算法分類

粒子群算法主要分為4個大的分支：

（1）標準粒子群算法的變形

?????? 在這個分支中，主要是對標準粒子群算法的慣性因子、收斂因子（約束因子）、“認知”部分的c1，“社會”部分的c2進行變化與調節，希望獲得好的效果。

???? ?慣性因子的原始版本是保持不變的，后來有人提出隨著算法迭代的進行，慣性因子需要逐漸減小的思想。算法開始階段，大的慣性因子可以是算法不容易陷入局部最優，到算法的后期，小的慣性因子可以使收斂速度加快，使收斂更加平穩，不至于出現振蕩現象。經過本人測試，動態的減小慣性因子w，的確可以使算法更加穩定，效果比較好。但是遞減慣性因子采用什么樣的方法呢？人們首先想到的是線型遞減，這種策略的確很好，但是是不是最優的呢？于是有人對遞減的策略作了研究，研究結果指出：線型函數的遞減優于凸函數的遞減策略，但是凹函數的遞減策略又優于線型的遞減，經過本人測試，實驗結果基本符合這個結論，但是效果不是很明顯。

?? 對于收斂因子，經過證明如果收斂因子取0.729,可以確保算法的收斂，但是不能保證算法收斂到全局最優，經過本人測試，取收斂因子為0.729效果較好。對于社會與認知的系數c2,c1也有人提出：c1先大后小，而c2先小后大的思想，因為在算法運行初期，每個鳥要有大的自己的認知部分而又比較小的社會部分，這個與我們自己一群人找東西的情形比較接近，因為在我們找東西的初期，我們基本依靠自己的知識取尋找，而后來，我們積累的經驗越來越豐富，于是大家開始逐漸達成共識（社會知識），這樣我們就開始依靠社會知識來尋找東西了。

??? 2007年希臘的兩位學者提出將收斂速度比較快的全局版本的粒子群算法與不容易陷入局部最優的局部版本的粒子群算法相結合的辦法，利用的公式是

??? v＝n*v(全局版本)＋（1－n）*v(局部版本)?????? 速度更新公式，v代表速度

?? w（k＋1）＝w(k)＋v???????????????????? 位置更新公式

該算法在文獻中討論了系數n取各種不同情況的情況，并且運行來了20000次來分析各種系數的結果。

（2）粒子群算法的混合

??????? 這個分支主要是將粒子群算法與各種算法相混合，有人將它與模擬退火算法相混合，有些人將它與單純形方法相混合。但是最多的是將它與遺傳算法的混合。根據遺傳算法的三種不同算子可以生成3中不同的混合算法。

?????? ?粒子群算法與選擇算子的結合，這里相混合的思想是：在原來的粒子群算法中，我們選擇粒子群群體的最優值作為p_g，但是相結合的版本是根據所有粒子的適應度的大小給每個粒子賦予一個被選中的概率，然后依據概率對這些粒子進行選擇，被選中的粒子作為p_g，其它的情況都不變。這樣的算法可以在算法運行過程中保持粒子群的多樣性，但是致命的缺點是收斂速度緩慢。

??????? 粒子群算法與雜交算子的結合，結合的思想與遺傳算法的基本一樣，在算法運行過程中根據適應度的大小，粒子之間可以兩兩雜交，比如用一個很簡單的公式

???? w（新）＝n×w1＋（1－n）×w2；

w1與w2就是這個新粒子的父輩粒子。這種算法可以在算法的運行過程中引入新的粒子，但是算法一旦陷入局部最優，那么粒子群算法將很難擺脫局部最優。

???? 粒子群算法與變異算子的結合，結合的思想：測試所有粒子與當前最優的距離，當距離小于一定的數值的時候，可以拿出所有粒子的一個百分比（如10％）的粒子進行隨機初始化，讓這些粒子重新尋找最優值。

???（3）二進制粒子群算法

????? 最初的PSO是從解決連續優化問題發展起來的.Eberhart等又提出了PSO的離散二進制版.用來解決工程實際中的組合優化問題。他們在提出的模型中將粒子的每一維及粒子本身的歷史最優、全局最優限制為1或0，而速度不作這種限制。用速度更新位置時，設定一個閾值，當速度高于該閾值時，粒子的位置取1，否則取0。二進制PSO與遺傳算法在形式上很相似，但實驗結果顯示，在大多數測試函數中，二進制PSO比遺傳算法速度快，尤其在問題的維數增加時

??（4）協同粒子群算法

?? 協同PSO，該方法將粒子的D維分到D個粒子群中，每個粒子群優化一維向量，評價適應度時將這些分量合并為一個完整的向量。例如第i個粒子群，除第i個分量外，其他D-1個分量都設為最優值，不斷用第i個粒子群中的粒子替換第i個分量，直到得到第i維的最優值，其他維相同。為將有聯系的分量劃分在一個群，可將D維向量分配到m個粒子群優化，則前D mod m個粒子群的維數是D/m的向上取整。后m－(D mod m)個粒子群的維數是D/m的向下取整。協同PSO在某些問題上有更快的收斂速度，但該算法容易被欺騙。

?? 基本的粒子群算法的分支就著4個，大部分的粒子群算法都圍繞著這4個分支在變化，其中粒子群算法的變形居多，從根本上來說，幾乎沒有什么新的思想的提出。

標準粒子群算法的實現?

????? 標準粒子群算法的實現思想基本按照粒子群算法(2)----標準的粒子群算法的講述實現。主要分為3個函數。第一個函數為粒子群初始化函數

InitSwarm(SwarmSize......AdaptFunc)其主要作用是初始化粒子群的粒子，并設定粒子的速度、位置在一定的范圍內。本函數所采用的數據結構如下所示：

表ParSwarm記錄的是粒子的位置、速度與當前的適應度值，我們用W來表示位置，用V來代表速度，用F來代表當前的適應度值。在這里我們假設粒子個數為N，每個粒子的維數為D。

W1,1	W1,2	...	W1,D	V1,1	V1,2	...	V1,D-1	V1,D	F1	第1個粒子
W2,1	W2,2	...	W2,D	V2,1	V2,2	...	V2,D-1	V2,D	F2	第2個粒子
...	...	...	...	...	...	...	...	...	...	.......
WN-1,1	WN-1,2	...	WN-1,D-1	VN-1,1	VN-1,2	...	VN-1,D-1	VN-1,D	FN-1	第N-1個粒子
WN,1	WN,2	...	WN,D	VN,1	VN,2	...	VN,D-1	VN,D	FN	第N個粒子

? 表OptSwarm記錄每個粒子的歷史最優解（粒子歷史最好的適應度）以及全部粒子搜索到的全局最優解。用Wg代表全局最優解，W.,1代表每個粒子的歷史最優解。粒子群初始化階段表OptSwarm的前N行與表ParSwarm中的相同，而Wg的值為表ParSwarm中適應度值的最大值對應的行。

Wj,1	Wj,2	...	Wj,D-1	Wj,D	第1個粒子的歷史最優解
Wk,1	Wk,2	...	Wk,D-1	Wk,D	第2個粒子的歷史最優解
...	...	...	...	...	...
Wl,1	Wl,2	...	Wl,D-1	Wl,D	第N-1個粒子的歷史最優解
Wm,1	Wm,2	...	Wm,D-1	Wm,D	第N個粒子的歷史最優解
Wg,1	Wg,2	...	Wg,D-1	Wg,D	全局粒子的歷史最優解

根據這樣的思想MATLAB代碼如下：

?

function?[ParSwarm,OptSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)
%功能描述：初始化粒子群，限定粒子群的位置以及速度在指定的范圍內
%[ParSwarm,OptSwarm,BadSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)
%
%輸入參數：SwarmSize:種群大小的個數
%輸入參數：ParticleSize：一個粒子的維數
%輸入參數：ParticleScope:一個粒子在運算中各維的范圍；
%?????????ParticleScope格式:
%???????????3維粒子的ParticleScope格式:
%???????????????????????????????????[x1Min,x1Max
%????????????????????????????????????x2Min,x2Max
%????????????????????????????????????x3Min,x3Max]
%
%輸入參數：AdaptFunc：適應度函數
%
%輸出：ParSwarm初始化的粒子群
%輸出：OptSwarm粒子群當前最優解與全局最優解
%
%用法[ParSwarm,OptSwarm,BadSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc);
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關的提示信息。
%
%編制人：XXX
%編制時間：2007.3.26
%參考文獻：無
%

%容錯控制
if?nargin~=4
????error('輸入的參數個數錯誤。')
end
if?nargout<2
????error('輸出的參數的個數太少，不能保證以后的運行。');
end

[row,colum]=size(ParticleSize);
if?row>1|colum>1
????error('輸入的粒子的維數錯誤，是一個1行1列的數據。');
end
[row,colum]=size(ParticleScope);
if?row~=ParticleSize|colum~=2
????error('輸入的粒子的維數范圍錯誤。');
end

%初始化粒子群矩陣

%初始化粒子群矩陣，全部設為[0-1]隨機數
%rand('state',0);
ParSwarm=rand(SwarmSize,2*ParticleSize+1);

%對粒子群中位置,速度的范圍進行調節
for?k=1:ParticleSize
????ParSwarm(:,k)=ParSwarm(:,k)*(ParticleScope(k,2)-ParticleScope(k,1))+ParticleScope(k,1);
????%調節速度，使速度與位置的范圍一致
????ParSwarm(:,ParticleSize+k)=ParSwarm(:,ParticleSize+k)*(ParticleScope(k,2)-ParticleScope(k,1))+ParticleScope(k,1);
end
????
%對每一個粒子計算其適應度函數的值

for?k=1:SwarmSize
????ParSwarm(k,2*ParticleSize+1)=AdaptFunc(ParSwarm(k,1:ParticleSize));
end

%初始化粒子群最優解矩陣
OptSwarm=zeros(SwarmSize+1,ParticleSize);
%粒子群最優解矩陣全部設為零
[maxValue,row]=max(ParSwarm(:,2*ParticleSize+1));
%尋找適應度函數值最大的解在矩陣中的位置(行數)
OptSwarm=ParSwarm(1:SwarmSize,1:ParticleSize);
OptSwarm(SwarmSize+1,:)=ParSwarm(row,1:ParticleSize);

下面的函數BaseStepPso實現了標準全局版粒子群算法的單步更新位置速度的功能

?

function?[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)
%功能描述：全局版本：基本的粒子群算法的單步更新位置,速度的算法
%
%[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)
%
%輸入參數：ParSwarm:粒子群矩陣，包含粒子的位置，速度與當前的目標函數值
%輸入參數：OptSwarm：包含粒子群個體最優解與全局最優解的矩陣
%輸入參數：ParticleScope:一個粒子在運算中各維的范圍；
%輸入參數：AdaptFunc：適應度函數
%輸入參數：LoopCount：迭代的總次數
%輸入參數：CurCount：當前迭代的次數
%
%返回值：含意同輸入的同名參數
%
%用法：[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關的提示信息。
%
%編制人：XXX
%編制時間：2007.3.26
%參考文獻：XXX
%參考文獻：XXX
%
%修改記錄
%----------------------------------------------------------------
%2007.3.27
%修改人：XXX
%?添加2*unifrnd(0,1).*SubTract1(row,:)中的unifrnd(0,1)隨機數，使性能大為提高
%參照基于MATLAB的粒子群優化算法程序設計
%
%?總體評價：使用這個版本的調節系數，效果比較好
%

%容錯控制
if?nargin~=8
????error('輸入的參數個數錯誤。')
end
if?nargout~=2
????error('輸出的個數太少，不能保證循環迭代。')
end

%開始單步更新的操作

%*********************************************
%*****更改下面的代碼，可以更改慣性因子的變化*****
%---------------------------------------------------------------------
%線形遞減策略
w=MaxW-CurCount*((MaxW-MinW)/LoopCount);
%---------------------------------------------------------------------
%w固定不變策略
%w=0.7;
%---------------------------------------------------------------------
%參考文獻：陳貴敏，賈建援，韓琪，粒子群優化算法的慣性權值遞減策略研究，西安交通大學學報，2006，1
%w非線形遞減，以凹函數遞減
%w=(MaxW-MinW)*(CurCount/LoopCount)^2+(MinW-MaxW)*(2*CurCount/LoopCount)+MaxW;
%---------------------------------------------------------------------
%w非線形遞減，以凹函數遞減
%w=MinW*(MaxW/MinW)^(1/(1+10*CurCount/LoopCount));
%*****更改上面的代碼，可以更改慣性因子的變化*****
%*********************************************

%得到粒子群群體大小以及一個粒子維數的信息
[ParRow,ParCol]=size(ParSwarm);
%得到粒子的維數
ParCol=(ParCol-1)/2;
SubTract1=OptSwarm(1:ParRow,:)-ParSwarm(:,1:ParCol);

%*********************************************
%*****更改下面的代碼，可以更改c1,c2的變化*****
c1=2;
c2=2;
%---------------------------------------------------------------------
%con=1;
%c1=4-exp(-con*abs(mean(ParSwarm(:,2*ParCol+1))-AdaptFunc(OptSwarm(ParRow+1,:))));
%c2=4-c1;
%----------------------------------------------------------------------
%*****更改上面的代碼，可以更改c1,c2的變化*****
%*********************************************
for?row=1:ParRow
???SubTract2=OptSwarm(ParRow+1,:)-ParSwarm(row,1:ParCol);
???TempV=w.*ParSwarm(row,ParCol+1:2*ParCol)+2*unifrnd(0,1).*SubTract1(row,:)+2*unifrnd(0,1).*SubTract2;
???%限制速度的代碼
???for?h=1:ParCol
???????if?TempV(:,h)>ParticleScope(h,2)
???????????TempV(:,h)=ParticleScope(h,2);
???????end
???????if?TempV(:,h)<-ParticleScope(h,2)
???????????TempV(:,h)=-ParticleScope(h,2)+1e-10;
???????????%加1e-10防止適應度函數被零除
???????end
???end??
???
???%更新速度
???ParSwarm(row,ParCol+1:2*ParCol)=TempV;
???
???%*********************************************
???%*****更改下面的代碼，可以更改約束因子的變化*****
???%---------------------------------------------------------------------
???%a=1;
???%---------------------------------------------------------------------
???a=0.729;
???%*****更改上面的代碼，可以更改約束因子的變化*****
???%*********************************************
???
???%限制位置的范圍
???TempPos=ParSwarm(row,1:ParCol)+a*TempV;
???for?h=1:ParCol
???????if?TempPos(:,h)>ParticleScope(h,2)
???????????TempPos(:,h)=ParticleScope(h,2);
???????end
???????if?TempPos(:,h)<=ParticleScope(h,1)
???????????TempPos(:,h)=ParticleScope(h,1)+1e-10;???????????
???????end
???end

???%更新位置?
???ParSwarm(row,1:ParCol)=TempPos;
???
???%計算每個粒子的新的適應度值
???ParSwarm(row,2*ParCol+1)=AdaptFunc(ParSwarm(row,1:ParCol));
???if?ParSwarm(row,2*ParCol+1)>AdaptFunc(OptSwarm(row,1:ParCol))
???????OptSwarm(row,1:ParCol)=ParSwarm(row,1:ParCol);
???end
end
%for循環結束

%尋找適應度函數值最大的解在矩陣中的位置(行數)，進行全局最優的改變?
[maxValue,row]=max(ParSwarm(:,2*ParCol+1));
if?AdaptFunc(ParSwarm(row,1:ParCol))>AdaptFunc(OptSwarm(ParRow+1,:))
????OptSwarm(ParRow+1,:)=ParSwarm(row,1:ParCol);
end

這兩個函數給出以后，需要一個函數來把這兩個函數組裝起來，以此實現一個完整的粒子群算法，這個函數就是PsoProcess

代碼如下：

function?[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot)
%功能描述：一個循環n次的PSO算法完整過程，返回這次運行的最小與最大的平均適應度,以及在線性能與離線性能
%[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot)
%輸入參數：SwarmSize:種群大小的個數
%輸入參數：ParticleSize：一個粒子的維數
%輸入參數：ParticleScope:一個粒子在運算中各維的范圍；
%?????????ParticleScope格式:
%???????????3維粒子的ParticleScope格式:
%???????????????????????????????????[x1Min,x1Max
%????????????????????????????????????x2Min,x2Max
%????????????????????????????????????x3Min,x3Max]
%
%輸入參數:InitFunc:初始化粒子群函數
%輸入參數:StepFindFunc:單步更新速度，位置函數
%輸入參數：AdaptFunc：適應度函數
%輸入參數：IsStep：是否每次迭代暫停；IsStep＝0，不暫停，否則暫停。缺省不暫停
%輸入參數：IsDraw：是否圖形化迭代過程；IsDraw＝0，不圖形化迭代過程，否則，圖形化表示。缺省不圖形化表示
%輸入參數：LoopCount：迭代的次數；缺省迭代100次
%輸入參數：IsPlot：控制是否繪制在線性能與離線性能的圖形表示；IsPlot=0,不顯示；
%?????????????????IsPlot=1；顯示圖形結果。缺省IsPlot=1
%
%返回值：Result為經過迭代后得到的最優解
%返回值：OnLine為在線性能的數據
%返回值：OffLine為離線性能的數據
%返回值：MinMaxMeanAdapt為本次完整迭代得到的最小與最大的平均適應度
%
%用法[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot);
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關的提示信息。
%
%編制人：XXX
%編制時間：2007.3.26
%參考文獻：XXXXX%

%修改記錄：
%添加MinMaxMeanAdapt，以得到性能評估數據
%修改人：XXX
%修改時間：2007.3.27
%參考文獻：XXX.

%容錯控制
if?nargin<4
????error('輸入的參數個數錯誤。')
end

[row,colum]=size(ParticleSize);
if?row>1|colum>1
????error('輸入的粒子的維數錯誤，是一個1行1列的數據。');
end
[row,colum]=size(ParticleScope);
if?row~=ParticleSize|colum~=2
????error('輸入的粒子的維數范圍錯誤。');
end

%設置缺省值
if?nargin<7
????IsPlot=1;
????LoopCount=100;
????IsStep=0;
????IsDraw=0;
end
if?nargin<8
????IsPlot=1;
????IsDraw=0;
????LoopCount=100;
end
if?nargin<9
????LoopCount=100;
????IsPlot=1;
end
if?nargin<10
????IsPlot=1;
end

%控制是否顯示2維以下粒子維數的尋找最優的過程
if?IsDraw~=0
????DrawObjGraphic(ParticleSize,ParticleScope,AdaptFunc);
end

%初始化種群???????
[ParSwarm,OptSwarm]=InitFunc(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)

%在測試函數圖形上繪制初始化群的位置
if?IsDraw~=0
????if?1==ParticleSize
????for?ParSwarmRow=1:SwarmSize
????????plot([ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,1)],[ParSwarm(ParSwarmRow,3),0],'r*-','markersize',8);
????????text(ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,3),num2str(ParSwarmRow));
????end
end

????if?2==ParticleSize
????????for?ParSwarmRow=1:SwarmSize
????????????stem3(ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,2),ParSwarm(ParSwarmRow,5),'r.','markersize',8);
????????end
????end
end
????
%暫停讓抓圖
if?IsStep~=0
????disp('開始迭代，按任意鍵：')
????pause
end

%開始更新算法的調用
for?k=1:LoopCount
????%顯示迭代的次數：
????disp('----------------------------------------------------------')
????TempStr=sprintf('第?%g?此迭代',k);
????disp(TempStr);
????disp('----------------------------------------------------------')
????
????%調用一步迭代的算法
????[ParSwarm,OptSwarm]=StepFindFunc(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,0.95,0.4,LoopCount,k)
????
????%在目標函數的圖形上繪制2維以下的粒子的新位置
????if?IsDraw~=0
????????if?1==ParticleSize
????????????for?ParSwarmRow=1:SwarmSize
????????????????plot([ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,1)],[ParSwarm(ParSwarmRow,3),0],'r*-','markersize',8);
????????????????text(ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,3),num2str(ParSwarmRow));
????????????end
????????end

????????if?2==ParticleSize
????????????for?ParSwarmRow=1:SwarmSize
????????????????stem3(ParSwarm(ParSwarmRow,1),ParSwarm(ParSwarmRow,2),ParSwarm(ParSwarmRow,5),'r.','markersize',8);
????????????end
????????end
????end
????
????XResult=OptSwarm(SwarmSize+1,1:ParticleSize);
????YResult=AdaptFunc(XResult);????
????if?IsStep~=0
????????XResult=OptSwarm(SwarmSize+1,1:ParticleSize);
????????YResult=AdaptFunc(XResult);????
????????str=sprintf('%g步迭代的最優目標函數值%g',k,YResult);
????????disp(str);
????????disp('下次迭代，按任意鍵繼續');
????????pause
????end
????
????%記錄每一步的平均適應度
????MeanAdapt(1,k)=mean(ParSwarm(:,2*ParticleSize+1));
end
%for循環結束標志

%記錄最小與最大的平均適應度
MinMaxMeanAdapt=[min(MeanAdapt),max(MeanAdapt)];
%計算離線與在線性能
for?k=1:LoopCount
????OnLine(1,k)=sum(MeanAdapt(1,1:k))/k;
????OffLine(1,k)=max(MeanAdapt(1,1:k));
end

for?k=1:LoopCount
????OffLine(1,k)=sum(OffLine(1,1:k))/k;
end

%繪制離線性能與在線性能曲線
if?1==IsPlot
????figure
????hold?on
????title('離線性能曲線圖')
????xlabel('迭代次數');
????ylabel('離線性能');
????grid?on
????plot(OffLine);

????figure
????hold?on
????title('在線性能曲線圖')
????xlabel('迭代次數');
????ylabel('在線性能');
????grid?on
????plot(OnLine);
end

%記錄本次迭代得到的最優結果
XResult=OptSwarm(SwarmSize+1,1:ParticleSize);
YResult=AdaptFunc(XResult);
Result=[XResult,YResult];

這里給出一個使用的例子代碼，并分別解釋各參數的含義：

%打開計時器
tic;
%
Scope=[-50?50
????-50?50
????-50?50
????-50?50
????-50?50
????-50?50
????-50?50
????-50?50
????-50?50
????-50?50];
[v,on,off,minmax]=PsoProcess(20,10,Scope,@initswarm,@BasestepPSO,@Griewank,0,0,4000,0);
toc

在上面的代碼中函數PsoProcess中的20代表粒子群的規模為20個，10代表每個粒子的維數為10，Scope是粒子的每一維的范圍，同時也是速度的范圍，@initswarm是初始化函數的句柄，@BasestepPSO是單步更新的函數句柄，@Griewank是適應度評價函數的句柄，4000代表真個算法循環4000次終止，其他參數參見說明文檔。

總結

以上是生活随笔為你收集整理的粒子群算法详解的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。