基于改進(jìn)蝙蝠算法的雙聚類算法設(shè)計(jì)與實(shí)現(xiàn)

崔衍，薛源

（南京郵電大學(xué)物聯(lián)網(wǎng)學(xué)院，南京 210003）

0 引言

各種生物的生命活動(dòng)都是由基因調(diào)控的，近年來(lái)DNA微陣列技術(shù)憑借其優(yōu)點(diǎn)成為研究基因的主要技術(shù)［1］。DNA微陣列技術(shù)產(chǎn)生的基因表達(dá)數(shù)據(jù)是一個(gè)矩陣，反映的是基因轉(zhuǎn)錄產(chǎn)物mRNA在細(xì)胞中的豐度，通過(guò)分析這些數(shù)據(jù)能夠得到當(dāng)前細(xì)胞的生理活動(dòng)狀態(tài)。該矩陣的行代表的是基因，列代表的是條件（實(shí)驗(yàn)環(huán)境）。對(duì)基因表達(dá)數(shù)據(jù)進(jìn)行分析的傳統(tǒng)方法是聚類，但是傳統(tǒng)的聚類只能對(duì)基因進(jìn)行聚類或者對(duì)條件進(jìn)行聚類，這樣得到的是部分基因在所有條件下相似的表達(dá)模式或者所有基因在部分條件下相似的表達(dá)模式［2］，并且得到的聚類元素集合沒(méi)有相交，即一個(gè)基因或者條件只能在一個(gè)聚類里，不能屬于多個(gè)聚類［3］，這并不能很好的應(yīng)用于對(duì)基因表達(dá)數(shù)據(jù)進(jìn)行分析，因?yàn)橐徊糠只蚩赡茉诙鄠€(gè)條件下調(diào)節(jié)同一生物功能［4］。Cheng和Church［5］在2000年首次將雙聚類應(yīng)用于基因表達(dá)數(shù)據(jù)分析。雙聚類顧名思義，就是同時(shí)對(duì)行（基因）和列（條件）進(jìn)行聚類，得到一個(gè)具有相似表達(dá)模式的子矩陣。雙聚類已經(jīng)被證明是一個(gè)NP-hard問(wèn)題，智能優(yōu)化算法在解決NP-hard問(wèn)題方面有著獨(dú)特的優(yōu)勢(shì)，因此將智能優(yōu)化算法應(yīng)用于雙聚類算法是一個(gè)很好的研究方向。

Yang［6］在2010年根據(jù)蝙蝠回聲定位提出了蝙蝠算法（bat algorithm，BA）。BA算法在解決復(fù)雜的優(yōu)化問(wèn)題方面有著很好的性能，引起了學(xué)者的極大關(guān)注，現(xiàn)在已經(jīng)應(yīng)用于多個(gè)領(lǐng)域來(lái)解決問(wèn)題。本文提出了一個(gè)基于改進(jìn)蝙蝠算法的雙聚類算法。改進(jìn)了原始蝙蝠算法的位置和速度更新公式，并引入變異算子來(lái)提高局部搜索能力。最后在酵母菌數(shù)據(jù)集上進(jìn)行實(shí)驗(yàn)，與多個(gè)雙聚類算法的實(shí)驗(yàn)結(jié)果進(jìn)行比較。結(jié)果表明本文提出的基于蝙蝠算法設(shè)計(jì)雙聚類算法是一個(gè)可行的方法，能夠挖掘出更優(yōu)的雙聚類。

1 雙聚類定義

Cheng和Church提出CC算法的時(shí)候給出了雙聚類的定義［5］，定義如下：

定義1設(shè)M行N列的基因表達(dá)矩陣為A，G為包含M個(gè)基因的基因集，C為包含N個(gè)條件的條件合，a i j為基因表達(dá)矩陣A的一個(gè)元素。雙聚類定義為B=(I,J)，其中I為G的一個(gè)子集，J為C的一個(gè)子集，b ij為子矩陣B的元素，子矩陣B的平均平方殘差為：

其中：a iJ、a Ij、a IJ分別為子矩陣B的第i行平均值、子矩陣B的第j列平均值和子矩陣B(I,J)的平均值。存在一個(gè)δ≥0，如果子矩陣B(I,J)滿足H(I,J)≤δ，則稱該子矩陣為一個(gè)δ-bicluster。

