閆群民，馬瑞卿，馬永翔，王俊杰

(1．西北工業(yè)大學(xué) 自動(dòng)化學(xué)院，陜西 西安 710072；2.陜西省工業(yè)自動(dòng)化重點(diǎn)實(shí)驗(yàn)室，陜西 漢中 723001 ；3.陜西理工大學(xué) 電氣工程學(xué)院，陜西 漢中 723001)

自然界中的生物均具有一定程度的群體行為。在計(jì)算機(jī)中建模分析這些行為是人工生命研究的一個(gè)重要方面，其背后隱藏的規(guī)律將對(duì)人類的生活產(chǎn)生很大的啟發(fā)。粒子群算法(Particle Swarm Optimization，PSO)[1]是KENNEDY博士和EBERHART教授通過對(duì)鳥群覓食模型進(jìn)行改進(jìn)所提出的。因?yàn)榱Ｗ尤核惴ㄐ柙O(shè)置的參數(shù)少且計(jì)算結(jié)構(gòu)簡單，目前已經(jīng)廣泛地應(yīng)用在實(shí)際工程中，如多目標(biāo)控制、神經(jīng)網(wǎng)絡(luò)研究、復(fù)雜多峰值優(yōu)化求解問題等。眾多的學(xué)者針對(duì)粒子群算法早熟和收斂速度慢等問題，做了大量的優(yōu)化研究，主要可以分為3個(gè)方向：(1) 最早出現(xiàn)的是針對(duì)粒子群算法自身迭代策略的改進(jìn)，將原本固定不變的學(xué)習(xí)因子和慣性權(quán)重系數(shù)改進(jìn)，為其設(shè)置隨著迭代次數(shù)線性或非線性動(dòng)態(tài)變化的調(diào)整策略[2-6]。在搜索的初期、中期及末期等階段有著不同的搜索重點(diǎn)，為達(dá)到平衡全局和局部搜索能力的要求，在保證搜索精度的情況下提升搜索速度。文獻(xiàn)[6]中針對(duì)不同的迭代策略進(jìn)行對(duì)比，經(jīng)過仿真證明在初始權(quán)值和最終權(quán)值相同的情況下，正弦函數(shù)遞減策略優(yōu)于線性策略，而線性策略優(yōu)于正切函數(shù)策略。(2) 將粒子群算法的優(yōu)勢和其他優(yōu)化算法進(jìn)行融合，結(jié)合其他算法的優(yōu)勢來改進(jìn)粒子群算法[7-8]。文獻(xiàn)[8]將牛頓最速下降因子引入粒子群算法中，提高了粒子群體的全局收斂性和收斂速度，克服了這種算法易“早熟”的缺點(diǎn)。(3) 針對(duì)粒子種群的拓?fù)浣Y(jié)構(gòu)進(jìn)行改進(jìn)[9-12]。文獻(xiàn)[11]將原本的D維空間轉(zhuǎn)化到混沌空間，根據(jù)混沌具有遍歷性等特點(diǎn)進(jìn)行迭代，提升了優(yōu)化效果；文獻(xiàn)[12]結(jié)合兩種鄰居拓?fù)浣Y(jié)構(gòu)，將種群分為3個(gè)子群，并采用不同的學(xué)習(xí)因子組合不同的搜索側(cè)重點(diǎn)，經(jīng)過一定迭代次數(shù)后對(duì)子種群的性能進(jìn)行評(píng)估，指導(dǎo)粒子在不同子群間遷移，提升了算法的綜合性能。筆者將前兩種類型進(jìn)行改進(jìn)，采用線性變化的策略控制社會(huì)及自我認(rèn)知因子，利用雙曲正切函數(shù)來控制慣性權(quán)重系數(shù)，在此基礎(chǔ)上融合模擬退火算法提升了算法跳出局部最優(yōu)解的能力。選出5種近幾年被提出的改進(jìn)粒子群算法的優(yōu)化算法[12-16]，與筆者提出的算法共同通過CEC2017中選出的7種測試函數(shù)進(jìn)行尋優(yōu)測試(其中包含1種單峰函數(shù)和6種典型復(fù)雜多峰函數(shù))，結(jié)果證明該優(yōu)化算法在尋優(yōu)速度、收斂精度及算法穩(wěn)定性等方面有顯著的提升。

