鼠疫傳染病優(yōu)化算法*

黃光球，陸秋琴

西安建筑科技大學(xué) 管理學(xué)院，西安 710055

1 引言

群智能優(yōu)化算法[1]求解具有大量局部最優(yōu)解優(yōu)化問題的全局最優(yōu)解具有獨(dú)到的優(yōu)勢(shì)。當(dāng)此類算法求解一個(gè)優(yōu)化問題時(shí)，同時(shí)啟動(dòng)很多試探解實(shí)施迭代計(jì)算，這種“群起而攻之”的方法可求解一些非常困難的多局部最優(yōu)解優(yōu)化問題。然而不幸的是，不存在一種群智能優(yōu)化算法可求解所有類型的優(yōu)化問題[2]。迄今為止，人們已提出很多已被廣泛應(yīng)用的群智能優(yōu)化算法，如遺傳算法[3]、蟻群算法[4-5]、粒子群算法[6-8]、生物地理學(xué)算法[9]、差分進(jìn)化算法[10-11]、蜂群算法[12]等。

在群智能優(yōu)化算法中，每個(gè)試探解被比喻成具有生物特征的個(gè)體。于是，一些特殊生物活動(dòng)場(chǎng)景被用來構(gòu)造群智能優(yōu)化算法[3-4,8,12]。然而，目前已提出的群智能優(yōu)化算法所源于的生物活動(dòng)場(chǎng)景都很簡單，且其數(shù)學(xué)基礎(chǔ)較弱，如遺傳算法、蟻群算法、生物地理學(xué)算法、蜂群算法等。有的群智能優(yōu)化算法甚至沒有生物活動(dòng)場(chǎng)景，如粒子群算法和差分進(jìn)化算法等。如果一個(gè)群智能優(yōu)化算法具有明確的較生動(dòng)的生物活動(dòng)場(chǎng)景，且該場(chǎng)景能夠用數(shù)學(xué)理論精確描述，那么該場(chǎng)景極有利于構(gòu)造群智能優(yōu)化算法，因?yàn)樵撊褐悄軆?yōu)化算法奠定在很好的數(shù)學(xué)基礎(chǔ)上，其性能易于被分析。本文提出的鼠疫傳染病優(yōu)化算法（plague infectious disease optimization，PIDO）正是基于這樣的生物活動(dòng)場(chǎng)景。

鼠疫是鼠疫耶爾森菌借鼠蚤傳播的烈性傳染病[13]，致死率極高，人類歷史上曾三次大流行，為廣泛流行于野生嚙齒動(dòng)物間的一種自然疫源性疾病。為了控制鼠疫的危害，全世界花費(fèi)了巨大的努力，構(gòu)建了大量數(shù)學(xué)模型來研究鼠疫的傳播機(jī)制。主要表現(xiàn)在如下四方面：

（1）運(yùn)用微分方程建立鼠間鼠疫和人間鼠疫的傳播動(dòng)力學(xué)模型，用于描述鼠疫的流行特征，分析模型中各參數(shù)與疫情發(fā)展變化的關(guān)系[14-17]。

（2）利用統(tǒng)計(jì)數(shù)學(xué)模型研究鼠疫疫源地的自然因素（如氣候、溫度、降雨量等）與鼠類密度、媒介密度、鼠疫感染率等參數(shù)的相關(guān)性及其滯后特征[18-19]。

（3）從生態(tài)適應(yīng)性、繁殖、生長發(fā)育、血液生化、冬眠和自帶病原譜等方面研究疫鼠的生物學(xué)特性，建立鼠疫動(dòng)物模型，從而為致病機(jī)理揭示、疫苗效果評(píng)價(jià)及鼠疫治療方案選擇提供依據(jù)[20-21]。

（4）篩選出與動(dòng)物鼠疫流行相關(guān)性最大的監(jiān)測(cè)指標(biāo)，構(gòu)建鼠疫疫源地鼠疫流行的預(yù)測(cè)預(yù)警模型，用于控制鼠疫的大規(guī)模爆發(fā)[22-23]。

在求解優(yōu)化問題時(shí)，群智能優(yōu)化算法采用的是啟發(fā)式搜索策略，同時(shí)啟動(dòng)大量個(gè)體進(jìn)行搜索。個(gè)體的行為以及個(gè)體之間的協(xié)同性變得十分重要。個(gè)體的行為可用算子來描述，個(gè)體之間的協(xié)同性就是群智能優(yōu)化算法的邏輯結(jié)構(gòu)。PIDO算法是基于具有脈沖預(yù)防接種的時(shí)滯鼠疫傳染病模型[13]而構(gòu)建的，其中個(gè)體的生物學(xué)含義就是某個(gè)生態(tài)系統(tǒng)中的人類個(gè)體，個(gè)體之間的協(xié)同性由鼠疫傳染病的傳播機(jī)制來體現(xiàn)，個(gè)體的行為由個(gè)體在易感、暴露、發(fā)病、治愈這4個(gè)狀態(tài)之間的隨機(jī)轉(zhuǎn)換來體現(xiàn)。