雙聚類的目的就是在基因表達(dá)數(shù)據(jù)矩陣中尋找滿足條件的子矩陣，使得子矩陣中基因集在對(duì)應(yīng)的條件集上具有相似的表達(dá)模式。

2 蝙蝠算法

蝙蝠算法是基于蝙蝠的回聲定位行為提出來(lái)的［7］。蝙蝠的回聲定位能夠幫助蝙蝠探測(cè)到目標(biāo)的方向和位置，還可以幫助蝙蝠躲避障礙物。蝙蝠算法中將蝙蝠的一些行為理想化，理想化的規(guī)則如下：

（1）所有的蝙蝠都利用回聲定位來(lái)感知自身與目標(biāo)的距離，并且以某種特殊的方式區(qū)別目標(biāo)與障礙物。

（2）每只蝙蝠都擁有相同的脈沖頻率f、不同的波長(zhǎng)λ和不同的響度A，在位置x以速度v隨機(jī)飛行。它們根據(jù)目標(biāo)的接近程度調(diào)整發(fā)射脈沖的波長(zhǎng)λ（或頻率f），并調(diào)整脈沖發(fā)射頻率r。

（3）假設(shè)響度從最大值A(chǔ)0變化到Amin。

該算法的主要思想是基于以上理想化規(guī)則，改變脈沖發(fā)射頻率、脈沖頻率、聲音響度，來(lái)尋找最優(yōu)解［8］。定義蝙蝠的頻率、速度、位置更新公式：

其中f i為第i只蝙蝠的脈沖頻率，fmax和fmin分別為脈沖頻率的最大值和最小值，β是一個(gè)服從均勻分布的隨機(jī)值且β∈[0,1]，和分別是第i只蝙蝠在第t代和t-1代的飛行速度；和分別為第t代和第t-1代的位置；x*為當(dāng)前群體中蝙蝠的最優(yōu)位置（最優(yōu)解）。

為了改善算法的局部搜索能力，使用如下公式更新進(jìn)行局部搜索：

其中，At為第t次迭代所有蝙蝠的平均響度，ε∈[-1,1]是一個(gè)隨機(jī)數(shù)。隨著迭代的進(jìn)行，響度A i和發(fā)射脈沖率r i都會(huì)進(jìn)行更新，更新公式如下：

其中，α和γ是常數(shù)，是第i只蝙蝠的初始脈沖發(fā)射率。

蝙蝠算法的偽代碼如下所示：

算法1原始蝙蝠算法

目標(biāo)函數(shù)：f(x),x=(x1,x2,…,x d)T;

1.初始化蝙蝠種群x i，速度v i，蝙蝠x i的頻率f i，脈沖發(fā)射率r i，響度A0，i=1,2,…,n

2.初始化參數(shù)蝙蝠種群數(shù)n,α,γ和最大迭代數(shù)ge nmax

3.根據(jù)目標(biāo)函數(shù)找到當(dāng)前全局最優(yōu)解

4.fori=1 togenmax:

5.根據(jù)公式（4），（5），（6）生成新解

6.ifrand>r i（rand是一個(gè)隨機(jī)數(shù)，且rand∈[0,1]）

7.根據(jù)公式（7）生成一個(gè)新解

8.end if

9.通過(guò)目標(biāo)函數(shù)評(píng)價(jià)新解

10.ifrand

11.將7生成的解作為新解；

12.根據(jù)公式（8）,（9）更新脈沖發(fā)射率r i和響度A i

13.end if

14.更新全局最優(yōu)解x*

15.end for

3 基于改進(jìn)蝙蝠算法的雙聚類算法

在這一部分我們提出基于改進(jìn)蝙蝠算法的雙聚類算法。原始的蝙蝠算法是適用于解決連續(xù)問(wèn)題，不適用于解決離散問(wèn)題，顯然原始蝙蝠算法并不適用于雙聚類算法，因此對(duì)蝙蝠算法做出如下改進(jìn)：