1 相關(guān)工作

1.1 粒子群算法

粒子群算法的尋優(yōu)過程通過粒子在搜索空間的飛行完成。粒子每次迭代中的移動(dòng)由3個(gè)部分組成，分別是對(duì)上一次速度的繼承、自身學(xué)習(xí)以及種群的信息交互，其表達(dá)式為

vi(k+1)=ωvi(k)+c1r1(Pbest.i(k)-xi(k))+c2r2(Gbest-xi(k)) ，

(1)

xi(k+1)=xi(k)+vi(k+1) ，

(2)

其中，ω為慣性權(quán)重系數(shù)；c1、c2分別表示自身認(rèn)知因子和社會(huì)認(rèn)知因子，是控制粒子群算法的迭代最重要的參數(shù)；xi(k)和vi(k)分別代表第i個(gè)粒子第k次迭代時(shí)的位置和速度；r1和r2為隨機(jī)數(shù)；Pbest.i為第i個(gè)粒子個(gè)體最優(yōu)位置；Gbest為種群最佳位置。飛行粒子的搜索能力不僅受到粒子之前的速度影響，還受到兩個(gè)學(xué)習(xí)因子指導(dǎo)。通過c1控制自我認(rèn)知部分，指導(dǎo)粒子向該粒子自身的歷史最優(yōu)解飛行；c2控制社會(huì)認(rèn)知部分，體現(xiàn)種群中各粒子的信息交換，指導(dǎo)粒子向種群最優(yōu)解的位置飛行。實(shí)際問題的尋優(yōu)的終止條件設(shè)置為全局最優(yōu)解滿足適應(yīng)度的要求，而實(shí)驗(yàn)中為比較不同算法的效果，設(shè)置迭代終止條件為達(dá)到最大迭代次數(shù)時(shí)終止。

每次粒子飛行前，先判斷vi(k)是否越過設(shè)置的速度范圍。若越過，則取速度邊界值替代當(dāng)前速度。飛行后再判斷xi(k)是否超出最大搜索空間。若超出，同樣取邊界值替代當(dāng)前值。根據(jù)相應(yīng)的適應(yīng)度值的變化來更新粒子群中Pbest.i和Gbest，更新方程為

(3)

(4)

1.2 模擬退火法

模擬退火算法(Simulated Annealing，SA )是由米特羅波利斯在20世紀(jì)50年代根據(jù)對(duì)固態(tài)物質(zhì)退火過程的研究而提出的一種算法，后在組合優(yōu)化領(lǐng)域中得到大規(guī)模的應(yīng)用推廣，如現(xiàn)代的玻爾茲曼機(jī)、圖像恢復(fù)處理、VLSI最優(yōu)設(shè)計(jì)、機(jī)器學(xué)習(xí)等科學(xué)領(lǐng)域中的組合優(yōu)化問題。該算法是一種基于蒙特卡羅思想設(shè)計(jì)的近似求解最優(yōu)化問題的方法，結(jié)構(gòu)簡單易實(shí)現(xiàn)，初始條件限制少，使用靈活且運(yùn)行效率高，與其他算法的結(jié)合能力較強(qiáng)，有著廣闊的前景。模擬退火算法的核心源自固體降溫?zé)崃W(xué)過程，退火指物體的能量隨著溫度的降低達(dá)到最低值的過程。米特羅波利斯準(zhǔn)則是模擬退火算法中的核心準(zhǔn)則，其特別之處在于溫度下降過程中不僅能夠接受優(yōu)值，也會(huì)根據(jù)溫度變量產(chǎn)生的概率接受較差的值，增大了算法在尋優(yōu)過程中跳出局部優(yōu)解的概率。將粒子群算法和模擬退火算法相結(jié)合，可以有效地減少粒子群算法陷入局部最優(yōu)解的概率。

