粒子群算法详解
粒子群算法簡介
一、粒子群算法的歷史
??????? 粒子群算法源于復(fù)雜適應(yīng)系統(tǒng)（Complex Adaptive System,CAS）。CAS理論于1994年正式提出，CAS中的成員稱為主體。比如研究鳥群系統(tǒng)，每個(gè)鳥在這個(gè)系統(tǒng)中就稱為主體。主體有適應(yīng)性，它能夠與環(huán)境及其他的主體進(jìn)行交流，并且根據(jù)交流的過程“學(xué)習(xí)”或“積累經(jīng)驗(yàn)”改變自身結(jié)構(gòu)與行為。整個(gè)系統(tǒng)的演變或進(jìn)化包括：新層次的產(chǎn)生（小鳥的出生）；分化和多樣性的出現(xiàn)（鳥群中的鳥分成許多小的群）；新的主題的出現(xiàn)（鳥尋找食物過程中，不斷發(fā)現(xiàn)新的食物）。
??? 所以CAS系統(tǒng)中的主體具有4個(gè)基本特點(diǎn)（這些特點(diǎn)是粒子群算法發(fā)展變化的依據(jù)）：
- 首先，主體是主動(dòng)的、活動(dòng)的。
- 主體與環(huán)境及其他主體是相互影響、相互作用的，這種影響是系統(tǒng)發(fā)展變化的主要?jiǎng)恿Α?
- 環(huán)境的影響是宏觀的，主體之間的影響是微觀的，宏觀與微觀要有機(jī)結(jié)合。
- 最后，整個(gè)系統(tǒng)可能還要受一些隨機(jī)因素的影響。
粒子群算法就是對(duì)一個(gè)CAS系統(tǒng)－－－鳥群社會(huì)系統(tǒng)的研究得出的。
???????粒子群算法（ Particle Swarm Optimization, PSO）最早是由Eberhart和Kennedy于1995年提出，它的基本概念源于對(duì)鳥群覓食行為的研究。設(shè)想這樣一個(gè)場(chǎng)景:一群鳥在隨機(jī)搜尋食物，在這個(gè)區(qū)域里只有一塊食物，所有的鳥都不知道食物在哪里，但是它們知道當(dāng)前的位置離食物還有多遠(yuǎn)。那么找到食物的最優(yōu)策略是什么呢?最簡單有效的就是搜尋目前離食物最近的鳥的周圍區(qū)域。???
????? PSO算法就從這種生物種群行為特性中得到啟發(fā)并用于求解優(yōu)化問題。在PSO中，每個(gè)優(yōu)化問題的潛在解都可以想象成d維搜索空間上的一個(gè)點(diǎn)，我們稱之為“粒子”（Particle），所有的粒子都有一個(gè)被目標(biāo)函數(shù)決定的適應(yīng)值(Fitness Value )，每個(gè)粒子還有一個(gè)速度決定他們飛翔的方向和距離，然后粒子們就追隨當(dāng)前的最優(yōu)粒子在解空間中搜索。Reynolds對(duì)鳥群飛行的研究發(fā)現(xiàn)。鳥僅僅是追蹤它有限數(shù)量的鄰居但最終的整體結(jié)果是整個(gè)鳥群好像在一個(gè)中心的控制之下.即復(fù)雜的全局行為是由簡單規(guī)則的相互作用引起的。
二、粒子群算法的具體表述
????? 上面羅嗦了半天，那些都是科研工作者寫論文的語氣，不過，PSO的歷史就像上面說的那樣。下面通俗的解釋PSO算法。
?????? PSO算法就是模擬一群鳥尋找食物的過程，每個(gè)鳥就是PSO中的粒子，也就是我們需要求解問題的可能解，這些鳥在尋找食物的過程中，不停改變自己在空中飛行的位置與速度。大家也可以觀察一下，鳥群在尋找食物的過程中，開始鳥群比較分散，逐漸這些鳥就會(huì)聚成一群，這個(gè)群忽高忽低、忽左忽右，直到最后找到食物。這個(gè)過程我們轉(zhuǎn)化為一個(gè)數(shù)學(xué)問題。尋找函數(shù)? y=1-cos(3*x)*exp(-x)的在[0,4]最大值。該函數(shù)的圖形如下：
當(dāng)x=0.9350-0.9450，達(dá)到最大值y=1.3706。為了得到該函數(shù)的最大值，我們?cè)赱0，4]之間隨機(jī)的灑一些點(diǎn)，為了演示，我們放置兩個(gè)點(diǎn)，并且計(jì)算這兩個(gè)點(diǎn)的函數(shù)值，同時(shí)給這兩個(gè)點(diǎn)設(shè)置在[0，4]之間的一個(gè)速度。下面這些點(diǎn)就會(huì)按照一定的公式更改自己的位置，到達(dá)新位置后，再計(jì)算這兩個(gè)點(diǎn)的值，然后再按照一定的公式更新自己的位置。直到最后在y=1.3706這個(gè)點(diǎn)停止自己的更新。這個(gè)過程與粒子群算法作為對(duì)照如下：
- 這兩個(gè)點(diǎn)就是粒子群算法中的粒子。
- 該函數(shù)的最大值就是鳥群中的食物?
- 計(jì)算兩個(gè)點(diǎn)函數(shù)值就是粒子群算法中的適應(yīng)值，計(jì)算用的函數(shù)就是粒子群算法中的適應(yīng)度函數(shù)。
- 更新自己位置的一定公式就是粒子群算法中的位置速度更新公式。
下面演示一下這個(gè)算法運(yùn)行一次的大概過程：
第一次初始化
第一次更新位置
第二次更新位置
第21次更新
最后的結(jié)果（30次迭代）
最后所有的點(diǎn)都集中在最大值的地方。
標(biāo)準(zhǔn)的粒子群算法
????????? 在上一節(jié)的敘述中，唯一沒有給大家介紹的就是函數(shù)的這些隨機(jī)的點(diǎn)（粒子）是如何運(yùn)動(dòng)的，只是說按照一定的公式更新。這個(gè)公式就是粒子群算法中的位置速度更新公式。下面就介紹這個(gè)公式是什么。在上一節(jié)中我們求取函數(shù)y=1-cos(3*x)*exp(-x)的在[0,4]最大值。并在[0,4]之間放置了兩個(gè)隨機(jī)的點(diǎn)，這些點(diǎn)的坐標(biāo)假設(shè)為x1=1.5； x2=2.5；這里的點(diǎn)是一個(gè)標(biāo)量，但是我們經(jīng)常遇到的問題可能是更一般的情況－－x為一個(gè)矢量的情況，比如二維的情況?z=2*x1+3*x22的情況。這個(gè)時(shí)候我們的每個(gè)粒子為二維，記粒子P1＝(x11,x12),P2=(x21,x22),P3=(x31,x32)，......Pn=(xn1,xn2)。這里n為粒子群群體的規(guī)模，也就是這個(gè)群中粒子的個(gè)數(shù)，每個(gè)粒子的維數(shù)為2。更一般的是粒子的維數(shù)為q，這樣在這個(gè)種群中有n個(gè)粒子，每個(gè)粒子為q 維。
?????? 由n個(gè)粒子組成的群體對(duì)Q維（就是每個(gè)粒子的維數(shù)）空間進(jìn)行搜索。每個(gè)粒子表示為：xi＝（xi1,xi2,xi3,...,xiQ），每個(gè)粒子對(duì)應(yīng)的速度可以表示為vi=(vi1,vi2,vi3,....,viQ)，每個(gè)粒子在搜索時(shí)要考慮兩個(gè)因素：
1。自己搜索到的歷史最優(yōu)值 pi ，pi=(pi1,pi2,....,piQ)，i=1,2,3,....,n。
2。全部粒子搜索到的最優(yōu)值pg，pg=(pg1,pg2,....,pgQ)，注意這里的pg只有一個(gè)。
下面給出粒子群算法的位置速度更新公式：
??????? 這里有幾個(gè)重要的參數(shù)需要大家記憶，因?yàn)樵谝院蟮闹v解中將會(huì)經(jīng)常用到：
它們是：
是保持原來速度的系數(shù)，所以叫做慣性權(quán)重。
是粒子跟蹤自己歷史最優(yōu)值的權(quán)重系數(shù)，它表示粒子自身的認(rèn)識(shí)，所以叫“認(rèn)知”。通常設(shè)置為2。
是粒子跟蹤群體最優(yōu)值的權(quán)重系數(shù)，它表示粒子對(duì)整個(gè)群體知識(shí)的認(rèn)識(shí)，所以叫做“社會(huì)知識(shí)”，經(jīng)常叫做“社會(huì)”。通常設(shè)置為2。
是[0,1]區(qū)間內(nèi)均勻分布的隨機(jī)數(shù)。
是對(duì)位置更新的時(shí)候，在速度前面加的一個(gè)系數(shù)，這個(gè)系數(shù)我們叫做約束因子。通常設(shè)置為1。
??? 這樣一個(gè)標(biāo)準(zhǔn)的粒子群算法就結(jié)束了。
下面對(duì)整個(gè)基本的粒子群的過程給一個(gè)簡單的圖形表示：
判斷終止條件可是設(shè)置適應(yīng)值到達(dá)一定的數(shù)值或者循環(huán)一定的次數(shù)。
???????? 注意：這里的粒子是同時(shí)跟蹤自己的歷史最優(yōu)值與全局（群體）最優(yōu)值來改變自己的位置預(yù)速度的，所以又叫做全局版本的標(biāo)準(zhǔn)粒子群優(yōu)化算法。
標(biāo)準(zhǔn)的粒子群算法(局部版本)
在全局版的標(biāo)準(zhǔn)粒子群算法中，每個(gè)粒子的速度的更新是根據(jù)兩個(gè)因素來變化的，這兩個(gè)因素是：1. 粒子自己歷史最優(yōu)值pi。2. ?粒子群體的全局最優(yōu)值pg。如果改變粒子速度更新公式，讓每個(gè)粒子的速度的更新根據(jù)以下兩個(gè)因素更新，A. 粒子自己歷史最優(yōu)值pi。B. 粒子鄰域內(nèi)粒子的最優(yōu)值pnk。其余保持跟全局版的標(biāo)準(zhǔn)粒子群算法一樣，這個(gè)算法就變?yōu)榫植堪娴牧Ｗ尤核惴ā?/span>
????? 一般一個(gè)粒子i 的鄰域隨著迭代次數(shù)的增加而逐漸增加，開始第一次迭代，它的鄰域?yàn)?，隨著迭代次數(shù)鄰域線性變大，最后鄰域擴(kuò)展到整個(gè)粒子群，這時(shí)就變成全局版本的粒子群算法了。經(jīng)過實(shí)踐證明：全局版本的粒子群算法收斂速度快，但是容易陷入局部最優(yōu)。局部版本的粒子群算法收斂速度慢，但是很難陷入局部最優(yōu)。現(xiàn)在的粒子群算法大都在收斂速度與擺脫局部最優(yōu)這兩個(gè)方面下功夫。其實(shí)這兩個(gè)方面是矛盾的。看如何更好的折中了。
???? 根據(jù)取鄰域的方式的不同，局部版本的粒子群算法有很多不同的實(shí)現(xiàn)方法。
第一種方法：按照粒子的編號(hào)取粒子的鄰域，取法有四種：1，環(huán)形取法 2，隨機(jī)環(huán)形取法 3，輪形取法 4，隨機(jī)輪形取法。
???? 1? 環(huán)形2 隨機(jī)環(huán)形
???? 3 輪形 4隨機(jī)輪形
因?yàn)楹竺嬗幸原h(huán)形取法實(shí)現(xiàn)的算法，對(duì)環(huán)形取法在這里做一點(diǎn)點(diǎn)說明：以粒子1為例，當(dāng)鄰域是0的時(shí)候，鄰域是它本身，當(dāng)鄰域是1時(shí)，鄰域?yàn)?，8；當(dāng)鄰域是2時(shí)，鄰域是2，3，7，8；......，以此類推，一直到鄰域?yàn)?，這個(gè)時(shí)候，鄰域擴(kuò)展到整個(gè)例子群體。據(jù)文獻(xiàn)介紹（國外的文獻(xiàn)），采用輪形拓?fù)浣Y(jié)構(gòu)，PSO的效果很好。
第二種方法：按照粒子的歐式距離取粒子的鄰域
????? ?在第一種方法中，按照粒子的編號(hào)來得到粒子的鄰域，但是這些粒子其實(shí)可能在實(shí)際位置上并不相鄰，于是Suganthan提出基于空間距離的劃分方案，在迭代中計(jì)算每一個(gè)粒子與群中其他粒子的距離。記錄任何2個(gè)粒子間的的最大距離為dm。對(duì)每一粒子按照||xa-xb||/dm計(jì)算一個(gè)比值。其中||xa-xb||是當(dāng)前粒子a到b的距離。而選擇閾值frac根據(jù)迭代次數(shù)而變化。當(dāng)另一粒子b滿足||xa-xb||/dm<frac時(shí)，認(rèn)為b成為當(dāng)前粒子的鄰域。
????? 這種辦法經(jīng)過實(shí)驗(yàn)，取得較好的應(yīng)用效果，但是由于要計(jì)算所有粒子之間的距離，計(jì)算量大，且需要很大的存儲(chǔ)空間，所以，該方法一般不經(jīng)常使用。
粒子群算法分類
粒子群算法主要分為4個(gè)大的分支：
（1）標(biāo)準(zhǔn)粒子群算法的變形
?????? 在這個(gè)分支中，主要是對(duì)標(biāo)準(zhǔn)粒子群算法的慣性因子、收斂因子（約束因子）、“認(rèn)知”部分的c1，“社會(huì)”部分的c2進(jìn)行變化與調(diào)節(jié)，希望獲得好的效果。
???? ?慣性因子的原始版本是保持不變的，后來有人提出隨著算法迭代的進(jìn)行，慣性因子需要逐漸減小的思想。算法開始階段，大的慣性因子可以是算法不容易陷入局部最優(yōu)，到算法的后期，小的慣性因子可以使收斂速度加快，使收斂更加平穩(wěn)，不至于出現(xiàn)振蕩現(xiàn)象。經(jīng)過本人測(cè)試，動(dòng)態(tài)的減小慣性因子w，的確可以使算法更加穩(wěn)定，效果比較好。但是遞減慣性因子采用什么樣的方法呢？人們首先想到的是線型遞減，這種策略的確很好，但是是不是最優(yōu)的呢？于是有人對(duì)遞減的策略作了研究，研究結(jié)果指出：線型函數(shù)的遞減優(yōu)于凸函數(shù)的遞減策略，但是凹函數(shù)的遞減策略又優(yōu)于線型的遞減，經(jīng)過本人測(cè)試，實(shí)驗(yàn)結(jié)果基本符合這個(gè)結(jié)論，但是效果不是很明顯。
?? 對(duì)于收斂因子，經(jīng)過證明如果收斂因子取0.729,可以確保算法的收斂，但是不能保證算法收斂到全局最優(yōu)，經(jīng)過本人測(cè)試，取收斂因子為0.729效果較好。對(duì)于社會(huì)與認(rèn)知的系數(shù)c2,c1也有人提出：c1先大后小，而c2先小后大的思想，因?yàn)樵谒惴ㄟ\(yùn)行初期，每個(gè)鳥要有大的自己的認(rèn)知部分而又比較小的社會(huì)部分，這個(gè)與我們自己一群人找東西的情形比較接近，因?yàn)樵谖覀冋覗|西的初期，我們基本依靠自己的知識(shí)取尋找，而后來，我們積累的經(jīng)驗(yàn)越來越豐富，于是大家開始逐漸達(dá)成共識(shí)（社會(huì)知識(shí)），這樣我們就開始依靠社會(huì)知識(shí)來尋找東西了。
??? 2007年希臘的兩位學(xué)者提出將收斂速度比較快的全局版本的粒子群算法與不容易陷入局部最優(yōu)的局部版本的粒子群算法相結(jié)合的辦法，利用的公式是
??? v＝n*v(全局版本)＋（1－n）*v(局部版本)?????? 速度更新公式，v代表速度
?? w（k＋1）＝w(k)＋v???????????????????? 位置更新公式
該算法在文獻(xiàn)中討論了系數(shù)n取各種不同情況的情況，并且運(yùn)行來了20000次來分析各種系數(shù)的結(jié)果。
（2）粒子群算法的混合
??????? 這個(gè)分支主要是將粒子群算法與各種算法相混合，有人將它與模擬退火算法相混合，有些人將它與單純形方法相混合。但是最多的是將它與遺傳算法的混合。根據(jù)遺傳算法的三種不同算子可以生成3中不同的混合算法。
?????? ?粒子群算法與選擇算子的結(jié)合，這里相混合的思想是：在原來的粒子群算法中，我們選擇粒子群群體的最優(yōu)值作為pg，但是相結(jié)合的版本是根據(jù)所有粒子的適應(yīng)度的大小給每個(gè)粒子賦予一個(gè)被選中的概率，然后依據(jù)概率對(duì)這些粒子進(jìn)行選擇，被選中的粒子作為pg，其它的情況都不變。這樣的算法可以在算法運(yùn)行過程中保持粒子群的多樣性，但是致命的缺點(diǎn)是收斂速度緩慢。
??????? 粒子群算法與雜交算子的結(jié)合，結(jié)合的思想與遺傳算法的基本一樣，在算法運(yùn)行過程中根據(jù)適應(yīng)度的大小，粒子之間可以兩兩雜交，比如用一個(gè)很簡單的公式
???? w（新）＝n×w1＋（1－n）×w2；
w1與w2就是這個(gè)新粒子的父輩粒子。這種算法可以在算法的運(yùn)行過程中引入新的粒子，但是算法一旦陷入局部最優(yōu)，那么粒子群算法將很難擺脫局部最優(yōu)。
???? 粒子群算法與變異算子的結(jié)合，結(jié)合的思想：測(cè)試所有粒子與當(dāng)前最優(yōu)的距離，當(dāng)距離小于一定的數(shù)值的時(shí)候，可以拿出所有粒子的一個(gè)百分比（如10％）的粒子進(jìn)行隨機(jī)初始化，讓這些粒子重新尋找最優(yōu)值。
???（3）二進(jìn)制粒子群算法
????? 最初的PSO是從解決連續(xù)優(yōu)化問題發(fā)展起來的.Eberhart等又提出了PSO的離散二進(jìn)制版.用來解決工程實(shí)際中的組合優(yōu)化問題。他們?cè)谔岢龅哪Ｐ椭袑⒘Ｗ拥拿恳痪S及粒子本身的歷史最優(yōu)、全局最優(yōu)限制為1或0，而速度不作這種限制。用速度更新位置時(shí)，設(shè)定一個(gè)閾值，當(dāng)速度高于該閾值時(shí)，粒子的位置取1，否則取0。二進(jìn)制PSO與遺傳算法在形式上很相似，但實(shí)驗(yàn)結(jié)果顯示，在大多數(shù)測(cè)試函數(shù)中，二進(jìn)制PSO比遺傳算法速度快，尤其在問題的維數(shù)增加時(shí)
??（4）協(xié)同粒子群算法
?? 協(xié)同PSO，該方法將粒子的D維分到D個(gè)粒子群中，每個(gè)粒子群優(yōu)化一維向量，評(píng)價(jià)適應(yīng)度時(shí)將這些分量合并為一個(gè)完整的向量。例如第i個(gè)粒子群，除第i個(gè)分量外，其他D-1個(gè)分量都設(shè)為最優(yōu)值，不斷用第i個(gè)粒子群中的粒子替換第i個(gè)分量，直到得到第i維的最優(yōu)值，其他維相同。為將有聯(lián)系的分量劃分在一個(gè)群，可將D維向量分配到m個(gè)粒子群優(yōu)化，則前D mod m個(gè)粒子群的維數(shù)是D/m的向上取整。后m－(D mod m)個(gè)粒子群的維數(shù)是D/m的向下取整。協(xié)同PSO在某些問題上有更快的收斂速度，但該算法容易被欺騙。
?? 基本的粒子群算法的分支就著4個(gè)，大部分的粒子群算法都圍繞著這4個(gè)分支在變化，其中粒子群算法的變形居多，從根本上來說，幾乎沒有什么新的思想的提出。
標(biāo)準(zhǔn)粒子群算法的實(shí)現(xiàn)?
????? 標(biāo)準(zhǔn)粒子群算法的實(shí)現(xiàn)思想基本按照粒子群算法(2)----標(biāo)準(zhǔn)的粒子群算法的講述實(shí)現(xiàn)。主要分為3個(gè)函數(shù)。第一個(gè)函數(shù)為粒子群初始化函數(shù)
InitSwarm(SwarmSize......AdaptFunc)其主要作用是初始化粒子群的粒子，并設(shè)定粒子的速度、位置在一定的范圍內(nèi)。本函數(shù)所采用的數(shù)據(jù)結(jié)構(gòu)如下所示：
表ParSwarm記錄的是粒子的位置、速度與當(dāng)前的適應(yīng)度值，我們用W來表示位置，用V來代表速度，用F來代表當(dāng)前的適應(yīng)度值。在這里我們假設(shè)粒子個(gè)數(shù)為N，每個(gè)粒子的維數(shù)為D。
| W1,1 | W1,2 | ... | W1,D | V1,1 | V1,2 | ... | V1,D-1 | V1,D | F1 | 第1個(gè)粒子 |
| W2,1 | W2,2 | ... | W2,D | V2,1 | V2,2 | ... | V2,D-1 | V2,D | F2 | 第2個(gè)粒子 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ....... |
| 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個(gè)粒子 |
| WN,1 | WN,2 | ... | WN,D | VN,1 | VN,2 | ... | VN,D-1 | VN,D | FN | 第N個(gè)粒子 |
? 表OptSwarm記錄每個(gè)粒子的歷史最優(yōu)解（粒子歷史最好的適應(yīng)度）以及全部粒子搜索到的全局最優(yōu)解。用Wg代表全局最優(yōu)解，W.,1代表每個(gè)粒子的歷史最優(yōu)解。粒子群初始化階段表OptSwarm的前N行與表ParSwarm中的相同，而Wg的值為表ParSwarm中適應(yīng)度值的最大值對(duì)應(yīng)的行。
| Wj,1 | Wj,2 | ... | Wj,D-1 | Wj,D | 第1個(gè)粒子的歷史最優(yōu)解 |
| Wk,1 | Wk,2 | ... | Wk,D-1 | Wk,D | 第2個(gè)粒子的歷史最優(yōu)解 |
| ... | ... | ... | ... | ... | ... |
| Wl,1 | Wl,2 | ... | Wl,D-1 | Wl,D | 第N-1個(gè)粒子的歷史最優(yōu)解 |
| Wm,1 | Wm,2 | ... | Wm,D-1 | Wm,D | 第N個(gè)粒子的歷史最優(yōu)解 |
| Wg,1 | Wg,2 | ... | Wg,D-1 | Wg,D | 全局粒子的歷史最優(yōu)解 |
根據(jù)這樣的思想MATLAB代碼如下：
?
function?[ParSwarm,OptSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)%功能描述：初始化粒子群，限定粒子群的位置以及速度在指定的范圍內(nèi)
%[ParSwarm,OptSwarm,BadSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)
%
%輸入?yún)?shù)：SwarmSize:種群大小的個(gè)數(shù)
%輸入?yún)?shù)：ParticleSize：一個(gè)粒子的維數(shù)
%輸入?yún)?shù)：ParticleScope:一個(gè)粒子在運(yùn)算中各維的范圍；
%?????????ParticleScope格式:
%???????????3維粒子的ParticleScope格式:
%???????????????????????????????????[x1Min,x1Max
%????????????????????????????????????x2Min,x2Max
%????????????????????????????????????x3Min,x3Max]
%
%輸入?yún)?shù)：AdaptFunc：適應(yīng)度函數(shù)
%
%輸出：ParSwarm初始化的粒子群
%輸出：OptSwarm粒子群當(dāng)前最優(yōu)解與全局最優(yōu)解
%
%用法[ParSwarm,OptSwarm,BadSwarm]=InitSwarm(SwarmSize,ParticleSize,ParticleScope,AdaptFunc);
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關(guān)的提示信息。
%
%編制人：XXX
%編制時(shí)間：2007.3.26
%參考文獻(xiàn)：無
%
%容錯(cuò)控制
if?nargin~=4
????error('輸入的參數(shù)個(gè)數(shù)錯(cuò)誤。')
end
if?nargout<2
????error('輸出的參數(shù)的個(gè)數(shù)太少，不能保證以后的運(yùn)行。');
end
[row,colum]=size(ParticleSize);
if?row>1|colum>1
????error('輸入的粒子的維數(shù)錯(cuò)誤，是一個(gè)1行1列的數(shù)據(jù)。');
end
[row,colum]=size(ParticleScope);
if?row~=ParticleSize|colum~=2
????error('輸入的粒子的維數(shù)范圍錯(cuò)誤。');
end
%初始化粒子群矩陣
%初始化粒子群矩陣，全部設(shè)為[0-1]隨機(jī)數(shù)
%rand('state',0);
ParSwarm=rand(SwarmSize,2*ParticleSize+1);
%對(duì)粒子群中位置,速度的范圍進(jìn)行調(diào)節(jié)
for?k=1:ParticleSize
????ParSwarm(:,k)=ParSwarm(:,k)*(ParticleScope(k,2)-ParticleScope(k,1))+ParticleScope(k,1);
????%調(diào)節(jié)速度，使速度與位置的范圍一致
????ParSwarm(:,ParticleSize+k)=ParSwarm(:,ParticleSize+k)*(ParticleScope(k,2)-ParticleScope(k,1))+ParticleScope(k,1);
end
????
%對(duì)每一個(gè)粒子計(jì)算其適應(yīng)度函數(shù)的值
for?k=1:SwarmSize
????ParSwarm(k,2*ParticleSize+1)=AdaptFunc(ParSwarm(k,1:ParticleSize));
end
%初始化粒子群最優(yōu)解矩陣
OptSwarm=zeros(SwarmSize+1,ParticleSize);
%粒子群最優(yōu)解矩陣全部設(shè)為零
[maxValue,row]=max(ParSwarm(:,2*ParticleSize+1));
%尋找適應(yīng)度函數(shù)值最大的解在矩陣中的位置(行數(shù))
OptSwarm=ParSwarm(1:SwarmSize,1:ParticleSize);
OptSwarm(SwarmSize+1,:)=ParSwarm(row,1:ParticleSize);
下面的函數(shù)BaseStepPso實(shí)現(xiàn)了標(biāo)準(zhǔn)全局版粒子群算法的單步更新位置速度的功能
?
function?[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)%功能描述：全局版本：基本的粒子群算法的單步更新位置,速度的算法
%
%[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)
%
%輸入?yún)?shù)：ParSwarm:粒子群矩陣，包含粒子的位置，速度與當(dāng)前的目標(biāo)函數(shù)值
%輸入?yún)?shù)：OptSwarm：包含粒子群個(gè)體最優(yōu)解與全局最優(yōu)解的矩陣
%輸入?yún)?shù)：ParticleScope:一個(gè)粒子在運(yùn)算中各維的范圍；
%輸入?yún)?shù)：AdaptFunc：適應(yīng)度函數(shù)
%輸入?yún)?shù)：LoopCount：迭代的總次數(shù)
%輸入?yún)?shù)：CurCount：當(dāng)前迭代的次數(shù)
%
%返回值：含意同輸入的同名參數(shù)
%
%用法：[ParSwarm,OptSwarm]=BaseStepPso(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,MaxW,MinW,LoopCount,CurCount)
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關(guān)的提示信息。
%
%編制人：XXX
%編制時(shí)間：2007.3.26
%參考文獻(xiàn)：XXX
%參考文獻(xiàn)：XXX
%
%修改記錄
%----------------------------------------------------------------
%2007.3.27
%修改人：XXX
%?添加2*unifrnd(0,1).*SubTract1(row,:)中的unifrnd(0,1)隨機(jī)數(shù)，使性能大為提高
%參照基于MATLAB的粒子群優(yōu)化算法程序設(shè)計(jì)
%
%?總體評(píng)價(jià)：使用這個(gè)版本的調(diào)節(jié)系數(shù)，效果比較好
%
%容錯(cuò)控制
if?nargin~=8
????error('輸入的參數(shù)個(gè)數(shù)錯(cuò)誤。')
end
if?nargout~=2
????error('輸出的個(gè)數(shù)太少，不能保證循環(huán)迭代。')
end
%開始單步更新的操作
%*********************************************
%*****更改下面的代碼，可以更改慣性因子的變化*****
%---------------------------------------------------------------------
%線形遞減策略
w=MaxW-CurCount*((MaxW-MinW)/LoopCount);
%---------------------------------------------------------------------
%w固定不變策略
%w=0.7;
%---------------------------------------------------------------------
%參考文獻(xiàn)：陳貴敏，賈建援，韓琪，粒子群優(yōu)化算法的慣性權(quán)值遞減策略研究，西安交通大學(xué)學(xué)報(bào)，2006，1
%w非線形遞減，以凹函數(shù)遞減
%w=(MaxW-MinW)*(CurCount/LoopCount)^2+(MinW-MaxW)*(2*CurCount/LoopCount)+MaxW;
%---------------------------------------------------------------------
%w非線形遞減，以凹函數(shù)遞減
%w=MinW*(MaxW/MinW)^(1/(1+10*CurCount/LoopCount));
%*****更改上面的代碼，可以更改慣性因子的變化*****
%*********************************************
%得到粒子群群體大小以及一個(gè)粒子維數(shù)的信息
[ParRow,ParCol]=size(ParSwarm);
%得到粒子的維數(shù)
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防止適應(yīng)度函數(shù)被零除
???????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;
???
???%計(jì)算每個(gè)粒子的新的適應(yīng)度值
???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循環(huán)結(jié)束
%尋找適應(yīng)度函數(shù)值最大的解在矩陣中的位置(行數(shù))，進(jìn)行全局最優(yōu)的改變?
[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
這兩個(gè)函數(shù)給出以后，需要一個(gè)函數(shù)來把這兩個(gè)函數(shù)組裝起來，以此實(shí)現(xiàn)一個(gè)完整的粒子群算法，這個(gè)函數(shù)就是PsoProcess
代碼如下：
function?[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot)%功能描述：一個(gè)循環(huán)n次的PSO算法完整過程，返回這次運(yùn)行的最小與最大的平均適應(yīng)度,以及在線性能與離線性能
%[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot)
%輸入?yún)?shù)：SwarmSize:種群大小的個(gè)數(shù)
%輸入?yún)?shù)：ParticleSize：一個(gè)粒子的維數(shù)
%輸入?yún)?shù)：ParticleScope:一個(gè)粒子在運(yùn)算中各維的范圍；
%?????????ParticleScope格式:
%???????????3維粒子的ParticleScope格式:
%???????????????????????????????????[x1Min,x1Max
%????????????????????????????????????x2Min,x2Max
%????????????????????????????????????x3Min,x3Max]
%
%輸入?yún)?shù):InitFunc:初始化粒子群函數(shù)
%輸入?yún)?shù):StepFindFunc:單步更新速度，位置函數(shù)
%輸入?yún)?shù)：AdaptFunc：適應(yīng)度函數(shù)
%輸入?yún)?shù)：IsStep：是否每次迭代暫停；IsStep＝0，不暫停，否則暫停。缺省不暫停
%輸入?yún)?shù)：IsDraw：是否圖形化迭代過程；IsDraw＝0，不圖形化迭代過程，否則，圖形化表示。缺省不圖形化表示
%輸入?yún)?shù)：LoopCount：迭代的次數(shù)；缺省迭代100次
%輸入?yún)?shù)：IsPlot：控制是否繪制在線性能與離線性能的圖形表示；IsPlot=0,不顯示；
%?????????????????IsPlot=1；顯示圖形結(jié)果。缺省IsPlot=1
%
%返回值：Result為經(jīng)過迭代后得到的最優(yōu)解
%返回值：OnLine為在線性能的數(shù)據(jù)
%返回值：OffLine為離線性能的數(shù)據(jù)
%返回值：MinMaxMeanAdapt為本次完整迭代得到的最小與最大的平均適應(yīng)度
%
%用法[Result,OnLine,OffLine,MinMaxMeanAdapt]=PsoProcess(SwarmSize,ParticleSize,ParticleScope,InitFunc,StepFindFunc,AdaptFunc,IsStep,IsDraw,LoopCount,IsPlot);
%
%異常：首先保證該文件在Matlab的搜索路徑中，然后查看相關(guān)的提示信息。
%
%編制人：XXX
%編制時(shí)間：2007.3.26
%參考文獻(xiàn)：XXXXX%
%修改記錄：
%添加MinMaxMeanAdapt，以得到性能評(píng)估數(shù)據(jù)
%修改人：XXX
%修改時(shí)間：2007.3.27
%參考文獻(xiàn)：XXX.
%容錯(cuò)控制
if?nargin<4
????error('輸入的參數(shù)個(gè)數(shù)錯(cuò)誤。')
end
[row,colum]=size(ParticleSize);
if?row>1|colum>1
????error('輸入的粒子的維數(shù)錯(cuò)誤，是一個(gè)1行1列的數(shù)據(jù)。');
end
[row,colum]=size(ParticleScope);
if?row~=ParticleSize|colum~=2
????error('輸入的粒子的維數(shù)范圍錯(cuò)誤。');
end
%設(shè)置缺省值
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維以下粒子維數(shù)的尋找最優(yōu)的過程
if?IsDraw~=0
????DrawObjGraphic(ParticleSize,ParticleScope,AdaptFunc);
end
%初始化種群???????
[ParSwarm,OptSwarm]=InitFunc(SwarmSize,ParticleSize,ParticleScope,AdaptFunc)
%在測(cè)試函數(shù)圖形上繪制初始化群的位置
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
%開始更新算法的調(diào)用
for?k=1:LoopCount
????%顯示迭代的次數(shù)：
????disp('----------------------------------------------------------')
????TempStr=sprintf('第?%g?此迭代',k);
????disp(TempStr);
????disp('----------------------------------------------------------')
????
????%調(diào)用一步迭代的算法
????[ParSwarm,OptSwarm]=StepFindFunc(ParSwarm,OptSwarm,AdaptFunc,ParticleScope,0.95,0.4,LoopCount,k)
????
????%在目標(biāo)函數(shù)的圖形上繪制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步迭代的最優(yōu)目標(biāo)函數(shù)值%g',k,YResult);
????????disp(str);
????????disp('下次迭代，按任意鍵繼續(xù)');
????????pause
????end
????
????%記錄每一步的平均適應(yīng)度
????MeanAdapt(1,k)=mean(ParSwarm(:,2*ParticleSize+1));
end
%for循環(huán)結(jié)束標(biāo)志
%記錄最小與最大的平均適應(yīng)度
MinMaxMeanAdapt=[min(MeanAdapt),max(MeanAdapt)];
%計(jì)算離線與在線性能
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('迭代次數(shù)');
????ylabel('離線性能');
????grid?on
????plot(OffLine);
????figure
????hold?on
????title('在線性能曲線圖')
????xlabel('迭代次數(shù)');
????ylabel('在線性能');
????grid?on
????plot(OnLine);
end
%記錄本次迭代得到的最優(yōu)結(jié)果
XResult=OptSwarm(SwarmSize+1,1:ParticleSize);
YResult=AdaptFunc(XResult);
Result=[XResult,YResult];
這里給出一個(gè)使用的例子代碼，并分別解釋各參數(shù)的含義：
%打開計(jì)時(shí)器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
在上面的代碼中函數(shù)PsoProcess中的20代表粒子群的規(guī)模為20個(gè)，10代表每個(gè)粒子的維數(shù)為10，Scope是粒子的每一維的范圍，同時(shí)也是速度的范圍，@initswarm是初始化函數(shù)的句柄，@BasestepPSO是單步更新的函數(shù)句柄，@Griewank是適應(yīng)度評(píng)價(jià)函數(shù)的句柄，4000代表真?zhèn)€算法循環(huán)4000次終止，其他參數(shù)參見說明文檔。
總結(jié)
- 上一篇： echart 表格_市政工程表格不会填？
- 下一篇： SwitchResX for Mac(屏