郭 佳，馬朝斌，張紹博，苗萌萌

1.北京交通大學 計算機與信息技術學院，北京100044

2.國家保密科技測評中心，北京100044

1 引言

人工蜂群算法（artificial bee colony algorithm，ABC）是Karaboga[1]于2005 年提出的一種全局優(yōu)化算法，具有計算簡潔、易于實現、控制參數少等優(yōu)點[2]，但該算法更新迭代過程中隨機選擇鄰域解，蜜源進化率較低。國內外學者對ABC 算法進行了諸多改進，主要有三方面：一是通過算法融合進行改進，如Panniem 等[3]提出Modified ABC，用螢火蟲算法產生新的位置來取代偵察蜂階段未更新的位置。耿璐等[4]利用差分進化具有良好的局部更新能力，將ABC 與DE（differential evolution）差分算子結合提出adaptive modified differential ABC。Tian 等[5]將元胞自動機與ABC 結合，利用元胞間的自動更新保持解的多樣性。這些算法將ABC 算法與其他進化算法進行整體融合，提高了搜索精度，但增加了搜索的復雜度。二是利用局部最優(yōu)信息引導算法，如受粒子群算法啟發(fā)，Zhu 等[6]提出GABC（gbeat ABC）算法，用全局最優(yōu)解gbest代替隨機生成的鄰域解引導搜索。周新宇等[7]在文獻[6]的基礎上使用正交實驗設計來產生新的食物源，確保偵察蜂在不同維度上保存放棄解和全局最優(yōu)解的有用信息。Jadon 等[8]在文獻[3]和文獻[6]的基礎上增加了選擇因子，以避免過快地收斂到局部最優(yōu)解。宋曉宇等[9]提出自適應搜索策略混合人工蜂群算法，利用目標函數值的反饋作用進行自適應調節(jié)，同時利用全局最優(yōu)解的引導信息來增強開發(fā)能力。這些方法均以解空間上相對最優(yōu)的某個值作為引導，逐步找到函數最優(yōu)解，但是這類方法容易使算法過分依賴這個值，從而導致陷入局部最優(yōu)解。三是針對解空間進行改進，如Karaboga[10]于2014年提出快速人工蜂群算法，引入周邊食物源的概念，根據周邊食物源參數來判別出周邊食物源種群大小，從而提高蜂蜜采集效率。Jadon 等[11]利用全局和局部鄰域的差分進化，更新食物源信息。Aydin 等[12]提出了增長式蜂群算法，通過迭代讓蜂群規(guī)模逐漸增長，以解決大規(guī)模持續(xù)優(yōu)化問題。Kumar 等[13]提出了一種兩步人工蜂群算法，用K-means 算法代替隨機初始化食物源。但這些思想僅改進了有限的局部解信息，沒有充分利用整個解空間的發(fā)展規(guī)律，對于多峰函數來說，較容易陷入局部最優(yōu)解。

馬爾可夫鏈在啟發(fā)式算法中已有一些工作，比如文獻[14]將馬爾可夫鏈與ABC 算法結合提出MABC（Markov ABC）算法，但該算法利用的二維解空間缺少馬爾可夫性的理論推導。Wang 等[15]將遺傳算法與馬爾可夫鏈結合，采用遺傳算法優(yōu)化初始參數，以馬爾可夫模型作為分類器對大型點焊機加工過程的運行狀態(tài)進行識別。Yang 等[16]將馬爾可夫鏈用于遺傳算法運行前的目標分類，并將改進算法應用于基于藍牙的交通檢測。這些算法均將Markov 鏈用于進化算法運行前或者后的數據處理階段，提高了進化算法的尋優(yōu)精度。

本文在文獻[14]的基礎上進行改進，將原來的二維解空間擴展至三維，從時間維度組成馬爾可夫鏈進行預測，提出IMABC（improved Markov ABC）算法。與文獻[15-16]中馬爾可夫鏈僅作用于進化算法部分數據和一段時間不同，IMABC 算法在空間上將馬爾可夫鏈作用于整個解空間，在時間上將馬爾可夫鏈作用于整個算法運行過程。

