陳 昆,張小駿

(西安電子科技大學計算機系,陜西西安710071)

AP聚類算法求解植入(l，d)模體識別問題

陳 昆,張小駿

(西安電子科技大學計算機系,陜西西安710071)

模體識別是運用計算機算法尋找一系列功能相近且形式相似的DNA序列片段，從而找出生物信息學中控制基因表達調(diào)控機制的轉(zhuǎn)錄因子結(jié)合位氛，將這種問題轉(zhuǎn)化為AP聚類算法可處理的模型，然后用AP聚類得到穩(wěn)定的候選模體聚類，最終利用貪心算法對問題進行求精，得出一組候選模體集，利用相對熵測度對候選模體集合進行評價并且擇優(yōu)輸出，從而構(gòu)造出一種新的模體識別算法.實驗結(jié)果分別從模擬數(shù)據(jù)和真實數(shù)據(jù)證明了所提算法的有效性.

基因轉(zhuǎn)錄;模體識別;AP聚類算法

0 引言

基因表達在生物上的表現(xiàn)是通過對基因轉(zhuǎn)錄、翻譯等過程形成蛋白質(zhì)分子來實現(xiàn).在DNA序列中可以啟動和控制基因轉(zhuǎn)錄、翻譯過程的相關(guān)結(jié)合位點,被稱為啟動子,這些啟動子序列中的特殊位點被稱作轉(zhuǎn)錄因子結(jié)合位點［1］.模體(motif)識別被作為計算手段來定位轉(zhuǎn)錄因子結(jié)合位點［2］.模體是一個短的DNA序列片段,通常長度是5～25個堿基對(base pairs).模體識別的難點在于模體并不是精確地出現(xiàn)在DNA序列中,而是以退化的形式出現(xiàn).模體識別通常被形式化地定義如下［3］:

植人(l,d)模體識別問題:給定t條長為n的字符表{A,C,G,T}上的DNA序列集合S= {s1,s2,…,st},以及非負整數(shù)l和d,滿足0≤d＜l＜n.模體發(fā)現(xiàn)問題就是找到一個l-mer(長為l的字符串)m,使得每條序列si中都存在一個與m至多有d個位置差異的l-mer mi.我們稱lmer m是一個模體,mi是一個模體實例.

如果植人的模體是充分保守的,那么(l,d)模體發(fā)現(xiàn)問題就可以簡化為計算l-mer在給定的序列S中出現(xiàn)的次數(shù).一般而言,模體識別都是基于序列中l(wèi)-mer的相似性比較.Buhler and Tompa采用了利用概率分析來估計(l,d)模體發(fā)現(xiàn)的求解難度［4］.

現(xiàn)有的模體識別算法主要有:一種是精確類算法,此種算法需遍歷整個樣本空間找到每條模體實例,算法的時間復雜度過高,使實際運用有著較大的困難,具有代表性的是W INNOWER算法［3］和PMSP算法［5］.另一種是非精確算法,此種算法構(gòu)造初始模型,通過模體特性多次迭代找到局部最優(yōu)解,因此時間復雜度較低,不足的是只能確定局部最優(yōu)解.具有代表性的算法是MEME算法［6］、基于模式驅(qū)動的PairMotif+算法［7］、基于Gibbs采樣DNA算法［8］、AP聚類等聚類算法.

筆者的思路是通過AP聚類［9］來解決模體識別問題.聚類算法按照特定標準,對各個樣本進行無監(jiān)督的分類,得出的各個聚類是一組相互相似的l-mer,它們是一組潛在的模體實例.

1 方法

1.1 AP聚類原理

