星型結構的多目標粒子群算法求解多模態多目標問題*

高海軍，潘大志

(西華師范大學數學與信息學院，四川 南充 637009)

1 引言

在實際生活中，存在很多需要同時優化2個或2個以上目標的問題，這類問題被稱為多目標優化問題[1]。在該類問題中，進行優化的目標之間常常會相互制約，因此要使多個目標同時達到最優十分困難。

目前，多目標優化問題的求解方法主要有2類：傳統的數學求解方法和進化算法。傳統的數學解析方法[2]主要是通過動態規劃求解問題的解，但當目標過多或決策變量過多時，數學解析方法求解非常困難。進化算法[3]是一種以種群進化為基礎的啟發式隨機搜索算法，有不斷演化逼近真實Pareto前沿的能力，目前出現了很多求解多目標優化問題的智能算法。Deb等人[4,5]提出NSGA-II(Non-dominated Sorting Genetic Algorithm Ⅱ)算法，使用帶精英策略的快速非支配排序，保證找到的最優解不被拋棄，提高了算法的全局搜索能力。林震等人[6]提出多策略差分進化算法，引入一種動態多策略差分進化模型，分析不同差分進化策略，依據每種策略對鄰域更新的貢獻度，動態地調整其子種群的大小，采用多策略相互協同進化，提高了算法性能。肖閃麗等人[7]提出不同維度的粒子向不同種群學習的策略，以增加種群的多樣性。夏星宇等人[8]在算法中引入均衡因子，提高種群的全局搜索能力。薛蒙蒙等人[9]針對粒子群算法容易陷入局部最優的缺點，提出了模擬退火粒子群算法，增強了算法全局搜索能力。以上算法中很少利用到多目標優化問題的特性：在m個目標的多目標優化問題中，目標空間中的Pareto前沿以及決策空間中的Pareto解集都可以形成一個m-1維的片段連續流型結構[10]。針對多目標優化問題的這種特性，Zhang等人[11]提出將自組織映射網絡應用于種群重組操作，有效地提高了多目標優化算法的性能。梁靜等人[12]提出了自組織映射網絡結構與多目標粒子群優化算法相結合的算法，為粒子結構構造新的鄰域關系，引導粒子在全局范圍內得到更優的解，平衡多目標粒子群算法的多樣性和收斂性。

上面提到的進化計算均是針對多目標優化問題的目標空間及性能分析，以算法獲得的Pareto解集與Pareto前沿的逼近程度進行評價，而沒有考慮Pareto解集在決策空間中的分布情況。決策空間的多樣性吸引了一些研究人員的注意。Shir等人[13]改變CMA-ES(Covariance Matrix Adaptation Evolution Strategies)[14]中的選擇算子和多樣性的度量，同時引入了空間聚合的概念，用于聚合空間中的多樣性維護，從而提高決策空間的多樣性。Tahernezhad等人[15]在優化系統中采用了基于聚類的創新方案，在解空間中獲得了更多樣化的非支配解集。Ulrich 等人[16]將決策空間多樣性納入超體積指標 ，以便同時優化這2個集合度量。這些算法旨在通過考慮決策空間中的多樣性來改善Pareto前沿的分布，但是沒有保留具有不同決策值的相同目標值的解。事實上，應該同時考慮決策空間中解的多樣性和收斂性。針對該問題，Liang等人[17,18]提出多模態多目標優化問題，即同一個目標值對應多個解現象的問題，在研究Pareto前沿的同時討論決策空間中Pareto解集的分布情況。

多模態多目標優化問題[17]由Liang教授提出，并展開相關研究的。Liang等人利用粒子群算法[18]、差分進化算法[19]對多模態多目標優化問題進行求解，取得了有效的成果。受Schütze等人[10]、Li[20]和Liang等人[17,18]的啟發，本文對Pareto解集的個體信息交換的鄰域結構進行改進，提出帶均勻計算方法的基于星型拓撲結構的多目標粒子群優化STMOPSONCMIU(Multi-Objective Particle Swarm Optimization algorithm using Star Topology and New Calculation Method of Individual Uniformity)算法，以提高個體之間信息交換的強度，增強算法的全局搜索能力。針對多模態多目標問題的決策空間，設計一種評價Perato解集個體分布均勻程度的新計算方法，以增強算法在求解多模態多目標優化問題時獲得的Pareto解集與實際的Pareto解集的逼近度。

