基于ACO的多序列間最長公共子序列查詢

王虹勝

（四川大學(xué)視覺合成圖形圖像技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室，成都 610065)

0 引言

多序列間的最長公共子序列（LCS）查詢廣泛地應(yīng)用于計(jì)算生物學(xué)、模式識別、文本比較和數(shù)據(jù)壓縮等。多序列間的最長公共子序列可以用來度量一組基因序列間的相似度，被廣泛應(yīng)用于基因測序和蛋白質(zhì)功能測序等領(lǐng)域。本文主要針對任意數(shù)量的生物數(shù)據(jù)查找最長公共子序列。

經(jīng)典的雙序列間LCS查詢算法是基于動態(tài)規(guī)劃（Dynamic Programming）的思想，可以精確查詢出雙序列間的所有最長公共子序列。但是，對于長度為M的N條序列，動態(tài)規(guī)劃算法的時(shí)間復(fù)雜度和空間復(fù)雜度均為O（MN），隨著數(shù)據(jù)規(guī)模的增大，其復(fù)雜度將以指數(shù)形式增加，這不利于多序列（N>10）或較長序列（M=500）的查詢。為了進(jìn)一步提升LCS問題的解決速度，人們提出了近似方法（Approximate Algorithm）和啟發(fā)式方法（Heuristic Algorithm）。這兩種方法均以降低求解精度為代價(jià)，獲得計(jì)算效率的大幅提升。近似方法是指在較短時(shí)間內(nèi)計(jì)算出最長公共子序列的近似最優(yōu)解。最早的LCS近似算法是Long Run方法，該方法在實(shí)際應(yīng)用中的利用率不高。2004年Huang等人提出的 BNMAS（Beat Next for Maximal Available Symbols）算法，但其時(shí)間復(fù)雜度嚴(yán)重依賴于序列集合的字符種類總數(shù)。啟發(fā)式方法是專為解決NP問題而提出的算法框架，常用于解決多序列間LCS查詢等無法在多項(xiàng)式時(shí)間內(nèi)獲得精確解的問題。2006年Hinkemeyer和Julstrom等人提出了基于遺傳算法（Genetic Algorithm）的LCS查詢算法。2008年Chiang等人對已有的GALCS算法做出部分改進(jìn)，進(jìn)一步降低了算法復(fù)雜度，該算法查詢質(zhì)量結(jié)果遠(yuǎn)高于Long Run算法。

針對現(xiàn)有算法的不足，本文主要利用啟發(fā)式算法中的蟻群優(yōu)化算法（ACO）求解LCS問題，ACO-LCS在查詢效率和結(jié)果近似度兩方面均高于GA-LCS算法；并對其中的啟發(fā)式因子計(jì)算部分做了改進(jìn)，進(jìn)一步提高算法的查詢結(jié)果質(zhì)量。

1 最長公共子序列問題模型

在計(jì)算機(jī)中，生物序列可表示為字符矩陣S={S0,S1,….,Sn-1}，矩陣每一行表示一條生物序列Si，且對于任意i(0 ≤i≤1)均有 |Si|=1。Si[j]表示序列Si中下標(biāo)為j的字符。子序列為一個(gè)嚴(yán)格遞增的下標(biāo)數(shù)組對應(yīng)的字符序列。例如，對于序列Si=CTATGA，下表數(shù)組p=[0 ,2,4,5]對應(yīng)的子序列為CAGA。同理，N條序列間的長度為D的任意公共子序列S，可以表示為一個(gè)D×N的下標(biāo)矩陣CS，即CS[i][j]表示S的第i個(gè)字符在原始序列Sj中的下標(biāo)。不同的公共子序列具有不同的下標(biāo)矩陣，因此，下標(biāo)矩陣可以作為公共子序列的唯一標(biāo)識。下標(biāo)矩陣CS具有如下兩個(gè)特點(diǎn)：

（1）對于給定的行號i，CS[i][j]和CS[i][k]對應(yīng)的字符相同，其中0≤i