AP聚類又稱近鄰傳播聚類,是2007年Frey等人提出的一種無監(jiān)督聚類算法［9］.通過對目標數(shù)據(jù)集建立相似度矩陣,引人Responsibility和Availability兩種信息傳遞,其中Responsibility表示從數(shù)據(jù)結(jié)點i信息傳向到候選聚類中心k,意味著k點作為i點聚類中心的可能性;Availability表示從候選聚類中心k信息傳向結(jié)點i,意味著i點歸屬于聚類中心k的可能性.初始目標數(shù)據(jù)集中每個元素作為聚類中心,圍繞這兩種信息關(guān)系進行迭代,每次迭代就是信息關(guān)系的進一步更新,直至穩(wěn)定的聚類出現(xiàn).

從理論上看AP聚類有這樣的特點:①不必要求設(shè)定初始的聚類中心;②聚類中心真實存在;③在相似度矩陣確定的情況下,聚類結(jié)果穩(wěn)定.所以較其他的聚類算法,此算法有著較穩(wěn)定的聚類效果,適合模式識別這類大規(guī)模數(shù)據(jù)處理.

1.2 用AP聚類求解模體識別問題

使用AP算法求解模體識別問題的步驟分為三步:第一步對給定數(shù)據(jù)集構(gòu)建圖結(jié)構(gòu);第二步用AP算法找到穩(wěn)定的聚類;第三步,對聚類結(jié)果求精得到植人(l,d)模體.

1.2.1 對目標數(shù)據(jù)集構(gòu)建圖結(jié)構(gòu)

對給定的t條DNA序列建立模體識別問題的圖模型,并定義圖中頂點的相似性,從而得到相似度矩陣S.對于輸人序列,每個長度為l的子序列對應于圖中的一個頂點,那么,模體識別問題就轉(zhuǎn)化成了在圖中尋找稠密子圖的問題.稠密子圖中的大多數(shù)頂點間的權(quán)值都比較大,為了保證聚類時的高效的時間性能,圖的規(guī)模不能太大.筆者采用如下的方式減小圖的規(guī)模:遍歷參考序列中的每個l-mer x,由x誘導出一個子圖G,并對每一個子圖G進行聚類.子圖G的性質(zhì)如下:

①頂點集合包含x,以及si(2≤i≤t)中所有的滿足與x小于等于2d個差異的l-mer.

②來自于同一條輸人序列中的兩個頂點間無邊相連.表示為負無窮大.

③來自于不同輸人序列中的兩個頂點間有邊相連.如果它們間的距離k滿足d＜k＜2d,邊的權(quán)重為l-k,否則邊的權(quán)重為10(l-k).計算邊權(quán)重的方法突顯了相似性高的l-mer在聚類中所起的作用.

1.2.2 用AP算法找到穩(wěn)定的聚類

利用兩個矩陣R和A迭代地計算來完成聚類,矩陣R代表著Responsibility信息傳遞關(guān)系,矩陣A代表著Availability信息傳遞關(guān)系.具體說來,矩陣R中的每個元素r(i,k)表示點xk作為點xi聚類中心的可能性,矩陣A中每個元素a(i,k)表示從點xi選擇xk當作聚類中心的可能性.步驟是:

(1)初始化相似度矩陣,初始設(shè)置a(i,k)= 0,對參考度S［k］［k］賦值,在模體識別應用中通常設(shè)置為中位數(shù).

(2)計算樣本點之間的Responsibility值.

r(i,k)=w(i,k)-max{a(i,k′)+w(i, k′)},1≤k′≤n,k′≠k,

(3)計算樣本點之間的Availability值.

(4)迭代計算.

rnew(i,k)=λrold(i,k)+(1-λ)r(i,k);

anew(i,k)=λaold(i,k)+(1-λ)a(i,k).

(5)輸出結(jié)果.

O=argmaxk{a(i,k)+r(i,k)}.

其中,第4個步驟引人的λ是阻尼系數(shù),λ取值空間為［0.5,1),下面取λ=0.8.

AP聚類算法實現(xiàn)偽代碼為:

A lgorithm 1 AP C luster

APCluster(S,n,k,m)

Input:matrix S

Output:a cluster