2 基本算法

2.1 馬爾可夫鏈

馬爾可夫（Markov）鏈是描述動態(tài)隨機現象的數學模型。它基于系統(tǒng)“狀態(tài)”和“狀態(tài)轉換”的概念[17-18]，根據概率分布，從一種狀態(tài)轉換為另一種狀態(tài)[19]，該方法已應用于眾多實際領域[20-21]，其定義如下：

在有限或可數個值的隨機過程{Xn,n=1,2,…}中，把所取可能值的全體記為狀態(tài)空間S，用{Xn=i}表示過程在時刻n處于狀態(tài)i，假設每當過程處于狀態(tài)i，則在下一個時刻將處于狀態(tài)j的概率是固定的pij，即對任意時刻n，P(Xn+1=j|Xn=i)=pij，若對任意狀態(tài)i1,i2,…,in，始終有P(Xn+1=j|Xn=in,Xn-1=in-1,…,X1=i1)=P(Xn+1=j|Xn=in)，這樣的離散過程被稱為馬爾可夫鏈。假設狀態(tài)空間是有限集，則稱為有限狀態(tài)馬爾可夫鏈。若對于任意的l，有P(Xt+l=j|Xl=i)=P(Xt=，即其轉移概率不依賴于l時,{Xt,t≥0}稱為齊次馬爾可夫鏈[16]，轉移概率為：

其中，Li表示狀態(tài)Si出現的總次數，表示狀態(tài)Si經n步轉移到狀態(tài)Sj的次數，根據轉移概率形成狀態(tài)轉移矩陣，并在下一個時刻進行預測：

其中，Δl是狀態(tài)區(qū)間的下限，Δu是狀態(tài)區(qū)間的上限，f(x)是滾動預測值。

2.2 人工蜂群算法

人工蜂群算法將蜜蜂分為雇傭蜂、跟隨蜂和偵察蜂。該算法的步驟如下：

（1）初始階段。確定蜜源大小SN；蜜蜂的數量與蜜源的數量一致，蜜源Xi=(xi1,xi2,…,xiD)代表候選解，其中i∈{1,2,…,SN}，D是蜜源的維度，蜜源范圍是[xmin,xmax]，首先按式（3）生成初始種群：

其中，rand是介于（0,1）之間的隨機數，j∈{1,2,…,D}，若xij超出了解空間范圍，則按式（4）變?yōu)檫吔缰担?/p>

求解后使用適應度函數來判斷蜜源質量：

（2）雇傭蜂階段。每次采蜜時在蜜源附近按式（6）隨機產生候選解：

其中，xij是搜索空間中的第i個蜜源，xkj是不同于xij的另一個蜜源，k∈{1,2,…,SN}，rand是介于（0，1）之間的隨機數。雇傭蜂通過貪婪選擇算法決定是否由υij替換xij。

（3）跟隨蜂階段。通過輪盤賭算法，按式（7）計算收益率Pi在整個種群中的概率，確定是否選擇該蜜源，選擇蜜源后，跟隨蜂轉變?yōu)楣蛡蚍洹?/p>

（4）偵察蜂階段。當一個蜜源達到開采的上限Limit時，將淘汰蜜源，并根據式（3）隨機選擇一個新的蜜源。

3 算法的改進

3.1 齊次馬爾可夫鏈性質的分析

（1）有限性論證。在初始化階段，解空間的范圍為[xmin,xmax]，每個解的維度D有限，蜂群數量為NP，生成時步為1，隨著時步的推移，生成時步為t，每個維度的參數都是有限的，故所有蜂群狀態(tài)空間有限。

（2）齊次馬爾可夫性論證。當雇傭蜂更新候選解的時候，t步概率表達式為：