模擬退火算法在迭代初始階段需要根據(jù)種群的初始狀態(tài)設(shè)置一個(gè)初始溫度。每次迭代對(duì)模擬固體內(nèi)部粒子在溫度下降情況下的移動(dòng)，根據(jù)米特羅波利斯準(zhǔn)則判斷是否由干擾產(chǎn)生的新解替代全局最優(yōu)解，其表達(dá)式如下：

(5)

其中，Ei(k)表示第i個(gè)粒子在第k次迭代時(shí)的內(nèi)能，即當(dāng)前粒子的適應(yīng)度值；Eg表示當(dāng)前種群最優(yōu)點(diǎn)的內(nèi)能；Ti表示當(dāng)前溫度。溫度每次迭代會(huì)以一定程度線性衰減，尋優(yōu)過程是不斷尋找新解和緩慢降溫的交替過程。Ei(k)完全決定了其下一次產(chǎn)生的新狀態(tài)Ei(k+1)，與之前的Ei(0)至Ei(k-1)無關(guān)，這個(gè)過程是一個(gè)馬爾可夫過程。使用馬爾可夫過程分析模擬退火步驟，經(jīng)過有限次的轉(zhuǎn)換在溫度Ti下的平衡態(tài)的分布如下：

(6)

(7)

其中，Smin是最優(yōu)值的搜索空間集合。當(dāng)溫度降到0時(shí)，pi的分布見式(7)。溫度下降的同時(shí)伴隨著大量的狀態(tài)轉(zhuǎn)移，達(dá)到熱平衡的狀態(tài)，則找到全局最優(yōu)的概率為1。模擬退火算法最大的優(yōu)點(diǎn)是跳出局部最優(yōu)點(diǎn)的能力突出。

2 自適應(yīng)模擬退火粒子群算法

以往應(yīng)用經(jīng)典的粒子群算法設(shè)置ω、c1和c2的值時(shí)，大多依賴經(jīng)驗(yàn)的判斷，或者根據(jù)大量的仿真實(shí)驗(yàn)來確定一個(gè)固定的值。但通過上述的分析知，如果這3個(gè)參數(shù)能夠隨著優(yōu)化的進(jìn)行不斷變化的話，粒子群算法將會(huì)有更加優(yōu)秀的效果。慣性權(quán)重系數(shù)ω直接影響著算法搜索能力的強(qiáng)弱，控制粒子繼承以往運(yùn)動(dòng)趨勢，即搜索飛行的慣性。若需要算法的全局搜索能力較強(qiáng)，則設(shè)置ω的值比較大，但受前一次迭代速度的影響較大，相應(yīng)的粒子飛行距離較大，不利于局部尋優(yōu)。若需要進(jìn)行局部精細(xì)解的搜索，提升局部區(qū)域的尋優(yōu)能力，則設(shè)置ω的值較小，但陷入局部最優(yōu)的概率會(huì)變大。ω值的選取與算法收斂速度和全局搜索能力成正比，與局部搜索能力成反比。筆者選取[-4，4]之間負(fù)的雙曲正切曲線來控制慣性權(quán)重系數(shù)的變化。雙曲正切曲線是一個(gè)非線性的控制策略，在搜索初期其遞減速度較慢，給粒子充分的時(shí)間進(jìn)行大范圍的全局搜索，減小陷入局部最優(yōu)的情況；中期近似線性遞減，逐漸加強(qiáng)局部搜索的能力；后期變化率再次減小，著重細(xì)致的局部搜索，精準(zhǔn)確定全局最優(yōu)解。慣性權(quán)重函數(shù)如下：

ω=(ωmax+ωmin)/2+tanh(-4+8×(kmax-k)/kmax)(ωmax-ωmin)/2 ，