（1）采用V型傳遞函數(shù)對(duì)位置更新公式作出改變，使蝙蝠算法適用于解決離散問(wèn)題［9］。

（2）采用動(dòng)態(tài)遞減慣性權(quán)重和速度調(diào)整因子，來(lái)增強(qiáng)全局搜索和局部搜索［10］。

（3）擴(kuò)展遺傳算法，采用一種變異算子進(jìn)行局部搜索［11］。

3.1 雙聚類編碼方式

對(duì)雙聚類進(jìn)行編碼，使之能夠表示蝙蝠算法的種群中的一個(gè)個(gè)體。編碼方式采用文獻(xiàn)［12］中的編碼方式。規(guī)則如下：用長(zhǎng)度為N+M的二進(jìn)制字符串進(jìn)行編碼，其中N和M分別代表基因表達(dá)式矩陣的行（基因）數(shù)和列（條件）數(shù)。二進(jìn)制字符串的前N位表示基因，后M位表示條件。若某位為1，則說(shuō)明該位對(duì)應(yīng)的基因或者條件屬于雙聚類，否則，不屬于雙聚類。如一個(gè)基因表達(dá)矩陣大小為7×6，包含7個(gè)基因，6個(gè)條件，該表達(dá)矩陣中的一個(gè)雙聚類編碼為1001001010010，表明該雙聚類含有3個(gè)基因和2個(gè)條件，大小為6，基因分別為基因表達(dá)矩陣中的第1、4、7個(gè)基因，條件分別為基因表達(dá)矩陣中的第2、5個(gè)條件。

3.2 初始化種群

蝙蝠算法是一種基于種群的進(jìn)化算法，因此需要對(duì)種群進(jìn)行初始化。采用文獻(xiàn)［13］的方法對(duì)種群進(jìn)行初始化。首先用K-means算法分別對(duì)基因和條件進(jìn)行聚類。共生成a個(gè)基因簇和b個(gè)條件簇，再將基因簇和條件簇進(jìn)行組合，得到a×b個(gè)雙聚類，對(duì)這些雙聚類進(jìn)行編碼，以此作為初始種群。

3.3 速度更新

我們使用動(dòng)態(tài)遞減慣性權(quán)重［10］來(lái)更新速度。在原始BA算法的速度更新公式中，上一代速度v t-1i的系數(shù)為1，不利于進(jìn)行全局搜索，容易陷入局部最優(yōu)。因此引入動(dòng)態(tài)遞減慣性權(quán)重，首先蝙蝠擁有較大的慣性權(quán)重，使蝙蝠的速度能夠在較大的范圍變化，有助于進(jìn)行全局搜索。隨著迭代的進(jìn)行，權(quán)重逐漸變小，使蝙蝠的速度在較小的范圍變化，有助于進(jìn)行局部搜索。公式如下：

其中，ωmax是最大權(quán)重，ωmin是最小權(quán)重，t代表當(dāng)前迭代次數(shù)，genmax代表最大迭代次數(shù)。arctan函數(shù)是單調(diào)遞增函數(shù)，所以ωi(t)是單調(diào)遞減函數(shù)，因此前期權(quán)重大，后期權(quán)重小。使算法得到前期速度快后期速度慢。前期快速全局搜索，后期慢速局部搜索，有效提高算法的速度和精度。

3.4 位置更新

在原始蝙蝠算法中，蝙蝠的速度負(fù)責(zé)更行蝙蝠的坐標(biāo)，且兩者都是連續(xù)值。由上文可知將蝙蝠算法應(yīng)用于雙聚類算法中時(shí)，每個(gè)蝙蝠的坐標(biāo)是由一串二進(jìn)制數(shù)字表示的，因此需要引入傳遞函數(shù)來(lái)負(fù)責(zé)由速度值變換到坐標(biāo)值。本文采用V型傳遞函數(shù)［9］來(lái)更新坐標(biāo)，如方程。

其中x ti(k)表示第i個(gè)蝙蝠在第t次迭代的第k維的坐標(biāo)值，v ti(k)表示第i個(gè)蝙蝠在第t次迭代的第k維的速度值，(x ti(k))-1是x ti(k)取反得到的。

3.5 局部搜索

