黃光球,陸秋琴
西安建筑科技大學(xué) 管理學(xué)院,西安710055
工程中存在一類非線性優(yōu)化問題,此類問題的數(shù)學(xué)表達(dá)式中存在不可導(dǎo)的數(shù)學(xué)表達(dá)式?;谔荻鹊膬?yōu)化方法無法求解此類問題。為了求解此類問題,人們提出了啟發(fā)式搜索方法[1],群智能優(yōu)化算法[2]就是其中的一種。因群智能優(yōu)化算法對工程優(yōu)化問題沒有特殊要求,故具有廣泛的適應(yīng)性。近幾年,群智能優(yōu)化算法得到快速發(fā)展。
群智能優(yōu)化算法是依據(jù)特定的生物進(jìn)化場景而構(gòu)建出來的,其算子依據(jù)個體間的相互作用而構(gòu)建,其邏輯結(jié)構(gòu)是依據(jù)生物進(jìn)化場景的內(nèi)涵而構(gòu)建[3]。早期的群智能算法典型代表有粒子群算法[4-7]、蟻群算法[8-9]、蜂群算法[10]、布谷鳥算法[11]、蝙蝠算法[12]、遺傳算法[13]等。目前,很多新型群智能算法先后被提出,如細(xì)菌覓食優(yōu)化、蛙跳算法、猴群算法、螢火蟲算法、鴿群算法、果蠅算法和頭腦風(fēng)暴算法等[14]。然而,目前所提出的群智能優(yōu)化算法的共同缺陷包括如下幾方面:(1)算法所依賴生物進(jìn)化場景均非常簡單,因而依據(jù)此類場景設(shè)計出來的群智能算法的邏輯結(jié)構(gòu)非常簡單;(2)所能開發(fā)出的算子非常少,從而導(dǎo)致個體之間的信息交換不充分,搜索易限于局部陷阱;(3)算法的全局收斂性難以保證。
隨著互聯(lián)網(wǎng)的高速發(fā)展,社會已進(jìn)入大數(shù)據(jù)時代,人們所遇到的問題變得越來越復(fù)雜,故需要提出高智能性的算法來解決這些復(fù)雜問題。對于群智能算法來說,如何增加群智能算法的智能特征,顯然是需要解決的關(guān)鍵問題。解決該問題的一種策略是:精心選擇出某種具有特殊性質(zhì)的生物進(jìn)化場景,巧妙解決目前群智能算法所存在問題,即可設(shè)計出高智能性的群智能算法。
若一類生物進(jìn)化場景具有豐富的演化內(nèi)涵,其中個體之間的相互作用關(guān)系十分豐富,則依據(jù)此類生物進(jìn)化場景就可以設(shè)計出具有算子眾多、邏輯結(jié)構(gòu)清晰的高智能性群智能優(yōu)化算法。能夠跨三個物種實(shí)施多級傳播的包蟲傳染病攻擊家犬、牛羊和人類的過程,就是這樣一個場景。
包蟲病是棘球蚴寄生在人體所致的一種人獸共患寄生蟲病。我國有囊型和泡型包蟲病兩種,分別由棘球蚴和泡球蚴寄生在人體組織器官中所致。囊型包蟲病呈世界性分布,畜牧業(yè)發(fā)達(dá)的國家和地區(qū)多見。我國包蟲病高發(fā)流行區(qū)主要集中在高山草甸地區(qū)及氣候寒冷、干旱少雨的牧區(qū)及半農(nóng)半牧區(qū),以新疆、青海、甘肅、寧夏、西藏、內(nèi)蒙古、陜西、河北、山西和四川等地較為嚴(yán)重。包蟲來源于動物的排泄物,沒有什么有效的治療方式,在西藏被視為西藏第一癌癥。
由于包蟲病可同時危害人類和家畜養(yǎng)殖業(yè),近年來全世界研究人員對該傳染病給予了高度關(guān)注,主要研究成果表現(xiàn)為如下幾方面:
(1)包蟲病流行病學(xué)分析。為了了解包蟲病的流行病學(xué)特征,通過對包蟲病病例進(jìn)行流行病學(xué)特征調(diào)查研究,探討病例可能的感染來源,為制定包蟲病防控措施提供依據(jù)[15]。
(2)包蟲病免疫學(xué)研究。通過對家畜和人的包蟲病免疫機(jī)理進(jìn)行研究,為包蟲病的防治和新型疫苗的生產(chǎn)提供科學(xué)參考[16]。
(3)包蟲病診斷方法研究。提出包蟲病的診斷方法,為包蟲病的治療提供科學(xué)、快速和準(zhǔn)確的依據(jù)[17]。
(4)包蟲病的治療方法研究。從手術(shù)和藥物等角度探討包蟲病的各種治療方法[18]。
包蟲病具有跨多物種傳播特征,能夠?qū)θ祟愒斐珊艽笸纯啵谶@悲慘的自然現(xiàn)象后面,卻隱藏著一個性能優(yōu)良的群智能優(yōu)化算法,本文稱其為包蟲傳染病優(yōu)化算法(hydatid disease optimization,HDO)。本文著重解決如下五個問題:
(1)如何將包蟲傳染病模型轉(zhuǎn)化為能求解復(fù)雜優(yōu)化問題的HDO 算法。
(2)如何使得HDO 算法中的算子能充分反映動物個體之間的相互作用關(guān)系,以便體現(xiàn)包蟲傳染病模型的思想。
(3)如何證明HDO 算法的全局收斂性。
(4)如何確定HDO 算法的最佳參數(shù)設(shè)置。
(5)如何用HDO 算法求解關(guān)聯(lián)區(qū)域VOCs(volatile organic compounds)聯(lián)防聯(lián)控最優(yōu)減排優(yōu)化問題。
設(shè)要求解的優(yōu)化問題為:
式中,Rn是n維歐氏空間;X=(x1,x2,…,xn)是一個n維決策向量;H為搜索空間;F(X)為目標(biāo)函數(shù);gi(X)≥0 為第i個不等式約束條件,I為不等式約束條件個數(shù)。目標(biāo)函數(shù)F(X)和約束條件gi(X)沒有限制。
在一個草原牧區(qū)生活有一群牧民,其數(shù)量有N3名,其編號分別是1,2,…,N3;牧民們以飼養(yǎng)藏羊為生,藏羊的數(shù)量有N2只,其編號分別是1,2,…,N2;與此同時,牧民們還飼養(yǎng)了一群狗用于控制和保護(hù)羊群,狗的數(shù)量有N1只,其編號分別是1,2,…,N1。一條狗、一只羊或一名牧民統(tǒng)稱為一個個體,每個個體均由n個特征(器官)來表征。令u表示個體類型,u=1表示狗類,u=2 表示羊類,u=3 表示牧民類,即對類型為u的個體i來說,其表征特征為。
該草原牧區(qū)有一種稱為棘球絳蟲的寄生蟲病在狗群中流行。棘球絳蟲寄生于狗的小腸中,蟲卵隨狗的糞便排出。由于狗隨羊群到處游動,其排泄的糞便廣泛地污染了牧區(qū)的水源、土壤和草場。棘球絳蟲蟲卵在自然界有非常強(qiáng)的適應(yīng)能力,它能在自然狀態(tài)下保持持續(xù)感染力。羊在吃草、喝水的過程中,就會把附著在草葉上或水中的蟲卵一同吃入體內(nèi)。生活在羊胃中的蟲卵會在胃酸的作用下大量繁殖,一部分蟲卵會在羊的體內(nèi)發(fā)育成包囊,另一部分蟲卵會隨羊的糞便又排泄到牧區(qū)的水源、土壤和草場中。
牧民不但與羊群有密切接觸,而且與狗也有更加密切的接觸。牧民與羊群密切接觸或吃羊肉、喝羊血時,就會把蟲卵吃入口中;狗有舔舐其肛門和身體的習(xí)慣,會把蟲卵從肛門帶到其身體上,牧民與狗密切接觸時,就會把蟲卵吃入口中;牧民喝被污染的水或食物時,更會將蟲卵吃入口中。
生活在人胃中的棘球絳蟲蟲卵在胃酸的作用下加速繁殖、數(shù)量大增。一部分蟲卵會隨血液鉆入人體的某些器官內(nèi),如肝、肺、腦中,然后發(fā)育成巨型包囊;另一部分蟲卵會隨人的糞便又排泄到牧區(qū)的水源、土壤和草場中。
棘球絳蟲蟲卵的傳播特征是糞口傳播,具有很強(qiáng)的傳染能力,既能夠在同物種類內(nèi)傳播,又能夠跨物種傳播。當(dāng)狗、羊和牧民傳染上包蟲病后,一部分能夠治愈,另一部分則死亡。棘球絳蟲蟲卵攻擊的是狗、羊和牧民的部分特征(器官)。
將上述場景映射到對優(yōu)化問題式(1)全局最優(yōu)解的搜索過程中,含義如下所述。
優(yōu)化問題式(1)的搜索空間與草原牧區(qū)相對應(yīng),該草原牧區(qū)中的一條狗、一只羊雞或一名牧民分別對應(yīng)于優(yōu)化問題式(1)的一個試探解,即對類型為u的Nu個個體所對應(yīng)的試探解集就是且類型為u的個體i的特征j與試探解的變量相對應(yīng)。對于優(yōu)化問題式(1),類型為u的個體i的適應(yīng)度Fit按下式計算:
在時期t,對于類型為u個體來說,共存在5 種狀態(tài):易感(未感染蟲卵)的個體,其在時期t的占比為Su(t),處于此狀態(tài)的個體用Su表示;已暴露(已感染蟲卵但未發(fā)?。┑膫€體,其在時期t的占比為Eu(t),處于此狀態(tài)的個體用Eu表示;已發(fā)病(感染蟲卵后已發(fā)?。┑膫€體,其在時期t的占比為Iu(t),處于此狀態(tài)的個體用Iu表示;已治愈的個體,其在時期t的占比為Ru(t),處于此狀態(tài)的個體用Ru表示;已死亡的個體,其在時期t的占比為Du(t),處于此狀態(tài)的用Du表示。
為了簡單起見,對上述場景進(jìn)行一些簡化:不考慮個體的自然死亡;當(dāng)一個個體因染上包蟲病而死亡后,立即就有一個新個體添加到該草原牧區(qū)中,從而確保各類個體的總數(shù)為常數(shù);棘球絳蟲蟲卵可以在同一物種內(nèi)傳播;不同物種間的傳播規(guī)律是:狗可以將蟲卵傳播給羊和牧民,羊只能將蟲卵傳播給牧民,牧民不會將蟲卵傳播給狗和羊。上述場景可采用傳染病傳播倉室模型[19-20]來描述,如圖1 所示。根據(jù)圖1,可以寫出其相應(yīng)的動力學(xué)方程組如下:
式中,βu表示包蟲病在類型為u的個體中的傳染率,0 <βu<1;λu表示類型為u的個體的復(fù)原率,0 <λu<1;γu表示類型為u的個體的發(fā)病率,0 <γu<1;δu表示類型為u的個體的治愈率,0 <δu<1;θu表示類型為u的個體的死亡率,0 <θu<1。
時期t,對于類型為u的任意一個個體來說,Su(t)、Eu(t)、Iu(t)、Ru(t)、Du(t)分別表示該個體處于狀態(tài)Su、Eu、Iu、Ru、Du的概率。必須指出,在同一時期,任意一個個體只能處在狀態(tài)Su、Eu、Iu、Ru、Du中的某一個。通常情況下模型式(3)中的參數(shù)βu、γu、δu、λu、θu的取值是隨時間變化的,可將式(3)應(yīng)用到類型為u的任意一個個體上,并將式(3)改寫成如式(4)所示離散遞推形式。
Fig.1 Flow chart of infectious disease transmission in interaction of dogs,sheep and herdsmen圖1 狗、羊與牧民相互作用的傳染病傳播流程圖
Fig.2 Legal types of state transition圖2 合法的狀態(tài)轉(zhuǎn)移類型
在時期t,采用模型式(4)計算類型為u的個體i的易感概率、暴露概率、發(fā)病概率、治愈概率和死亡概率;類型為u的個體i在時期t處于狀態(tài)Su、Eu、Iu、Ru、Du中的哪個狀態(tài),由它們之中的最大者確定。依據(jù)圖1,可以識別出所有合法的狀態(tài)轉(zhuǎn)移類型,如圖2 所示,圖2 中所描述的狀態(tài)轉(zhuǎn)移類型的含義及其所對應(yīng)的算子如表1 所示。
Table 1 Legal state transitions表1 合法狀態(tài)轉(zhuǎn)換
從表1 可知,圖2 所確定的算子有3×9=27 個。較多的算子可提升HDO 算法的智能性。
2.4.1 特征集合生成方法
設(shè)當(dāng)前個體類型u,個體編號為i,個體狀態(tài)s∈{Su,Eu,Iu,Ru,Du},則:
2.4.2 狀態(tài)轉(zhuǎn)移算子設(shè)計方法
對圖2 進(jìn)行分解,可得如圖3 所示的三種狀態(tài)轉(zhuǎn)移,存在下列三種情況:
(1)個體從時期t的狀態(tài)A轉(zhuǎn)移到時期t+1 的狀態(tài)B,如圖3(a)所示,其中A、B∈{Su,Eu,Iu,Ru,Du|u=1,2,3}但A≠B。有如下兩種情形:
Fig.3 Three situations of individual state transition圖3 個體三種狀態(tài)轉(zhuǎn)移情形
情形1當(dāng)A≠Du,B≠Du時,大量的個體狀態(tài)轉(zhuǎn)移都屬于這種情況。例如Su→Eu、Eu→Iu等。為了實(shí)現(xiàn)圖3(a)所示的狀態(tài)轉(zhuǎn)移,可將已處于狀態(tài)B的若干個體的某些特征的狀態(tài)值經(jīng)加工處理后傳給處于狀態(tài)A的當(dāng)前個體的對應(yīng)特征,從而使其具有處于狀態(tài)B的個體的特征。例如,對于狀態(tài)轉(zhuǎn)移Su→Eu,將已處于暴露狀態(tài)(Eu)的若干個體的某些特征的狀態(tài)值經(jīng)加工處理后傳給處于易感狀態(tài)(Su)的當(dāng)前個體,即可使其感染上包蟲病。此相當(dāng)于已感染棘球絳蟲蟲卵但未發(fā)病的暴露個體將其帶有蟲卵的東西傳給當(dāng)前易感者,使其也感染上蟲卵。其他轉(zhuǎn)移的含義可類似解釋。
情形2當(dāng)A=Iu,B=Du時,當(dāng)前個體從時期t的發(fā)病狀態(tài)Iu轉(zhuǎn)移到時期t+1 的死亡狀態(tài)Du,此時可將適應(yīng)度極低的虛弱個體用適應(yīng)度較高的強(qiáng)壯個體來取代,從而實(shí)現(xiàn)虛弱個體的死亡。此情形可隨機(jī)淘汰極差個體。
(2)當(dāng)前個體在時期t處于某個狀態(tài)A時,A∈{Su,Eu,Iu,Ru|u=1,2,3},在時期t+1 沒有發(fā)生狀態(tài)轉(zhuǎn)移,即相當(dāng)于A→A,如圖3(b)所示。圖2 中的每個節(jié)點(diǎn)包含了圖3(b)所示的情形。例如,Su→Su、Eu→Eu、Iu→Iu和Ru→Ru。注意,不存在Du→Du。
生物個體在生存競爭過程中總是盡量使自身強(qiáng)壯,以便更好地生存發(fā)展。為達(dá)到此目的,可將處于相同狀態(tài)的若干個強(qiáng)壯個體的某些特征的狀態(tài)值經(jīng)加工處理后傳給當(dāng)前個體的對應(yīng)特征,使得當(dāng)前個體獲得強(qiáng)壯個體的特征,從而向更好的方向發(fā)展。
(3)在時期t,當(dāng)前個體在處于狀態(tài)C的個體的作用下從狀態(tài)A轉(zhuǎn)移到狀態(tài)B,如圖3(c)所示,其中A∈{S2,S3},B∈{E2,E3},C∈{E1∪I1,E2∪I2}。例如,對于圖3(c)中的節(jié)點(diǎn)S2,有。
設(shè)當(dāng)前個體類型u∈{1,2,3},個體編號為i,則HDO 算法中各個算子的數(shù)學(xué)表達(dá)式如表2 所示。
對于優(yōu)化問題式(1),其生長算子可以描述為:
(1)初始化:①令時期t=0,按第4 章介紹的方法設(shè)置本算法中的參數(shù):演化時期數(shù)G,全局最優(yōu)解誤差ε,個體特征被包蟲病攻擊的最大概率E0、L、N1、N2、N3;②分別初始化個體,u∈{1,2,3}。
(3)計算類型為u的個體i的Su、Eu、Iu、Ru、Du狀態(tài),u∈{1,2,3}。函數(shù)GetSEIRD()用于確定類型為u的個體i將處于Su、Eu、Iu、Ru、Du狀態(tài)中的哪一個狀態(tài)。
Table 2 Mathematical expressions of operators in HDO algorithm表2 HDO 算法中各個算子的數(shù)學(xué)表達(dá)式
(4)執(zhí)行下列操作:
(5)結(jié)束。
函數(shù)GetSEIRD(pS,pE,pI,pR,pD)的定義如下:
FunctionGetSEIRD(pS,pE,pI,pR,pD)//pS、pE、pI、pR、pD分別表示易感狀態(tài)S、暴露狀態(tài)E、發(fā)病狀態(tài)I、治愈狀態(tài)R、死亡狀態(tài)D 的概率
(1)HDO 算法演化過程具有Markov 特性。從Su-Su、Su-Eu、Eu-Eu、Eu-Iu、Iu-Iu、Iu-Ru、Iu-Du、Ru-Ru、Ru-Su的定義知,當(dāng)i、k為種群編號,j=1~n,u∈{1,2,3},時,這些算子均有如下特征:
從式(6)和式(7)可知,時期t+1(新一代)的試探解的產(chǎn)生只與該試探解在時期t(當(dāng)前代)的狀態(tài)有關(guān),而與該試探解在時期t以前是如何演變到時期t所對應(yīng)的狀態(tài)的歷程無關(guān)。因此,依據(jù)Markov 特性的定義[21],HDO 算法演化過程具有Markov 特性。
(2)HDO 算法演化過程具有“步步不差”特性。從生長算子的定義可知,因總有,故新一代試探解的適應(yīng)度不會劣于其老一代試探解的適應(yīng)度。
(3)HDO 算法性能優(yōu)良。由于HDO 算法擁有27個算子,是迄今為止擁有算子最多的群智能算法,由文獻(xiàn)[22]已經(jīng)證明,群智能算法的性能優(yōu)良程度與算子個數(shù)成正比,因此HDO 算法具有性能優(yōu)良的特征。
(4)HDO 算法收斂速度快。HDO 算法利用包蟲病病毒只攻擊個體的極少部分特征這一優(yōu)勢而獲得每次只需要處理n/1 000~n/100 個變量這一能力。因被處理變量大幅減少,故當(dāng)求解復(fù)雜優(yōu)化問題,特別是高維優(yōu)化問題時,能夠顯著提升收斂速度。
令N=(N1+N2+N3)/3,HDO 算法的時間復(fù)雜度計算過程如表3 所示。
由HDO 算法知,草原牧區(qū)與搜索空間H等價,將N1只狗只羊和N3個牧民重新排列成M=N1+N2+N3個個體,形成新的個體序列為。每個個體是優(yōu)化問題式(1)的一個試探解,其目標(biāo)函數(shù)值按式(1)計算為F(Xti),所有個體的狀態(tài)所形成的集合為:
Table 3 Time complexity table of HDO algorithm表3 HDO 算法的時間復(fù)雜度計算表
進(jìn)一步令:
不失一般性,令F1即為所求的全局最優(yōu)解。由式(8)的下標(biāo)形成的集合為:
U={1,2,…,M}
?i∈U,i就是個體i執(zhí)行隨機(jī)搜索時可能處的狀態(tài)。假設(shè)在時期t個體i搜索到的最好目標(biāo)函數(shù)值為Fi,其對應(yīng)的狀態(tài)為i。顯然,由式(8)知,在時期t+1 個體i進(jìn)行搜索時,若向更優(yōu)的狀態(tài)k轉(zhuǎn)移,則應(yīng)滿足k<i;相反,若向更差的狀態(tài)k轉(zhuǎn)移,則應(yīng)滿足k>i,如圖4 所示。
?X∈H,有F1≤F(X)≤FM,將H劃分為如下非空子集:
另外,顯然有:
Fig.4 State transition in random search of HDO algorithm圖4 HDO 算法隨機(jī)搜索時的狀態(tài)轉(zhuǎn)移
引理1在HDO 算法中,i=1,2,…,M,滿足:
證明(1)引理式(10)的證明:設(shè)狀態(tài)i為時期t個體i的狀態(tài),其在H 中對應(yīng)的位置為Xt,由2.6 節(jié)知,HDO 算法的隨機(jī)搜索過程具有步步不差的特性,故在時期t+1,個體i不會轉(zhuǎn)移到任何更差的狀態(tài)上去,如圖4 所示,故有:
(2)引理式(11)的證明:設(shè)狀態(tài)i為時期t個體i的狀態(tài),在時期t+1,個體i隨機(jī)選擇算子Su-Su、Su-Eu、Eu-Eu、Eu-Iu、Iu-Iu、Iu-Ru、Iu-Du、Ru-Ru、Ru-Su進(jìn)行演化以便轉(zhuǎn)移到更好的狀態(tài)k上。此時,存在有如下兩種情況:
①若狀態(tài)i就是全局最優(yōu)狀態(tài),也即i=1,則下一步轉(zhuǎn)移仍留在原狀態(tài),即k=i=1,這是因為由2.6節(jié)知,個體i不會轉(zhuǎn)移到比原狀態(tài)i更差的其他狀態(tài)上去,故必以概率p1,1=1 留在原狀態(tài)i上。因p1,1=1 >0,命題得證。
②若狀態(tài)i不是全局最優(yōu)狀態(tài),則在狀態(tài)1 和當(dāng)前狀態(tài)i之間必至少存在一個中間狀態(tài)k,如圖4 所示,使得F1≤Fk<Fi,也就是1 ≤k<i,此時個體i可以從當(dāng)前狀態(tài)i轉(zhuǎn)移到更優(yōu)的新狀態(tài)k上去,也即pi,k>0。命題得證。
綜上所述,可得?k<i,pi,k>0。 □
定理1[23]設(shè)P′是一n階可歸約隨機(jī)矩陣,即通過相同的行和列變換后可得到,其中C是m階本原隨機(jī)矩陣,且T≠0,R≠0,則有:
上述矩陣是一個穩(wěn)定隨機(jī)矩陣,P′∞=1′P′∞,P′∞=P′0P′∞唯一確定且與初始分布無關(guān),P′∞滿足條件:
定理2HDO 算法具有全局收斂性。
證明由2.6 節(jié)知,HDO 算法的搜索過程具有Markov 特性。對于每個,i=1,2,…,M可看作是有限Markov 鏈上的一個狀態(tài),根據(jù)引理中式(10)的結(jié)論可得該Markov 鏈的狀態(tài)轉(zhuǎn)移矩陣為:
由式(9)知,P′中每行的概率之和等于1。又根據(jù)引理中式(11)的結(jié)論可得:
由以上可知,P′是一個M階可歸約隨機(jī)矩陣,即Markov 狀態(tài)轉(zhuǎn)移矩陣,滿足定理1 的條件,故有:
因C∞=C=(1),T∞=0,故必有R∞=(1,1,…,1)T,這是因為P′中任意一行的概率之和總為1,故有:
上式表明,當(dāng)k→∞時,pi,1=1,i=1,2,…,M,即無論各個體初始狀態(tài)如何,最后都能以概率1 轉(zhuǎn)移到全局最優(yōu)狀態(tài)1 上去。于是有:
因此,HDO 算法具有全局收斂性。 □
HDO 算法參數(shù)包括兩部分:一部分是包蟲病傳染病模型參數(shù),該部分參數(shù)為算法內(nèi)置參數(shù),無需用戶再進(jìn)行設(shè)置;另一部分是算法運(yùn)行控制參數(shù),此類參數(shù)需要用戶根據(jù)情況進(jìn)行設(shè)置。
(1)包蟲病傳染病模型參數(shù)確定方法。包蟲病傳染病模型參數(shù)的選擇依據(jù)是確保具有足夠的隨機(jī)性,u∈{1,2,3}。依據(jù)文獻(xiàn)[19]介紹的參數(shù)取值方法并經(jīng)隨機(jī)化后,可得βu=Rand(0.42,0.93),γu=Rand(0.12,0.24),δu=Rand(0.31,0.57),λu=Rand(0.32,0.53),θu=Rand(0.18,0.43)。應(yīng)用此取值策略,任取測試情況如圖5 所示。從圖5 可知,具有極好的隨機(jī)性。參數(shù)βu、γu、δu、λu、θu的取值方法可作為算法內(nèi)置參數(shù)進(jìn)行設(shè)置,無須用戶干預(yù)。
Fig.5 Randomicity of state E2圖5 狀態(tài)E2出現(xiàn)的隨機(jī)性
(2)算法運(yùn)行控制參數(shù)設(shè)置方法。HDO 算法的運(yùn)行控制參數(shù)有G、ε、E0、L、N1、N2、N3。G和ε是兩個互補(bǔ)參數(shù),只要滿足一個即可,ε由所求解的工程問題決定,通??扇ˇ?10-5~10-8即可;G可取充分大,不妨設(shè)G=108。HDO 算法關(guān)鍵參數(shù)只有E0、L、N1、N2和N3,可令N=N1=N2=N3。下面主要討論關(guān)鍵參數(shù)E0、L、N的取值方法。由于BUMP 優(yōu)化問題與本文所要求解的工程問題形態(tài)相似,且BUMP優(yōu)化問題極難求解,故下面以BUMP 優(yōu)化問題為例來探明E0、L、N的取值方法。BUMP 優(yōu)化問題如下所示。
設(shè)平均最優(yōu)目標(biāo)函數(shù)值用Y表示,平均CPU 時間用T表示。當(dāng)L取不同值時,采用HDO 算法求解BUMP 優(yōu)化問題,令n=50,G=108,N=100,E0=0.01,運(yùn)行100 次,表4 描述了L與Y和T之間的關(guān)系。結(jié)果表明,當(dāng)L=3~7 時,Y的精度達(dá)到最高,而T增加較低。由此可見,L=3~7 為L的最佳取值區(qū)間。
Table 4 Relationship of L with Y and T表4 L 與Y 和T 之間的關(guān)系
令n=50,G=108,N=100,E0=0.01,L=3,HDO 算法運(yùn)行100 次。表5 描述了參數(shù)E0與Y和T之間的關(guān)系。結(jié)果表明,當(dāng)E0=0.006~0.100 時,Y的精度較高,且T較少;當(dāng)E0>0.100 時,T增加很大,且Y的精度也大大降低;特別是當(dāng)E0=1.000 時,無法獲得最佳解。由此可見,E0的最佳取值區(qū)間為E0=0.006~0.100。
令G=108,N=100,L=3,HDO 算法運(yùn)行100 次。表6 描述了Y、n、N和T之間的關(guān)系。從表6 可以看到:
(1)當(dāng)n增加時,T大大增加;
(2)對于給定的n,如果N增加,T大大增加;
Table 5 Rrelationship of E0 with Y and T表5 E0與Y 和T 之間的關(guān)系
(3)對于給定的n和N,如果E0增加,Y的精度會降低,但T增加。
因此,如果n>500,N=50~100 就足夠了;如果n<500,N=100 就足夠了。
某個關(guān)聯(lián)區(qū)域由n個區(qū)域組成,在氣象因素作用下,每個區(qū)域所排放的VOCs(揮發(fā)性有機(jī)廢氣)一部分留存在本區(qū)域,另一部分會擴(kuò)散到其他區(qū)域。關(guān)聯(lián)區(qū)域VOCs 聯(lián)防聯(lián)控減排方案的含義是,如何控制關(guān)聯(lián)區(qū)域內(nèi)每個區(qū)域的VOCs 排放量,才能使關(guān)聯(lián)區(qū)域內(nèi)VOCs 對大氣環(huán)境的影響降到最低。假設(shè)VOCs減排工作已進(jìn)行了t-1 期,在當(dāng)前時期t,VOCs 在n個區(qū)域的預(yù)測總產(chǎn)生量分別為Q1(t),Q2(t),…,Qn(t);對區(qū)域i來說,其他區(qū)域遷入到區(qū)域i的VOCs 分別為R1i(t),R2i(t),…,Rni(t),Rji(t)≥0,j=1~n,i=1~n;而從區(qū)域i遷出的VOCs 分別為Si1(t),Si2(t),…,Sin(t),Sij(t)≥0,j=1~n,i=1~n。在時期t,n個區(qū)域的減排量分別為x1(t),x2(t),…,xn(t),它們是待求的變量,{x1(t),x2(t),…,xn(t)}的一種取值方案構(gòu)成一個時期t的減排方案。由于該方案既考慮到了前t-1 期的減排情況,又考慮了n個區(qū)域間的相互影響,因此該減排方案具有跨時間和跨空間的特征,是一種跨時空協(xié)同減排方案。關(guān)聯(lián)區(qū)域VOCs聯(lián)防聯(lián)控最優(yōu)減排模型如式(12)所示。
Table 6 Relationship of E0,n,N with Y and T表6 E0、n、N 與Y 和T 之間的關(guān)系
式中,F(xiàn)(X(t))為減排方案的總效用;減排總量控制目標(biāo);減排任務(wù)控制目標(biāo);政府補(bǔ)貼效用控制目標(biāo);總成本控制目標(biāo)f4i(xi(t))=-ci(t)xi(t);X(t)=(x1(t),x2(t),…,xn(t));ci(t)為區(qū)域i的單位VOCs 減排成本;Pi(s)為時期s區(qū)域i的VOCs 實(shí)際減排量;Vi為受區(qū)域i影響的其他區(qū)域集合;V0為受區(qū)域間相互影響的最大允許值,簡稱區(qū)域影響系數(shù);wi為區(qū)域i的補(bǔ)貼系數(shù),區(qū)域越重要,補(bǔ)貼系數(shù)越大;O1~O4為4 個目標(biāo)函數(shù)的優(yōu)先級,優(yōu)先級次序要求滿足O1>O2>O3>O4,因此可設(shè)O1=1 000,O2=100,O3=10,O4=1;qt為前t-1 期未完成的減排量在時期t所分?jǐn)偟谋壤?,簡稱分?jǐn)偙壤?/p>
計算時,減排方案X(t)也稱為試探解;若減排方案X(t) 不滿足約束條件,則令f(X(t))=-∞,f(X(t))∈{f1(X(t)),f2(X(t)),f3(X(t)),f4(X(t))}。優(yōu)化模型式(12)是一個非線性優(yōu)化問題,傳統(tǒng)的基于函數(shù)連續(xù)性和可導(dǎo)性的數(shù)學(xué)優(yōu)化方法無法求解該優(yōu)化模型。本文采用HDO 算法對其求解。
本文以西安市為例來說明本文所提出方法的使用過程。結(jié)合西安市各區(qū)縣2018 年10 月至12 月VOCs 排放情況來制定2019 年1 月份的VOCs 最佳 減排方案,以此來說明算法HDO 算法的使用方法。表7 給 出了該市2018 年10 月至12 月VOCs 實(shí)際減 排量,西安市區(qū)縣數(shù)n=13。
表8為該市2018年10月至12月VOCs產(chǎn)生量。利用EViews8軟件可以預(yù)測出該市各區(qū)縣在2019年1月份的VOCs產(chǎn)生量,如表8 的最后一列所示。
依據(jù)該市的氣象規(guī)律,每年10 月至12 月的氣象特性較為相似。因此,依據(jù)該市各區(qū)縣2018年10月至12 月的氣象資料,采用HYSPLIT4 模式軟件可以計算出2018 年10 月至12 月各區(qū)縣VOCs 遷入和遷出量,如表9 所示。表9 最后一列是對2019 年1 月各區(qū)縣VOCs 遷入和遷出百分比的預(yù)測值。從表9 也可以確定出不同時期受某個區(qū)縣影響的其他區(qū)縣的集合。
Table 7 Actual emission reduction of VOCs in various counties from October to December 2018表7 2018 年10 月至12 月各區(qū)縣VOCs實(shí)際減排情況 億m3
Table 8 VOCs production amount in various counties from October to December 2018表8 2018年10月至12月各區(qū)縣VOCs產(chǎn)生量 億m3
表10 給出了該市各區(qū)縣VOCs 減排成本和減排補(bǔ)貼系數(shù)。前期未完成的減排量在當(dāng)前時期所分?jǐn)偟谋壤畹蜑閝t=7%;區(qū)域間相互影響的最大允許值V0=5%。
選擇7 種優(yōu)化算法與HDO 算法進(jìn)行比較,這些算法包括NP-PSO(non-parametric particle swarm optimization)[7]、DASA(differential ant-stigmergy algor-ithm)[8]、ABC(artificial bee colony)[10]、RCGA(realcoded genetic algorithm)[13]、MBBO(metropolis biogeography-based optimization)[24]、DE(differential evolution)[25]和SaDE(self-adaptive differential evolution)[26]。RCGA 是一種實(shí)數(shù)編碼型遺傳算法,搜索效率高;DASA 是一種仿蟻群算法設(shè)計思路的新型算法,收斂速度很快;NP-PSO 是一種新型粒子群算法,不需要參數(shù)設(shè)置,收斂速度很快;MBBO 是BBO 算法的改進(jìn)版,強(qiáng)化了島嶼的都市特征,搜索效率極佳;DE 通過在基本DE 算法中引入了局部誘導(dǎo)遺傳算子,大幅提升了其收斂速度;SaDE 是在傳統(tǒng)自適應(yīng)差分算法中引入了高效自調(diào)整參數(shù)演化算子,使得該算法的搜索效率大幅提升;ABC 是一種基于基因組合參數(shù)控制策略的蜂群算法,收斂速度很快。
Table 9 Emigration percentage of VOCs in counties from October 2018 to January 2019表9 2018 年10 月至2019 年1 月各區(qū)縣VOCs遷出百分比
Table 10 Cost of VOCs emission reduction and subsidies coefficient for emission reduction in counties表10 各區(qū)縣VOCs減排成本和減排補(bǔ)貼系數(shù)
計算時,7 種優(yōu)化算法的參數(shù)按表11 進(jìn)行初始化;HDO 算法的參數(shù)設(shè)置為G=108,N=100,L=3,E0=0.01。
采用HDO 算法和7 個比較算法進(jìn)行求解,每個算法共運(yùn)行50 次,得到2019 年1 月西安市各區(qū)縣VOCs最佳減排方案平均值如表12 所示。
從表12 可以看出,HDO 算法所獲得的目標(biāo)函數(shù)值要優(yōu)于其他7 種算法。對新城區(qū)來說,HDO 算法所獲得的VOCs 最佳減排方案最好,SaDE 和DE 算法所獲得的VOCs 最佳減排方案略遜于HDO 算法,但好于RCGA、DASA 和ABC 算法,而NP-PSO、MBBO所獲得的VOCs 最佳減排方案最差;對其他區(qū)縣來說,各算法所獲得的VOCs 最佳減排方案對比可作類似的分析。與其他算法相比,HDO 算法發(fā)現(xiàn)VOCs最佳減排方案的平均計算時間只有43 s,略優(yōu)于RCGA(66 s)、ABC(74 s)、SaDE(79 s)和DE(87 s),明顯優(yōu)于MBBO(102 s)、NP-PSO(243 s)和DASA(253 s);HDO 算法發(fā)現(xiàn)VOCs 最佳減排方案的平均適應(yīng)度計算次數(shù)只有58 473 次,明顯優(yōu)于RCGA(202 193 次)和ABC(243 192 次),大幅優(yōu)于SaDE(338 794 次)和DE(357 597 次),極優(yōu)于MBBO(464 381 次)、NP-PSO(523 675 次)和DASA(572 247 次)。圖6 給出了各算法求解優(yōu)化模型式(12)時的樣本收斂曲線,從該圖也可以看出,HDO 算法的收斂過程要明顯快于其他7種參與比較的算法。
Table 11 Parameters of 7 optimization algorithms表11 7 種優(yōu)化算法的參數(shù)
Table 12 Average of VOCs optimum emission reduction schemes for counties in January 2019表12 2019 年1 月各區(qū)縣VOCs最佳減排方案平均值 億m3
Fig.6 Sample convergence curve圖6 樣本收斂曲線
與其他算法相比,HDO 算法具有如下優(yōu)點(diǎn):
(1)HDO 算法中包括形態(tài)為Su-Su、Su-Eu、Eu-Eu、Eu-Iu、Iu-Iu、Iu-Ru、Iu-Du、Ru-Ru、Ru-Su的27 個算子,擁有3 個不同物種的個體,可顯著地提升算法的搜索能力。
(2)在HDO 算法中,Su-Su、Eu-Eu、Iu-Iu、Ru-Ru算子可利用強(qiáng)壯個體的特征來改善虛弱個體的特征,從而提升算法的求精(exploitation)能力;Su-Eu、Eu-Iu、Iu-Ru、Ru-Su算子可改良個體的適應(yīng)度分布特征,從而提升算法的探索(exploration)能力;Iu-Du算子可使極虛弱個體得到有效清除,從而降低算法陷入局部陷阱的概率。
(3)HDO 算法進(jìn)行演化計算時,每次只需要處理極少部分?jǐn)?shù)量,具有自動降維特性,具有求解高維優(yōu)化問題能力。
(4)HDO 算法具有全局收斂性。
HDO 算法今后的改進(jìn)方向:
(1)已經(jīng)發(fā)現(xiàn),某些傳染病能夠跨不低于4 個物種,可以利用HDO 算法的設(shè)計思路,提出跨多物種的傳染病優(yōu)化算法。
(2)將HDO 算法的狀態(tài)數(shù)從當(dāng)前的S(易感)、E(暴露)、I(染?。(治愈)、D(死亡)等5個狀態(tài)擴(kuò)展到S(易感)、E(暴露)、I(發(fā)病)、V(免疫)、R(治愈)、D(死亡)等6個狀態(tài),從而使HDO算法擁有更多的算子。
(3)深入分析Su-Su、Su-Eu、Eu-Eu、Eu-Iu、Iu-Iu、Iu-Ru、Iu-Du、Ru-Ru、Ru-Su的性能。
(4)在HDO 算法中納入DNA 機(jī)制、免疫機(jī)制,可使HDO 算法的研究更加深入。