(8)

圖1 慣性系數(shù)自適應(yīng)變化圖

其中，ωmax，ωmin是慣性權(quán)重系數(shù)的最大值和最小值，文中取ωmax=0.95，ωmin=0.40。經(jīng)大量實(shí)驗(yàn)證明，采用如上取值時(shí)算法性能會(huì)大幅度提升[2]。k是當(dāng)前迭代次數(shù)，kmax是最大迭代次數(shù)，其圖像見圖1。

當(dāng)c1>c2時(shí)，粒子的運(yùn)動(dòng)更偏向個(gè)體最優(yōu)的方向；反之，則更偏向群體最優(yōu)方向。筆者所提的優(yōu)化算法初始階段注重全局搜索，著重突出粒子的自我認(rèn)知能力，注重粒子運(yùn)動(dòng)的遍歷性，減小陷入局部最優(yōu)解的概率。隨著迭代的進(jìn)行，加強(qiáng)粒子間的交流，使種群最優(yōu)解的位置對(duì)每個(gè)粒子的運(yùn)行起到更大的影響作用，著重對(duì)Gbest的附近進(jìn)行局部搜索。對(duì)應(yīng)著這樣的搜索策略，隨著迭代次數(shù)的增加，ω不斷減小，c1逐漸減小，c2逐漸增大。筆者采取如下參數(shù)變化策略：

c1=c1max-k(c1max-c1min)/kmax，

(9)

c2=c2min-k(c2min-c2max)/kmax，

(10)

其中，c1max，c1min分別是自我學(xué)習(xí)因子的最大值和最小值；c2max，c2min分別是社會(huì)學(xué)習(xí)因子的最大值和最小值。參照文獻(xiàn)[16]并進(jìn)行大量實(shí)驗(yàn)后，取c1max=2.50，c2max=1.25，c1min=1.25，c2min=2.50。隨著迭代次數(shù)k的增加，c1由2.50線性減小至1.25，而c2則由1.25逐步線性增加到2.50，這樣設(shè)置就滿足了初期注重粒子在空間上的遍歷性要求，增強(qiáng)了全局搜索的能力。在迭代次數(shù)過半時(shí)，c1

將模擬退火算法中的米特羅波利斯準(zhǔn)則引入迭代中。根據(jù)最初的粒子最優(yōu)值設(shè)置初始溫度，并且每次迭代后以一定的降溫系數(shù)μ衰減。具體操作如下：

(11)

其中，T為初始溫度。每次迭代后，計(jì)算更新后位置的內(nèi)能(適應(yīng)度)與種群最優(yōu)點(diǎn)內(nèi)能的差距，根據(jù)式(5)計(jì)算得出的概率與rand()進(jìn)行對(duì)比，判斷是否接受較差的解。取降溫系數(shù)μ=0.95。筆者所提出的自適應(yīng)模擬退火粒子群算法對(duì)粒子群算法的3個(gè)重要參數(shù)進(jìn)行了自適應(yīng)變化且加入了模擬退火操作，增加了尋優(yōu)的精度和速度。具體步驟如下：

步驟1 設(shè)置搜索空間及搜索速度的邊界值，設(shè)置種群規(guī)模及最大迭代次數(shù)kmax。

步驟2 隨機(jī)產(chǎn)生種群中所有粒子的初始位置和初始速度。

步驟3 評(píng)價(jià)全局粒子的適應(yīng)度值并記錄Gbest，根據(jù)式(11)設(shè)置模擬退火的初始溫度。

步驟4 根據(jù)式(8)～式(10)自適應(yīng)地改變?chǔ)?、c1和c2。

步驟5 根據(jù)式(1)改變粒子速度，根據(jù)式(2)進(jìn)行一次迭代尋優(yōu)。

步驟6 計(jì)算移動(dòng)后粒子的適應(yīng)度。

步驟7 根據(jù)式(3)更新粒子的自身的歷史最優(yōu)位置。

