基于Kriging模型的多目標代理優(yōu)化算法及其收斂性評估

張建俠，宋明順，方興華，鄧鈺佳

(中國計量大學 經濟與管理學院，浙江 杭州 310018)

0 引言

實踐中的優(yōu)化設計問題通常是多目標優(yōu)化問題(Multiobjective Optimization Problem, MOP)。由于優(yōu)化目標之間的沖突性，MOP一般不存在單個絕對的最優(yōu)解，其目標是尋求一組Pareto最優(yōu)意義下的折衷解[1](Pareto Set, PS)。智能優(yōu)化算法(如多目標進化算法)是求解MOP的有效方法，但由于在優(yōu)化進程中需要大量調用函數(shù)，不能滿足昂貴仿真情形下的優(yōu)化設計需求[2]。為提高復雜系統(tǒng)的優(yōu)化設計效率，基于近似模型的代理優(yōu)化方法應運而生[3-4]。

Kriging模型可以同時提供任意試驗點上的預測值及預測誤差，易于構造引導優(yōu)化過程快速收斂的自適應優(yōu)化策略[5]。以有效全局優(yōu)化(Efficient Global Optimization, EGO)算法[6]為代表的自適應代理優(yōu)化算法，合理利用設計的加點準則序貫選取能收斂到全局最優(yōu)解的新試驗點，優(yōu)化效率比傳統(tǒng)進化類方法高出1～2個數(shù)量級[7]。近年來，自適應優(yōu)化思想逐漸推廣至MOP。由于容易理解和執(zhí)行，學者們常采用直接優(yōu)化Pareto解集質量指標(performance indicator)的方式優(yōu)化原MOP[8]。

文獻[9]基于超體積(HyperVolume, HV)質量指標，提出了期望超體積改進(Expected Hypervolume Improvement, EHVI)準則及其精確計算方法，并證明了EHVI準則的單調性質;文獻[10]對ParEGO等4種多目標代理優(yōu)化算法進行了比較，指出基于EHVI準則的代理優(yōu)化算法是唯一一種在求解人工構造問題和實際工程問題時，都表現(xiàn)優(yōu)異的算法?？尚行愿怕?Probability of Feasibility, PoF)準則能估計試驗點落入可行域的概率，在代理優(yōu)化研究中常被用于設計黑箱約束的應對策略[11]。文獻[12]進一步指出，基于PoF準則的約束應對策略不需要擬合拉格朗日乘子函數(shù)和選擇罰因子，與基于增廣拉格朗日乘子法的約束優(yōu)化策略[13]相比，更易于執(zhí)行且優(yōu)化效果好;文獻[14]將PoF準則分別與Kriging模型預測標準差和EHVI準則結合，提出了探索MOP可行域和改進PS質量的加點策略。當優(yōu)化問題包含多個非連通可行子區(qū)域時，其可行域探索準則不能確保找出所有子區(qū)域，且近似PS改進準則逼近可行域邊界上最優(yōu)解的效率不高[15]。

評估算法的收斂性是代理優(yōu)化研究的另一課題。文獻[16]歸類分析了多目標進化算法研究中的收斂準則，并強調了綜合應用多個評估指標(如Epsilon、R2等)的重要性。對多目標代理優(yōu)化算法，代理模型的近似誤差將會影響求得的近似PS，因此PS的不確定性反映了算法的收斂性。條件模擬(conditional simulation)與Kriging模型有相同理論基礎，但更側重于對區(qū)域變量空間分布的隨機模擬，被廣泛應用于地質、土壤、水文、生態(tài)等領域的不確定性評價[17]。針對固定翼飛機采集地面火災數(shù)據(jù)存在采樣不足的問題，文獻[18]提出了利用Kriging模型和條件模擬擬合火輻射能量密度的方法;文獻[19]針對不考慮約束影響的優(yōu)化問題，分析了條件模擬方法評估算法收斂性的可能性。