對(duì)于局部搜索部分，原始蝙蝠算法使用隨機(jī)游走公式（4）為每個(gè)蝙蝠生成一個(gè)新的解。然而，本文提出的算法的搜索空間是離散空間，因此需要使用另一種方法來(lái)進(jìn)行局部搜索。本文使用位翻轉(zhuǎn)變異算子［11］。變異算子如下：

其中rand是一個(gè)隨機(jī)值且rand∈[0,1]，Pm是一個(gè)變異概率，x*全局最優(yōu)解。

3.6 適應(yīng)度函數(shù)

本文提出的算法的主要目標(biāo)是找到具有較低均方殘差的同時(shí)還要兼顧具有較大體積的雙聚類，適應(yīng)度函數(shù)如下：

其中ω1和ω2分別為雙聚類的均方殘差和雙聚類的體積的權(quán)值。均方殘差residue(B)的計(jì)算如公式（1），s ize(A)為基因表達(dá)矩陣的大小，計(jì)算公式如下：

size(B)為雙聚類B的體積，計(jì)算公式如下：

3.7 基于蝙蝠算法的雙聚類算法實(shí)現(xiàn)

基于蝙蝠算法的雙聚類算法的實(shí)現(xiàn)過(guò)程如下：首先，用3.2節(jié)描述的種群初始化方法初始化種群，再依據(jù)3.1節(jié)描述的編碼規(guī)則對(duì)每個(gè)蝙蝠個(gè)體進(jìn)行編碼；其次，以3.4節(jié)、3.5節(jié)和3.6節(jié)描述的速度更新、位置更新、局部搜索綜合而成的搜索方法對(duì)解空間進(jìn)行搜索尋找質(zhì)量更高的雙聚類；最后重復(fù)前兩個(gè)步驟，直到算法完成規(guī)定的迭代次數(shù)，輸出雙聚類。該算法的偽代碼如算法2所示。

算法2 BIBA

輸入：基因表達(dá)矩陣，蝙蝠種群數(shù)n，MSR閾值δ，最大迭代次數(shù)genmax，雙聚類數(shù)量BiNum,變異概率Pm

輸出：雙聚類集合

1.fornum=1 toBi Num:

2.初始化響度A i，脈沖發(fā)射率r i，頻率f i，速度v i，i=1,2,…,n

3.通過(guò)K-means算法分別對(duì)矩陣行和列進(jìn)行聚類，得到行聚類和列聚類

4.將行聚類和列聚類相乘，得到一定個(gè)數(shù)的雙聚類，選取M SR值最小的前n個(gè)雙聚類進(jìn)行編碼，將其作為初始種群x i

5.for種群迭代次數(shù)iter=1 to max_gen

6.fori=1 ton

7.根據(jù)公式（4），（11），（10），（12），（13）分別更新頻率，速度權(quán)重，速度，坐標(biāo)

8.ifrand>r i（rand是一個(gè)隨機(jī)數(shù)，且rand∈[0,1]）

10.end if

12.ifrand

13.接受9產(chǎn)生的新解，然后按照公式（8）,（9）更新響度A i和脈沖發(fā)射率r i

14.end if

15.end for

16.更新全局最優(yōu)解x*

17.end for

18.將迭代得到的雙聚類存入雙聚類結(jié)果集中

19.end for

20.輸出結(jié)果

4 結(jié)果與分析

4.1 數(shù)據(jù)集

實(shí)驗(yàn)所使用的數(shù)據(jù)集是釀酒酵母菌數(shù)據(jù)集［14］。數(shù)據(jù)集由2884個(gè)基因和17個(gè)條件組成，所有基因在所有條件下的表達(dá)值都位于0~600之間，其中有34個(gè)缺失值，用0~800之間的隨機(jī)數(shù)進(jìn)行填充。

4.2 參數(shù)設(shè)置

初始化參數(shù)設(shè)置如下：蝙蝠數(shù)量n=100，響度A=0.6，脈沖發(fā)射率r=0.5,α=γ=0.9，頻率fmin=0,fmax=1，最大種群迭代次數(shù)t=2500，變異算子，遞減慣性權(quán)重ωmax=0.9,ωmin=0.42,ω1=0.95,ω2=0.05。