步驟8 根據(jù)式(5)計(jì)算接受新解的概率pi(k)。

步驟9 以米特羅波利斯準(zhǔn)則為依據(jù)，對(duì)比概率pi(k)與rand()，判斷是否由產(chǎn)生的新解替代全局最優(yōu)解進(jìn)行退火操作，更新溫度。

步驟10 判斷是否達(dá)到最大迭代次數(shù)kmax，若未達(dá)到則返回步驟4。

步驟11 輸出當(dāng)前最優(yōu)粒子，即尋優(yōu)結(jié)果，算法終止。

3 實(shí)驗(yàn)及結(jié)果分析

3.1 測試函數(shù)及參數(shù)設(shè)置

為了驗(yàn)證筆者提出的模擬退火粒子群優(yōu)化算法的效果，現(xiàn)從CEC2017測試集中選取出1種單峰函數(shù)f1及6種多峰函數(shù)f2-f7來對(duì)筆者所提出的算法進(jìn)行測試，根據(jù)要求設(shè)置其搜索空間。具體的函數(shù)表達(dá)式見表1。

同時(shí)在不同類型的粒子群算法的優(yōu)化算法中選取5種，分別是MSMPSO[12](Multi-population based Self-adaptive Migration PSO)，RPPSO[13](PSO with Random Parameters)，PSO-DAC[14](PSO baesd on Dynamic Acceleration Coefficients)，NDPSO[15](PSO with inertia weight decay by Normal Distribution)和APSO-DA[16](Adaptive PSO via Disturbing Acceleration coefficents)進(jìn)行仿真。通過對(duì)比筆者所提算法與這5種算法在搜索空間維數(shù)分別是D=10和D=30時(shí)，粒子群包含的粒子數(shù)相同，最大迭代次數(shù)相同的情況下的運(yùn)行效果，來體現(xiàn)筆者所改進(jìn)的算法的優(yōu)越性。

表1 測試函數(shù)信息表

3.2 算法精度對(duì)比

使用MATLAB 2016b環(huán)境進(jìn)行仿真驗(yàn)證。所選取的5種對(duì)比算法的慣性權(quán)重系數(shù)和認(rèn)知因子等參數(shù)都按照文獻(xiàn)[12-16]中的設(shè)置，粒子種群規(guī)模統(tǒng)一設(shè)置為30，最大迭代次數(shù)kmax=1 000。每個(gè)測試函數(shù)的維數(shù)在D=10和D=30的情況下分別運(yùn)行30次，結(jié)果記錄如表2及表3所示。從表2～表3中根據(jù)30次尋優(yōu)結(jié)果的平均值和標(biāo)準(zhǔn)差可以看出，無論是針對(duì)單峰函數(shù)還是多峰函數(shù)，本文算法針對(duì)多種測試函數(shù)的搜索精度都優(yōu)于以往的算法。f2～f5這4個(gè)復(fù)雜的多峰函數(shù)的尋優(yōu)精度，優(yōu)于以往優(yōu)化算法數(shù)十倍。f1～f5的尋優(yōu)結(jié)果的標(biāo)準(zhǔn)差最小，有極強(qiáng)的搜索穩(wěn)定性。

表2 D=10時(shí)6種優(yōu)化算法的運(yùn)行結(jié)果對(duì)比

表3 D=30時(shí)6種優(yōu)化算法的運(yùn)行結(jié)果對(duì)比

搜索空間維度的增長，表示著每次迭代所需要處理的數(shù)據(jù)成倍提升，全局最優(yōu)點(diǎn)位置信息也更加復(fù)雜，尋優(yōu)任務(wù)變得十分困難。在實(shí)際工況下，面對(duì)高維問題首先要做的是降低維數(shù)，減少計(jì)算量。D=30的這類高維數(shù)情況較少，但為證明優(yōu)化算法面對(duì)高維問題時(shí)的性能，進(jìn)行了針對(duì)D=30時(shí)的尋優(yōu)實(shí)驗(yàn)。實(shí)驗(yàn)中面對(duì)搜索空間維數(shù)從10增長到30，大部分算法的收斂精度都大幅度下降，而筆者所提出的算法依然有著精準(zhǔn)的尋優(yōu)精度，針對(duì)高維數(shù)的復(fù)雜多峰測試函數(shù)有很強(qiáng)的尋優(yōu)結(jié)果。