PIDO算法采用與現(xiàn)有群智能算法不同的設(shè)計(jì)思路，提出了將脈沖預(yù)防接種的時(shí)滯鼠疫傳染病模型轉(zhuǎn)化為能求解一類大量局部最優(yōu)解優(yōu)化問題的方法，構(gòu)造出的算子可以充分反映脈沖預(yù)防接種的時(shí)滯鼠疫傳染病模型的相互作用關(guān)系，從而體現(xiàn)出鼠疫傳染病動(dòng)力學(xué)理論的基本思想，該算法具有全局收斂性。

2 鼠疫傳染病優(yōu)化算法設(shè)計(jì)方法

設(shè)要求解的優(yōu)化問題為：

式中，Rn是n維歐氏空間，n為該優(yōu)化問題所包含的實(shí)變量總數(shù)；X=(x1,x2,…,xn)是一個(gè)n維解向量；f(X)為目標(biāo)函數(shù)，f(X)的數(shù)學(xué)表達(dá)式?jīng)]有限制條件；H為搜索空間。

2.1 PIDO算法場(chǎng)景設(shè)計(jì)

假設(shè)在一個(gè)村莊Z中生活有N個(gè)村民，每個(gè)村民用編號(hào)表示為1,2,…,N；每個(gè)村民均由n個(gè)特征來表征，一個(gè)特征等價(jià)于村民的某個(gè)器官，即對(duì)第i個(gè)村民來說，其用特征表示為(xi,1,xi,2,…,xi,n)，i=1～N；該村莊Z中流行鼠疫，村民通過與帶鼠疫病毒的老鼠進(jìn)行有效接觸，如被老鼠咬、吃老鼠肉或誤食被老鼠糞便污染的食物，就會(huì)傳染上鼠疫，鼠疫會(huì)在村民中傳播；鼠疫病毒攻擊的是村民的某些器官；該村莊中未染上鼠疫的村民稱為易感者，他們?nèi)旧鲜笠吆?，不?huì)立即發(fā)病，其體內(nèi)的病毒進(jìn)入潛伏期；體內(nèi)鼠疫病毒處于潛伏期的村民稱為暴露者，他們也會(huì)將病毒傳給與其有效接觸的其他村民；潛伏期過后的暴露者會(huì)發(fā)病，這些村民稱為發(fā)病者，他們更會(huì)將其體內(nèi)的病毒傳給與其有效接觸的其他村民；發(fā)病者和暴露者可通過醫(yī)學(xué)治療被治愈；發(fā)病者和暴露者被治愈后稱為治愈者；為了防止鼠疫病毒對(duì)村民們的危害，村民每隔一段時(shí)間接種一次疫苗，接種后的村民不會(huì)100%地成功獲得免疫；成功獲得免疫的村民在一段時(shí)間內(nèi)自身不會(huì)染上鼠疫病毒；沒有成功獲得免疫的村民仍然是易感者；成功獲得免疫的村民的免疫力在一段時(shí)間后會(huì)自動(dòng)消失而喪失免疫力；沒有實(shí)施免疫或喪失免疫力的村民會(huì)再次染上鼠疫病毒。在鼠疫病毒影響之下，該村莊中的每個(gè)村民的生長狀態(tài)將在易感、暴露、發(fā)病、治愈這4個(gè)狀態(tài)之間轉(zhuǎn)換。這種轉(zhuǎn)換映射到優(yōu)化問題式（1）的H，相當(dāng)于每個(gè)試探解在H中從一個(gè)位置跳轉(zhuǎn)到另一個(gè)位置，于是達(dá)到了對(duì)H隨機(jī)搜索的目的。每個(gè)村民的體質(zhì)強(qiáng)壯程度由該村民的特征決定，體質(zhì)強(qiáng)壯的村民能繼續(xù)生長，體質(zhì)虛弱的村民則停止生長。

把上述場(chǎng)景投射到對(duì)優(yōu)化問題式（1）的全局最優(yōu)解的搜索中，其含義如下：