4.3 實(shí)驗(yàn)結(jié)果

應(yīng)用改進(jìn)蝙蝠算法共產(chǎn)生100個(gè)雙聚類，隨機(jī)挑選其中四個(gè)雙聚類來(lái)繪制它們的表達(dá)值，每個(gè)雙聚類選擇了50個(gè)基因進(jìn)行繪制，如圖1所示，橫軸代表的是條件，縱軸代表的是基因表達(dá)值，從這些圖中看出每個(gè)雙聚類的基因表達(dá)值變化趨勢(shì)具有一定的相似性。在表1中分別列出了這四個(gè)雙聚類的MSR值和體積SIZE。圖2顯示了該算法的收斂曲線，橫軸代表的是迭代次數(shù)，縱軸代表的是適應(yīng)度函數(shù)值，可以發(fā)現(xiàn)該算法前期全局搜索收斂較快，后期局部搜索較慢，能夠找到更優(yōu)的結(jié)果。

圖1 4個(gè)雙聚類的基因表達(dá)值

圖2 算法收斂曲線

表1 4個(gè)雙聚類的M SR、SIZ E和C I

4.4 蝙蝠算法改進(jìn)前后對(duì)比

表2是基于原始二進(jìn)制蝙蝠算法的雙聚類算法（bi-clustering based on original binary bat algo?rithm，BOBBA）和基于改進(jìn)蝙蝠算法的雙聚類算法（BIBA）比較結(jié)果。

表2 BOBBA與BIBA的比較

其中，AMSR是所有雙聚類的平均均方殘差，MMSR是所有雙聚類中最小MSR值，ASIZE是所有雙聚類的平均體積，ACON是所有雙聚類的條件數(shù)平均值，AGENE是所有雙聚類基因平均值。并采用CI［15］對(duì)雙聚類算法結(jié)果進(jìn)行評(píng)估：

CI越小，說(shuō)明雙聚類的體積越大，均方殘差越小，則質(zhì)量越好。

由表2可以看出，BIBA算法無(wú)論是在A M S R方面還是在A S IZ E方面，結(jié)果都是要優(yōu)于BOBBA的，CI指標(biāo)也是低于BOBBA的。

4.5 不同算法的對(duì)比分析

表3列出了七種算法在酵母菌數(shù)據(jù)集下的運(yùn)行結(jié)果。評(píng)價(jià)標(biāo)準(zhǔn)和4.4的評(píng)價(jià)標(biāo)準(zhǔn)一致。

表3 7種算法結(jié)果比較

由表3可以看出，Single-object GA算法得到的AM S R最小，但是ASI ZE也是最小，因此得到的雙聚類綜合質(zhì)量并不是最好，其它算法同樣存在著相同的問(wèn)題。而B(niǎo)IBA算法得到的雙聚類AMSR較小的同時(shí)，ASIZE也較大，C I為0.21，與其他算法相比值更低，因此得到的雙聚類的綜合質(zhì)量也更優(yōu)。

5 結(jié)語(yǔ)

基因表達(dá)數(shù)據(jù)隱藏著豐富的生物信息，通過(guò)雙聚類對(duì)基因表達(dá)數(shù)據(jù)矩陣進(jìn)行挖掘可以發(fā)現(xiàn)這些生物信息。本文提出了一種新的基于改進(jìn)蝙蝠算法的雙聚類算法。它將較新的、性能較好的智能優(yōu)化算法——蝙蝠算法應(yīng)用到尋找雙聚類里。通過(guò)在酵母菌數(shù)據(jù)集上進(jìn)行實(shí)驗(yàn)，并與原始的蝙蝠算法和其他6種雙聚類算法的實(shí)驗(yàn)結(jié)果進(jìn)行比較，該方法表現(xiàn)出了較為理想的效果。但是，該算法也存在許多不足之處：計(jì)算消耗時(shí)間較長(zhǎng)，尋優(yōu)仍需繼續(xù)改進(jìn)，評(píng)價(jià)函數(shù)也需要進(jìn)一步調(diào)整。對(duì)于未來(lái)的工作，將繼續(xù)修改該算法，提升算法的性能，尋找質(zhì)量更高的雙聚類。