準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化定量構(gòu)效特征選擇方法

李小林，王 靜，張?jiān)?，黃世國

(福建農(nóng)林大學(xué) 計(jì)算機(jī)與信息學(xué)院，福州 350002) E-mail：fjhsg25@126.com

1 引 言

帝王蝶算法(Monarch Butterfly Optimization，MBO)是一種新的群體智能優(yōu)化算法，通過模擬帝王蝶遷徙的行為實(shí)現(xiàn)尋優(yōu)過程[1].為了避免算法過早收斂，通過結(jié)合和聲搜索操作[2]或結(jié)合自適應(yīng)交叉算子和貪婪策略[3]對該算法進(jìn)行了改進(jìn).馮艷紅等[4]在帝王蝶算法上結(jié)合差分變異策略來生成新個(gè)體，提高了獲得折扣{0，1}KP問題全局最優(yōu)解的概率.針對TSP問題，Wang等[5]提出了一種二進(jìn)制帝王蝶優(yōu)化算法.

定量構(gòu)效關(guān)系(Quantitative Structure-Activity Relationship，QSAR)的目標(biāo)是建立化合物的結(jié)構(gòu)性質(zhì)和相應(yīng)的生物活性之間的關(guān)系[6].化合物的結(jié)構(gòu)由各種特征(分子描述符)進(jìn)行編碼來描述.然而，并不是所有的特征都與其生物活性相關(guān).特征太多可能導(dǎo)致模型過擬合，特征太少則可能會導(dǎo)致欠擬合[7].只包含相關(guān)和無冗余特征的最優(yōu)特征子集能夠提高QSAR模型預(yù)測的精度和可解釋性[8].因此，特征選擇成為QSAR研究的關(guān)鍵預(yù)處理步驟.

元啟發(fā)式算法因其較強(qiáng)的全局搜索能力和潛力已成功應(yīng)用于QSAR特征選擇.這些代表性的方法包括遺傳算法[9，10]，粒子群優(yōu)化[11，12]，差分進(jìn)化算法等.

針對QSAR特征選擇需在優(yōu)先滿足最大化準(zhǔn)確度前提下同時(shí)最小化特征數(shù)的多目標(biāo)要求，本研究在MBO算法的粒子排序步驟中采用了非支配排序算法，并修改了調(diào)整算子；采用準(zhǔn)確度優(yōu)先的策略，改善了傳統(tǒng)多目標(biāo)算法在低準(zhǔn)確度區(qū)域耗費(fèi)大量計(jì)算資源的缺陷；采用分小組使用不同突變策略的方法解決了算法早熟的問題.實(shí)驗(yàn)結(jié)果表明，這些策略的綜合應(yīng)用提升了MBO算法的性能，實(shí)現(xiàn)了QSAR建模提高準(zhǔn)確度并減少特征數(shù)的目的.

2 定量構(gòu)效建模

為了建立化合物結(jié)構(gòu)特征與其活性之間的關(guān)系，采用機(jī)器學(xué)習(xí)方法來建立二者的關(guān)聯(lián).其中偏最小二乘法(Partial Least Squares，PLS)常用于預(yù)測化合物結(jié)構(gòu)特征與其相應(yīng)的生物活性或化學(xué)性質(zhì)之間的關(guān)系.偏最小二乘法通過具有因子分析的線性多元模型將兩個(gè)數(shù)據(jù)矩陣x和y聯(lián)系起來[13].在QSAR建模中，x表示所有化合物的特征向量，y表示相應(yīng)化合物的活性.

用Q2值衡量QSAR建模的精度，定義如式(1)所示：

(1)

3 帝王蝶優(yōu)化算法

帝王蝶優(yōu)化算法模擬了帝王蝶的遷移行為[1].其行為簡單描述為：1)帝王蝶的整個(gè)種群被分為兩個(gè)子種群，稱為子種群1和子種群2；2)每個(gè)子種群中更優(yōu)的子代個(gè)體將分別由遷移操作算子或調(diào)整算子產(chǎn)生，并取代其父代個(gè)體；3)當(dāng)前最優(yōu)的個(gè)體自動存活到下一代.