STMOPSONCMIU在解決多模態多目標優化問題時具有2個較好的能力，一是盡可能多地找到Pareto最優解的能力，二是能保持對應目標空間中同一點的Pareto最優解的搜索能力。在STMOPSONCMIU算法中，借助于粒子的星型拓撲結構，使得每個粒子可以與與之相鄰的4個粒子交換信息，構造小生物環境，從而提高種群的多樣性，搜索出更多的Pareto最優解。另外，針對多模態問題，如圖1所示，解決方案A1和A2均對應目標A，它們在目標空間中的擁擠距離為零，但這2個解在決策空間中的距離很大，本文通過設計新的均勻度距離計算方法來選擇決策空間中的粒子，使得粒子在決策空間中分布更加均勻。

Figure 1 Illustration of multi-modal multi-objective optimization problem圖1 多模態多目標優化問題示意圖

1.1 多目標優化問題

多目標優化問題的數學模型(以最小化為例)[1]可表示為：

minF(X)=(f1(X),f2(X),…,fm(X))

其中,X=(x1,x2,…,xn)∈Ω是決策空間中的n維決策變量；Ω=[ai,bi]是搜索空間的可行域；m是待優化的目標函數個數；F:Ω→Rm是m個待優化的目標函數由決策空間Ω到目標空間Rm的映射關系；gi(X)(i=1,2,…,k)為不等式約束條件；hj(X)(j=1,2,…,l)為等式約束條件。一些基礎的定義[12]如下所示：

定義1(Pareto支配性[1]) 決策變量X支配決策變量Y(記為XY)當且僅當滿足以下條件：