（2）下標(biāo)矩陣每列存儲的下標(biāo)值嚴(yán)格遞增。即CS[i][j]

對 于 給 定 的 序 列 集 合 S={S0：ATCACG，S1：TTCTAG，S2：TCCAGA}，其長度為 4的公共子序列s=TCAG可由一個(gè)4×3的下標(biāo)矩陣表示。該矩陣存儲了公共子序列s中的四個(gè)字符在所有原始序列中的下標(biāo)。例如，矩陣的第0行[1 0 0]表示s的第0個(gè)字符T在原始序列S0中的下標(biāo)為1，在原始序列S1和S2中的下標(biāo)為0。矩陣的第0列則表示公共子序列s中的四個(gè)字符分別對應(yīng)原始序列 S0中的下標(biāo)為 1、2、3、5 的字符。

當(dāng)下標(biāo)矩陣CS的行數(shù)D達(dá)到最大值時(shí)，即為最長公共子序列的字符下標(biāo)矩陣，遍歷矩陣的任意一列即可計(jì)算出序列集合的最長公共子序列。因此，LCS問題可以轉(zhuǎn)換為計(jì)算最大行數(shù)的字符下標(biāo)矩陣的問題。在本文中，我們將利用蟻群優(yōu)化算法來解決該問題。

（3）蟻群優(yōu)化算法在LCS問題中的應(yīng)用

對于給定的序列矩陣S={S0,S1,….,Sn-1}，其最長公共子序列是所有原始序列的子序列，因此通過遍歷任意一條原始序列均可得出序列間的LCS下標(biāo)矩陣。蟻群優(yōu)化算法模擬了一個(gè)規(guī)模為M的人工蟻群A={A0,A1,…,AM-1}，每只螞蟻Ai隨機(jī)地分配一條序列Sj，以遍歷該序列為基準(zhǔn)來查詢序列間的最長公共子序列。

螞蟻Ai把原始序列Sj看作自然環(huán)境中的地面，把序列中每個(gè)字符s及其下標(biāo)p組成的組對(s ,p)視作地面上的一個(gè)位置點(diǎn)。所有滿足下標(biāo)矩陣特點(diǎn)的未訪問位置點(diǎn)為Ai的可轉(zhuǎn)向位置點(diǎn)，即若s為所有原始序列的公共字符，且s在每條原始序列中的下標(biāo)p大于該序列已被訪問位置點(diǎn)的下標(biāo)值，則(s ,p)為可轉(zhuǎn)向位置點(diǎn)。螞蟻Ai模擬自然界中蟻群覓食時(shí)的尋路方式，評估所有位置點(diǎn)的可轉(zhuǎn)向性，若(s ,p)滿足下標(biāo)矩陣的特點(diǎn)，即將字符s在所有原始序列中的下標(biāo)作為新的一行加入螞蟻Ai構(gòu)建的LCS下標(biāo)矩陣中。螞蟻Ai以(s ,p)作為當(dāng)前位置繼續(xù)向前訪問，不斷地將滿足下標(biāo)矩陣特點(diǎn)的位置點(diǎn)對應(yīng)的下標(biāo)添加進(jìn)LCS下標(biāo)矩陣中，直到遍歷完序列Sj的最后一個(gè)位置點(diǎn)，便可得到該螞蟻構(gòu)建的近似LCS下標(biāo)矩陣。圖1為螞蟻Ai的一次構(gòu)建過程。

圖1 LCS下標(biāo)矩陣構(gòu)建過程