優(yōu)化問題的H與村莊Z相對(duì)應(yīng)，該村莊中一個(gè)村民對(duì)應(yīng)于優(yōu)化問題式（1）的一個(gè)試探解，N個(gè)村民所對(duì)應(yīng)的試探解集合就是X={X1,X2,…,XN}，Xi=(xi,1,xi,2…,xi,n)，i=1～N。村民i的一個(gè)特征j對(duì)應(yīng)于試探解Xi的一個(gè)變量xi,j。因此，村民i與試探解Xi是等價(jià)的。村民的體質(zhì)強(qiáng)弱用指數(shù)HHI（human health index，HHI）來表示，村民i的HHI指數(shù)按下式計(jì)算：

2.2 具有脈沖預(yù)防接種的時(shí)滯鼠疫傳染病模型

為論述方便，把在一個(gè)村莊生活的村民分成4類：

（1）易感者（S 類）：即所有未染鼠疫的村民，若S類村民有效接觸帶鼠疫病毒的村民，就可能染上鼠疫。

（2）暴露者（E 類）：即已染上鼠疫病毒但還未發(fā)病的村民，他們是潛在發(fā)病者。

（3）發(fā)病者（I類）：即已染上鼠疫病毒且已發(fā)病的村民，I類村民有效接觸S類村民，就可能把病毒傳給S類村民。

（4）治愈者（R 類）：即E 類或I 類村民的治愈者，這群村民暫時(shí)尚不會(huì)得病，但經(jīng)過一段時(shí)間后，若有效接觸E類或I類病村民，還會(huì)重新染病。

具有脈沖預(yù)防接種的時(shí)滯鼠疫傳染病模型[13]如式（3）所示。

式中，t為時(shí)期；R(t)、I(t)、E(t)、S(t)分別為時(shí)期t屬于R類、I類、E類、S類人群的比例，R(t)≥0，I(t)≥0，E(t)≥0，S(t)≥0，R(t)+I(t)+E(t)+S(t)=1；其他參數(shù)如表1[13]所示。表中，Rand(A,B)表示在區(qū)間[A,B]產(chǎn)生一個(gè)均勻分布隨機(jī)數(shù)；INT(w)表示將實(shí)數(shù)w按四舍五入取整。

Table 1 Meanings of parameters and method of taking value表1 參數(shù)含義與取值方法

時(shí)期t，因任意一個(gè)人只可處于R類、I類、E類、S類中的某個(gè)類，故R(t)、I(t)、E(t)、S(t)實(shí)際上是一個(gè)人屬于R 類、I 類、E 類、S 類的概率；當(dāng)一個(gè)人屬于R類、I 類、E 類、S 類時(shí)，就表示一個(gè)人處于R 狀態(tài)、I 狀態(tài)、E狀態(tài)或S狀態(tài)。

因此，可以將式（3）應(yīng)用于人群的任何一個(gè)人，即：

式（4）用于計(jì)算時(shí)期t每個(gè)人處于R狀態(tài)、I狀態(tài)、E狀態(tài)和S狀態(tài)的概率。

記時(shí)期t參數(shù)b、r、p、q、v、β、μ的取值分別為bt、rt、pt、qt、vt、βt、μt；為計(jì)算方便，將式（4）改為遞推離散表達(dá)式，即：

式（5）和式（6）中，Ri(t)、Ii(t)、Ei(t)、Si(t)分別為時(shí)期t個(gè)體i屬于R 類、I 類、E 類、S 類人群的概率，Ri(t)≥0，Ii(t)≥0，Ei(t)≥0，Si(t)≥0，Ri(t)+Ii(t)+Ei(t)+Si(t)=1；參數(shù)bt、rt、pt、qt、vt、βt、μt，τi的取值方法如表1所示。

時(shí)期t，隨機(jī)產(chǎn)生村民的bt、rt、pt、qt、vt、βt、μt、τi，采用式（5）和式（6）計(jì)算村民i的Ri(t)、Ii(t)、Ei(t)和Si(t)。村民i在時(shí)期t處于R 狀態(tài)、I 狀態(tài)、E 狀態(tài)和S 狀態(tài)4 個(gè)狀態(tài)中的哪一個(gè)狀態(tài)，由Ri(t)、Ii(t)、Ei(t)和Si(t)所構(gòu)成的概率分布決定，即Ri(t)、Ii(t)、Ei(t)和Si(t)中的哪個(gè)值越大，其所對(duì)應(yīng)的狀態(tài)被選中的概率也越大。表2 給出了鼠疫病毒在村民中傳播情形。

Table 2 Transformation of legal states of plague virus transmission表2 鼠疫病毒傳播的合法狀態(tài)轉(zhuǎn)換

表2 所示的9 種合法狀態(tài)轉(zhuǎn)換可用圖1 來表示，圖中的每條弧對(duì)應(yīng)一個(gè)算子。表2 所示的9 種合法狀態(tài)轉(zhuǎn)換可用9 個(gè)算子描述：S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S。