3.3 收斂速度對(duì)比

為了更加直觀地對(duì)比不同算法在精度上的差別以及迭代速度的優(yōu)劣性，選出部分測試函數(shù)的優(yōu)化曲線圖。圖2～圖5分別顯示在高維數(shù)D=30時(shí)，筆者所提算法與5種算法的對(duì)比。圖2顯示在面對(duì)簡單的單峰函數(shù)f1時(shí)，文中算法的速度優(yōu)勢并不突出，個(gè)別算法也具有較快的速度和精度。由圖3～圖5可知，在面對(duì)f3、f5和f7復(fù)雜多峰函數(shù)時(shí)，筆者所提算法在收斂速度上的優(yōu)勢得到顯著的體現(xiàn)，尋優(yōu)速度快的同時(shí)也保證了收斂的精度。觀察圖5，面對(duì)搜索空間極大、峰值極多且分布極不規(guī)則的Schwefel測試函數(shù)，并且這個(gè)函數(shù)的各個(gè)局部優(yōu)解的差值極小，尋優(yōu)難度較大?？梢钥闯鲈诘街衅?500次左右)時(shí)大部分算法已經(jīng)陷入局部優(yōu)解，因?yàn)檫@些算法的尋優(yōu)過程在這個(gè)階段，種群中的大部分粒子已經(jīng)聚集在局部最優(yōu)處。隨著越來越多粒子受種群最優(yōu)的影響被引向局部最優(yōu)解處，所以跳出局部最優(yōu)解的能力逐步減弱至喪失；但筆者所提的算法雖然在初始迭代的前400次，其精度優(yōu)勢并不突出，但其優(yōu)勢突出體現(xiàn)在其因?yàn)槿诤狭四M退火操作，極大程度地增強(qiáng)了跳出局部最優(yōu)的能力，在迭代800多次的后期，仍然可以不斷地更新全局最優(yōu)值，保證跳出局部最優(yōu)的能力，提升搜索結(jié)果的精準(zhǔn)性。綜上所述可知，筆者所提的算法在收斂速度上也具有一定程度的優(yōu)勢。

圖2 f1函數(shù)的6種方法優(yōu)化曲線對(duì)比圖

圖3 f3函數(shù)的6種方法優(yōu)化曲線對(duì)比圖

圖4 f5函數(shù)的6種方法優(yōu)化曲線對(duì)比圖

圖5 f7函數(shù)的6種方法優(yōu)化曲線對(duì)比圖

4 結(jié)束語

采用雙曲正切函數(shù)控制慣性權(quán)重系數(shù)隨著迭代進(jìn)行非線性變化，用線性變化策略控制認(rèn)知因子c1、c2，從而保證了3個(gè)控制搜索性能的參數(shù)都能隨著迭代進(jìn)行自適應(yīng)變化；同時(shí)，在粒子群算法尋優(yōu)過程中加入模擬退火操作，通過米特羅波利斯準(zhǔn)則引導(dǎo)種群以一定的概率接受差解，彌補(bǔ)了粒子群算法易陷入局部最優(yōu)的缺陷，極大地提升了粒子群算法的搜索能力，滿足了算法在不同階段的搜索要求。面對(duì)多種高維多峰的測試函數(shù)，與以往的粒子群優(yōu)化算法的尋優(yōu)效果進(jìn)行對(duì)比，其尋優(yōu)精度提升數(shù)十倍，能精準(zhǔn)地收斂在全局最優(yōu)處，以較快的收斂速度完成尋優(yōu)目標(biāo)。