初 始 時(shí) 下 標(biāo) 矩 陣 為 空 ，位 置 點(diǎn)(A ,0)，(T ,1)，(C ,2)，(A ,4)均為可轉(zhuǎn)向位置點(diǎn)。假設(shè)Ai選擇轉(zhuǎn)向位置點(diǎn)(A ,0)，則搜索(A ,0)在所有原始序列中的下標(biāo)[0 ,3,1]，并將其作為第0行加入下標(biāo)矩陣中。Ai繼續(xù)評估下一個(gè)位置點(diǎn)(T ,1)，分別以3和1為當(dāng)前下標(biāo)在原始序列S1和S2中搜索字符T，未搜尋到則表示(T ,1)不是當(dāng)前的可轉(zhuǎn)向點(diǎn)。Ai繼續(xù)以[0 ,3,1]為當(dāng)前下標(biāo)評估下一個(gè)位置點(diǎn)(C ,2)，并將其對應(yīng)的下標(biāo)數(shù)組[2 ,4,2]作為新的一行加入LCS下標(biāo)矩陣中。直到S0遍歷結(jié)束螞蟻便可構(gòu)建出路徑(A ,0)→(C ,2)對應(yīng)的行數(shù)為2的下標(biāo)矩陣，其對應(yīng)的公共子序列為AC。然而，若要構(gòu)建出行數(shù)較大的下標(biāo)矩陣，還需要螞蟻有效地選擇“較好”的可轉(zhuǎn)向點(diǎn)。

螞蟻在遍歷原始序列時(shí)，所有滿足下標(biāo)矩陣要求的未訪問位置點(diǎn)均為其候選轉(zhuǎn)向的位置點(diǎn)，螞蟻通過評估轉(zhuǎn)向每個(gè)位置點(diǎn)可獲得更長公共子序列的可能性，來選擇出“較好”的位置點(diǎn)。在圖1中，若螞蟻首次選擇時(shí)放棄位置點(diǎn)(A ,0)，選擇轉(zhuǎn)向位置點(diǎn)(T ,1)，則可構(gòu)建出路徑(T ,1)→(C ,2)→(A ,4)對應(yīng)的行數(shù)為3的下標(biāo)矩陣，其對應(yīng)的公共子序列TCA。因此，字符T是比A更好的候選轉(zhuǎn)向。在實(shí)際LCS問題中，螞蟻通常無法直觀地判斷那個(gè)字符“較好”，只能通過位置點(diǎn)(s ,p)的特性計(jì)算其能獲取更長公共子序列的可能性來覺得是否選擇。這種可能性被稱作字符的轉(zhuǎn)換概率。轉(zhuǎn)換概率的決定因素主要是位置點(diǎn)(s ,p)的信息素因子和啟發(fā)式因子。

位置點(diǎn)(s ,p)的信息素因子實(shí)現(xiàn)了蟻群系統(tǒng)的正反饋性，使每只螞蟻有了經(jīng)驗(yàn)學(xué)習(xí)的能力，從而有效地選擇“較好”字符。螞蟻每次選擇一個(gè)位置點(diǎn)(s ,p)后，便會在改點(diǎn)及其在所有原始序列中的對應(yīng)位置點(diǎn)上排放一定量的信息素。隨著蟻群不斷地選擇。位置點(diǎn)(s ,p)上信息素的積累量，表示為τ(s,p)。τ(s,p)越大表示該字符在蟻群歷史構(gòu)建過程中被選擇的次數(shù)越多，是“較好”字符的可能性也就越大。位置點(diǎn)(s ,p)的啟發(fā)式因子是指局部范圍內(nèi)，轉(zhuǎn)向位置點(diǎn)(s ,p)可獲得更長公共子序列的可能性，表示為η(s,p)。η(s,p)值反映了每只螞蟻的局部搜索能力。信息素因子和啟發(fā)式因子共同決定了螞蟻構(gòu)建LCS過程中的每次字符選擇。啟發(fā)式因子具體計(jì)算公式如式（1）所示:

蟻群系統(tǒng)的交流和優(yōu)化主要由信息素機(jī)制實(shí)現(xiàn)。螞蟻通過每個(gè)位置點(diǎn)的信息素累積量，獲取整個(gè)蟻群對于該位置點(diǎn)的歷史選擇經(jīng)驗(yàn)，從而實(shí)現(xiàn)蟻群系統(tǒng)的學(xué)習(xí)和協(xié)作。每輪構(gòu)建結(jié)束后，蟻群優(yōu)化算法會在所有螞蟻構(gòu)建的LCS下標(biāo)矩陣中，選出行數(shù)最大的一個(gè)作為本輪迭代的最優(yōu)解。同時(shí)維護(hù)一個(gè)歷史最優(yōu)解，若當(dāng)前輪產(chǎn)生的本輪最優(yōu)解行數(shù)大于當(dāng)前歷史最優(yōu)解，則更新歷史最優(yōu)解。在LCS問題中，信息素表示為一個(gè)與序列矩陣S同等規(guī)模的矩陣τ[N][D]，τ[N][D]存儲了字符Si[j]位置的信息素累積量。蟻群在每輪構(gòu)建LCS下標(biāo)矩陣完成之后便會對全局信息素矩陣進(jìn)行更新操作，以反饋蟻群的當(dāng)前構(gòu)建情況，為下一輪構(gòu)建提供經(jīng)驗(yàn)學(xué)習(xí)。信息素的更新主要包括信息素隨時(shí)間揮發(fā)的減少量，以及優(yōu)質(zhì)的LCS下標(biāo)矩陣儲存位置點(diǎn)的排放量。信息素更新公式如（2）所示，其中ρ為每輪構(gòu)建過程中信息素的揮發(fā)率。Δτ[i][j]為本輪構(gòu)建過程中Si[j]位置信息素的總排放量。φ為信息素排放系數(shù)，且φ∈(0,1)，以防信息素排放量過高，使得某條路徑信息素積累量明顯高于其他路徑，造成蟻群系統(tǒng)過早停滯。

2 實(shí)驗(yàn)結(jié)果分析

蟻群優(yōu)化算法是一個(gè)群體智能算法，群體的規(guī)模對算法構(gòu)建的質(zhì)量和效率都有一定影響。蟻群優(yōu)化算法中螞蟻數(shù)量越多時(shí)，構(gòu)建的公共子序列種類也就越多，可查詢到更長的公共子序列的可能性也就越大。但是，當(dāng)構(gòu)建的公共子序列種類數(shù)過多時(shí)，大部分位置點(diǎn)的信息素便會趨于平均，這極大地抑制了蟻群系統(tǒng)的正反饋性，導(dǎo)致算法收斂速度變慢，影響最長公共子序列的查詢效率。而若蟻群規(guī)模太小，算法很容易出現(xiàn)停止現(xiàn)象，不利于最優(yōu)解的構(gòu)建。本實(shí)驗(yàn)以不同規(guī)模的DNA序列集合為例，研究不同蟻群規(guī)模對算法查詢效率和查詢結(jié)果長度的結(jié)果。如表1所示，隨著蟻群規(guī)模不斷增大，算法查詢的公共子序列長度隨之變化，并且查詢時(shí)間消耗成倍增加，但螞蟻構(gòu)建出的下標(biāo)矩陣多樣性增加，查詢到的公共子序列長度也逐漸增加。但當(dāng)蟻群規(guī)模增加到一定程度時(shí)，其查詢出的公共子序列長度增加量逐漸減小，查詢長度提升逐漸變得不明顯。

表1 蟻群規(guī)模對LCS查詢性能影響

3 結(jié)語

本文分析了蟻群優(yōu)化算法查詢最長公共子序列的基本思想和大致流程，采用基于蟻群優(yōu)化算法實(shí)現(xiàn)了多序列間的最長公共子序列查詢。本文通過大量的實(shí)驗(yàn)研究并分析了不同蟻群規(guī)模對查詢算法的性能影響，確保蟻群優(yōu)化算法可有效解決LCS問題。

[1]Dorigo M,Birattari M,Stutzle T.Ant colony optimization[J].IEEE Computational Intelligence Magazine,2006,1（4）:28-39.

[2]Bullnheimer B,Hartl R,Strauss C.A New Rank Based Version of the Ant System.A Computational Study[J].1997,7（1）:25-38.

[3]Stutzle T,Hoos H.MAX-MIN Ant System[J].Future Generation Computer Systems,2000,16（8）:889-914.

[4]Dorigo M,Maniezzo V,Colorni A.Positive Feedback as a Search Strategy[J].,1991.