Fig.1 9 legal state transitions圖1 9種合法狀態(tài)轉(zhuǎn)換

2.3 特征集合的生成方法

（1）優(yōu)勢(shì)村民集PSs的產(chǎn)生方法：從處于狀態(tài)s的村民中隨機(jī)選出L個(gè)村民，這些村民的HHI 指數(shù)比當(dāng)前村民i高，組成優(yōu)勢(shì)村民集PSs，s∈{R,I,E,S}；L稱為施加影響村民數(shù)。

（2）類別村民集CSs的產(chǎn)生方法：從處于狀態(tài)s的村民中隨機(jī)選出L個(gè)村民，組成類別村民集CSs，s∈{R,I,E,S}。

2.4 PIDO演化算子設(shè)計(jì)方法

（1）S_S算子。S_S算子描述的是：在時(shí)期t處于S 狀態(tài)的村民，在時(shí)期t+1 仍未染上傳染病的情形。對(duì)于處于易感狀態(tài)S的村民i，有：

式中，vi,j(t+1)、xs,j(t)分別為時(shí)期t+1 和時(shí)期t村民i的特征j的狀態(tài)值；αs為常數(shù)，αs=Rand(0.3,0.7)。

（2）S_E算子。S_E算子描述的是：在時(shí)期t處于S 狀態(tài)的村民，在時(shí)期t+1 通過與已處于E 狀態(tài)或I狀態(tài)的村民有效接觸后染上鼠疫病毒的情形。讓L個(gè)已處于E狀態(tài)或I狀態(tài)的村民的隨機(jī)選擇的特征j的狀態(tài)值經(jīng)處理后傳給處于S狀態(tài)的村民i的對(duì)應(yīng)特征j，使其染上鼠疫病毒，即成為暴露者。即對(duì)處于S狀態(tài)的村民i來說，有：

（3）E_E算子。E_E算子描述的是：在時(shí)期t已處于E狀態(tài)的村民，在時(shí)期t+1因潛伏期未滿而還處于潛伏期的情形。讓L個(gè)處于E 狀態(tài)但其HHI 指數(shù)高于當(dāng)前村民i的村民的特征j的狀態(tài)值經(jīng)處理后傳給處于E狀態(tài)的村民i的對(duì)應(yīng)特征j，使其HHI指數(shù)提升。即處于E狀態(tài)的村民i來說，有：

（4）E_I 算子。E_I 算子描述的是：在時(shí)期t已處于E狀態(tài)的村民，在時(shí)期t+1因潛伏期已滿鼠疫病毒開始發(fā)作的情形。L個(gè)處于I狀態(tài)的村民的特征j的狀態(tài)值經(jīng)處理后傳給已處于E 狀態(tài)的村民i的對(duì)應(yīng)特征j，使其發(fā)病。即對(duì)處于E狀態(tài)的村民i來說，有：

（5）E_R 算子。E_R 算子描述的是：在時(shí)期t已處于E狀態(tài)的村民，在時(shí)期t+1通過治療使其痊愈的情形。讓L個(gè)處于R 狀態(tài)的村民的特征j的狀態(tài)值經(jīng)處理后傳給已處于E 狀態(tài)的村民i的對(duì)應(yīng)特征j，使其痊愈。即對(duì)處于E狀態(tài)的村民i來說，有：

（6）I_I算子。I_I算子描述的是：在時(shí)期t已處于I狀態(tài)的村民，在時(shí)期t+1 仍處于I狀態(tài)的情形。讓L個(gè)處于I狀態(tài)但其HHI指數(shù)高于當(dāng)前村民i的村民的特征j的狀態(tài)值經(jīng)處理后傳給處于I狀態(tài)的村民i的對(duì)應(yīng)特征j，使其HHI指數(shù)升高。即對(duì)處于I狀態(tài)的村民i來說，有：

（7）I_R算子。I_R算子描述的是：在時(shí)期t處于I狀態(tài)的村民，在時(shí)期t+1 通過治療使其痊愈的情形。讓L個(gè)處于R 狀態(tài)的村民的特征j的狀態(tài)值經(jīng)處理后傳給當(dāng)前村民i的對(duì)應(yīng)特征j，使其痊愈。即對(duì)處于I狀態(tài)的村民i來說，有：

（8）R_R 算子。R_R 算子描述的是：在時(shí)期t已處于R 狀態(tài)的村民，在時(shí)期t+1 仍處于R 狀態(tài)的情形。讓L個(gè)處于R狀態(tài)但其HHI指數(shù)高于當(dāng)前村民i的村民的特征j的狀態(tài)值經(jīng)處理后傳給處于R 狀態(tài)的村民i的對(duì)應(yīng)特征j，使其HHI 指數(shù)升高。即對(duì)處于R狀態(tài)的村民i來說，有：