(1)iter←0

(2)λ←0.8

(3)while iter＜m and exemp lar changed do

(4) for i←0 to n

(5) for k←0 to n

(6) for j←0 to n

(7) R(i,k)←R(i,k)-max{A(i, j)+S(i,j)}

(8) for i←0 to n

(9) for k←0 to n

(10) for j←0 to n

(11) A(i,k)←min{0,R(k,k)+sum (max{0,R(j,k)})}

(12) for i←0 to n

(13) for k←0 to n

(14) R(i,k)←λR(i,k)+(1-λ)R(i -1,k)

(15) A(i,k)←λA(i,k)+(1-λ)A(i -1,k)

(16) for i←0 to n

(17) put A(i,k)+R(i,k)into cluster

(18) iter←iter+1

(19)return cluster

1.2.3 對聚類結(jié)果求精得到植入(l，d)模體

對得到的聚類結(jié)果即O(它使得所有數(shù)據(jù)點到最近的類代表點的相似度之和最大)中的各個聚類進一步求精來得到最終的模體.

從子圖G中通過AP算法可以得到多個聚類(cluster).得到的cluster雖然包含了一組相互相似的l-mer,但并不一定能夠涵蓋所有的模體實例,需要對聚類進一步求精,得出其他輸人序列中的模體實例.通過算法1對得到的聚類進行貪心地擴展,使得聚類中的元素能夠遍布在每條輸人序列之中.

A lgorithm 2 C luster Refinem ent

Inpu t:a cluster C

Outpu t:a motif

(1)Q←the sequences that contain l-mers in C

(2)while｜Q｜＜t do

(3) select a random sequence siin S-Q

(4) find an l-mer x in sisuch that the distance from x to the l-mers in C is smallest

(5) add x to C

(6) add sito Q

(7)get a motif m by aligning the l-mers in C

(8)return m

每個得到的聚類經(jīng)過擴展之后,對齊聚類中的l-mer,之后便可以得到一個模體.

1.3 整體算法

模體識別的整體算法如下所示.圖1為整體算法的流程圖.

A lgorithm 3 M otif D iscovery A lgorithm

Inpu t:l,d,S={s1,s2,…,st}

Output:(l,d)motifs

(1)M←Φ

(2)select a random reference sequence sifrom S

(3)for each l-mer x in sido

(4) construct a graph induced by x

(5)call Algorithm 1 and get the cluster of the graph using AP algorithm

(6)for each obtained cluster C do

(7) if｜C｜＞10 then

(8)refine C using Algorithm 2 and get a motif m

(9)add m to M

(10)sort themotifs in M

(11)return motifs with high score用相對嫡對模體進行打分排序［10］:

式中:pjk是輸人序列中所有模體出現(xiàn)的第j列上堿基rk∈{A,C,G,T}的觀測頻率;bk是堿基rk的背景頻率.相對嫡通過用于測量模體的保守性,以及模體與背景分布差異程度.

圖1 整體算法的流程圖Fig.1 Flow chart of the whole algorithm

2 實驗與討論

2.1 數(shù)據(jù)和評價方法

分別使用模擬數(shù)據(jù)和真實數(shù)據(jù)對所提算法進行測試.對于模擬數(shù)據(jù),采用文獻［3］中提供的標準生成方法:隨機生成t條獨立同分布的DNA序列和一個長為l的模體m;在每條DNA序列中的一個隨機位置上植人一個模體m的實例m′,使得m′和m的最大距離為d.

對于真實數(shù)據(jù),選用了preproinsulin、DHFR、c-fos、metallothionein和Yeast ECB5等5組生物基因序列.

用C++實現(xiàn)了所提算法,并與比較算法在相同的實驗環(huán)境下進行了測試:惠普工作站,CPU為英特爾至強(W3505)雙核2.53 GHz,內(nèi)存4GB.模擬數(shù)據(jù)的所得結(jié)果是在五組隨機數(shù)據(jù)上計算結(jié)果的平均值.采用核苷酸水平的性能系數(shù)(NPC)來評估識別準確率［11］:

式中:nTP表示預測的位點中所包含的已知位點的數(shù)量;nFP表示預測的位點中不是已知位點的數(shù)量;nFN表示已知的位點中沒有被預測出的位點數(shù)量;nPC在準確性測量中可以同時體現(xiàn)特異性和敏感性.NPC的取值范圍是0到1,取值越高,準確率越高.

2.2 模擬數(shù)據(jù)集上的結(jié)果

選取了3個典型的算法與所提算法進行比較,即A lignACE、MotifCut和VINE.分別是基于Gibbs采樣的統(tǒng)計方法、基于圖聚類的算法和啟發(fā)式方法.

首先,通過植人不同的(l,d)實例來對算法進行比較.表1為各個算法的預測準確率.在4個比較算法中,相比AlignACE和VINE,所提算法和MotifCut都利用了全局信息,但采用了不同的聚類求精方法.在7組數(shù)據(jù)集中,所提算法有4組結(jié)果優(yōu)于MotifCut,因此,在現(xiàn)實中所提算法的結(jié)果可以與MotifCut的結(jié)果進行相互補充.

表1 不同(l，d)問題實例上的預測準確率Tab.1 Identification accuracy on different(l，d) problem instances

其次,固定(l,d)問題實例為(15,4),表2為不同序列長度的預測準確率.可以發(fā)現(xiàn),所提算法的預測準確率隨著序列長度的增加呈現(xiàn)出一個較緩的降低趨勢,表明所提算法對不同的序列長度具有較好的適應性.

2.3 真實數(shù)據(jù)集上的結(jié)果

本節(jié)驗證本文算法在真實數(shù)據(jù)集中識別模體的有效性.一般而言,真實數(shù)據(jù)集與模擬數(shù)據(jù)集有較大的區(qū)別.真實數(shù)據(jù)集沒有模擬數(shù)據(jù)集規(guī)整,模體的保守性沒有模擬數(shù)據(jù)集強,從而進一步增加了模體識別的難度.取輸出結(jié)果中得分最高的模體為檢測出的模體,如果檢測出的模體所對應的位點能夠與已知的位點有部分重疊,則認為是有效的檢測;也即如果預測準確率(NPC值)大于0,則認為是有效的檢測.

表2 不同序列長度上的預測準確率Tab.2 Identification accuracy on different sequence lengths

表3 真實數(shù)據(jù)集上的結(jié)果Tab.3 Results on real data sets

圖2為檢測出模體的序列l(wèi)ogo圖,它形象地展示出了模體的序列保守性.可以發(fā)現(xiàn),檢測出的模體與公布的模體有較好的相似性.更為具體的,表3給出了檢測每組數(shù)據(jù)集所使用的(l,d)、得到的模體以及預測準確率.其中,檢測出的模體中下劃線部分表示了與公布模體重疊的部分.

圖2 模體的序列l(wèi)ogo圖Fig.2 Sequence logos of detected motifs

從表3可以看出,對于Yest ECB數(shù)據(jù),識別準確率達到了最高,主要原因是其中模體的保守性很強.對于c-fos數(shù)據(jù),識別準確率最低,主要原因是所采用的(9,2)實例對應的序列保守性低,受到背景序列的干擾強.對于DHFH數(shù)據(jù),雖然檢測出的模體與公布的模體能夠完全重疊,但是識別準確率仍然較低,原因是沒有識別出所有的結(jié)合位點.

3 結(jié)論