其中，pxij→vij(t)為經t步由狀態(tài)xij轉換到vij的概率，L為通過貪婪算法對解的選擇結果，若fit(xij)≤fit(vij)，L=1，否則等于0。根據式（6），為xij至xkj空間中任意一點的概率。由式（8）可以看出，除了考慮空間的上下界外，解從狀態(tài)xij轉換到vij的概率只與xij和隨機選擇的值xkj有關，也就是說，個體的狀態(tài)轉移只與前一時刻有關。因此，人工蜂群算法的種群狀態(tài)序列為有限齊次Markov 鏈。

3.2 優(yōu)化的改進算法

IMABC 算法將采蜜過程分為兩個階段，分割參數記為ω。第一階段為原始ABC 算法，運行ω次后進入第二階段，從初始狀態(tài)到第一階段結束，過程中所有歷史解記為Markov 解空間，即由xij擴充為xijt，記為Xt=(xij1,xij2,…,xijω)，將Xt按照當前解生成的時步順序ω進行排序，最小時步的生成解在最前面，依次形成按時步遞進的Markov 鏈解空間，再將Markov 鏈解空間劃分為N個狀態(tài)子空間，邊界值Qn如式（9）所示：

其中，n∈{1,2,…,N+1}，在確定D和SN的前提下，Xmin為Markov 三維解空間在時間維度的最小值，Xmax為Markov 三維解空間在時間維度的最大值，t∈{1,2,…,ω}。將狀態(tài)子空間記為Sn，按式（10）所示劃分：

將Markov解空間內的解Xt分別置于N個子空間中。按式（1）計算m步狀態(tài)轉移概率，m∈{1,2,…,M}，形成M個狀態(tài)轉移概率矩陣P(m)，選取離預測值最近的M個時步的值進行預測，以離預測時步的遠近為順序，時步小在前，在狀態(tài)轉移矩陣中，選取與m時步相對應的起始狀態(tài)所在的行向量，組成預測矩陣。計算狀態(tài)子空間Sn轉移到Sn′的總概率，如式（11）所示：

Table 1 State transition prediction表1 狀態(tài)轉移預測

根據表1 計算得出：

pbest所對應的狀態(tài)空間即為預測解所在的最優(yōu)狀態(tài)子空間。獲得最優(yōu)狀態(tài)子空間后，按式（2）計算新解，代替適應度最差的解。隨著t的增加，Xt中的元素逐漸增多，優(yōu)化過程變慢，需對Xt進行處理，過濾掉初始階段對預測沒有意義的解，形成新的Markov 解空間V，如式（12）所示：

其中，t為迭代次數，λ為丟棄參數，介于(1,lnt)之間，用于控制丟棄比例。在早期階段，基于初始隨機解生成的解空間不能充分代表解空間的發(fā)展趨勢，應被丟棄，故在時步t較小時，前期解丟棄的比例越大對預測進行的干擾越少。隨著算法運行時步t的增加，解空間具有兩方面的特點：一是范圍迅速擴大；二是隨著適用度值逐漸增加，越來越趨于某一值。故如果按前期的丟棄數量進行丟棄，解空間范圍沒有明顯縮小，算法運行效率會越來越低；如果仍按前期的丟棄比例進行丟棄，會加大收斂速度，但在局部最優(yōu)解較多的函數尋優(yōu)過程中，會減少算法的開拓能力，導致算法尋優(yōu)結果為局部最優(yōu)解。故在Markov 解空間的變化過程中應保證：隨著時步t的增加，丟棄解的數量增加，但占總時步比例減小。按式（12），取λ=1，計算對應于不同循環(huán)時步的新Markov 解空間范圍如表2 所示。

Table 2 Range of solution space changing with time表2 隨時步增加而變化的解空間范圍

從表2 可以看出，在循環(huán)次數較小時，早期解丟棄的數量小，但占總時步比例大；隨著循環(huán)次數的增加，解丟棄的數量大，但占總時步比例小。在算法運行效率和結果準確性之間達到了一定平衡。