（9）R_S算子。R_S算子描述的是：在時(shí)期t處于R狀態(tài)的村民，在時(shí)期t+1 因免疫力消失而轉(zhuǎn)為S狀態(tài)的情形。讓L個(gè)處于S 狀態(tài)的村民的特征j的狀態(tài)值經(jīng)處理后傳給當(dāng)前村民i的對(duì)應(yīng)特征j，使其轉(zhuǎn)為S狀態(tài)。即對(duì)處于R狀態(tài)的村民i來說，有：

（10）生長算子。該算子描述的是村民的生長，即：

式中，Xi(t)=(xi,1(t),xi,2(t),…,xi,n(t))，Vi(t+1)=(vi,1(t+1),vi,2(t+1),…,vi,n(t+1))。HHI(Xi(t))、HHI(Vi(t+1))按式（2）進(jìn)行計(jì)算。

2.5 PIDO算法的構(gòu)造

步驟1初始化：（1）令t=0。演化時(shí)期數(shù)G=8 000～60 000，誤差要求ε=10-5～10-10，N=50～500，村民接種疫苗成功的概率Q0=0.5～0.9，村民受影響概率E0=1/1 000～1/100，L=3～6，T=3～10。（2）在H中隨機(jī)選擇N個(gè)村民所對(duì)應(yīng)的試探解[X1(0),X2(0),…,XN(0)]。（3）令V(i)=0，i=1～N。V(i)=0 表示村民i疫苗接種未成功或未接種；V(i)>0 表示村民i疫苗接種成功。（4）隨機(jī)選擇全局最優(yōu)解Y*初始值。

步驟2Ii(0)，i=1～N。

步驟3計(jì)算村民i的RIES狀態(tài)，RIESi(0)=RIES(Si(0),Ei(0),Ii(0),Ri(0))，i=1～N。其中RIESi(0)為時(shí)期0 村民i所處的狀態(tài)；函數(shù)RIES(Si(0),Ei(0),Ii(0),Ri(0))用于確定村民i將處于哪種狀態(tài)。

步驟4結(jié)束。

函數(shù)RIES(pS，pE，pI，pR)的定義如下：

2.6 PIDO算法的特點(diǎn)

（1）時(shí)間復(fù)雜度。PIDO 算法的時(shí)間復(fù)雜度計(jì)算過程如表3所示。

Table 3 Time complexity table for PIDO algorithm表3 PIDO算法的時(shí)間復(fù)雜度計(jì)算表

（2）Markov 特性。從S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S等算子的定義的式（7）～式（15）知，任意一個(gè)在時(shí)期t+1 的新試探解X(t+1)的計(jì)算生成只與其在時(shí)期t的狀態(tài)X(t)有關(guān)，而與其以前是如何演變到當(dāng)前狀態(tài)的歷程無關(guān)，表明PIDO算法的演進(jìn)過程具有Markov特性。

（3）從生長算子的定義式（16）知，時(shí)期t+1 任一村民的HHI 指數(shù)永遠(yuǎn)不會(huì)低于其在時(shí)期t時(shí)的HHI指數(shù)，表明PIDO 算法的演進(jìn)過程具有“步步不差”特性。

3 PIDO算法全局收斂性證明

由PIDO算法知，村莊Z與搜索空間H等價(jià)，將村莊內(nèi)的N個(gè)村民{X1(t),X2(t),…,XN(t)}排列成序列村民數(shù)為N，一個(gè)村民等價(jià)于一個(gè)試探解，其目標(biāo)函數(shù)值為，所有村民的狀態(tài)所形成的集合為：

進(jìn)一步令：

不失一般性，令F1即為所求的全局最優(yōu)解。由式（17）的下標(biāo)可得集合：

U={1,2,…,N}

?i∈U，i就是村民i執(zhí)行隨機(jī)搜索時(shí)可能處的狀態(tài)。假設(shè)在時(shí)期t村民i搜索到的最好目標(biāo)函數(shù)值為Fi，其對(duì)應(yīng)的狀態(tài)為i。顯然，由式（17）知，在時(shí)期t+1村民i進(jìn)行搜索時(shí)，若向更優(yōu)的狀態(tài)k轉(zhuǎn)移，則應(yīng)滿足ki，如圖2所示。

Fig.2 State transition of PIDO algorithm in random search圖2 PIDO算法隨機(jī)搜索時(shí)的狀態(tài)轉(zhuǎn)移