為提高黑箱系統(tǒng)多目標優(yōu)化設計的效率，綜合考慮Kriging模型的預測不確定性、試驗點的PoF和EHVI、試驗點之間的距離，并利用隨機集理論，提出一種改進的多目標代理優(yōu)化算法及收斂性分析方法。該算法對可行域非連通的問題也有效；近似PS改進準則兼具刻畫可行域邊界的能力；收斂性評估從分析近似PS不確定性的角度評估算法的收斂性。通過兩個算例與已有算法作比較，計算結果驗證了所提算法的高效性。

1 多目標優(yōu)化和Kriging模型

1.1 問題描述

一個典型的包含m個目標函數(shù)和r個約束條件的MOP可以表示為

s.t.

gi(x)≤0,i=1,…,r;

x∈X?d。

(1)

其中：m>1;x為d維設計變量;y(x)為目標函數(shù)向量;gi(x)為第i個約束條件。由于優(yōu)化目標之間的沖突性，MOP通常不存在絕對的最優(yōu)解，而是尋求一組Pareto意義下的折衷解[1]。設G為問題(1)的可行域，x,x′∈G是兩個可行解，稱x支配x′(記作xx′)當且僅當：?i∈{1,…,m},yi(x)≤yi(x′)∧y(x)≠y(x′)。Pareto最優(yōu)解集(PS)是所有Pareto最優(yōu)解構成的集合：PS={x∈G|x*∈G:x*x}。PS在目標空間上的投影稱為Pareto最優(yōu)前沿(Pareto Front, PF)：PF={y(x)|x∈PS}。

1.2 Kriging模型

Kriging模型將仿真模型的響應y(x)看作是隨機過程的實現(xiàn)，即

(2)

(3)

1.3 其他相關概念

辨識優(yōu)化問題可行域、改進PS質量是多目標代理優(yōu)化算法的基本功能。為此，首先介紹以下概念。

(1)期望超體積改進[9]

超體積指的是被PF支配的空間的體積，其取值越大越好。超體積改進(HV Improvement，HVI)則被用來度量新試驗數(shù)據(jù)(x,y(x))帶來的超體積的增量。

HVI(x)=HVI(y(x),PF,r)=

(4)

(5)

(2)可行性概率[11]

考慮約束條件的影響，令G={x∈X|g(x)≤0}表示可行域，則新試驗點x帶來的可行超體積改進為

I(x)=I(y(x),PF,r)=1g(x)≤0·

HVI(y(x),PF,r)

(6)

其中1g(x)≤0為示性函數(shù)。若目標函數(shù)y獨立于約束條件g，則新試驗點x帶來的期望的可行超體積改進可表示為

HVI(y(x),PF,r))

(7)

(3)條件模擬

(8)

2 基于Kriging模型的多目標代理優(yōu)化算法

本文提出的多目標代理優(yōu)化算法(Surrogate-based Multiobjective Optimization Algorithm, SMOA)的基本流程如圖1所示,包括辨識可行域、改進解集質量及評估解集不確定性3個主要部分。

2.1 辨識可行域和改進近似PF

(1)辨識可行域

若MOP的約束條件較難滿足(如可行域較小)，初始樣本點集不含可行試驗點，則先用以下加點準則選取一個可行試驗點：

(9)

(10)

(2)改進近似PF

約束優(yōu)化問題的最優(yōu)解通常位于可行域邊界上。為此，提出如下兼顧目標改進和可行域邊界刻畫的近似PF改進準則：

Cpf(x)=maxx∈Ppf(EHVI(x)·PoF(x))。

(11)

式中Ppf表示備選樣本集(同時優(yōu)化EHVI和PoF時求得的Pareto集),

Ppf=maxx∈X(EHVI(x),PoF(x))。

(12)

式中EHVI(x)=EHVI(x,PF,r)表示試驗點x的期望超體積改進。從具有Pareto最優(yōu)性的備選樣本集中選取新試驗點，不僅縮小了新試驗點選擇范圍、避免在無效區(qū)域選取新試驗點，還使新試驗點兼具了刻畫可行域邊界的能力，進而提高Cpf(x)準則選取新試驗點和改進近似PF的針對性，提高優(yōu)化效率和穩(wěn)健性。