IMABC 的偽代碼如下，其中SN為蜜源數量，MaxCycles為最大迭代次數，D為解的維度；ub為解空間的上限，lb為解空間的下限，limit為單個蜜源的開采次數上限，ω為分割參數，N為劃分的子空間數量。

IMABC 主函數偽代碼如下：

3.3 算法性能分析

3.3.1 算法收斂效率分析

IMABC 算法的迭代分為初期、中期和后期三部分。前期為迭代次數小于分割參數階段，因IMABC算法依賴于已經產生的解空間，故需先運行ABC 算法，產生用于迭代的初始解空間，此時算法的收斂效率和運行時間與原始的ABC 算法一致，按式（6），算法復雜度為TABC=sn×(d+1)，其中搜索的復雜度是sn×d，計算每個食物源適應度值的復雜度是sn。中期為根據初始解空間進行迭代更新尋找最優(yōu)值的階段，此時解空間的維度是sn×d×t，其中t從ω開始。在已有初始解空間的基礎上，IMABC 算法依時步的增加不斷更新解空間，既保留了已有解的有效信息，又避免了依賴某一最優(yōu)解導致的算法早熟，可以高效準確地找到最優(yōu)點。但同時，隨著迭代次數t的增加，解空間也隨之增加，按式（12）進行處理后解空間維度變?yōu)?，此時算法復雜度為。后期為t較大的階段，此時IMABC 算法的搜索時間明顯變長，故IMABC算法雖然可以顯著提高算法的精度，卻要犧牲一定的時間成本。但同時也要考慮，在實際尋優(yōu)過程中，往往會在算法達到一定精度后觸發(fā)停止條件，因此IMABC 算法會比ABC 提前達到收斂所需精度，減少IMABC 算法的復雜度和運行時間。

3.3.2 算法尋優(yōu)能力分析

ABC 算法本身具有的特點就是重勘探輕開采，即具有良好的全局尋優(yōu)能力，但是局部尋優(yōu)能力差，IMABC 算法在初期采用原始的ABC 算法，且通過分割參數控制原始ABC 算法在整體算法中的比重，可以有效地保留算法的全局尋優(yōu)能力。

當逐漸趨于最優(yōu)值的時候，需要盡量多地保留解空間的優(yōu)質解信息，但在ABC 算法中，雇傭蜂階段按式（6）找到新蜜源，由于xkj的隨機性，蜂群的先驗信息丟失，適應度高的蜜源被放棄，造成算法震蕩。但是若按以文獻[6]為代表的最優(yōu)值引導收斂算法，會導致過分依賴gbest，在有多個局部極值點的情況下算法易趨向于局部最小值。IMABC 最大限度地保留了解空間的優(yōu)質解信息，并根據適應度值預測出最佳子解空間，進一步更新迭代，避免了搜尋過程中的盲目和過分引導，提高了局部尋優(yōu)能力。

4 實驗驗證

本章包括4 個實驗，第1、2 個實驗分別從結果的準確性、收斂效率和運行時間上比較了IMABC 算法、GABC 算法和ABC 算法，第3 個實驗比較了不同分割參數對IMABC 算法的影響，第4 個實驗比較了不同維度對IMABC 算法的影響。

選擇9 種經典測試函數進行實驗，不同函數的公共參數設置相同，如表3 所示。

這9 個函數[22-24]中f1、f4、f5、f7和f8為單峰函數，僅有一個極值點，其中f5會在給定的間隔上出現不同的階躍現象，并伴隨產生大量局部極值，f2、f3、f6和f9是典型的非線性多模態(tài)函數，在解空間中存在大量的局部極小值點，算法對此類函數進行尋優(yōu)極易陷入局部最優(yōu)解。

4.1 結果準確性比較