?X∈H，有F1≤F(X)≤FN，將H劃分成子集：

引理1在PIDO 算法中，i=1,2,…,，滿足：

（1）引理式（19）的證明：設(shè)狀態(tài)i為時(shí)期t村民i的狀態(tài)，其在H中對(duì)應(yīng)的位置為Xt，由2.6 節(jié)知，PIDO算法的隨機(jī)搜索過程具有步步不差的特性。故在時(shí)期t+1，村民i不會(huì)轉(zhuǎn)移到任何更差的狀態(tài)上去，如圖2所示，故有：

（2）引理式（20）的證明：設(shè)狀態(tài)i為時(shí)期t村民i的狀態(tài)，在時(shí)期t+1，村民i隨機(jī)選擇S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S 等算子進(jìn)行演化以便轉(zhuǎn)移到更好的狀態(tài)k上。此時(shí)，存在有如下兩種情況：

①若狀態(tài)i就是全局最優(yōu)狀態(tài)，也即i=1，則下一步轉(zhuǎn)移仍留在原狀態(tài)，即k=i=1，這是因?yàn)橛?.6節(jié)知，村民i不會(huì)轉(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之間必至少存在一個(gè)中間狀態(tài)k，如圖2 所示，使得F1≤Fk0。命題得證。

綜上所述，可得?k0。 □

定理1[24]設(shè)P′是一n階可歸約隨機(jī)矩陣，即通過相同的行和列變換后可得到，其中C是m階本原隨機(jī)矩陣，且T≠0，R≠0，則有：

上述矩陣是一個(gè)穩(wěn)定隨機(jī)矩陣，P′∞=1′P′∞，P′∞=P′0P′∞唯一確定且與初始分布無關(guān)，P′∞滿足條件：

定理1的證明過程可參見文獻(xiàn)[24]。

定理2PIDO算法具有全局收斂性。

證明由2.6 節(jié)知，PIDO 算法的搜索過程具有Markov 特性。對(duì)于每個(gè),i=1,2,…,2N可看作是有限Markov 鏈上的一個(gè)狀態(tài)，根據(jù)引理中式（19）的結(jié)論可得該Markov鏈的狀態(tài)轉(zhuǎn)移矩陣為：

由式（18）知，P′中每行的概率之和等于1。又根據(jù)引理中式（20）的結(jié)論可得：

由以上可知，P′是一個(gè)N階可歸約隨機(jī)矩陣，即Markov狀態(tài)轉(zhuǎn)移矩陣，滿足定理1的條件，故有：

因C∞=C=(1)，T∞=0，故必有R∞=(1,1,…,1)T，這是因?yàn)镻′中任意一行的概率之和總為1，故有：

上式表明，當(dāng)k→∞時(shí)，pi,1=1，i=1,2,…,N，即無論各村民初始狀態(tài)如何，最后都能以概率1轉(zhuǎn)移到全局最優(yōu)狀態(tài)1上去。于是有：

因此，PIDO算法具有全局收斂性。 □

4 PIDO算法性能分析

4.1 村民數(shù)量對(duì)PIDO算法性能的影響分析

下面以著名的Bump 優(yōu)化問題和Michalewicz 優(yōu)化問題為例對(duì)PIDO 算法在不同村民數(shù)量設(shè)置下的性能進(jìn)行分析，這兩個(gè)優(yōu)化問題求解難度很大。

Bump優(yōu)化問題：

Michalewicz優(yōu)化問題：

利用PIDO 算法求解Bump 優(yōu)化問題和Michalewicz 優(yōu)化問題，計(jì)算結(jié)果如表4 所示。計(jì)算時(shí)，G=107，ε=10-7，Q0=0.9，L=3，T=3。

從表4可以看到：

（1）當(dāng)N增加時(shí)，消耗的CPU時(shí)間大大增加。

（2）對(duì)于給定的n，如果N增加，CPU消耗時(shí)間也大大增加，但平均最佳目標(biāo)函數(shù)值的精度同時(shí)增加。

（3）對(duì)于給定的n和N，如果E0增加，平均最佳目標(biāo)函數(shù)值的精度也增加，但是消耗的CPU 時(shí)間可能增加或減少。

因此，如果n＞500，N=100～200 就足夠了；如果n＜500，N=200就足夠了。

Table 4 Comparisons of calculation results under different parameter settings表4 不同參數(shù)設(shè)置下的計(jì)算結(jié)果對(duì)比

4.2 與其他群智能算法的比較