根據(jù)算法原理和實驗數(shù)據(jù)可以發(fā)現(xiàn),所提AP聚類算法應用于模體識別對模體比較保守且生物背景對模體干擾較小時會取得比較好的效果.筆者雖然提出了用AP聚類算法解決模體識別的一種方法,但實際運用中AP聚類算法時間復雜度較高為O(mkN3),而且只完全適用于OOPS類型的DNA序列,對于ZOOPS類型部分適用.對于TCM類型基本不適用.若使之也可以對后兩種類型適用,則m、k、N的值將更大,時間復雜度更高,所以,所提算法的改進方向是著重對初始的相似度矩陣降維修剪.

［1］ HAESELEER P D.How does DNA sequencemotif discovery work［J］.Nature Biotechnology,2006,24 (8),68-74.

［2］ ZAMBELLI F,PESOLE G.Motif discovery and transcription factor binding sites before and after the nextgeneration sequencing era［J］.Briefings in Bioinformatics,2013,14(2):225-237.

［3］ PEVZNER PA,SZE SH.Combinatorial approaches to finding subtle signals in DNA sequences［C］//Proceedings of the Eighth International Conference on Intelligent Systems for Molecular Biology.California: Spring-verlag,2000:269-278.

［4］ BUHLER J,TOMPA M.Finding motifs using random projections［J］.Journal of Computational Biology, 2002,9(2):225-242.

［5］ DAVILA J,BALLA S,RAJASEKARAN S.Space and time efficient algorithms for planted motif search［C］//Proceedings of the Second International W orkshop on Bioinformatics Research and App lications.［s.l.］:CRC Press Inc,UK.2006:822-829.

［6］ BAILEY T L,ELKAN C.Fitting a mixture model by expectation maximization to discover motifs in biopolymers［C］//Proceedings of the 2nd International Conference on Intelligent Systems for Molecular Biology. Berlin:Spring-verlag,1994:28-36.

［7］ YU Qiang,HUO Hong-wei,ZHANG Yi,et al.Pair-Motif+:a fast and effective algorithm for De Novo motif discovery in DNA sequences［J］.International Journal of Biological Sciences,2013,9(4):412 -424.

［8］ LAWRENCE C E,ALTSCHUL S F,BOGUSKIM S, et al.Detecting subtle sequence signals:a Gibb's sampling strategy formultiple alignment［J］.Science, 1993,262:208-214.

［9］ SBRENDAN J,FERY,DELBERT D.Clustering by passingmessages between data points［J］.SCIENCE, 2007,2(3):315-328

［10］GIULIO P,GIANCARLO M,GRAZIANO P.An algorithm for finding signals of unknown lengthinDNAsequences［J］.Bioinformatics,2001,4(3):207-214.

［11］TOMPA M,LI N,BAILEY T L,et al.Assessing computational tools for the discovery of transcription factor binding sites［J］.Nat Biotechnol,2005,23 (1):137-149.

AP Clustering Algorithm Solving Planted(L,d)Motif Identification

CHEN Kun,ZHANG Xiao-jun
(School of Computer Science,Xidian University,Xi′an 710071,China)

Transcription factors can be combined with the special DNA sequence that can control gene transcription process.The special DNA sequence is called the motifs.The motif identification is to find a set of DNA fragments with both similar functions and similar forms.It plays a crucial role in the research on the structure and function of genes.The problem was converted to themodelwhich can be processed by AP clustering algorithm.Then we get steady candidate motifs by using AP clustering.Finally we use the greedy algorithm to refine the clustering results.We can get a group of candidatemotifs set,evaluate candidatemotifs set by information content and output the optimalmotif set.Thereby the new algorithm is designed for the problem.The experimental results on both simulated data and real data demonstrate the validity of the proposed algorithm.

gene transcription;motif identification;AP clustering algorithm

TP39

A

10.3969/j.issn.1671-6833.2015.03.024

1671-6833(2015)03-0110-05

2015-01-10;

2015-03-21

中央高校基本科研項目(K50513100011)

陳昆(1975-),河南鄭州人,西安電子科技大學博士研究生,高級工程師,研究方向為生物信息學,E-mail:553058474@qq.com.