該實驗從結果的準確性上比較IMABC 算法、GABC 算法和ABC 算法，設置維度D=30，分割參數ω=1 000，后續(xù)實驗將擴展比較D和ω。不同文獻給出的SN參考值不同，本文采用文獻[1]中的取值30。Limit的大小直接決定IMABC 算法的性能，其值越小，食物源被放棄的可能性越大，算法執(zhí)行率越高。IMABC 算法的Limit參數有3 種典型值100、200和SN×D，本文采用Limit的最小值100，以實現函數的最高執(zhí)行率。設置最大循環(huán)數MaxCycle=10 000。進行30 次實驗，記錄每個值并計算平均值、標準差、最佳值和最差值。測試結果如表4 所示。

從表4 可以看出，在本節(jié)設定的實驗參數前提下，ABC 算法、GABC 算法和IMABC 算法均求出了函數f5的最優(yōu)解，尋優(yōu)精度一致。在其他函數尋優(yōu)過程中，改進的IMABC 算法比原始ABC 算法和GABC 算法的結果精度有明顯提高。

4.2 收斂速度比較

該實驗從算法收斂效率和運行時間上比較了IMABC 算法、GABC 算法和ABC 算法，設置維度D=30，MaxCycle=500，分割參數ω=100，SN=30 。首先，為每個測試功能設置預設精度，如果在最大循環(huán)時間內達到預設精度，則認為該操作成功，否則為失敗。進行10 次實驗，測試結果列于表5。其中，“最小運行次數”為函數重復尋優(yōu)過程10 次中達到預設精度值時的最小運行次數，“最大運行次數”為函數重復尋優(yōu)過程10 次中達到預設精度值時的最大運行次數，“算法成功率”為函數重復尋優(yōu)過程10 次中達到預設精度的比率，“平均運行次數”為10 次函數重復尋優(yōu)過程所有運行次數的平均值，其中沒有達到預設精度的運行次數不算在內，“/”表示函數重復尋優(yōu)過程中沒有一次達到預設精度，“運行時間”為10 次函數尋優(yōu)過程運行時間的平均數。

Table 3 Standard test function表3 測試函數

Table 4 Comparison of optimization accuracy results of different algorithms表4 不同算法的優(yōu)化精度結果比較

Table 5 Comparison of algorithm success rate under given accuracy表5 在給定精度下算法成功率的比較

根據表5 可以看出，在指定預設精度的函數尋優(yōu)過程中，對于函數f2、f3、f6和f8，算法ABC 和GABC 的成功率為0，而IMABC 達到了100%。對于函數f1和f4，算法ABC 的成功率為0，GABC 有部分達到了收斂精度，而IMABC 的收斂成功率達到了100%。對于函數f7，算法IMABC 比ABC 和GABC算法的成功率高。對于函數f5和f9，3 個算法均在500 次循環(huán)內找到了函數的最優(yōu)值，算法成功率達到了100%，但平均運行次數差不多。因為f5為階躍函數，函數值計算的過程中采用取底運算，f9的函數值計算過程中采用絕對值運算，故這兩個函數均存在多個不同的解對應于一個函數值的情況，也就是不同的解對應相同的適應度值，導致IMABC 算法無法有效地判斷Markov 解空間內子空間的優(yōu)劣，故在迭代次數較少時IMABC算法的收斂效率與ABC相差不多。

表5 最后一列為10 次函數重復尋優(yōu)運行時間的平均數，從函數f9所用時間可以看出，IMABC 算法的運行時間約為ABC 和GABC 算法的4 倍，可見IMABC 算法以時間換取精度。但因在預設精度的尋優(yōu)過程中，IMABC 算法達到停止條件所用迭代次數少，所用尋優(yōu)時間明顯降低，f3在IMABC 算法下的尋優(yōu)時間甚至還要少于ABC 及GABC 算法。

圖1 為本次實驗中達到預設精度時的函數收斂曲線。