本文采用國際上通用的智能優(yōu)化算法測(cè)試包CEC2013[25]所提供的6 個(gè)基準(zhǔn)函數(shù)優(yōu)化問題來測(cè)試PIDO 算法的性能，如表5 所示。該測(cè)試包括有28 個(gè)經(jīng)過精心設(shè)計(jì)的基準(zhǔn)函數(shù)優(yōu)化問題。

表5 中O是可以任意設(shè)定的理論全局最優(yōu)解。F14、F25、F28的理論全局最優(yōu)解目前尚未發(fā)現(xiàn)，表5中的全局最優(yōu)目標(biāo)函數(shù)值是本文算法發(fā)現(xiàn)的迄今為止最好的目標(biāo)函數(shù)值。本文用PIDO算法求解表5中的6 個(gè)函數(shù)優(yōu)化問題，其參數(shù)是N=200，n=50，G=107，ε=10-7，Q0=0.9，E0=0.01，L=3，T=3。與PIDO算法進(jìn)行比較的7種智能優(yōu)化算法均選自國際著名期刊近期刊登的算法，這些算法如表6 所示，即RCGA（real-coded genetic algorithm）、DASA（differential ant-stigmergy algorithm）、NP-PSO（non-parametric particle swarm optimization）、MpBBO（metropolis biogeography-based optimization）、MDE-LiGO（modifieddifferential evolution with locality induced genetic operators）、SLADE（symmetric Latin-based adaptive differential evolution）、ABC（artificial bee colony algorithm）。這7 種算法的終止運(yùn)行條件為：進(jìn)化代數(shù)G=107或者最優(yōu)解誤差ε=10-7。

Table 5 6 benchmark function optimization problems表5 6個(gè)基準(zhǔn)函數(shù)優(yōu)化問題

Table 6 Parameters setting of compared intelligent optimization algorithms表6 參與比較的智能優(yōu)化算法的參數(shù)設(shè)置

表6中各參數(shù)的含義可參見相關(guān)文獻(xiàn)。針對(duì)表5中列出的6個(gè)優(yōu)化問題，用PIDO算法和表6列出的7種算法進(jìn)行求解，每種算法均獨(dú)立求解51 次。表7給出了8種算法所求得的最優(yōu)目標(biāo)函數(shù)值的平均值、標(biāo)準(zhǔn)差、中值、最小值、最大值、適應(yīng)度評(píng)價(jià)次數(shù)。

總排名1 是這些算法基于平均最優(yōu)目標(biāo)函數(shù)值進(jìn)行的排名，總排名2是這些算法基于平均最優(yōu)目標(biāo)函數(shù)值和適應(yīng)度評(píng)價(jià)次數(shù)進(jìn)行的排名；最終總排名1和最終總排名2分別是基于總排名1和總排名2所進(jìn)行的排名。

非參數(shù)Friedman檢驗(yàn)[26-27]是基于PIDO算法所得的結(jié)果與7 種被比較算法所獲得的結(jié)果之間進(jìn)行的非參數(shù)檢驗(yàn)。Friedman檢驗(yàn)的結(jié)果顯示在表8中，其中顯著性＝1表示PIDO算法的性能與被比較算法具有99%的統(tǒng)計(jì)學(xué)差異，顯著性＝0 表示沒有顯著差異。在表8 中，顯著性＝1 的案例數(shù)和顯著性＝0 的案例數(shù)分別表示PIDO算法與7種被比較算法顯著不同和沒有顯著差異地求解基準(zhǔn)函數(shù)優(yōu)化問題的數(shù)目。

從表7 可以看出PIDO、RCGA、DASA、NP-PSO、MpBBO、MDE-LiGO、SLADE、ABC按平均最優(yōu)目標(biāo)函數(shù)值的精度進(jìn)行排序所得的結(jié)果如下：

按平均最優(yōu)目標(biāo)函數(shù)值和適應(yīng)度評(píng)價(jià)次數(shù)進(jìn)行排序所得的結(jié)果如下：

從表8可以知道，PIDO算法求解6個(gè)基準(zhǔn)函數(shù)優(yōu)化問題的顯著性案例總數(shù)為39，明顯大于不顯著性案例總數(shù)3，表明PIDO算法的性能明顯優(yōu)于7種被比較算法。

圖3（a）～圖3（f）說明PIDO、RCGA、DASA、NPPSO、MpBBO、MDE-LiGO、SLADE、ABC 算法求解表5 所示的6 個(gè)基準(zhǔn)函數(shù)優(yōu)化問題時(shí)的樣本收斂曲線，其中水平和垂直軸采用對(duì)數(shù)刻度。從圖3 和表7均可以看出，當(dāng)PIDO 算法求解6 個(gè)基準(zhǔn)函數(shù)優(yōu)化問題時(shí)，均能發(fā)現(xiàn)質(zhì)量最好的全局最優(yōu)解。綜合看來，PIDO 算法的綜合性能要優(yōu)于7 種被比較算法，表明其求解精度高且計(jì)算速度快。