2.2 度量近似PF的不確定性

工程實踐中MOP的真實PF是未知的，使優(yōu)化算法的收斂性評估遇到困難(Epsilon、世代距離等解集質量指標需要用到真實PF)。為此，提出一種利用條件模擬度量近似PF不確定性，并將其作為算法收斂性評價依據(jù)的間接方法。

(1)條件Pareto前沿(Conditional Pareto Fronts，CPFs)

條件模擬方法生成N個CPFs的步驟見算法1，這些CPFs反映了近似PF的不確定。

算法1生成N個CPFs。

步驟1選取模擬點。

(1)隨機生成n個模擬點{e1,…,en}?X?d。

(2)利用約束條件的Kriging模型，選出“可行”的模擬點{e1′,…,en′}。

步驟2利用條件模擬方法預測目標函數(shù)響應值：

步驟3利用“可行”的模擬點及響應數(shù)據(jù)，生成一個條件PF及條件PS。

步驟4重復以上步驟N次，得到N個CPFs。

文獻[22]證明了CPFs在分布上與被其支配的隨機集等價。因此，CPFs的不確定性又可借助Vorob’ev期望和方差的概念加以描述[23]。

(2)Vorob’ev期望和方差

定義1收斂函數(shù)和上水平集[22-23]。設Y是拓撲空間B?m中的隨機閉集，定義Y的收斂函數(shù)(coverage function)為pY:y∈BP(y∈Y)。利用收斂函數(shù)，定義Y的上水平集(upper level set)為Qβ={z∈B,pY(z)≥β}，并稱β為分位數(shù)。

定義2Vorob’ev期望和方差[22-23]。令m表示m上的勒貝格測度且E(m(Y))≤+∞，若存在β*滿足E(m(Y))=m(Qβ*),則定義Y的Vorob’ev期望為其上水平集Qβ*;否則，先用式子{m(Qβ)≤E(m(Y))≤m(Qβ*),?β>β*}求分位數(shù)β*,再令對應的水平集Qβ*為Y的Vorob’ev期望。定義Y與其Vorob’ev期望Qβ*的對稱差的期望E(μ(Qβ*ΔY))為Y的Vorob’ev方差。

(3)近似PF的不確定性

算法2評估近似PF的不確定性(求Vorob’ev期望和方差)。

步驟1利用算法1生成近似PF的CPFs。

步驟4二分法求Vorob’ev期望的β*分位數(shù)

(1)取a=0,b=1。

步驟5利用收斂函數(shù)和β*，確定Y的Vorob’ev期望Qβ*(上水平集)。

3 算例與結果分析

本章通過兩個典型算例，將提出的SMOA算法與文獻[14]的KEMOCO代理優(yōu)化算法及U-NSGA-Ⅲ[24]進化算法作比較，以驗證SMOA算法的有效性、高效性及條件模擬方法評估算法收斂性的可行性。

3.1 數(shù)值算例

算例1

f1(x)=(x1-10)2+(x2-15)2+10,

(13)

其中x=(x1,x2)∈[-5,10]×[0,15]。算例1為自構造函數(shù)。

算例2[25]

f1(x)=-(25(x1-2)2+(x2-2)2+

(x3-1)2+(x4-4)2+(x5-1)2),

g1=2-x1-x2≤0,

g2=x1+x2-6≤0,

g3=x2-x1-2≤0,

g4=x1-3x2-2≤0

g5=(x3-3)2+x4-4≤0,

g6=4-(x5-3)2-x6≤0。

(14)

其中：x1,x2,x6∈[0,10],x3,x5∈[1,5],x4∈[0,6]。

算例1的可行域由3個非連通子區(qū)域組成且可行域占設計變量空間的比例小，可用于驗證SMOA算法可行域探索準則的有效性。算例2包含6個變量和6個約束條件，可用于驗證SMOA算法在求解較復雜MOP時的高效性和通用性。

3.2 試驗設定