假設(shè)問題對應(yīng)于一個(gè)d維解空間，初始化產(chǎn)生一個(gè)大小為NP的種群P=(X1，X2，…，XN)，第i個(gè)個(gè)體在解空間的位置可以表示為一個(gè)d維向量Xi=(xi，1，xi，2，…，xi，d)T，顯然，每個(gè)個(gè)體的位置都對應(yīng)著目標(biāo)問題的一個(gè)潛在解.其遷移算子和調(diào)整算子描述如下.

設(shè)NP1和NP2分別為兩個(gè)子種群中個(gè)體數(shù)量，p為子種群1個(gè)體的占比，即NP1=ceil(p×NP)，NP2=NP-NP1.其中，函數(shù)ceil(c)表示向上取整.遷移算子的數(shù)學(xué)模型描述如式(2)所示：

(2)

其中，i=1，2，…，N1，t表示當(dāng)前迭代次數(shù)，xr1和xr2分別表示帝王蝶個(gè)體r1和r2的位置，這兩只帝王蝶分別隨機(jī)取自兩個(gè)子種群.標(biāo)量r=rand×peri，rand是[0，1]間的一個(gè)隨機(jī)數(shù)；peri代表遷移周期.

在調(diào)整算子中，若rand≤p，則子種群2中個(gè)體按式(3)更新：

(3)

其中，i=1，2，…，N2.xbest表示全局最優(yōu)個(gè)體的位置.

否則，個(gè)體將按式(4)更新：

(4)

其中，xr3表示從子種群2中隨機(jī)選擇的一個(gè)帝王蝶的位置.

(5)

(6)

其中，α是權(quán)重因子，由α=Smax/t2計(jì)算所得，Smax表示帝王蝶的最大步長.dx表示帝王蝶的步長，由Levy函數(shù)計(jì)算.BAR表示調(diào)整率.以最小值優(yōu)化為例，MBO算法的步驟如算法1所示.

算法1.MBO算法偽代碼

初始化參數(shù)：種群大小NP，子種群1中個(gè)體的比例p，調(diào)整率BAR，遷移周期peri，權(quán)重因子α，最大步長Smax，t=1，最大評價(jià)次數(shù)T；

隨機(jī)初始化大小為NP的種群；

計(jì)算種群中每個(gè)個(gè)體的適應(yīng)度值；

Whilet≤Tdo

根據(jù)適應(yīng)度值對所有個(gè)體進(jìn)行排序；

適應(yīng)度值較小的前p×NP個(gè)個(gè)體組成子種群1，其余個(gè)體組成子種群2；

For子種群1中的每個(gè)個(gè)體

根據(jù)公式(2)進(jìn)行個(gè)體的更新；

End

For子種群2中的每個(gè)個(gè)體

根據(jù)公式(3-6)進(jìn)行個(gè)體的更新；

End

將兩個(gè)更新后的子種群合并成一個(gè)種群；

計(jì)算種群中每個(gè)個(gè)體的適應(yīng)度值；

t=t+1；

End

輸出最優(yōu)個(gè)體的位置xbest及其對應(yīng)的適應(yīng)度值f(xbest).

4 準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化QSAR特征選擇

4.1 二進(jìn)制帝王蝶特征選擇算法

傳統(tǒng)的帝王蝶優(yōu)化算法的目標(biāo)是在連續(xù)數(shù)據(jù)空間中尋找最優(yōu)解.但特征選擇是在二進(jìn)制空間中搜索最優(yōu)解.為了解決QSAR特征選擇的問題，每個(gè)個(gè)體均表示為一個(gè)二進(jìn)制字符串.字符串的長度等于特征的數(shù)量.字符串中的“1”表示選擇了對應(yīng)的特征，“0”表示未選擇對應(yīng)的特征.同時(shí)，偏最小二乘法用于計(jì)算每個(gè)個(gè)體的適應(yīng)度值.

二進(jìn)制帝王蝶優(yōu)化(Binary Monarch Butterfly Optimization，BMBO)算法的過程與標(biāo)準(zhǔn)帝王蝶優(yōu)化算法相似，只在兩個(gè)地方做了修改：