(fi(X)≤fi(Y))∧(fj(X)

其中，i,j=1,2,…,m0。

定義2(Pareto解) 決策變量X是Pareto解[1]:?Y∈Ω:YX，即在可行域內，不存在任何可行解支配X。

定義3(Pareto解集PS(Pareto optional Set)[1])PS={X∈Ω|?Y∈Ω：YX}。

定義4(Pareto前沿PF(Pareto Front)[1])PF={F(X)|X∈PS}

定義5(多模態多目標優化問題MMO(Multimodal Multi-Objective optimization problems)[2]) 對于多目標優化問題，當目標空間的同一個區域對應在決策空間中的解有2個或2個以上時，該問題被稱為多模態多目標優化問題。

圖1簡潔地說明了同一個Pareto前沿對應2個Pareto解集的情況。

1.2 粒子群優化算法

粒子群優化PSO(Particle Swarm Optimization)[21 - 23]算法是一種基于種群群體演化的算法，算法思想來源于鳥類群體性的社會活動，引導鳥群飛向食物。為求解多目標優化問題，將PSO算法擴展為多目標粒子群優化MOPSO(Multi-Objective Particle Swarm Optimization)[21]算法。在該算法中，一個粒子經歷過的歷史最優位置標記為pbest，其鄰域內所有粒子經歷的歷史最優位置標記為nbest。種群中的每個粒子都由pbest和nbest引導，在可行域中飛行。設當前種群中有n個粒子，第i個粒子的第t次迭代的位置為Xi(t)，速度為Vi(t)，粒子位置與速度更新公式可表示為：

Xi(t)=Xi(t-1)+Vi(t)

Vi(t)=ωVi(t-1)+c1r1(Xpbesti-Xi(t))+

c2r2(Xnbesti-Xi(t))

其中,ω為慣性權重，c1與c2為學習因子，r1與r2是在[0,1]內均勻生成的隨機數。

2 星型拓撲結構的多目標粒子群算法

針對多目標問題的Pareto前沿對應決策空間中的多個非支配個體的問題，Liang等人[17]對NSGAII[4]算法進行擴展提出DN-NSGAII(Decision space based Niching NSGAII)算法，該算法旨在定位更多的PS。Yue等人[18]利用環形拓撲結構構建Pareto解集的鄰域關系，同時設計了一種特殊的擁擠度距離計算方法SCD(Special Crowding Distance)，對Pareto解集進行排序，提出環形拓撲結構的多目標優化MO_Ring_PSO_SCD(Multi-Objective Particle Swarm Optimization using Ring topology and Special Crowding Distance)算法，有效地解決了多模態多目標優化問題。受此啟發，為了進一步加強粒子間信息的交換強度，本文修改Pareto解集中個體的鄰域結構，增加鄰域中的個體數，將其變為星型結構，提出基于星型拓撲結構的多目標粒子群優化STMOPSO(Multi-Objective Particle Swarm Optimization using Star Topology)算法。

2.1 Pareto解集的星型拓撲結構

一般的環形拓撲結構的粒子關系如圖2所示。

圖2和圖3中黑色的點表示粒子，從圖2中可知，粒子i只與粒子i-1和粒子i+1交換信息，形成線性結構。在這種結構中，鄰域內個體少，粒子在位置更新過程中無法充分使用其他個體的信息。為增強粒子之間的信息共享度，提高算法的全局搜索能力，對粒子鄰域結構進行擴展，擴充鄰域半徑，擴充鄰域內的個體，由原來的3個個體擴大為5個個體，其結構如圖3所示。由圖3可知，粒子i同時與粒子i-2、粒子i-1、粒子i+1、粒子i+2共4個粒子共享信息，擴大了信息交換程度。分析粒子i的鄰域結構可知，其鄰域結構為星型結構，由線性結構連接變成星型結構，當前粒子的信息交換對象由原來的前后2個變成4個，這使得粒子之間的信息交流更充分，全局搜索能力更強。

Figure 3 Star topology with each particle圖3 粒子的星型拓撲結構

2.2 PS的均勻度計算方法

文獻[18]中提出的特殊擁擠距離計算法(SCD)沒有很好地控制PS中個體分布的均勻度，因而針對該問題提出一種評價PS中個體均勻度的新計算方法NCMIU(New Calculation Method of Individual Uniformity)，通過個體均勻度的值更新PS，使求得的Pareto解集中的個體分布更加均勻，更好地逼近真實的PS。

圖4給出了在二維情況下粒子i的均勻度新計算方法所涉及到的相關個體。個體i的均勻度表示為個體i與其鄰域的個體i-1,i-2,i+1和i+2在各維度上的距離比值。設k表示決策空間的第k維；xi,k表示第i個個體的第k(k=1,2,…,n)維分量；cdi,k表示決策空間中第i個個體在第k維的均勻度；gcdi表示第i個個體的綜合均勻度。因此，在決策空間上個體的均勻度計算公式如式(1)和式(2)所示：

i=3,4,…,m-2

(1)

(2)

Figure 4 Illustration of new calculation method of individual uniformity圖4 個體均勻度的新計算方法示意圖

當個體位于決策空間邊緣時，其第k維的均勻度無法按照式(1)和式(2)進行計算，需單獨計算。第1個個體第k維的均勻度計算公式為：

第2個個體第k維的均勻度計算公式為：

第m-1個個體第k維的均勻度計算公式為：

第m個個體第k維的均勻度計算公式為:

通過上述計算公式計算個體的均勻度，按值升序排序。將該計算方法應用到非支配解集的排序函數中，在決策空間上得到一組分布相對均勻的外部存儲集。以下描述新的均勻度計算方法與粒子群算法結合求解多模態多目標優化問題(STMOPSONCMIU)的算法過程。

2.3 STMOPSONCMIU算法過程

在STMOPSONCMIU算法中，pbestA表示個體最優外部存儲集；nbestA表示鄰近最優外部存儲集；pop表示整個種群；popi(t)表示第t代的第i個粒子；pbestA{i}表示前i個粒子發現的最優位置，nbestA{i}表示第i個粒子鄰域內最佳位置。每個鄰域中有5個粒子，每個粒子與其直接鄰近的左右側粒子相互作用，使得每個粒子與其鄰域最佳位置上的粒子進行信息交互，以避免粒子群體收斂到某個局部最優點。nbestA的引入限制了群體之間的整體信息交流，因此在搜索期間可以形成多個穩定的外部存儲集，偽代碼如算法1所示。

算法1STMOPSONCMIU

// 初始化種群pop(0)

Evaluation(pop(0));

// 初始化pbestA和nbestA

fori=1→mdo

pbestA(i)?pop(0);

nbestA(i)?pop(0);

endfor

fori=1→Max_Gendo

// 更新pbestA

forj=1 →mdo

ifj==1then

temp_nbestA? [pbestA{m-1,:};pbestA{m,:}];

temp_nbestA?[temp_nbestA;pbestA{1,:};pbestA{2,:};pbestA{3,:}];

elseifj==2then

temp_nbestA?[pbestA{m,:};pbestA{1,:}];

temp_nbestA? [temp_nbestA;pbestA{2,:};pbestA{3,:};pbestA{4,:}];

elseifj==m-1then

temp_nbestA?[pbestA{m-3,:};pbestA{m-2,:}];

temp_nbestA?[temp_nbestA;pbestA{m-1,:}];

temp_nbestA?[temp_nbestA;pbestA{m,:};pbestA{1,:}];

elseifj==mthen

temp_nbestA?[pbestA{m-2,:};pbestA{m-1,:}];

temp_nbestA?[temp_nbestA;pbestA{m,:};pbestA{1,:};pbestA{2,:}];

else

temp_nbestA?[pbestA{j-2,:};pbestA{j-1,:};pbestA{j,:}];

temp_nbestA?[temp_nbestA;pbestA{j+1,:};pbestA{j+2,:}];

endif

// 個體均勻度計算

nbestA{j}?NCMIU(temp_NBA(:,1:k+n),k,n);

endfor

forj=1→mdo

Sort particles inpbestAandnbestA;

Selectpbestandnbest;

Upadatepop;

Evalutionpop;

UpadatepbestA;

endfor

endfor

returnthe non_dominated particles innbestA;

2.4 STMOPSONCMIU性能分析

STMOPSONCMIU是基于MO_Ring_PSO_SCD改進的，主要做了以下的改進:一是將環形的2個粒子擴展到星型的5個粒子，增加粒子之間選擇區域和信息交流；二是針對粒子之間的擁擠度距離提出新的計算方法，根據粒子之間的鄰域關系和均勻度的排序選擇全局最優和局部最優。文獻[18]用多模態多目標測試函數4(MMF4)對MO_Ring_PSO_SCD、Omni-optimizer和DN-NSGAII進行收斂性分析，本文對STMOPSONCMIU與MO_Ring_PSO_SCD、Omni-optimizer和DN-NSGAII進行對比分析。

為了更好地對比測試，本文選擇測試函數MMF4，設置種群數量為800，最大迭代次數為100，累計運行30次取平均值，所有的其他參數與第3節實驗參數一致。測試函數MMF4有4個可行解，將MMF4的可行域劃分為4個子區域，每個區域是1個PS解集，分別為Region1{x1∈[-1,0],x2∈[1,2]},Region2{x1∈[0,1],x2∈[1,2]},Region3 {x1∈[-1,0],x2∈(0,1]}和Region4 {x1∈(0,1],x2∈[0,1)}，在測試函數MMF4上，算法的收斂性體現為每次迭代在每個區域中得到解的比例情況，越趨近于25%說明該算法的收斂性越好。理想情況下，算法在4個區域中解的比例各為25%。

如圖5所示為測試函數MMF4從1到100次迭代得到的每個區域解的占比。針對4種算法，STMOPSONCMIU算法得到的解的可行區域占比差距(即縱坐標到0.25的距離)為0.37%，MO_Ring_PSO_SCD算法差距為0.45%，STMOPSONCMIU算法優于MO_Ring_PSO_SCD算法。整體看來， STMOPSONCMIU算法迭代20次后在24.6%～25.3%的小范圍內波動，而MO_Ring_PSO_SCD算法則在50代之后才能達到24.6%～25.3%，說明STMOPSONCMIU算法在收斂時間上占優。但是，Omni-optimizer算法和DN-NSGAII算法波動范圍一直很大，沒有均勻收斂到25%的小范圍內(25%±0.5%)。因此， STMOPSONCMIU算法的收斂性相對優于其他3種算法的。

通過以上分析可知，STMOPSONCMIU主要通過擴大粒子信息交流的范圍，來增強粒子之間的信息轉移；新的擁擠度計算方法可更好地引導粒子收斂到全局最優和局部最優，能有效地解決多模態多目標優化問題。由圖5可知，STMOPSONCMIU的收斂速度比算法MO_Ring_PSO_SCD、Omni-optimizer、DN-NSGAII的快。

Figure 5 Converging behavior of four algorithms on MMF4圖5 MMF4上4個算法的收斂性

3 實驗結論分析

為了評估算法獲得的近似PF的收斂性和均勻性，通過Pareto解集的近似度PSP(Pareto Sets Proximity)[18]和反轉世代距離IGD(Inverted Generational Distance)[12]性能指標進行評比。設P為真實 Pareto最優解集，P*是得到的Pareto最優解集，反轉世代距離IGD的計算如下所示:

Pareto解集的近似度PSP的計算公式為：

其中,IGDx表示決策空間中的反轉世代距離，CR表示Pareto前沿最大差值的覆蓋率。PSP可以很好地反映算法求得的Pareto解集與真實的PS之間的相似度，PSP值越大，說明Pareto解集相似度越高，同時滿足Pareto前沿的分布。其中，IGDx和覆蓋率CR的計算方法參考文獻[18]。同時，超體積值hv(hyper-volume)也可以作為評估算法獲得的近似PF的收斂性和均勻性的性能指標。實驗分析中，用IGDx和hv的值來衡量種群的收斂性和均勻度，IGDx值越小(hv值越大)，得到的PS越好，越接近真實的PS。在STMOPSONCMIU中，ω=0.7298，c1=c2=2.05，其他參數設置參考文獻[4,5,17,19,24,25]。

為了說明本文算法的有效性，通過文獻[18]中的測試函數進行求解，所有的測試函數均獨立運行30次，再與算法MOPSONCMIU、STMOPSO、MOPSO和算法MO_Ring_PSO_SCD[18]求得的PSP值進行比較。表1列出了以上5種算法在11個多目標測試函數上的PSP性能指標。表1中的數值表示同一個算法在同一個測試函數上獨立運行30次的平均值，粗體數值表示所有對比算法在每一個測試函數上最大的PSP值。表2給出了測試函數上性能指標hv、IGDx、IGDf、CR和PSP的值。IGDf表示在目標空間中的反轉世代距離，IGDf的值越小，表示得到的PF越好，越接近真實的PF，且能最大程度覆蓋真實的PF。另外，TTS表示各算法在求解11個測試函數時，對應性能參數獲得最佳的次數。

由表1可知， 除了在函數MMF1和Omni_test 上STMOPSONCMIU算法的PSP平均值比MO_Ring_PSO_SCD算法的低之外，其他9個測試函數上STMOPSONCMIU算法的PSP平均值均為最優。因此，STMOPSONCMIU算法對求解多模態多目標優化問題是有效的。

Table 1 PSPs (mean) of different algorithms 表1 不同算法的PSP值(均值)

通過表2可以看出，由TTS值可知，STMOPSONCMIU算法計算得到的性能指標hv、IGDx、PSP優于其他算法的；測試函數SYM_PART_simple上，算法STMOPSO的CR優于其他算法的，其他的測試函數均是STMOPSONCMIU算法的計算結果占優；在測試函數Omni_test上，STMOPSO算法的性能指標IGDf優于其他算法的，其他的測試函數上，STMOPSONCMIU算法的性能指標占優。為了更好地對比各個算法的性能，將各個測試函數的IGDf、IGDx以及PSP指標取對數分別作折線圖，如圖6所示。由圖6a可知，測試函數MMF1、MMF2上是MOPSO算法占優，其他9個測試函數上都是STMOPSONCMIU算法占優；由圖6b可知，測試函數Omni_test上的IGDf值是STMOPSO算法占優；而由圖7和圖8可知，所有的測試函數中，改進的算法STMOPSONCMIU的性能指標IGDx的PSP值均占優。由此可以說明，STMOPSONCMIU算法可以兼顧PF，從而找到所有的Pareto解集。為更好地體現STMOPSONCMIU算法求解多模態多目標問題的效果，圖9給出了其某一次計算11個測試函數獲得PS與真實PS之間的逼近程度圖，圖9中橫坐標為x1，縱坐標為x2。由圖9可以看出，STMOPSONCMIU算法與MO_Ring_PSO_SCD差距比較小，說明STMOPSONCMIU有效可行，性能上是占優的。

Table 2 Relevant performance index (mean) of each algorithm on each test function表2 各算法在各測試函數上相關性能指標(均值)

Figure 6 Each test function corresponds to IGDf performance value of different algorithms and takes logarithmic polyline graph圖6 各個測試函數對應不同算法的IGDf性能值取對數折線圖

Figure 7 Logarithmic polyline graph was taken for each test function corresponding to IGDx performance values of different algorithms圖7 各個測試函數對應不同算法的IGDx性能值取對數折線圖

Figure 8 Logarithmic polyline graph was taken for each test function corresponding to PSP performance values of different algorithms圖8 各個測試函數對應不同算法的PSP性能值折線圖

Figure 9 PSs obtained by STMOPSONCMIU圖9 STMOPSONCMIU算法在計算各測試函數時得到的PS

4 結束語

本文通過多目標優化問題的Pareto解集多樣性和粒子之間的拓撲結構關系，對其鄰域結構進行擴展，提出一種基于星型拓撲結構的多目標粒子群優化算法STMOPSONCMIU。為更好地確定粒子群算法中鄰域內的最優個體，設計出一種個體均勻度計算方法計算PS中個體的均勻度，通過均勻度來選擇最優個體。將該算法用于求解多模態多目標優化問題，得到的PS分布均勻且與真實的Pareto前沿的逼近程度達到95%以上，表明該算法能較好地解決多模態多目標優化問題。