馮艷紅 王改革 李明亮 李 晰
組合優(yōu)化問題(Combinatorial Optimization Pro-blem, COP)的求解一直是運籌學領域重要的研究主題之一.典型的組合優(yōu)化問題包括:作業(yè)車間調(diào)度問題(Job Shop Scheduling Problem, JSSP)[1-2],圖著色問題(Graph Coloring Problem, GCP)[3],車輛路徑規(guī)劃問題(Vehicle Routing Problem, VRP)[4],旅行商問題(Travelling Salesman Problem, TSP)[5],多維背包問題(Multidimensional Knapsack Problem, MKP)[6].
多需求多維背包問題(Multidemand MKP, MDMKP)[7]是MKP的一個重要擴展,也被證明是一類復雜的NP-hard問題.相比MKP,MDMKP包含一組與背包約束沖突的需求約束,是一類更難求解的多約束0-1規(guī)劃問題.現(xiàn)實生活中的很多問題可以建模為MDMKP,包括設備配置、選址方案設計、投資預算、證券選擇等,但是,MDMKP并未引起廣大學者更多的關(guān)注.
目前,針對MDMKP問題的求解分為兩類:確切性算法(Exact Algorithm)和啟發(fā)式算法(Heuristic Algorithm).Cappanera等[7]提出嵌套禁忌搜索啟發(fā)式算法,在標準禁忌搜索算法中增加擾動策略.Gortázar等[8]介紹一種黑盒分散搜索算法,特別提出一種靜態(tài)罰函數(shù)法,處理MDMKP問題的約束.Arntzen等[9]提出交替控制樹搜索框架,交替使用確切性算法和啟發(fā)式算法求解問題.Lai等[10]提出高效的基于兩階段解的禁忌搜索算法.Al-Shihabi[11]提出基于核思想的優(yōu)化框架,更新28個已知最優(yōu)解.因此,尋找有效的方法解決MDMKP問題是有意義的.
群智能算法[12-13]是利用群體智能高效求解優(yōu)化問題的方法,具有不依賴梯度信息、對待求解問題無需連續(xù)、可導等特點,已被廣泛應用于求解各類優(yōu)化問題,特別是組合優(yōu)化問題.受蛾子趨光性和萊維飛行的啟發(fā),Wang[14]提出蛾子搜索算法(Moth Search, MS),基于萊維飛行算子(Lévy Flight Ope-rator, LFO)和直接飛行算子(Straight Flight Ope- rator, SFO)實現(xiàn)種群中個體的更新.目前,MS已在數(shù)值優(yōu)化、約束優(yōu)化、多目標優(yōu)化、組合優(yōu)化等諸多領域被廣泛應用[15-16].因此,將MS應用到MDMKP問題,無疑為該問題的求解提供一種新的方法.
因此,本文提出改進蛾子搜索算法(Modified Moth Search, MMS),用于求解MDMKP問題.主要思想是:引入突變概率,設計突變直接飛行算子(Mutation SFO, MSFO),增強MMS的搜索能力;擴大參數(shù)α的調(diào)節(jié)范圍,設計自適應萊維飛行算子(Adaptive LFO, ALFO),幫助算法跳出局部最優(yōu);增加均勻變異算子(Uniform Mutation Operator, UMO),增強種群多樣性.
本文方法的創(chuàng)新之處在于:1)采用自適應ε約束處理法解決MDMKP中兩組相互沖突的約束條件,從而易于采用群智能算法求解;2)對自適應ε約束處理法的參數(shù)使用正交實驗方法進行驗證,得到一組針對MDMKP問題的最適合參數(shù);3)基于空間轉(zhuǎn)換思想,設計應用連續(xù)型群智能算法求解離散優(yōu)化問題的框架.
給定具有n個物品的集合V={1,2,…,n},MDMKP問題[10]可描述為:從集合V中選取一組滿足所有問題約束的物品組合并使所選物品總價值最大.數(shù)學描述如下:
xj∈{0,1}, ?j∈{1,2,…,n}.
(1)
其中,cj表示物品j的價值;xj表示物品是否被選中,xj=1表示物品被選擇,xj=0表示物品未被選擇;ai,j表示物品j對第i種資源的消耗量;bi表示第i種資源總量.特別地,式(1)的小于等于約束稱為背包約束,大于等于約束稱為需求約束.
不失一般性,有如下假設:
在MS中,假設待求解問題對應d維解空間,隨機初始化具有N個個體的種群
P={X1,X2,…,XN}
分布于此解空間,第i個蛾子個體在解空間的位置向量
Xi=[xi,1,xi,2,…,xi,d]T.
顯然,每個個體位置對應目標問題的一個潛在解.MS按照適應度值升序排序后將種群等分為子種群1(SP1)和子種群2(SP2).子種群1中的個體依據(jù)LFO進行位置更新,子種群2中的個體依據(jù)SFO進行位置更新.
1)LFO.SP1中的個體i位置更新公式如下:
2)SFO.SP2中的個體i位置更新公式如下:
考慮最小優(yōu)化問題,MS流程如算法1所示.
算法1MS
輸入初始化所有參數(shù),適應度函數(shù)f(X),
隨機初始化N個個體,確定問題維度D,
最大迭代次數(shù)T
輸出最優(yōu)蛾子個體的位置,對應適應度函數(shù)值
step 1 WHILE未達到最大迭代次數(shù)或終止條件DO.
step 2 計算個體適應度值并進行升序排列,記錄最優(yōu)個體信息.
step 3 種群劃分.適應度值較小的N/2個體組成SP1,余下個體組成SP2.
step 4 應用LFO更新SP1.
step 5 應用SFO更新SP2.
step 6 合并兩個新生成子種群為一個種群.
step 7 判斷終止條件是否滿足.若未滿足,t=t+1,轉(zhuǎn)step 2.
一般而言,MS能夠求得全局最優(yōu)解.但是,面對一些復雜的約束優(yōu)化問題,特別是問題的可行解介于可行域和不可行域邊界之時,MS可能會陷入局部最優(yōu).
MDMKP問題是復雜的離散約束優(yōu)化問題,問題的可行域是離散的,如何使搜索能穿梭于各個可行區(qū)域間,是該問題的主要挑戰(zhàn).多種群可能是解決該問題的方案之一.MS是基于兩個子種群進化的群智能算法,為此,提出如下改進策略.
1)自適應萊維飛行算子(ALFO).在LFO中,參數(shù)α控制蛾子飛行的步長,該步長只與當前迭代次數(shù)有關(guān),忽略總迭代次數(shù).受文獻[17]啟發(fā),將參數(shù)α修改為
2)突變直接飛行算子(MSFO).在SFO中,參數(shù)λ是調(diào)節(jié)因子,控制算法的收斂速度和種群多樣性.在MS中,參數(shù)λ∈U(0,1),為了使蛾子朝著當前最優(yōu)個體飛去,使用一個與最大迭代次數(shù)和當前迭代次數(shù)相關(guān)的突變率λ對蛾子個體的局部搜索進行擾動.突變率
顯然:過高的突變率增加算法在搜索空間中搜索更多區(qū)域的概率,防止種群收斂到任何最優(yōu)解;過小的突變率可能導致種群早熟收斂到局部最優(yōu).參數(shù)λ從0.9線性遞減到0,突變率增加算法的搜索能力.
3)均勻變異算子(UMO).為了增加種群多樣性,對種群中每個個體的每個基因,以較小的概率pm進行變異,即更新當前基因為隨機生成的(0,1)內(nèi)的隨機實數(shù).變異公式
xi,j=rand(),rand()≤pm,
其中,xi,j表示個體i的第j個分量,rand()是在(0,1)內(nèi)服從均勻分布的實數(shù).
將上述3種策略應用到MS框架中,得到MMS.
MS的提出是用來求解連續(xù)優(yōu)化問題,LFO和SFO在連續(xù)空間上的計算是封閉的,而MDMKP問題
的每個潛在解對應一個二進制向量,因此,直接應用MS求解MDMKP問題是不可行的.
本文基于空間轉(zhuǎn)換機制[18],采用映射的方法,將搜索空間的實值向量轉(zhuǎn)換為問題空間的二進制向量.假設
Ω={(x1,x2,…,xn)|Li≤xi≤Ui,1≤i≤n}
為n維搜索空間,
Φ={(y1,y2,…,yn)|yi∈{0,1},1≤i≤n}
為由所有可行解和不可行解組成的問題空間,則可構(gòu)建如下映射:
m∶Ω→Φ.
(2)
本文選取一種簡單有效的映射函數(shù)m(x)=x,映射方法
種群中的任何個體表示為一個二元組〈X,Y〉,其中,X構(gòu)成問題的搜索空間,Y形成問題的候選解空間.在本文中,X∈[-5,5]n,則Y∈{0,1}n.
下面以n=10為例,具體實值向量映射為二進制向量的映射過程如表1所示.
表1 實值向量映射為二進制向量Table 1 Mapping process from real valued vectors to binary vectors
由式(1)可知,MDMKP包含m個小于等于不等式約束,q個大于等于不等式約束,顯然兩組約束是對立的.采用傳統(tǒng)的約束處理技術(shù)如罰函數(shù)法,基于特定問題的不可行解修復法等很難對該約束進行處理.
因此,本文采用自適應ε約束處理法[19-21].該方法的核心思想是:基于人為設定的ε值,對個體的違約程度進行排序,從而將個體劃分為不同的區(qū)域.在不同的區(qū)域中,可行解和不可行解采用不同的評價準則.此方法有效利用不可行域中目標函數(shù)值較好的不可行解信息,具有更好的收斂性能.ε水平控制方法如下:
(3)
其中,Xθ表示種群中個體按照約束違反程度升序排序后的第θ個個體,Tc表示代數(shù)控制參數(shù),cp在區(qū)間[2,10]中取值.
設兩個個體X1、X2的目標函數(shù)值分別為f1、f2,約束違約度分別為φ1、φ2,則基于ε水平的個體比較規(guī)則如下:
具體算法框架如算法2所示.
算法2自適應約束處理法
輸入初始化所有參數(shù),初始種群,當前代數(shù)t
輸出第t代的ε值
step 1 對初始種群中所有個體的違約目標進行降序排列.
step 2 計算違約目標值第θ·N大的值,記錄ε(0).
step 3 如果t≥Tc返回0值.
step 4 否則
step 4.1 計算當前代所有個體違約目標值.
step 4.2 如果違約目標值為0的個體個數(shù)大于等于ap·N, 返回0值.
step 4.3 否則基于式(3)的第1項計算ε(t).
step 5 輸出ε(t).
融合2.1節(jié)~2.4節(jié)的內(nèi)容,給出MMS求解MDMKP問題的流程,如算法3所示.
算法3MMS求解MDMKP框架
輸入初始化所有算法參數(shù),最大迭代次數(shù)T,
隨機初始化N個個體,適應度函數(shù)f(X),
問題實例I
輸出最優(yōu)蛾子個體的位置,
對應的適應度函數(shù)值
step 1 WHILE未達到最大迭代次數(shù)或終止條件DO.
step 2 應用式(2)對個體編碼為二元組〈X,Y〉.
step 3 基于自適應ε約束處理法對個體進行排序,記錄最優(yōu)個體信息.
step 4 種群劃分.較優(yōu)的N/2個體組成SP1,余下個體組成SP2.
step 5 應用ALFO更新SP1.
step 6 應用MSFO更新SP2.
step 7 應用UMO更新SP1和SP2.
step 8 合并兩個新生成子種群為一個種群.
step 9 判斷終止條件是否滿足.若未滿足,t=t+1,轉(zhuǎn)step 2.
應用MMS求解MDMKP問題主要由種群初始化、自適應ε約束處理、ALFO、MSFO、UMO組成.設種群規(guī)模為N,測試用例中物品個數(shù)為D,則種群初始化時間復雜度為O(ND).
在迭代過程中,對種群中所有個體按照自適應ε約束處理,時間復雜度為O(N),SP1應用ALFO進行位置更新,時間復雜度為O(DN/2),SP2應用MSFO進行位置更新,時間復雜度為O(ND/2),UMO的時間復雜度為O(ND).
因此,每次迭代過程中算法的時間復雜度為
與MS一致.
本節(jié)通過實驗評價MMS的綜合性能.實驗包括如下3部分:1)考察自適應ε約束處理法參數(shù)對算法性能的影響;2)考察3種策略對算法性能的影響;3)對比MMS與MS,驗證MMS的優(yōu)越性.
為了驗證MMS的優(yōu)化性能,選取文獻[10]中的兩組測試用例與MS進行對比.每組包含48個測試用例,表示為n-m-q-g-id.例如100-5-2-0-0表示該測試用例包含100個物品、5個背包約束、2個需求約束、第0組的第0個實例.第1組中n=100,第2組中n=250,兩組實例其余參數(shù)設置如下:
m=5,10,q=2,5,10,g=0,1,id=0,1,2,3,4,5.
與g=0組不同的是,g=1的組中項目對資源的消耗值可以為負值.具體實驗數(shù)據(jù)可在https://grafo.etsii.urjc.es/optsicom/binaryss.html#instances下載.
為了確保實驗的公平性,依據(jù)文獻[14],MS和MMS在所有測試用例上設置相同參數(shù).種群規(guī)模為50,迭代次數(shù)為1 000.算法參數(shù)設置如下:φ=0.618,MAXSTEP=1.0,λ取值為在(0,1)區(qū)間上服從均勻分布的隨機數(shù).此外,為了消除隨機性對算法性能的評估,2種算法在所有96個測試用例上均進行30次獨立實驗.
本文采用最好目標值(Best)、平均目標值(Avg.)、最差目標值(Worst)這3個基本性能指標.同時,為了對比最好目標值與已知最優(yōu)值(BKV)之間的相對偏差,引入Gap值考察算法逼近BKV的程度:
本文所有實驗的硬件環(huán)境為:Intel(R) Core(TM) i5-2415M CPU 2.30 GHz、內(nèi)存4.00 GB,操作系統(tǒng)為Microsoft Windows 7.編程語言為Python,解釋器為Python3.9.
自適應ε約束處理法包含3個主要參數(shù):θ,cp,Tc.一般而言,算法性能對參數(shù)設置具有敏感性.因此,基于正交實驗設計(Orthogonal Experimental Design, OED)[22]尋找一組較優(yōu)的參數(shù)組合.實驗包含3個因子,每個因子包含3個水平,設計具有9種組合的正交表.
選取100-5-2-0-0作為測試用例,每組參數(shù)組合獨立運行30次,取平均值作為實驗對比值,則不同參數(shù)值的組合如表2所示.正交實驗參數(shù)及實驗結(jié)果如表3所示.因子分析如表4所示,表中黑體數(shù)字表示最優(yōu)值.
表2 不同參數(shù)值組合Table 2 Combination of different parameter values
表3 正交實驗參數(shù)及實驗結(jié)果Table 3 Parameters of orthogonal experiment and experimental results
表4 基于正交設計方法的因子分析Table 4 Factor analysis based on orthogonal design method
從表4的標準差可以發(fā)現(xiàn):Tc是3個參數(shù)中最重要的參數(shù),需要合理選擇.基于上述OED,分析最適合的參數(shù)組合為θ=0.2,cp=5,Tc=500.下面實驗將采用這組參數(shù).
為了驗證3種策略對MMS性能的影響,對比研究MS,包括MS應用ALFO(簡記為MS-I)、MS應用MSFO(簡記為MS-II)、MS應用UMO(簡記為MS-III)、融合3種策略的MMS.選取100-5-2-0-id和100-5-2-1-id形式的12個測試用例,每種算法對每個測試用例均獨自運行30次,實驗結(jié)果如表5所示,表中黑體數(shù)字表示5種算法求解同個測試用例的最優(yōu)值,total欄表示4種改進算法與MS相比取得更優(yōu)值的個數(shù).
表5 不同改進策略的實驗結(jié)果Table 5 Experimental results of different improvement strategies
由表5可知,ALFO和UMO均能在可行域復雜的情況下,充分挖掘搜索空間,對所有的12個測試用例的求解質(zhì)量均有所提高.MSFO對9個測試用例的最優(yōu)值有更高的精度,說明MSFO可以有效增加種群多樣性.同時,融合3種策略的MMS的最優(yōu)值和Gap值都優(yōu)于單一改進策略的MS,尋優(yōu)性能顯著,表明MMS能夠充分發(fā)揮不同改進策略的優(yōu)點.
為了從統(tǒng)計意義上說明4種改進算法與MS是否有顯著性差異,采用Wilcoxon秩和檢驗驗證MS的12個測試用例的Best值在α= 0.05的顯著性水平下與其它算法的顯著性差異.其中:p-value>0.05,表明接受原假設H0,即認為兩種優(yōu)化算法之間無顯著性差異,算法尋優(yōu)能力相當;p-value<0.05,拒絕原假設H0,即兩種算法存在顯著性差異.由此可見,MS-I、MS-III、MMS與MS存在顯著性差異,MS-II與MS無顯著性差異.由表5還可以看出,1)MMS在幾乎所有的12個測試用例中,取得的Best值均超過其它四種算法.說明三種策略的組合有益于提高算法求解性能.2)MS-I、MS-III在所有的12個測試用例上均優(yōu)于MS,MS-II在絕大多數(shù)測試用例上優(yōu)于MS,說明三種改進策略均能夠提高MS的優(yōu)化性能.3)相比MS-I和MS-III,MS-II的改進效果不是很明顯,可能的原因是過高的突變率增加在搜索空間搜索更多可行域的概率,延緩算法收斂速度,太小的突變率導致算法早熟收斂.
為了驗證MMS對于MDMKP的求解能力,與MS對96個測試用例進行對比,實驗結(jié)果如表6~表9所示,表中黑體數(shù)字表示最優(yōu)值.
表6 MS和MMS求解100-5-x-x-x測試用例的實驗結(jié)果Table 6 Experimental results of MS and MMS on test instances(100-5-x-x-x)
表7 MS和MMS求解100-10-x-x-x測試用例的實驗結(jié)果Table 7 Experimental results of MS and MMS on test instances(100-10-x-x-x)
表8 MS和MMS求解250-5-x-x-x測試用例的實驗結(jié)果Table 8 Experimental results of MS and MMS on test instances(250-5-x-x-x)
續(xù)表8Table 8(Continued)
表9 MS和MMS求解250-10-x-x-x測試用例的實驗結(jié)果Table 9 Experimental results of MS and MMS on test instances (250-10-x-x-x)
由表6可以看出:1)就求得最優(yōu)結(jié)果的總個數(shù)而言,MMS的Best指標、Avg.指標、Worst指標分別有20、24、20個優(yōu)于MS;2)就4種評價指標平均值而言,MMS均優(yōu)于MS.
由表7可以看出,除了Worst指標以外,MMS劣于MS.由此可以推斷,MMS可能對于測試用例的屬性具有敏感性.
由表8可以看出:1)就求得最優(yōu)結(jié)果的總個數(shù)而言,MMS的Best指標、Avg.指標、Worst指標分別有20、8、6個優(yōu)于MS;2)就4種評價指標平均值而言,MMS均優(yōu)于MS.
由表9可以看出,MMS的Avg.指標有12個優(yōu)于MS,Worst指標有10個優(yōu)于MS,Best指標有4個劣于MS.
從表6~表9的對比結(jié)果可以看出,隨著問題維數(shù)的增大,MMS優(yōu)于MS的求解結(jié)果更加明顯,說明改進策略在搜索復雜的可行域時優(yōu)勢突顯.
為了圖形化顯示MS和MMS求解n=100,250兩組實例的Gap指標分布情況,圖1給出小提琴圖,圖中紅點表示中位數(shù),黑色盒型表示最集中的50%區(qū)域.從圖1可以看出:1)兩種算法均存在較明顯的離散值(上側(cè)的須和下側(cè)的須較長);2)MS求解兩組測試用例的Gap值分布相對集中,MMS的Gap值分布不太均勻.說明相比MS,MMS能夠求得更優(yōu)解,但算法的穩(wěn)定性需要進一步提高.
(a)n=100
(b)n=250圖1 2種算法求解兩組測試用例的Gap值分布Fig.1 Gap value distribution of 2 algorithms on 2 sets of instances
在MDMKP的每一組測試用例中,n-m-q-g-0測試用例中背包中物品價值系數(shù)全部為正整數(shù),n-m-q-g-1測試用例中背包中物品價值系數(shù)包含負整數(shù).為了說明算法的求解性能是否與測試用例的特征具有相關(guān)性,對MMS求得的兩組測試用例(n=100,250)的Gap值繪制散點圖,如圖2所示.
(a)n=100
(b)n=250圖2 96個測試用例的Gap值分布散點圖Fig.2 Scatter plot of Gap value distribution of 96 instances
從圖2可以得出如下結(jié)論.
1)總體而言,Gap值的分布與測試用例的分組基本一致,即同組的6個測試用例Gap值可以作為一組,其中100-10-5-1、250-5-2-0、250-5-5-0、250-10- 10-1四組測試用例尤其明顯,說明測試用例的屬性對算法的求解性能具有一定影響.
2)對于背包約束個數(shù)和需求約束個數(shù)相同的兩組實例,即g=0,1,組別為0的正系數(shù)測試用例求解效果明顯優(yōu)于組別為1的兼具正負系數(shù)的測試用例,前者整體處于后者的底端,說明收益可為正負值,增加求解難度.
3)背包約束個數(shù)的增加使可行域的搜索變得困難,解的質(zhì)量變差,Gap值相對較大.然而,需求約束個數(shù)增加,求解效果影響不明顯.可能的原因是約束條件比較寬泛.
MDMKP問題包含兩組對立的約束條件,使得該問題的求解非常困難.如果采用合適的方法,處理好約束,那么有可能降低問題的難度.本文提出基于ε約束處理法的改進蛾子搜索算法(MMS),用于求解該問題,基于96個標準測試函數(shù)的實驗表明,MMS比MS具有更高的尋優(yōu)精度.實驗分析表明3種改進策略對算法性能具有各自的影響:自適應萊維飛行算子(ALFO)擴大蛾子飛行步長的調(diào)節(jié)范圍,對個體進行局部擾動,以更大概率跳出局部最優(yōu);變異直接飛行算子(MSFO)引入合適的突變率,防止種群過早收斂,增加在搜索空間中搜索更多區(qū)域的概率;均勻變異算子(UMO)以較小的概率更新種群中部分個體的基因,增加種群多樣性.同時,本文還分析ε約束處理法的3個主要參數(shù)對算法性能的影響.
本研究工作豐富現(xiàn)存的求解MDMKP問題的技術(shù),下一步將考慮應用不同約束處理方法,如可行性法則、多目標優(yōu)化法等.此外,該框架也可應用于求解其它的多約束組合優(yōu)化問題,如凸可分背包問題、動態(tài)背包問題、多目標背包問題等.盡管MMS優(yōu)于MS,但解的質(zhì)量還需要進一步提高.因此,提出更好的策略在復雜的搜索空間中尋求高質(zhì)量解將是下一步的研究工作.