1)在初始化階段，所有個(gè)體的位置均隨機(jī)生成，然后根據(jù)給定閾值轉(zhuǎn)換成二進(jìn)制字符串.

2)因在調(diào)整算子中獲得的值是連續(xù)值，使用Tanh()傳遞函數(shù)[14]將其映射為二進(jìn)制值.故當(dāng)rand>BAR時(shí)，公式(5)修改為：

(7)

需要注意的是，為了適應(yīng)最小值優(yōu)化，本研究用準(zhǔn)確度的負(fù)值即-Q2作為適應(yīng)度函數(shù).

4.2 快速非支配排序策略

為了同時(shí)滿足最小化準(zhǔn)確度負(fù)值和最小化特征數(shù)的要求，引入了非支配排序算法.

首先，快速非支配排序涉及到的概念有如下兩點(diǎn)：

1)對于兩個(gè)向量x和y，當(dāng)且僅當(dāng)：

(8)

我們稱x支配y，表示為xy，否則，這兩個(gè)向量互不支配.

2)對于一個(gè)解x∈X，當(dāng)且僅當(dāng)：

{?x∈X|xjxi，i≠j}

(9)

則稱該解x為帕累托最優(yōu).所有帕累托最優(yōu)解的集合稱為帕累托解集，定義如式(10)所示：

Ps={xi，xj∈X|?xjxi，i≠j}

(10)

一個(gè)包含帕累托解集目標(biāo)函數(shù)值的集合稱為帕累托前沿，表示為式(11)：

Pf={f(xi)|xi∈Ps}

(11)

因此，快速非支配排序算法對每個(gè)解x均需計(jì)算兩個(gè)參數(shù)：1)種群中支配x的個(gè)體的數(shù)目nx；2)種群中被x支配的個(gè)體的集合Sx.

快速非支配排序的整個(gè)過程描述如下：

1)找到nx=0的個(gè)體，然后將這些個(gè)體存儲到當(dāng)前集合Fi中，其中i代表級別；

2)找出在Fi中由x支配的個(gè)體的集合Sx.對于Sx中的每個(gè)個(gè)體x，用nx-1更新nx.如果nx=0，則將x存儲在集合H中.

3)F1中的個(gè)體為第一非支配水平的個(gè)體.H代表當(dāng)前集合.重復(fù)上述這些步驟，直到所有的個(gè)體被分配到不同級別的集合為止.

同時(shí)，還對調(diào)整算子進(jìn)行了修改.與傳統(tǒng)的非支配排序算法不同，我們只考慮第一層帕累托前沿F1，并且只計(jì)算F1中個(gè)體的擁擠度.當(dāng)rand≤p時(shí)，調(diào)整算子中的個(gè)體按照式(12)進(jìn)行更新：

(12)

4.3 準(zhǔn)確度優(yōu)先策略

該策略的思想是當(dāng)子代個(gè)體擁有更好的精度時(shí)，對應(yīng)的父代個(gè)體將會被其取代.增加該策略的原因是：在QSAR建模中，模型的準(zhǔn)確度指標(biāo)優(yōu)先于特征數(shù)指標(biāo).采用準(zhǔn)確度優(yōu)先策略可以避免多目標(biāo)算法耗費(fèi)大量計(jì)算資源來搜索低準(zhǔn)確度的空間，即圖1中矩形框內(nèi)的區(qū)域.

圖1 低準(zhǔn)確度的區(qū)域示意圖Fig.1 Illustration of region with the low accuracy

4.4 突變策略

為了解決早熟的問題，整個(gè)種群分為子組1、子組2和子組3，每個(gè)子組的個(gè)體數(shù)量均為NP/3.子組1中的突變率設(shè)置為隨著迭代次數(shù)的增加而降低，子組2的突變率則設(shè)置為隨機(jī)數(shù)，子組3中的個(gè)體不進(jìn)行突變.

4.5 準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化QSAR特征選擇的基本框架

綜合上述3個(gè)策略，準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化(Accuracy Prior Multi-Objective Monarch Butterfly Optimization，APMOMBO)QSAR特征選擇的步驟如算法2所示.該算法最終輸出的最優(yōu)粒子是準(zhǔn)確度最高的粒子.