從圖1 可以看出除函數f5和f9外，其余7 個函數通過IMABC 算法的收斂效率明顯高于ABC 算法和GABC 算法。從圖1 還可以看出IMABC 算法的收斂曲線相對平滑，ABC 和GABC 算法的收斂曲線呈階梯狀較明顯。因為IMABC 算法的尋優(yōu)過程是基于整個解空間的迭代更新，函數迭代過程中基本每一次均會更新解，而ABC 算法迭代過程中新解的更新依賴于隨機選取的鄰域解，存在多次迭代仍未更新的現象，故在函數尋優(yōu)過程中呈階梯下降狀態(tài)。同樣，GABC 利用最優(yōu)解引導尋優(yōu)過程，會出現一次明顯的函數值變化后長期處于不更新的狀態(tài)。

在本次實驗中，函數的尋優(yōu)過程如圖2所示，因圖形較多，僅選取單峰函數f1、f7和多峰函數f2。

Fig.1 Curve of function convergence result圖1 函數收斂結果曲線

根據圖2 的（a）至（c）可以看出，因本次實驗分割參數設為100，故在迭代次數小于100 時IMABC 與ABC、GABC算法收斂速度基本一致，但是隨著迭代次數的增加，根據圖2的（d）至（f）可以看出，在150到200階段，IMABC 算法已能顯出較好的收斂速度，隨著迭代次數的進一步增加，就達到如圖1所示的收斂結果。

4.3 分割參數的擴展比較

分割參數值ω的不同直接影響函數優(yōu)化結果，本實驗驗證不同分割參數對IMABC 算法的影響。設維度D=30，MaxCycle=1 000，SN=10。本實驗從ω=100 開始，每增加100 進行一次運算。當ω=1 000 時，表示未運行IMABC 算法，等效于ABC 算法。分別選擇單峰測試函數f1，多峰測試函數f2和f3，進行30次實驗，記錄每個值并計算平均值、標準差、最優(yōu)值和最差值，測試結果如表6 所示。

從表6 可看出，ω值小于1 000 時的實驗結果均優(yōu)于ω等于1 000 時的性能，進一步證明了IMABC算法的性能優(yōu)于ABC 算法。

Fig.2 Curve of function convergence process圖2 函數收斂過程曲線

Table 6 Comparison of optimization results of different segmentation parameters表6 不同分割參數的優(yōu)化結果比較

同時可以看出，在函數f1的尋優(yōu)過程中，最優(yōu)值出現在ω=500 時，而在f2和f3的尋優(yōu)過程中，最優(yōu)值出現在ω=300 時。因為f2和f3均是多峰函數，有多個局部極值點，尋找極值點較困難，故越早使用IMABC 算法越有助于函數尋優(yōu)。

4.4 不同維度的擴展比較

本實驗驗證不同維度對IMABC 算法性能的影響。設分割參數ω=100，MaxCycle=5 000，SN=10，維度D從10 開始，每增加20 進行一次運算，增至90。分別選擇單峰測試函數f1，多峰測試函數f2和f3，進行30 次實驗，記錄每個值并計算平均值、標準差、最優(yōu)值和最差值，測試結果如表7 所示。

Table 7 Comparison of optimization results of different dimensions表7 不同維度的優(yōu)化結果比較

如表7 所示，隨著維度的增加，尋優(yōu)難度增大，精度逐漸降低。

5 結束語

本文將馬爾可夫鏈與人工蜂群算法融合，提出IMABC 算法。IMABC 算法分為兩個階段，第一階段運行ABC 算法得出初始解空間，第二階段利用馬爾可夫鏈對第一階段產生的解空間進行重構，并進一步預測新解。通過實驗驗證得出IMABC算法具有較高的精度和收斂速度，但算法運行時間較長。在后續(xù)研究中，將進一步完善IMABC 算法，研究最優(yōu)解狀態(tài)空間的動態(tài)變化，提取有用信息，逐漸縮小解空間，減少馬爾可夫鏈預測過程中消耗的運行步驟，從而進一步減小算法的運行時間。