為方便分析，首先對算例的變量作歸一化處理，使x∈[0,1]d。其次，假定算例的目標和約束條件相互獨立，并分別構建他們的Kriging模型(取初始設計的樣本容量為k=5×d，其中d為設計變量的維數(shù))。這是由于：優(yōu)化問題函數(shù)之間的相關性信息通常是未知的，且構建多函數(shù)聯(lián)合Kriging模型的過程復雜、易出現(xiàn)數(shù)值不穩(wěn)定，在實踐中并未表現(xiàn)出更好的預測精度[26]。

SMOA算法包含兩個需要預先設定的參數(shù)q和tfea。其中，距離項(D(x))q中的參數(shù)q對可行域探索準則Cfea(x)的探索能力有很大影響，經過反復測試，建議對有多個非連通可行子區(qū)域的優(yōu)化問題(如算例1)取q=2，對可行域單連通的優(yōu)化問題(如算例2)取q=1。鑒于Cfea(x)準則中的距離項能迫使新試驗點彼此遠離，進而在可行域中分布均勻，此處將用于辨識可行域的可行試驗點的數(shù)量設定為tfea=10(足以刻畫出可行域的輪廓)。此外，本文選取有20個隨機初始點的擬牛頓法BFGS算法優(yōu)化Cfea(x)準則，以保證新試驗點選取的準確性；在改進近似PF時，選取種群大小為100、進化代數(shù)為200的U-NSGA-Ⅲ算法生成Cpf(x)準則的備選試驗點集。

在比較SMOA、KEMOCO以及U-NSGA-Ⅲ算法的優(yōu)化性能時，本文以試驗的總次數(shù)T作為SMOA和KEMOCO算法的終止條件，取U-NSGA-Ⅲ算法的種群大小為100、進化代數(shù)為200，并選取相對超體積(Relative HV，RHV)、Epsilon和R2指標評估算法的優(yōu)化性能，用這些指標分別描述近似PF與真實PF之間的相對誤差、最大最小距離和平均距離[27]。

3.3 優(yōu)化結果分析

KEMOCO算法直接以EHVI與PoF的乘積作為近似PF的改進準則。本文SMOA算法首先以同時最大化EHVI和PoF為目標生成具有Pareto最優(yōu)性的備選點集Ppf，再從中選取新試驗點，不僅提高了新試驗點選取的目的性，還使其兼具了改進最優(yōu)解和刻畫可行域邊界的能力，保證了算法的優(yōu)化效率和精度。通過兩個數(shù)值算例，在45組隨機初始試驗設計下，用SMOA和KEMOCO算法的PF改進準則改進可行域探索階段得到的近似PF，并與U-NSGA-Ⅲ算法的優(yōu)化結果作比較，得到的解集質量指標數(shù)據(jù)如表1所示。表中，算例1的真實PF是通過比較細密網(wǎng)格點上目標和約束函數(shù)的取值得到的；算例2的自變量維度高，難以用網(wǎng)格點法確定其真實PF，為此將由U-NSGA-Ⅲ算法(種群大小取200，進化代數(shù)取500)求得的“高精度”PF當作真實的PF。

表1 不同優(yōu)化算法求得的近似PF的質量指標數(shù)據(jù)(45組)

由表1的指標數(shù)據(jù)可以看出：用SMOA算法求解算例1和算例2，得到的解集質量指標數(shù)據(jù)明顯優(yōu)于KEMOCO算法的相應數(shù)據(jù)(差異性檢驗結果如表2)，且SMOA算法優(yōu)化結果的標準差更小。說明本文所提近似PF改進準則(式(11))，兼具了改進目標函數(shù)和刻畫可行域邊界的功能，具有更高的優(yōu)化精確和穩(wěn)定性。算例1的優(yōu)化結果數(shù)據(jù)還顯示SMOA算法和U-NSGA-Ⅲ算法的優(yōu)化精度在同一數(shù)量級，證明SMOA算法具有類似于多目標進化算法的優(yōu)化能力但優(yōu)化效率更高(SMOA算法的試驗次數(shù)是T=80，遠小于U-NSGA-Ⅲ算法的試驗次數(shù)T=20 000)。