算法2.準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化QSAR特征選擇算法偽代碼

初始化參數(shù)：種群大小NP，子種群1的個(gè)體占比p，調(diào)整率BAR，遷移周期peri，最大步長Smax，t=1，最大評價(jià)次數(shù)T；

隨機(jī)初始化大小為NP的種群；

通過偏最小二乘法計(jì)算種群中每個(gè)個(gè)體的適應(yīng)度值，并統(tǒng)計(jì)每個(gè)個(gè)體所選擇的特征數(shù)量；

通過非支配排序策略對種群中所有個(gè)體進(jìn)行排序；

Whilet≤Tdo

根據(jù)p將種群分成子種群1和子種群2；

For子種群1中的每個(gè)個(gè)體

根據(jù)公式(2)進(jìn)行個(gè)體的更新；

End

For子種群2中的每個(gè)個(gè)體

根據(jù)公式(7)和公式(12)進(jìn)行個(gè)體的更新；

End

將兩個(gè)更新后的子種群合并成一個(gè)種群，并將該種群分成子組1、子組2和子組3；

根據(jù)突變策略對前兩個(gè)子組中的個(gè)體進(jìn)行變異操作；

通過偏最小二乘法計(jì)算種群中每個(gè)個(gè)體的適應(yīng)度值，并統(tǒng)計(jì)每個(gè)個(gè)體所選擇的特征數(shù)量；

根據(jù)準(zhǔn)確度優(yōu)先策略更新子組1和子組2中的父代個(gè)體；

將所有子組合并為一個(gè)種群；

通過非支配排序策略對種群中所有個(gè)體進(jìn)行排序；

t=t+1；

End

輸出最優(yōu)個(gè)體的位置xbest及其對應(yīng)的適應(yīng)度值f(xbest).

5 算法性能測試

5.1 實(shí)驗(yàn)設(shè)置及性能指標(biāo)

為了測試我們提出的算法的特征選擇性能，本文選擇了文獻(xiàn)[11]中的3個(gè)基準(zhǔn)數(shù)據(jù)集，如表1所示.

表1 測試數(shù)據(jù)集Table 1 List of datasets used in the experiments

同時(shí)，將提出的APMOMBO算法和BMBO、PSO和WS-PSO進(jìn)行比較.在本文中BMBO和APMOMBO算法的相關(guān)參數(shù)設(shè)置為：NP=150，p=0.5，BAR=0.5，peri=1.2，Smax=1.0，T=20000.除種群數(shù)量和評價(jià)次數(shù)分別為150和20000外，PSO和WS-PSO算法的其它參數(shù)設(shè)置與文獻(xiàn)[11]相同，即PSO的學(xué)習(xí)率a=0.5，WS-PSO的學(xué)習(xí)率a=0.5，b=0.8，突變率c=0.05以及加權(quán)抽樣比例為α=0.5.偏最小二乘法中潛在化合物的最大數(shù)量被設(shè)置為20，并且通過5折交叉驗(yàn)證來確定潛在化合物的最佳數(shù)量.為了消除隨機(jī)性，所有算法在3個(gè)數(shù)據(jù)集上均獨(dú)立運(yùn)行100次.本文所有實(shí)驗(yàn)在硬件配置為Intel(R)Core(TM)i5-7500 CPU @3.40GHz、內(nèi)存為16.0GB的計(jì)算機(jī)上進(jìn)行，操作系統(tǒng)為Microsoft Windows 10.編程語言為MATLAB，編譯環(huán)境為MATLAB R2018b.我們使用文獻(xiàn)[11]中描述的5折交叉驗(yàn)證方法來評估QSAR特征選擇的性能.本文的性能指標(biāo)包括準(zhǔn)確度Q2、均方根誤差[15](表示為RMSECV)和特征數(shù)量(表示為NSD).準(zhǔn)確度的值越大越好，其它兩個(gè)指標(biāo)則越小越好.同時(shí)，我們還計(jì)算了各指標(biāo)的標(biāo)準(zhǔn)差.

5.2 性能比較與分析