4.3 局部尋優(yōu)能力和全局尋優(yōu)能力分析

本文采用文獻(xiàn)[28]提出的方法來對(duì)PIDO算法的局部和全局尋優(yōu)能力的平衡性進(jìn)行分析。其中穿透行為用來說明PIDO算法的全局尋優(yōu)能力，而膨脹行為用來說明PIDO算法的局部尋優(yōu)能力。

圖4顯示了當(dāng)PIDO算法求解基準(zhǔn)優(yōu)化問題F3[25]時(shí)，穿透和膨脹行為在時(shí)間間隔[0，694]和[0，966]內(nèi)的協(xié)調(diào)關(guān)系。從圖4 可以看出，穿透行為在0 s 到120 s占優(yōu)，而膨脹行為在之后120 s占優(yōu)。

Table 7 Optimal solutions obtained by 8 compared algorithms表7 8種參與比較的算法所獲得的最優(yōu)解

續(xù)表7

Table 8 Results comparison of Friedman test(α=0.01)表8 Friedman檢驗(yàn)結(jié)果比較(α=0.01)

圖5分別說明了當(dāng)PIDO算法求解F3時(shí)，在時(shí)間間隔[0，694]和[0，966]內(nèi)膨脹和粘滯行為之間的協(xié)調(diào)性。從圖5 中可以看出，膨脹行為總是占優(yōu)，粘滯行為根本不占優(yōu)，但是膨脹行為的占優(yōu)隨著時(shí)間的流逝而減少；相反，粘滯行為出現(xiàn)的次數(shù)隨時(shí)間而增加。

5 結(jié)束語

PIDO具有如下特點(diǎn)：

（1）PIDO算法包括S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S等9個(gè)算子，這些算子可明顯地增加該算法的搜索能力。

（2）采用隨機(jī)方法確定PIDO算法中的鼠疫傳染病模型中的參數(shù)和各算子中的相關(guān)參數(shù)，既可使模型更能表達(dá)實(shí)際情況，又大幅減少了輸入的參數(shù)個(gè)數(shù)。

Fig.3 Convergence curves of samples圖3 樣本收斂曲線

Fig.4 Coordination of penetration and expansion behavior圖4 穿透行為和膨脹行為之間的協(xié)調(diào)性分析

（3）PIDO 算法中的算子S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S是利用鼠疫傳染病模型來構(gòu)造的，不需要與求解的問題相關(guān)，因此PIDO 算法具有通用性。

Fig.5 Coordination of expansion and stickiness behavior圖5 膨脹行為和粘滯行為之間的協(xié)調(diào)性分析

（4）在PIDO算法中，算子S_S、E_E、I_I和R_R能使HHI指數(shù)高的村民向HHI指數(shù)低的村民傳遞強(qiáng)壯特征信息，從而使后者能向好的方向發(fā)展；算子S_E、E_I、E_R、I_R 和R_S 既能使處于不同狀態(tài)的村民之間交換信息，又能使村民獲得其他村民經(jīng)處理后的特征信息，降低了村民落入局部最優(yōu)陷阱的幾率，從而提升了PIDO 算法的局部求精(exploitation)能力；脈沖式預(yù)防接種擁有使村民躍出局部最優(yōu)陷阱的特性。因此，PIDO 算法能從多個(gè)角度實(shí)現(xiàn)村民之間的信息交換，這對(duì)提升PIDO算法的搜尋(exploration)能力有意義。

（5）在PIDO 算法進(jìn)行演化計(jì)算時(shí)，由于只有村民的極少部分特征被鼠疫病毒攻擊，因此每次村民之間交換信息，只有極少部分特征參與計(jì)算，從而當(dāng)求解維數(shù)較高的優(yōu)化問題時(shí)，PIDO 算法收斂速度可顯著提升，PIDO算法適于求解維數(shù)較高的優(yōu)化問題。

（6）PIDO 算法的演化過程具有Markov 特性和“步步不差”特性，從而確保該算法具有全局收斂性。

PIDO算法未來的改進(jìn)方向如下：

（1）利用鼠疫傳染病動(dòng)力學(xué)模型優(yōu)化PIDO算法的相關(guān)參數(shù)，使得這些參數(shù)設(shè)置更合理。

（2）深入研究S_S、S_E、E_E、E_I、E_R、I_I、I_R、R_R、R_S等算子的動(dòng)態(tài)特征。

（3）深入研究PIDO算法求解過程中村民們的動(dòng)態(tài)特征。