表2 SMOA算法與KEMOCO算法優(yōu)化效果差異性檢驗(配對t檢驗，差異性水平0.05)

3.4 算法收斂性評估

當真實PF未知時，解集質量指標難以評價代理優(yōu)化算法的收斂性。本文用條件模擬方法分析近似PF的不確定性，并將此不確定性作為算法收斂性的評價依據(jù)。用SMOA算法求解算例1，并用條件模擬方法度量近似PF的不確定性，結果如圖3～圖5所示。圖5中圓點對應真實的PF，菱形表示隨機集Vorob’ev期望的PF，上三角形表示SMOA算法在優(yōu)化過程中添加的可行試驗數(shù)據(jù)(它們的連線對應于求得的近似PF)；陰影區(qū)域則反映了CPFs的不確定性(y屬于對稱差Qβ*ΔY的概率)。

由圖3可知：在可行域探索階段選取的前10個試驗點中，有4個是可行的且它們成功地找到了算例1的全部3個可行子區(qū)域；但此時CPFs的不確定程度高(圖中陰影區(qū)域大)、SMOA算法求得的近似PF與隨機集Vorob’ev期望的PF相差大，算法未收斂。用加點準則Cfea(x)繼續(xù)選取10個新試驗點(如圖4)，此時CPFs的不確定性已較小(圖中陰影區(qū)域較小)，Vorob’ev期望的PF接近于真實的PF，說明當前試驗數(shù)據(jù)包含了PF的較多信息；但SMOA算法求得的近似PF與Vorob’ev期望的PF差異仍然較大，算法尚未收斂。用SMOA算法繼續(xù)添加新試驗點直到總試驗次數(shù)T=80(如圖5)，這時圖中的陰影區(qū)域已經非常小，且真實的PF、Vorob’ev期望的PF、近似PF三者幾乎重合，SMOA算法收斂。

表征近似PF不確定性的Vorob’ev期望和方差的超體積(分別記為VE和VD)，也可由條件模擬方法(算法2)得到。受變異系數(shù)概念的啟發(fā)，在此用VD/VE作為定量分析近似PF不確定性的指標，并與RHV、Epsilon、R2等作比較，結果如圖6所示。可知：①VD/VE與RHV、Epsilon、R2等解集質量指標的總體變化趨勢一致，表明近似PF的不確定性確實反映了代理優(yōu)化算法的收斂性；②VD/VE指標的靈敏度略有不足：VE和VD數(shù)據(jù)的獲取涉及隨機抽樣(用以生成CPFs)，由此產生的隨機誤差導致VD/VE指標在優(yōu)化進程的后期不能迅速趨于零；③用條件模擬方法間接分析代理優(yōu)化算法的收斂性，即便真實PF未知，也能從定性、定量兩方面評估優(yōu)化進程，是一種新穎、可行的算法收斂性評估方法。

4 結束語

針對包含黑箱約束的多目標優(yōu)化設計問題，本文提出一種有效辨識可行域和提高優(yōu)化效率的SMOA算法。該算法的可行域探索準則包含了考慮試驗點之間距離的項，辨識可行域的能力強；近似PF改進準則從具有Pareto最優(yōu)性的點集中選取新試驗點，選點目的性強且兼具刻畫可行域邊界，保證了算法的優(yōu)化效率和精度；此外，用隨機模擬方法間接評估算法的收斂性，避免了真實PF未知情況下解集質量指標難以評估算法收斂性的問題。最后，通過與KEMOCO及U-NSGA-Ⅲ算法作比較，驗證了SMOA算法的有效性和高效性。由于條件模擬涉及大量模擬點的選取與比較，其計算量較大。如何進一步提高算法收斂性評估中條件模擬的效率以使其適用于高維問題，以及如何在算法的優(yōu)化迭代中同時選取多個新試驗點以充分利用仿真模型的并行計算能力，是需要進一步研究的問題。