由于BMBO未在QSAR建模特征選擇中得到應(yīng)用，我們給出了BMBO和PSO和WS-PSO的比較結(jié)果，如表2所示，其中粗體為較優(yōu)結(jié)果.在準(zhǔn)確度Q2和均方根誤差RMSECV方面，BMBO的性能均優(yōu)于PSO及WS-PSO算法，且標(biāo)準(zhǔn)差相對更低.這表明將BMBO算法用于QSAR建模特征選擇，可以保持較高的準(zhǔn)確度和穩(wěn)定性.但BMBO算法在Artemisinin數(shù)據(jù)集和Selwood數(shù)據(jù)集的降維效果較差，這說明需要進(jìn)一步改進(jìn)BMBO的降維性能.

為了評估不同策略與BMBO算法結(jié)合的性能和效率，我們將3個(gè)策略逐個(gè)添加到BMBO算法中，得到3個(gè)不同的算法，即僅采用快速非支配排序策略的二進(jìn)制帝王蝶優(yōu)化算法(稱為IBMBO-1)，同時(shí)采用快速非支配排序策略和準(zhǔn)確度優(yōu)先策略的二進(jìn)制帝王蝶優(yōu)化算法(稱為IBMBO2)，以及同時(shí)采用快速非支配排序策略、準(zhǔn)確度優(yōu)先策略以及突變策略的二進(jìn)制帝王蝶優(yōu)化算法即APMOMBO.表3給出了BMBO和上述不同策略結(jié)合的實(shí)驗(yàn)結(jié)果，其中粗體為較優(yōu)結(jié)果.

從表3可以看出，IBMBO-2和APMOMBO在準(zhǔn)確度和均方根誤差方面比IBMBO-1表現(xiàn)更好，而在特征數(shù)量方面比IBMBO-1表現(xiàn)得差.這表明，準(zhǔn)確度優(yōu)先策略和準(zhǔn)確度優(yōu)先+突變策略有利于提高特征選擇的準(zhǔn)確性和穩(wěn)定性，僅使用快速非支配排序算法很難搜索到所有有用的特征來保持較高的準(zhǔn)確性.與IBMBO-2相比，APMOMBO除了在Selwood數(shù)據(jù)集的特征數(shù)量之外的所有性能指標(biāo)上都有提升.結(jié)果表明，采用突變策略避免了早熟，提高了精度，減少了均方根誤差和特征數(shù).

對比表2和表3中的數(shù)據(jù)可以看出，APMOMBO在3個(gè)性能指標(biāo)上都比BMBO、WS-PSO和PSO表現(xiàn)得更好.這說明改進(jìn)的BMBO算法通過使用快速非支配排序、準(zhǔn)確度優(yōu)先和突變策略，提高了精度，減少了均方根誤差和特征數(shù)量.

表2 3種算法用于QASR建模特征選擇的性能指標(biāo)對比Table 2 Comparison of the algorithms for feature selection in QSAR

表3 BMBO算法結(jié)合3種策略的實(shí)驗(yàn)結(jié)果Table 3 Experimental results of BMBO algorithms with three strategies

6 總 結(jié)

在本文中，我們應(yīng)用二進(jìn)制帝王蝶優(yōu)化來選擇分子描述符用于QSAR建模，并且結(jié)合快速非支配排序、準(zhǔn)確度優(yōu)先和突變3種策略來提高BMBO的性能.

對BMBO和準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化QSAR特征選擇方法在3個(gè)基準(zhǔn)數(shù)據(jù)集上進(jìn)行了測試.實(shí)驗(yàn)結(jié)果表明BMBO算法可以用于QSAR特征選擇，其在準(zhǔn)確度和均方根誤差上的性能均優(yōu)于PSO和WS-PSO，但在降維方面略有不足.而準(zhǔn)確度優(yōu)先的多目標(biāo)帝王蝶優(yōu)化QSAR特征選擇方法則在所有指標(biāo)上均優(yōu)于PSO、WS-PSO和BMBO算法.這表明快速非支配排序、準(zhǔn)確度優(yōu)先和突變策略提高了特征選擇的準(zhǔn)確度，降低了均方根誤差，并減少了特征的數(shù)量.該算法可有效實(shí)現(xiàn)QSAR建模中的特征選擇.