• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    SegHMC:一種基于Segmental HMM模型的順式調(diào)控模塊識(shí)別算法

    2016-12-17 08:23:58郭海濤霍紅衛(wèi)于強(qiáng)
    自動(dòng)化學(xué)報(bào) 2016年11期
    關(guān)鍵詞:模體集上調(diào)控

    郭海濤 霍紅衛(wèi) 于強(qiáng)

    SegHMC:一種基于Segmental HMM模型的順式調(diào)控模塊識(shí)別算法

    郭海濤1霍紅衛(wèi)1于強(qiáng)1

    順式調(diào)控模塊(Cis-regulatory module,CRM)在真核生物基因的轉(zhuǎn)錄調(diào)控中起著重要作用,識(shí)別順式調(diào)控模塊是當(dāng)前計(jì)算生物學(xué)的一個(gè)重要課題.雖然當(dāng)前有許多計(jì)算方法用于識(shí)別順式調(diào)控模塊,但識(shí)別準(zhǔn)確率仍有待進(jìn)一步提高.將順式調(diào)控模塊的多種特征信息結(jié)合在一起,有助于提高識(shí)別順式調(diào)控模塊的準(zhǔn)確率.基于此,本文提出了一種識(shí)別順式調(diào)控模塊的算法SegHMC(Segmental HMM model for discovery of cis-regulatory module).該算法建立了一種關(guān)于順式調(diào)控模塊識(shí)別問題的Segmental HMM模型,進(jìn)一步擴(kuò)展了順式調(diào)控模塊調(diào)控結(jié)構(gòu)(或調(diào)控語法)的表示,不僅將順式調(diào)控模塊表示為模體(Motif)的組合,還進(jìn)一步將模體共同出現(xiàn)的頻率、模體順序偏好以及順式調(diào)控模塊中相鄰模體間的距離分布等特征引入到順式調(diào)控模塊的調(diào)控語法中.在模擬數(shù)據(jù)集和真實(shí)生物數(shù)據(jù)集上的實(shí)驗(yàn)結(jié)果表明,本文方法識(shí)別順式調(diào)控模塊的準(zhǔn)確率顯著優(yōu)于當(dāng)前的主要方法.

    基因的轉(zhuǎn)錄調(diào)控,模體,Segmental HMM,順式調(diào)控模塊識(shí)別

    在基因的表達(dá)調(diào)控系統(tǒng)中,轉(zhuǎn)錄因子(Transcription factor,TF)通過與所調(diào)控基因附近被稱為轉(zhuǎn)錄因子結(jié)合位點(diǎn)(Transcription factor binding site,TFBS)或模體(Motif)的特定DNA序列片段相結(jié)合,來啟動(dòng)基因的轉(zhuǎn)錄調(diào)控[1?2].在真核生物中,多個(gè)轉(zhuǎn)錄因子對(duì)基因的轉(zhuǎn)錄調(diào)控,并不是孤立進(jìn)行的,而是轉(zhuǎn)錄因子之間或者各個(gè)轉(zhuǎn)錄因子與它們的模體之間通過一系列的時(shí)空交互來實(shí)施更復(fù)雜、更精確的轉(zhuǎn)錄調(diào)控.在被調(diào)控基因的調(diào)控區(qū)(Transcriptional regulatory region)中,模體非均勻地聚集為一系列的稱為順式調(diào)控模塊(Cis-regulatory module,CRM)的離散區(qū)域,如啟動(dòng)子(Promoter)、增強(qiáng)子(Enhancer)、沉寂子(Silencer)、絕緣子(Insulator)等.結(jié)合這些順式調(diào)控模塊的轉(zhuǎn)錄因子通過相互協(xié)作、相互競爭,激活或抑制所調(diào)控基因的轉(zhuǎn)錄表達(dá).一個(gè)順式調(diào)控模塊包含單個(gè)或多個(gè)轉(zhuǎn)錄因子的多個(gè)模體實(shí)例,其長度通常約為幾百到幾千個(gè)堿基對(duì)(Base pair,bp).一個(gè)真核生物基因調(diào)控區(qū)序列可能的調(diào)控結(jié)構(gòu)的簡單例子如圖1所示.

    圖1 順式調(diào)控模塊結(jié)構(gòu)示意圖(順式調(diào)控模塊是包含多個(gè)轉(zhuǎn)錄因子相應(yīng)模體的序列區(qū);模體的方向、模體間的間隔距離、模體間的相互關(guān)系可能包含了給定順式調(diào)控模塊的重要性質(zhì).)Fig.1 The structure discription of cis-regulatory modules(A cis-regulatory module is a sequence region that contains multiple motifs of multiple transcription factors;motif orientation,the interval distance between motifs and their cooperation relationship may imply the important regulatory properties of the cis-regulatory module.)

    識(shí)別順式調(diào)控模塊是理解基因轉(zhuǎn)錄調(diào)控分子機(jī)制的基礎(chǔ),同時(shí)也是構(gòu)建基因調(diào)控網(wǎng)絡(luò)[3?4]的關(guān)鍵步驟.此外,識(shí)別具有特定調(diào)控功能的順式調(diào)控模塊對(duì)疾病機(jī)理的研究也有重要的意義.許多疾病的發(fā)生都與基因的異常表達(dá)有關(guān),調(diào)控基因表達(dá)的順式調(diào)控模塊發(fā)生變異是造成基因異常表達(dá)的主因.有證據(jù)表明特定順式調(diào)控模塊中協(xié)作調(diào)控元素的破壞,可以導(dǎo)致畸形和疾病;例如,Kleinjan等[5]發(fā)現(xiàn)PAX6的任何遠(yuǎn)端調(diào)控元素的缺失都會(huì)改變其表達(dá)水平,從而造成先天性眼球畸形、無虹膜以及大腦缺陷等疾病.

    通過生物實(shí)驗(yàn),例如高通量檢測與順式調(diào)控模塊相關(guān)的表觀遺傳標(biāo)記特征[6],可以識(shí)別順式調(diào)控模塊,但這種方法費(fèi)時(shí)費(fèi)力代價(jià)較大,并且許多時(shí)候受限于實(shí)驗(yàn)條件而很難實(shí)施.因此,使用計(jì)算方法直接從DNA序列中識(shí)別順式調(diào)控模塊已成為一個(gè)非常有吸引力的手段.

    然而,使用計(jì)算方法識(shí)別順式調(diào)控模塊也存在著許多挑戰(zhàn):1)同一個(gè)順式調(diào)控模塊不同實(shí)例內(nèi)的模體排列順序并不完全相同,但也并非完全無序的.此外,模塊內(nèi)的模體之間的距離也不確定,即同一模塊不同實(shí)例內(nèi)相同的兩相鄰模體間的距離也都不相同.因此,很難確定性地刻畫這種結(jié)構(gòu).2)真核生物的調(diào)控區(qū)通常很長,構(gòu)成順式調(diào)控模塊的模體通常較短且存在退化,根據(jù)已知的模體或者借助于現(xiàn)有的模體庫(如TRANSFAC[7]、JASPAR[8])直接搜索,會(huì)找出大量的假陽性匹配,所以很難通過直接搜索相關(guān)模體的方式來識(shí)別包含這些模體的順式調(diào)控模塊.

    至今,已有多種用于識(shí)別順式調(diào)控模塊的模型和方法[6,9?27].為了識(shí)別真核生物基因的順式調(diào)控模塊,不同方法利用順式調(diào)控模塊的不同特征(如模體的聚集和物種間的保守性),使用不同搜索策略.

    其中一類方法基于窗聚集,利用模體傾向于聚集的特性來搜索順式調(diào)控模塊.這類方法用最簡單的方式表示順式調(diào)控模塊,通過概率統(tǒng)計(jì)度量給定長度窗口內(nèi)的模體組合的統(tǒng)計(jì)顯著性,相應(yīng)方法如MSCAN[22]和MCAST[28]等,或者使用組合方法,在指定窗口大小范圍內(nèi)搜索在多個(gè)序列共同出現(xiàn)一組模體實(shí)例的最小區(qū)域,將其作為候選順式調(diào)控模塊,如CMStalker[29]等.這類方法本質(zhì)上假定了序列窗內(nèi)的模體之間獨(dú)立同分布.這類方法雖然簡單直接,但需要合理確定窗大小以及度量統(tǒng)計(jì)顯著性的打分閾值,而這些參數(shù)在實(shí)際應(yīng)用中通常很難確定;此外,這類方法也忽略了順式調(diào)控模塊可能的調(diào)控結(jié)構(gòu)(或調(diào)控語法),如模體間的順序和距離.

    另一類方法基于概率模型,通過對(duì)序列或順式調(diào)控模塊建立概率模型,進(jìn)而找出待搜索的目標(biāo)序列中的順式調(diào)控模塊.基于概率模型的方法,除了少數(shù)采用判別模型的方法,如HexDiff[30]、Regulatory Potential[31]等外,大部分的方法使用生成模型,主要是隱馬爾科夫模型(Hidder Markov model, HMM).HMM模型的主要優(yōu)勢在于,它可以對(duì)順式調(diào)控模塊的出現(xiàn)進(jìn)行可靠的統(tǒng)計(jì)度量,并能刻畫順式調(diào)控模塊的調(diào)控語法.此外,HMM模型所使用的期望最大化的參數(shù)估計(jì)算法,可以自動(dòng)調(diào)節(jié)大量的參數(shù),避免了手動(dòng)設(shè)置的麻煩.基于HMM的順式調(diào)控模塊識(shí)別方法通常將順式調(diào)控模塊看作由一組過表達(dá)的模體和背景組合生成的序列片段.與窗聚集方法相比,它們不僅考慮了構(gòu)成順式調(diào)控模塊的模體組合,也同時(shí)考慮了構(gòu)成順式調(diào)控模塊的模體之間的距離.最初的一些方法,如CisModule[20],使用HMM間接捕捉順式調(diào)控模塊內(nèi)部以及順式調(diào)控模塊之間背景的概率分布.但這種方法僅使用了一般的順式調(diào)控模塊內(nèi)部背景,并未推斷順式調(diào)控模塊內(nèi)的模體之間的任何順序.后續(xù)的方法,進(jìn)一步擴(kuò)展這種模型的表示,如Stubb[32]方法,創(chuàng)建了一個(gè)僅包含模體和背景兩個(gè)狀態(tài)的HMM模型,使用統(tǒng)計(jì)方法度量模體對(duì)間共同出現(xiàn)的顯著性,進(jìn)而決定是否引入相應(yīng)的轉(zhuǎn)移概率.通過使用定義在指定長度窗口內(nèi)序列上的打分函數(shù)度量窗口內(nèi)模體聚集顯著性來預(yù)測順式調(diào)控模塊.該方法僅利用了HMM模型的轉(zhuǎn)移概率,沒有使用任何其他HMM模型特性.后來的方法CORECLUST[33]和TSHAS[14]建立了更復(fù)雜的HMM模型,引入了順式調(diào)控模塊內(nèi)部背景狀態(tài),加入了模體到模體的概率轉(zhuǎn)移,更細(xì)致地刻畫了順式調(diào)控模塊的調(diào)控結(jié)構(gòu).訓(xùn)練和解碼算法使用標(biāo)準(zhǔn)的Baum-Welch算法[34]和Veterbi算法[34].這類方法使用了HMM的特性,但局限于HMM的表達(dá)能力,建立的調(diào)控模型并不直觀,添加了大量的輔助狀態(tài).另一類HMM相關(guān)的順式調(diào)控模塊識(shí)別方法,使用加強(qiáng)的HMM來刻畫順式調(diào)控模塊的調(diào)控結(jié)構(gòu);如BayCis算法[35],使用貝葉斯層次HMM模型,對(duì)順式調(diào)控模塊和包含順式調(diào)控模塊調(diào)控序列進(jìn)行建模.利用了層次HMM的特性建立了模體之間的概率轉(zhuǎn)移.該模型將HMM狀態(tài)轉(zhuǎn)移參數(shù)看作隨機(jī)變量,并引入貝葉斯先驗(yàn).該模型雖然結(jié)構(gòu)直觀,表達(dá)能力強(qiáng),但模型的訓(xùn)練和解碼需要大量的計(jì)算.Lemnian等提出的基于Extended sunflower HMM的順式調(diào)控模塊識(shí)別方法[19],使用Extended sunflower HMM,對(duì)模型的刻畫深入到模體內(nèi)部,不僅刻畫了模體間的依賴,還刻畫了模體內(nèi)部堿基之間的依賴關(guān)系,但該方法僅用于同型順式調(diào)控模塊的識(shí)別.

    此外,還有一類方法利用相近物種的進(jìn)化保守性來識(shí)別順式調(diào)控模塊,如 MorphMS[26]、MultiModule[36]和ReLA[27]等.這類方法首先通過雙序列或多序列比對(duì)同源基因的調(diào)控區(qū)找出其中的保守區(qū)域,然后使用其他方法在保守區(qū)域中搜索順式調(diào)控模塊.由于大多數(shù)基因的調(diào)控區(qū)中存在著大量重復(fù)(Duplication)和改組(Shuffling)的序列片段,很難進(jìn)行序列比對(duì),所以這類方法并不總是有效.

    將順式調(diào)控模塊的多種特征信息結(jié)合在一起,有助于提高識(shí)別順式調(diào)控模塊的準(zhǔn)確率.順式調(diào)控模塊雖然結(jié)構(gòu)具有不確定性,很難刻畫,但作為頻繁出現(xiàn)在多個(gè)被調(diào)控基因調(diào)控區(qū)中的調(diào)控功能單位,很可能包含某些保守成分.MOPAT[37]從一組同源基因的調(diào)控區(qū)中搜索保守模體模塊(這里的保守模體模塊被定義為在一組同源基因調(diào)控區(qū)中頻繁出現(xiàn)且具有一定距離約束的相鄰模體對(duì),不考慮兩個(gè)模體的相對(duì)順序)出發(fā),查找包含多個(gè)這種保守模體模塊的區(qū)域,將其作為候選順式調(diào)控模塊,找出了一些保守的順式調(diào)控模塊.這一事實(shí)間接說明了順式調(diào)控模塊中存在保守成分.這里,我們同樣利用保守模體模塊這種順式調(diào)控模塊的保守成分,但限定保守模體模塊內(nèi)的模體具有特定的次序.然后,進(jìn)一步將順式調(diào)控模塊保守性假設(shè)從同源基因(不同物種)推廣到共調(diào)控基因(同一物種).最終,我們將順式調(diào)控模塊表示為由單模體和保守模體模塊混合組成的,具有部分保守特征的調(diào)控結(jié)構(gòu)(也稱調(diào)控語法),從而將順式調(diào)控模塊結(jié)構(gòu)保守性和其內(nèi)部模體傾向于聚集的特征結(jié)合起來.為了刻畫這種復(fù)雜的順式調(diào)控模塊調(diào)控結(jié)構(gòu),我們使用一種被稱為Segmental HMM[38]的增強(qiáng)HMM模型來表達(dá).

    基于此,本文提出了一種識(shí)別順式調(diào)控模塊的概率模型方法SegHMC(Segmental HMM model for discovery of cis-regulatory module).該方法使用Segmental HMM在給定候選模體集上構(gòu)建同源或共調(diào)控基因調(diào)控區(qū)序列和順式調(diào)控模塊的調(diào)控語法結(jié)構(gòu).同一般的識(shí)別順式調(diào)控模塊的HMM模型相比,我們不僅將順式調(diào)控模塊表示為模體的組合,還將模體共同出現(xiàn)的頻率、模體順序偏好以及順式調(diào)控模塊中的相鄰模體之間距離分布等特征引入到順式調(diào)控模塊的調(diào)控語法當(dāng)中,這些特征可以有效提高順式調(diào)控模塊的識(shí)別精度.此外,為了處理真核生物基因長的調(diào)控區(qū),我們對(duì)模型進(jìn)行了降低搜索空間的優(yōu)化.這種優(yōu)化通過提前進(jìn)行片段分割,顯式建立Segmental HMM狀態(tài)轉(zhuǎn)換圖,去除了大量不必要的搜索路徑,降低了搜索空間,同時(shí)又不失精度.得到的模型可用于待搜索目標(biāo)基因調(diào)控區(qū)中甚至整個(gè)基因組中的相似順式調(diào)控模塊識(shí)別.我們分別在一個(gè)模擬數(shù)據(jù)集和兩個(gè)真實(shí)生物數(shù)據(jù)集: Muscle數(shù)據(jù)集和果蠅早期發(fā)育數(shù)據(jù)集上對(duì)我們的方法進(jìn)行測試,并選取當(dāng)前主要方法進(jìn)行比較,所有方法識(shí)別順式調(diào)控模塊的準(zhǔn)確率使用通用評(píng)價(jià)指標(biāo)相關(guān)系數(shù)(Correlation coefficient,CC)和F1-score來度量.實(shí)驗(yàn)結(jié)果表明,我們的方法識(shí)別順式調(diào)控模塊的準(zhǔn)確率顯著優(yōu)于當(dāng)前的主要方法.

    1 SegHMC算法

    1.1 SegHMC的Segmental HMM模型

    Segmental HMM[38]是HMM的一個(gè)擴(kuò)展,也稱Generalized HMM.與一般HMM的每個(gè)狀態(tài)僅能發(fā)射一個(gè)堿基相比,Segmental HMM的每個(gè)狀態(tài)可以發(fā)射可變長度的堿基序列片段;狀態(tài)所發(fā)射的堿基序列可由一個(gè)片段模型來表示.該片段模型,給出了生成長度為u的觀察序列o=o1o2···ou的聯(lián)合概率,可由下式表示.

    因此,片段模型由兩個(gè)分布組成,一個(gè)是描述片段長度似然的片段長度分布ds(u),另一個(gè)為表示不同長度觀察序列發(fā)射概率的發(fā)射模型es(o).因此,在順式調(diào)控模塊的識(shí)別模型中,可以根據(jù)對(duì)順式調(diào)控模塊和調(diào)控序列結(jié)構(gòu)的抽象,對(duì)這兩個(gè)分布給出具體的定義.

    本文使用Segmental HMM,在片段層次上對(duì)順式調(diào)控模塊和調(diào)控序列的調(diào)控結(jié)構(gòu)進(jìn)行建模,具有更強(qiáng)的表達(dá)能力;例如,可以對(duì)片段之間的依賴進(jìn)行建模.本節(jié)將詳細(xì)闡述Segmental HMM模型,該模型主要包括:模型的構(gòu)建、狀態(tài)轉(zhuǎn)移概率、片段長度分布和生成狀態(tài)的發(fā)射模型.

    1.1.1 Segmental HMM模型構(gòu)建

    我們將轉(zhuǎn)錄調(diào)控序列的調(diào)控結(jié)構(gòu)定義如下.轉(zhuǎn)錄調(diào)控序列由一系列的順式調(diào)控模塊和順式調(diào)控模塊之間的背景(稱為全局背景)構(gòu)成,而每個(gè)順式調(diào)控模塊又由一組具有特定次序的模體和模體之間的背景(稱為局部背景)構(gòu)成,這種抽象具有明顯的層次性.基于這種結(jié)構(gòu)定義,給定的轉(zhuǎn)錄調(diào)控序列可由下列過程生成:

    1)定位給定候選模體集中的模體在目標(biāo)轉(zhuǎn)錄調(diào)控序列中所有可能出現(xiàn)實(shí)例;

    2)以這些被定位的模體實(shí)例為錨點(diǎn),使用兩模體實(shí)例之間的背景(全局背景或局部背景,具體類別待定)序列連接這些模體實(shí)例,從而生成整個(gè)調(diào)控序列.

    上述過程中,我們?cè)试S模體實(shí)例在空間上存在重疊,模體之間可能通過多種類型的背景序列相連;因此,存在許多平行的生成路徑.從這些路徑中,找出最可能的生成路徑,即可得出該轉(zhuǎn)錄調(diào)控序列最可能的調(diào)控結(jié)構(gòu),從而找出相應(yīng)的順式調(diào)控模塊.

    將每個(gè)具體片段(模體、全局背景和局部背景)表示為Segmental HMM的一個(gè)狀態(tài),片段之間的連接對(duì)應(yīng)了兩個(gè)狀態(tài)之間的轉(zhuǎn)移,根據(jù)上述生成過程,我們顯式構(gòu)造Segmental HMM的狀態(tài)轉(zhuǎn)換圖.顯式構(gòu)建狀態(tài)轉(zhuǎn)換圖一方面移除了不必要的狀態(tài)路徑,減小算法的搜索空間;另一方面,更便于構(gòu)建模體的二元語法,模型順式調(diào)控模塊內(nèi)的相鄰模體間的一階依賴關(guān)系.為了標(biāo)識(shí)順式調(diào)控模塊,我們?cè)黾酉鄳?yīng)的輔助狀態(tài):順式調(diào)控模塊開始狀態(tài)和順式調(diào)控模塊結(jié)束狀態(tài).圖2給出了表示一個(gè)調(diào)控序列調(diào)控語法結(jié)構(gòu)的Segmental HMM狀態(tài)轉(zhuǎn)換圖的具體例子,整個(gè)模型所包含的狀態(tài)如下:

    2)模體狀態(tài)M={m1,m2,···,mK};

    4)順式調(diào)控模塊狀態(tài)C,又由順式調(diào)控模塊開始狀態(tài)Cs和順式調(diào)控模塊結(jié)束狀態(tài)Ce構(gòu)成,即C=Cs∪Ce=

    Segmental HMM狀態(tài)轉(zhuǎn)換圖的具體構(gòu)造過程如算法1所示.

    算法1.Segmental HMM狀態(tài)轉(zhuǎn)換圖的構(gòu)造

    輸入:一組Motif的PWM集PWMS,用于搜索Motif的p-value閾值和一個(gè)調(diào)控序列

    輸出:狀態(tài)轉(zhuǎn)換圖的狀態(tài)集Q和這些狀態(tài)之間的連接集T

    1)創(chuàng)建該模型的初始狀態(tài) S和終止?fàn)顟B(tài) E

    2)對(duì)PWMS中每個(gè)PWM,在所給調(diào)控序列中找出小于給定p-value閾值所有Motif匹配

    3)根據(jù)找出的Motif匹配對(duì)所給調(diào)控序列進(jìn)行分割,標(biāo)記狀態(tài)類型,創(chuàng)建相應(yīng)的狀態(tài)集M,Bg和Cs

    5)以狀態(tài)在序列中的位置為關(guān)鍵字對(duì)狀態(tài)集Q進(jìn)行排序

    6)Ce←?

    7)Bc←?

    8)for每個(gè)狀態(tài)qi∈Q do

    9)if qi為模型的初始狀態(tài)

    10)從Q中順序取出下一狀態(tài)qi+1

    11)T←T∪{qi→qi+1}

    12)else if qi為一個(gè)全局背景狀態(tài)then

    13)從Q中找出qi位置之后的第一個(gè)全局背景狀態(tài),記為qj

    14)T←T∪{qi→qj}

    15)從Q中找出qi位置之后的第一個(gè)CRM初始狀態(tài),記為qj

    16)T←T∪{qj→qi}

    17)for qi的每個(gè)前端模體狀態(tài)m do

    18)創(chuàng)建一個(gè)CRM終止?fàn)顟B(tài)ce

    19)Ce←Ce∪{ce}

    20)T←T∪{m→ce}

    21)T←T∪{ce→qi}

    22)else if qi是一個(gè)CRM初始狀態(tài)then

    23)T←T∪{qi→mi}

    24)else if qi是一個(gè)Motif狀態(tài)then

    25)從Q中找出qi位置之后且不與它重疊的下一Motif狀態(tài),記為qj

    26)創(chuàng)建一個(gè)局部背景狀態(tài)bc

    27)Bc←Bc∪{bc}

    28)T←T∪{qi→bc}

    29)T←T∪{bc→qj}

    30)Q←Q∪Ce∪Bc

    31)以狀態(tài)在序列中的位置為關(guān)鍵字對(duì)狀態(tài)集Q進(jìn)行重新排序

    32)return Q和T

    圖2 Segmental HMM狀態(tài)轉(zhuǎn)移圖Fig.2 The state transition diagram of segmental HMM

    在我們的模型中使用位置權(quán)重矩陣(Position weight matrice,PWM)[39]表示相應(yīng)的模體,在算法中提前給定待搜索順式調(diào)控模塊所包含的可能模體集,所對(duì)應(yīng)的PWM集表示為PWMS.算法所要求的其他輸入包括:搜索模體的p-value,以及待建模的轉(zhuǎn)錄調(diào)控序列.

    在上述Segmental HMM狀態(tài)轉(zhuǎn)換圖的構(gòu)造算法中,第1行創(chuàng)建模型的初始狀態(tài)和結(jié)束狀態(tài);第2行,根據(jù)所給p-value找出給定模體集中模體及其反向互補(bǔ)的模體在序列中所有的出現(xiàn)實(shí)例.第3行,根據(jù)所找出的模體實(shí)例,分別構(gòu)建模體狀態(tài)、全局背景狀態(tài)和順式調(diào)控模塊開始狀態(tài),對(duì)應(yīng)的狀態(tài)集分別為M、Bg和Cs.第5行對(duì)所有狀態(tài)的集合進(jìn)行排序.第6~7行,分別初始化順式調(diào)控模塊結(jié)束狀態(tài)集Ce和模體間局部背景狀態(tài)集Bc.第8~31行,確定各狀態(tài)間的轉(zhuǎn)移,具體為:對(duì)于模型的初始狀態(tài)(對(duì)應(yīng)于第9~11行),只需連接下一任意有效狀態(tài);對(duì)于全局背景狀態(tài),需要與下一全局背景狀態(tài)(對(duì)應(yīng)于第13~14行)、相鄰順式調(diào)控模塊初始狀態(tài)(對(duì)應(yīng)于第15~16行)和順式調(diào)控模塊終止?fàn)顟B(tài)(對(duì)應(yīng)于第17~21行)相連;對(duì)于順式調(diào)控模塊初始狀態(tài),則只需連接到對(duì)應(yīng)的模體狀態(tài),這對(duì)應(yīng)于第22~23行;對(duì)于模體狀態(tài),對(duì)應(yīng)于第24~29行,找出后續(xù)與當(dāng)前狀態(tài)不重疊的模體狀態(tài),然后創(chuàng)建相應(yīng)的局部背景狀態(tài),并依次連接這些狀態(tài).第30~31行,將新創(chuàng)建的順式調(diào)控模塊結(jié)束狀態(tài)和模體間局部背景狀態(tài)加入到總狀態(tài)集Q,并對(duì)Q重新排序.

    關(guān)于上述構(gòu)造算法的幾點(diǎn)說明:

    1)在第3行中,對(duì)每個(gè)模體,分別創(chuàng)建了位于模體前后的兩個(gè)可能的全局背景狀態(tài).所創(chuàng)建的全局背景狀態(tài),僅有一端位置是確定的(即模體的前端或模體的后端),在后面第13~14行的操作中,會(huì)進(jìn)一步將這些半連接的全局背景片段連接起來,形成一個(gè)大的全局背景;

    2)在第5行和第31行中,按照兩個(gè)關(guān)鍵字(位置、狀態(tài)的類型)將前面分別生成的、無次序的各種狀態(tài)按照生成轉(zhuǎn)錄調(diào)控序列的時(shí)空順序進(jìn)行排序,以便于后面確定各狀態(tài)的連接轉(zhuǎn)移的操作;

    3)在算法中,我們顯式地創(chuàng)建順式調(diào)控模塊的開始狀態(tài)和結(jié)束狀態(tài),一方面,可以使結(jié)構(gòu)更清晰,另一方面,也便于在推斷時(shí)確定順式調(diào)控模塊的邊界;

    4)對(duì)于第27行,為了簡化順式調(diào)控模塊的結(jié)構(gòu)表示模型,我們假設(shè)順式調(diào)控模塊內(nèi)相鄰的模體間是非重疊的;

    5)假定所給序列的長度為T,在整個(gè)算法中,耗時(shí)的操作主要集中在:模體的查找,最壞時(shí)間復(fù)雜度為O(KT),其中K為PWM的個(gè)數(shù);查找后繼狀態(tài)的操作,最壞時(shí)間復(fù)雜度為O(T2);對(duì)狀態(tài)集的排序,最壞時(shí)間復(fù)雜度為O(T2).因此,整個(gè)算法的時(shí)間復(fù)雜度為O(T2).

    1.1.2 狀態(tài)轉(zhuǎn)移概率

    在Segmental HMM狀態(tài)轉(zhuǎn)換圖中,每個(gè)狀態(tài)對(duì)應(yīng)模體、全局背景或局部背景這些類型的一個(gè)實(shí)例,狀態(tài)之間的轉(zhuǎn)移概率即為相應(yīng)狀態(tài)類型之間轉(zhuǎn)移的概率.

    對(duì)于順式調(diào)控模塊內(nèi)的模體狀態(tài)mi和模體狀態(tài)mj之間的轉(zhuǎn)移概率,可由下式估計(jì)得到:

    其中,t(m)表示模體狀態(tài)(對(duì)應(yīng)于模體實(shí)例)m所對(duì)應(yīng)的模體類型,A為模體狀態(tài)間的轉(zhuǎn)移計(jì)數(shù).

    對(duì)于全局背景類型狀態(tài)bg,其反映了在長的序列區(qū)域中出現(xiàn)順式調(diào)控模塊的概率.由狀態(tài)bg到順式調(diào)控模塊狀態(tài)的轉(zhuǎn)移概率,或順式調(diào)控模塊到bg狀態(tài)的轉(zhuǎn)移概率,由于即使在很長的調(diào)控序列中順式調(diào)控模塊的數(shù)目也相對(duì)較少,所能獲得的數(shù)據(jù)難以訓(xùn)練出可靠的模型參數(shù).為避免過擬合,可由經(jīng)驗(yàn)估計(jì)得出,作為常量參數(shù),在系統(tǒng)運(yùn)行時(shí)設(shè)定.

    1.1.3 片段長度分布

    全局背景長度和局部背景長度分別表示了順式調(diào)控模塊之間以及順式調(diào)控模塊內(nèi)的模體之間的空白區(qū)域的長度分布.對(duì)于全局背景狀態(tài)bg和局部背景狀態(tài)bc,我們假定其序列長度分別滿足期望為wg和wc的幾何分布;這種假定一方面反映了我們對(duì)順式調(diào)控模塊結(jié)構(gòu)的不確定性,另一方面又為模型順式調(diào)控模塊內(nèi)的模體的二元語法特征提供足夠的適應(yīng)性.在該假設(shè)下,背景序列長度為d的概率為:

    這里,b表示bg或bc,wi表示wg或wc.

    對(duì)于模體狀態(tài)m,由于模體所對(duì)應(yīng)的位置權(quán)重矩陣[39]是直接從數(shù)據(jù)庫中獲取的,其長度wm及其特定位置堿基的概率都是已知的,所以模體狀態(tài)上的序列長度d′的概率分布是特定的,即:

    1.1.4 生成狀態(tài)的發(fā)射模型

    在本文模型中,只有模體狀態(tài)、全局背景狀態(tài)和局部背景狀態(tài)為生成狀態(tài).每種生成狀態(tài)發(fā)射長度滿足特定分布的堿基序列片段.

    對(duì)于全局背景和局部背景狀態(tài),我們分別使用k階Markov模型和m階的局部Markov模型.在局部Markov模型中,位置t處堿基的條件概率僅由以位置t為中心長度為2D的窗口內(nèi)的序列片段來估計(jì).這可采用記筆記的方式預(yù)先計(jì)算出每個(gè)位置堿基的條件概率,并存儲(chǔ)計(jì)算的結(jié)果,需要時(shí)直接查表即得.

    對(duì)于模體的生成概率,在本文模型中,使用經(jīng)典的PM 模型[40].假定模體實(shí)例為O,該模體的PWM Θ=[θ1,θ2,···,θL],其中θi(1≤i≤L)為堿基頻率的列向量,則模體狀態(tài)所對(duì)應(yīng)的堿基序列片段為O的概率為:

    這里oi為模體實(shí)例O中第i位置的堿基.

    1.2 解碼和訓(xùn)練算法

    在我們的模型中,將輸入的序列分為訓(xùn)練集和測試集.在訓(xùn)練集上訓(xùn)練出模型參數(shù)后,使用已訓(xùn)練的模型識(shí)別給定測試集中所有序列的順式調(diào)控模塊,這一過程表現(xiàn)為解碼出模型的最優(yōu)狀態(tài)路徑過程.在Segmental HMM模型中,最優(yōu)狀態(tài)路徑可形式化定義為:

    給定長度為T的轉(zhuǎn)錄調(diào)控序列(即觀測序列) O=o1o2···oT,設(shè)其對(duì)應(yīng)的狀態(tài)序列為Π=(π1,···,πT),則該轉(zhuǎn)錄調(diào)控序列所對(duì)應(yīng)的最優(yōu)狀態(tài)序列可表示為:

    這里 Psi(di)表示狀態(tài) si的序列片段長度的概率分布,asi,si+1為狀態(tài)si到狀態(tài)si+1的轉(zhuǎn)移概率,e(oti+1···ti+1|si)表示狀態(tài)si生成觀測序列片段oti+1···ti+1的概率.

    由于缺少足夠的標(biāo)注數(shù)據(jù),本文模型使用無監(jiān)督的Baum-Welch算法[34]直接從訓(xùn)練集中訓(xùn)練系統(tǒng)的模型參數(shù).對(duì)于模型的初始狀態(tài)概率,由于它只確定了輸入序列的第一個(gè)位置的初始功能狀態(tài),在沿序列的后續(xù)操作中,其影響完全可以忽略,在本文模型中簡單地由均勻分布隨機(jī)生成.

    我們已提前標(biāo)出對(duì)應(yīng)片段的可能狀態(tài),創(chuàng)建了Segmental HMM的狀態(tài)轉(zhuǎn)換圖.因此,在解碼時(shí),不再需要通過使用像最大似然之類的方法去推斷最可能的片段分割位置,可直接使用解碼算法找出最優(yōu)路徑.為了求式(6)所對(duì)應(yīng)的最優(yōu)狀態(tài)路徑,本文使用基于動(dòng)態(tài)規(guī)劃的Veterbi算法[34],記為SegHMC Veterbi,并把它作為模型的缺省設(shè)置.此外,為了提供足夠的彈性,代替求最優(yōu)狀態(tài)路徑,本文還給出了類似于MAP(Maximum a posteriori probability)算法[34]基于閾值的后驗(yàn)解碼算法,該算法給出了最可能的狀態(tài)路徑,記為SegHMC threshold.與MAP算法輸出每個(gè)后驗(yàn)概率最大的序列區(qū)域相比,SegHMC threshold輸出后驗(yàn)概率大于指定閾值包含順式調(diào)控模塊的序列區(qū)域.在SegHMC threshold算法中,本文搜索后驗(yàn)概率大于給定閾值且至少包含兩個(gè)模體的連續(xù)區(qū)域作為候選順式調(diào)控模塊順式調(diào)控模塊.候選順式調(diào)控模塊區(qū)域的邊界定義為首個(gè)模體的起始位置,和最后一個(gè)模體的結(jié)束位置.在本文模型中,選擇的閾值范圍為[0.45,0.70],在模型的后驗(yàn)推斷中該范圍內(nèi)的閾值通常能給出好的性能.相對(duì)于完全輸出后驗(yàn)概率最大的MAP輸出來講,能通過合理地選取相應(yīng)的閾值在精度和召回率之間達(dá)到一個(gè)平衡.

    1.3 SegHMC算法的整體框架

    為了找出目標(biāo)序列中的順式調(diào)控模塊,本文算法需要給定相應(yīng)的輸入,它主要包括:訓(xùn)練集、測試集、常量參數(shù)、候選模體的PWM集、篩選順式調(diào)控模塊的閾值.算法通過執(zhí)行如下的過程,輸出所有測試集中大于給定閾值的順式調(diào)控模塊:

    1)在訓(xùn)練集上,使用第1.1.1節(jié)中的Segmental HMM構(gòu)造算法構(gòu)造相應(yīng)的Segmental HMM模型,使用Baum-Welch算法訓(xùn)練出Segmental HMM的狀態(tài)轉(zhuǎn)移概率;

    2)在測試集上,使用第1.1.1節(jié)中的Segmental HMM構(gòu)造算法構(gòu)造相應(yīng)的Segmental HMM模型,利用第1)步訓(xùn)練得出的模型參數(shù),使用Veterbi算法解碼或基于閾值的后驗(yàn)解碼算法找出最優(yōu)或最可能的狀態(tài)路徑,預(yù)測出測試集中所有序列的順式調(diào)控模塊,進(jìn)一步過濾得出大于給定閾值的順式調(diào)控模塊集.

    模型的常量參數(shù)主要包括:搜索候選位點(diǎn)的pvalue;開始一個(gè)順式調(diào)控模塊的概率參數(shù)ps和結(jié)束一個(gè)順式調(diào)控模塊的概率參數(shù)pe;相鄰位點(diǎn)距離的幾何分布參數(shù)ml和相鄰順式調(diào)控模塊距離的幾何參數(shù)mg;以及順式調(diào)控模塊的權(quán)重閾值w.這些參數(shù)的取值首先結(jié)合真實(shí)生物數(shù)據(jù)特征的先驗(yàn)知識(shí)來選取;例如,一個(gè)序列中通常的順式調(diào)控模塊的含量、順式調(diào)控模塊的平均長度、順式調(diào)控模塊內(nèi)模體間的平均距離等.其次,通過在模擬數(shù)據(jù)和搜集到的順式調(diào)控模塊數(shù)據(jù)集上,試驗(yàn)一個(gè)取值范圍內(nèi)不同的取值,選擇在多數(shù)情況下能產(chǎn)生好的結(jié)果的值作為模型的參數(shù).

    算法所輸出的順式調(diào)控模塊的具體信息主要包括:序列中所有找到的順式調(diào)控模塊的位置及其分值,以及構(gòu)成順式調(diào)控模塊的模體、相應(yīng)模體的位置和分值;其中,順式調(diào)控模塊分值和順式調(diào)控模塊內(nèi)的各個(gè)模體的分值,分別標(biāo)識(shí)了所找到的順式調(diào)控模塊和相應(yīng)模體的統(tǒng)計(jì)顯著性.我們將順式調(diào)控模塊的分值定義為對(duì)應(yīng)的序列片段由順式調(diào)控模塊狀態(tài)生成的后驗(yàn)概率和該序列片段由全局背景狀態(tài)生成的后驗(yàn)概率的log似然比值;類似地,模體的分值定義為該模體所對(duì)應(yīng)的序列片段由模體狀態(tài)生成的后驗(yàn)概率和該序列片段由全局背景狀態(tài)生成的后驗(yàn)概率的log似然比值.

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

    2.1 實(shí)驗(yàn)數(shù)據(jù)

    我們?cè)谝粋€(gè)模擬數(shù)據(jù)集和兩個(gè)真實(shí)生物數(shù)據(jù)集上測試本文方法:脊椎動(dòng)物的Muscle特定表達(dá)系統(tǒng)和果蠅早期胚胎發(fā)育系統(tǒng),以下簡稱Muscle數(shù)據(jù)集和果蠅早期發(fā)育數(shù)據(jù)集.這兩個(gè)真實(shí)生物數(shù)據(jù)集是評(píng)價(jià)當(dāng)前順式調(diào)控模塊方法性能的標(biāo)準(zhǔn)數(shù)據(jù)集[9?10],分別代表了兩類不同的序列數(shù)據(jù);其中Muscle數(shù)據(jù)上的序列為共調(diào)控序列,果蠅早期發(fā)育數(shù)據(jù)集中的序列為同源序列.

    模擬數(shù)據(jù)集由30條序列和一個(gè)包含52個(gè)模體的模體集構(gòu)成,這52個(gè)模體從TRANSFAC數(shù)據(jù)庫中隨機(jī)提取;其中,每個(gè)序列的長度均為20kbp.序列集中的每條序列包含0~3個(gè)順式調(diào)控模塊;這些順式調(diào)控模塊的長度隨機(jī)為200~1500bp,模體位點(diǎn)間的空白的平均長度為50bp;每個(gè)順式調(diào)控模塊包含2~6種不同的模體,大約有15個(gè)模體位點(diǎn)實(shí)例;此外,為了模擬某些模體傾向于共同出現(xiàn)的偏好,設(shè)定每個(gè)順式調(diào)控模塊內(nèi)包含40%的有序模體對(duì);順式調(diào)控模塊內(nèi)部和順式調(diào)控模塊間的背景序列均采用3階Markov模型,模型參數(shù)通過掃描D.melanogaster基因組間的序列獲得.

    Muscle數(shù)據(jù)集最初由Wassermen和Fickett收集編錄[41];后來,Klepper等選取該數(shù)據(jù)集的一個(gè)子集并進(jìn)行了擴(kuò)展,作為文獻(xiàn)[9]中的基準(zhǔn)數(shù)據(jù).我們所使用的Muscle數(shù)據(jù)集即為文獻(xiàn)[9]中的基準(zhǔn)數(shù)據(jù),該數(shù)據(jù)集共包含5個(gè)模體和來自不同物種24條共調(diào)控序列.所包含的模體具體為Mef2、Myf、Sp1、SRF和Tef,這些模體在肌肉的轉(zhuǎn)錄調(diào)控中起重要作用.所包含的24條序列,平均長度為850bp,分別來自老鼠、人類、雞等物種.24條序列共包含5個(gè)模體的84個(gè)實(shí)例;每條序列均包含1個(gè)順式調(diào)控模塊,這些順式調(diào)控模塊的平均長度為120bp,其中最短的為14bp,最長的為294bp.

    本文模型可以在一組具有相似調(diào)控結(jié)構(gòu)基因的調(diào)控區(qū)上進(jìn)行訓(xùn)練,然后用于搜索其他基因甚至整個(gè)基因組中的相似順式調(diào)控模塊.但限于收集完整標(biāo)注全基因組數(shù)據(jù)的困難,我們僅在已標(biāo)注的調(diào)控果蠅早期前后軸發(fā)育系統(tǒng)的一組基因上進(jìn)行了測試.相較于Muscle數(shù)據(jù)集,果蠅早期發(fā)育數(shù)據(jù)集中的基因有相對(duì)較長的順式調(diào)控模塊.在該數(shù)據(jù)集中,我們選取參與果蠅早期胚胎發(fā)育轉(zhuǎn)錄調(diào)控的7個(gè)模體:Bcd、Hb、Cad、Kr、Kni、Tll和Gt,并從iDMMPWM數(shù)據(jù)庫[42]中下載這些模體的PWM.我們從在果蠅早期胚胎發(fā)育過程中起重要作用的基因中選取一個(gè)包含7個(gè)基因的子集.這7個(gè)基因具體為:kni、kr、hb、tll、btd、eve和h,這些基因主要在果蠅早期胚胎的前后軸發(fā)育中調(diào)控果蠅幼體的體軸和體節(jié)發(fā)育.對(duì)于每個(gè)基因,我們選擇D.melanogaster和其12個(gè)同源的對(duì)應(yīng)基因序列作為搜索集.基因的染色體坐標(biāo)和同源信息從FlyBase數(shù)據(jù)庫[43]中獲取.典型的真核生物基因的調(diào)控區(qū)的長度通常為包含轉(zhuǎn)錄起始位點(diǎn)上游15kbps和下游5kbps的區(qū)域,但為了盡可能地覆蓋包含順式調(diào)控模塊的區(qū)域,我們將搜索區(qū)域定義為基因的轉(zhuǎn)錄起始位點(diǎn)附近的上游區(qū)20kbp和下游區(qū)20kbp,總共40kbp的序列區(qū)域,使用RepeatMasker軟件(http://www.repeatmasker.org)掩住(即替換為N字符串)搜索區(qū)內(nèi)的重復(fù)子序列.所選的7個(gè)模體和7個(gè)基因上的約91條序列構(gòu)成了我們實(shí)驗(yàn)中所使用的果蠅早期發(fā)育數(shù)據(jù)集.我們從REDfly數(shù)據(jù)庫[44]中收集這些基因的順式調(diào)控模塊,并對(duì)存在重疊的順式調(diào)控模塊進(jìn)行合并,將整理后的集合作為真實(shí)注釋的順式調(diào)控模塊基準(zhǔn)集.

    2.2 評(píng)價(jià)

    為了對(duì)本文方法在上述數(shù)據(jù)集上的預(yù)測結(jié)果進(jìn)行客觀評(píng)價(jià),并且使得評(píng)價(jià)結(jié)果不傾向于某一單個(gè)指標(biāo),我們使用評(píng)價(jià)指標(biāo)相關(guān)系數(shù)(Correlation coefficient,CC)[45]和定義在精度(Precision)和召回率(Recall)上的F1-score[46]作為堿基水平上的主要評(píng)價(jià)指標(biāo)對(duì)方法的預(yù)測準(zhǔn)確率進(jìn)行評(píng)價(jià),這些指標(biāo)被大多數(shù)順式調(diào)控識(shí)別方法所使用.此外,為了評(píng)價(jià)本文方法在位點(diǎn)水平上的精度/召回率,我們畫出了SegHMC threshold在不同閾值下的P/R曲線,以及模型缺省設(shè)置SegHMC Veterbi下的Precision/Recall值.這些評(píng)價(jià)指標(biāo)的具體定義為:

    這里TP、FP、TN和FN分別表示真陽性、假陽性、真陰性和假陰性;式(9)中的精度(Precision)和召回率(Recall)的具體定義為:

    1)精度(Precison):Precision=TP/(TP+ FP),度量了預(yù)測的正確率,表示預(yù)測結(jié)果中被正確識(shí)別的順式調(diào)控模塊的比率;

    2)召回率(Recall):Recall=TP/(TP+FN),度量了順式調(diào)控模塊的識(shí)別率,表示基準(zhǔn)集中被正確識(shí)別的順式調(diào)控模塊的比率.

    CC度量了預(yù)測集和基準(zhǔn)集中堿基位置的相關(guān)性,同時(shí)考慮了預(yù)測結(jié)果的真陽性率和真陰性率. CC取值范圍在–1和+1之間,+1表示預(yù)測結(jié)果與基準(zhǔn)集完全一致,–1表示預(yù)測結(jié)果與基準(zhǔn)集完全相反;當(dāng)預(yù)測結(jié)果接近隨機(jī)時(shí),CC趨向于0.

    從精度和召回率的定義,可以看到,它們分別衡量了預(yù)測算法的兩個(gè)方面:查全率和查準(zhǔn)率.雖然僅從定義上看這兩個(gè)指標(biāo)之間并沒有必然的聯(lián)系,但事實(shí)上它們?cè)诙鄶?shù)情況下是相互制約、相互矛盾的.這種制約和矛盾主要是由搜索策略的不完善造成的.通常情況下,如果希望獲得較高的召回率,則不得不降低搜索策略標(biāo)準(zhǔn),以盡可能地覆蓋基準(zhǔn)集中結(jié)果,從而引入了大量不相關(guān)的結(jié)果,降低了精度;如果想獲得較高的精度,則需要提高搜索策略標(biāo)準(zhǔn),這不可避免地過濾掉了相關(guān)的結(jié)果,造成了召回率的下降.因此,通常需要在這兩個(gè)指標(biāo)之間找到一個(gè)合理的平衡點(diǎn).F1-score是精度和召回率的調(diào)和平均值,綜合了方法的精度和召回率,度量了方法在精度和召回率之間達(dá)到平衡的能力;其值在0和1之間,值越高,通常表示算法的性能越好(在真陽性率方面).

    我們使用 Cython 語言,將本文方法SegHMC開發(fā)成相應(yīng)的工具,并在 Intel Xeon E5640@2.67GHz處理器,4GB內(nèi)存的Windows 7 64位系統(tǒng)的平臺(tái)上進(jìn)行測試.

    我們選取當(dāng)前的主要方法:BayCis[35]、Stubb[32]、MSCAN[22]、MotEvo[25]、ReLA[27]和CMStalker[29],在 Muscle 數(shù)據(jù)集上與本文方法SegHMC進(jìn)行預(yù)測性能的比較.

    在這6個(gè)被比較的方法中,BayCis和Stubb屬于概率模型方法,它們構(gòu)造相應(yīng)的HMM模型,訓(xùn)練模型參數(shù),使用訓(xùn)練的HMM解碼出給定序列中潛在的順式調(diào)控模塊.MSCAN是一個(gè)窗聚集方法,它計(jì)算給定序列滑動(dòng)的窗內(nèi)出現(xiàn)的所有模體聚集的統(tǒng)計(jì)顯著性,將統(tǒng)計(jì)顯著性超過給定閾值的序列區(qū)域判定為候選順式調(diào)控模塊.MotEvo是窗聚集和概率模型的結(jié)合,它使用貝葉斯模型對(duì)滑動(dòng)窗口內(nèi)的模體聚集區(qū)域進(jìn)行打分,將分值大于給定閾值的區(qū)域作為候選順式調(diào)控模塊輸出.CMStalker是一個(gè)組合方法,它綜合了約束滿足規(guī)劃與參數(shù)松馳技術(shù),能有效遍歷模體組合的解空間,在給定的序列中找出可能的順式調(diào)控模塊.ReLA是一個(gè)基于局部比對(duì)的順式調(diào)控模塊搜索方法,它首先通過使用第三方模體搜索工具找出序列集中的所有模體,然后使用修改的Smith-Waterman算法,找出參考序列和被比較序列中得分最高的局部比對(duì)作為候選順式調(diào)控模塊.

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

    在模擬數(shù)據(jù)集上,模型固定參數(shù)具體設(shè)置為: ps=0.001,pe=0.1,p-value=0.01,和 w=100. mg和ml分別被設(shè)置為500bp和50bp,表示一個(gè)序列內(nèi)的順式調(diào)控模塊間的平均距離為500bp,順式調(diào)控模塊內(nèi)相鄰模體間的平均距離為50bp.對(duì)于其他被比較的方法,我們使用其默認(rèn)設(shè)置.

    圖3(a)給出了所有方法的CC值和F1-score值.從圖中可以看出,我們的方法在這兩個(gè)評(píng)價(jià)指標(biāo)上明顯優(yōu)于其他方法.圖3(b)給出了所有方法的P/R值.從圖中可以看出,除MSCAN外的其他方法都有很高的召回率(Recall值),即所有方法的預(yù)測結(jié)果大都覆蓋了所植入的順式調(diào)控模塊,但預(yù)測的精度(Precision值)表現(xiàn)出很大的差異.但我們方法的兩個(gè)版本在保持高召回率的同時(shí),精度超過絕大部分的其他方法,在精度和召回率之間達(dá)到一個(gè)很好的平衡.

    圖3 SegHMC threshold和SegHMC Veterbi在模擬數(shù)據(jù)集上的預(yù)測性能Fig.3 The prediction performances of SegHMC threshold and SegHMC Veterbi on the synthetic dataset

    2.4 Muscle數(shù)據(jù)集上的結(jié)果

    考慮到該數(shù)據(jù)集上的序列所包含的順式調(diào)控模塊,大多具有模體高度聚集但總體長度較短的特征,我們?cè)O(shè)置模型參數(shù)ml為20bp,同時(shí)保持其他參數(shù)值不變.

    圖4(a)給出了所有方法在Muscle數(shù)據(jù)集上的CC和F1-score值.從圖中可以看出,本文方法SegHMC在這兩個(gè)指標(biāo)上均優(yōu)于現(xiàn)存的方法(這里為模型缺省設(shè)置SegHMC Veterbi下的實(shí)驗(yàn)結(jié)果).

    圖4 SegHMC threshold和SegHMC Veterbi在Muscle集上的預(yù)測性能Fig.4 The prediction performances of SegHMC threshold and SegHMC Veterbi on the muscle dataset

    為了考查SegHMC在該數(shù)據(jù)集上的具體表現(xiàn),我們根據(jù)不同的閾值給出SegHMC threshold在位點(diǎn)水平上預(yù)測性能的P/R曲線.對(duì)于模型缺省設(shè)置SegHMC Veterbi和其他被比較的方法,我們根據(jù)它們的缺省輸出在圖中給出相應(yīng)P/R值的點(diǎn),如圖4(b)所示.在圖4(b)中,SegHMC threshold給出了一個(gè)平衡范圍,SegHMC Veterbi輸出位于曲線中間的一個(gè)點(diǎn).

    從圖4(b)中可以看出,本文方法的兩個(gè)版本在該數(shù)據(jù)集上較其他方法在預(yù)測精度和召回率之間均達(dá)到一個(gè)更好的平衡,具體為:在預(yù)測精度(Precision值)上高于其他大多數(shù)被比較的方法,召回率(Recall值)超過一半的方法.對(duì)于其他方法,分析圖中的數(shù)據(jù),還可以看到,在該數(shù)據(jù)集上不同的方法達(dá)到不同的P/R平衡,表現(xiàn)出不同的性能趨向.對(duì)于MotEvo和Stubb方法,它們傾向于給出較高召回率的預(yù)測,但整體上并未達(dá)到合理的平衡點(diǎn),表現(xiàn)為預(yù)測精度(Precision值)明顯低于所有方法的平均值,僅約為本文方法的一半.很顯然,這種高的召回率是通過放寬順式調(diào)控模塊的篩選條件,盡可能覆蓋可能存在順式調(diào)控模塊的序列區(qū),以給出盡可能多的預(yù)測結(jié)果,犧牲預(yù)測精度換來的.而對(duì)于BayCis和CMStalker方法,為確保有較高的預(yù)測精度,它傾向于給出保守的預(yù)測結(jié)果,具體表現(xiàn)為其召回率(Recall值)明顯低于其他方法.對(duì)于MSCAN方法,它并未完全傾向于某一單個(gè)指標(biāo),而是在保證召回率的同時(shí),給出了較高的預(yù)測精度,表現(xiàn)出良好的預(yù)測性能.在所有的方法中,ReLA并不傾向任何單個(gè)指標(biāo),給出了最保守的預(yù)測.

    2.5 果蠅早期發(fā)育數(shù)據(jù)集上的結(jié)果

    我們將這些方法進(jìn)一步在果蠅早期發(fā)育數(shù)據(jù)集上進(jìn)行測試,模型的所有常量參數(shù)設(shè)置與模擬數(shù)據(jù)集上相同.其中,對(duì)于Stubb,我們選擇其多物種版StubbMS(Stubb中的一個(gè)模塊,在以下的描述中我們?nèi)杂洖镾tubb).StubbMS利用相近物種間的保守性來提高方法的預(yù)測性能,具體為:通過雙序列比對(duì),找出兩物種間的保守區(qū)域,對(duì)保守區(qū)內(nèi)的片段進(jìn)行打分,將分值高于給定閾值的保守片段作為候選順式調(diào)控模塊.所有方法在該數(shù)據(jù)集上的整個(gè)測試方式如下:

    1)對(duì)于 BayCis和 SegHMC,使用和D.melanogaster同源的所有其他果蠅物種的基因作為訓(xùn)練集,將D.melanogaster的對(duì)應(yīng)基因作為測試集(REDfly數(shù)據(jù)庫僅搜集了D.melanogaster基因的順式調(diào)控模塊標(biāo)注信息,其他同源物種的順式調(diào)控模塊暫時(shí)無法搜集);

    2)對(duì)于需要雙序列比對(duì)的 Stubb方法,選取其原文中所使用的物種 (D.melanogaster和D.pseudoobscura),在相應(yīng)基因上進(jìn)行測試;

    3)對(duì)于僅在單序列上預(yù)測順式調(diào)控模塊的方法MSCAN和MotEvo,僅在D.melanogaster相應(yīng)的基因上進(jìn)行測試;

    4)對(duì)于ReLA,選取D.melanogaster相應(yīng)的基因作為參考序列,其余同源基因序列作為被比對(duì)的序列;

    5)對(duì)于CMStalker,在所有同源基因上進(jìn)行測試.

    圖5(a)給出了所有方法在該數(shù)據(jù)集堿基水平上的CC和F1-score.從圖中可以看出,與其他方法相比,本文方法SegHMC在整個(gè)數(shù)據(jù)集的CC和F1-score值(這里為模型缺省設(shè)置SegHMC Veterbi下的實(shí)驗(yàn)結(jié)果)上仍高于其他方法,表現(xiàn)出穩(wěn)定的預(yù)測性能.此外,通過對(duì)比可以發(fā)現(xiàn),與Muscle數(shù)據(jù)集上的結(jié)果相比,這些指標(biāo)的分值均有較大幅度的下降.考查這兩個(gè)數(shù)據(jù)集上數(shù)據(jù)的特征,發(fā)現(xiàn)這主要是由于搜索區(qū)長度的增加,而造成的信噪比(順式調(diào)控模塊相對(duì)于背景)的降低造成的.與Muscle數(shù)據(jù)集相比(序列平均長度為850bp),果蠅早期發(fā)育數(shù)據(jù)集中的序列更長(長度均為40kbp),各方法在識(shí)別單個(gè)模體時(shí),不可避免地給出大量假陽性的預(yù)測,進(jìn)而影響到順式調(diào)控模塊的預(yù)測結(jié)果,造成了預(yù)測結(jié)果中的大量假陽性預(yù)測.

    圖5(b)給出了所有方法位點(diǎn)水平上的P/R曲線.從圖中可以看出,本文方法的兩個(gè)版本都達(dá)到了很好的平衡.并且閾值方法在某些閾值上較缺省設(shè)置達(dá)到更好的平衡,盡管如此,我們也應(yīng)該明白:在實(shí)際應(yīng)用中,我們所能達(dá)到的平衡完全基于數(shù)據(jù)的輸入,而P/R是完全未知的,我們很難選擇合適的閾值參數(shù),得到好的預(yù)測結(jié)果.因此,方法自動(dòng)選擇合適的參數(shù)達(dá)到適當(dāng)?shù)腜/R平衡就顯得更為重要.從圖4(b)中,可以發(fā)現(xiàn):在Muslce數(shù)據(jù)集上表現(xiàn)良好的MSCAN方法,在該數(shù)據(jù)集上表現(xiàn)得很差,其精度(Precision值)僅處于中間水平;召回率(Recall值)明顯低于其他方法.在所有方法中MSCAN是唯一的基于窗聚集的方法,對(duì)于包含長度多變的順式調(diào)控模塊的序列,窗聚集方法很難確定一個(gè)合理的窗大小和打分閾值,從而給出好的預(yù)測結(jié)果(這里使用其默認(rèn)輸出);另一方面,MotEvo也是基于窗聚集的方法,但與純粹的基于窗聚集的MSCAN方法相比,它結(jié)合了貝葉斯概率模型對(duì)滑動(dòng)窗內(nèi)的模體進(jìn)行打分,從而有更好的性能,這也在某種程度上說明了,概率模型方法對(duì)不同類型的數(shù)據(jù)有更好的適應(yīng)性.對(duì)于引入系統(tǒng)發(fā)生的方法Stubb和ReLA,通過序列比對(duì)利用物種的保守性確實(shí)提高了預(yù)測精度(具有僅次于本文方法的Precision值),但完全基于物種間的保守性很難保證有較高的召回率(其Recall值最低).雖然我們的方法未直接通過序列比對(duì)利用物種間的保守性,但通過概率模型刻畫共調(diào)控或同源序列間結(jié)構(gòu)的相似性,在某種程度上也利用了物種間的系統(tǒng)發(fā)生關(guān)系,這也是本文方法優(yōu)于其他方法的一個(gè)原因.對(duì)于基于概率模型的方法BayCis,雖然在該數(shù)據(jù)集上給出了較高召回率(Recall值)的預(yù)測,但其預(yù)測的精度是最低的(Precision值).對(duì)于純粹組合搜索方法CMStalker,給出了很保守的預(yù)測,與在Muscle數(shù)據(jù)集上的結(jié)果相比,預(yù)測精度和召回率都下降很多.考慮到該數(shù)據(jù)集上的順式調(diào)控模塊的長度特點(diǎn),我們推測可能超出了CMStalker的所能處理的空間范圍.

    圖5 SegHMC threshold和SegHMC Veterbi在果蠅早期發(fā)育數(shù)據(jù)集上的預(yù)測性能Fig.5 The prediction performances of SegHMC threshold and SegHMC Veterbi on the early drosophila development dataset

    為了考查不同的方法在該數(shù)據(jù)集中單個(gè)數(shù)據(jù)上的預(yù)測性能的差異,我們畫出了所有方法在該數(shù)據(jù)集中每個(gè)基因上的預(yù)測性能的箱線圖,如圖6所示.

    圖6給出了每個(gè)方法CC值的中位數(shù)、下四分位數(shù)和上四分位數(shù)以及最小值和最大值.從圖中可以看出:所有方法在不同基因上的預(yù)測性能變化很大,這種變化要比方法總體性能之間的差異要大得多;這也說明了即使表現(xiàn)最好的方法在許多情況下也可能給出錯(cuò)誤的預(yù)測.本文方法總體上較其他方法有更高的CC值,表現(xiàn)出更好的預(yù)測性能.雖然CMStalker在所有基因上的CC值變化最小,表現(xiàn)得最穩(wěn)定,但其CC值總體上明顯要低于大多數(shù)的方法.

    圖6 所有方法在果蠅早期發(fā)育數(shù)據(jù)集中各基因上CC的變化Fig.6 A boxplot describing variation for all methods in CC across the genes in the early drosophila development dataset

    2.6 引入的保守結(jié)構(gòu)信息帶來的性能提升

    為了考查所引入的結(jié)構(gòu)信息(相鄰位點(diǎn)的相關(guān)性和位點(diǎn)間的距離分布)對(duì)方法性能的影響,我們回退模型到不考慮這些結(jié)構(gòu)信息一般模型,記為SegHMC-simple,在所有3個(gè)數(shù)據(jù)集上進(jìn)行實(shí)驗(yàn)對(duì)比,如圖7所示.

    從圖中可以看出,加入結(jié)構(gòu)信息的模型SegHMC相對(duì)未加入相關(guān)信息的SegHMC-simple總體上具有明顯的性能提升,盡管對(duì)于不同的數(shù)據(jù)集性能提升不同.總體上,SegHMC在模擬數(shù)據(jù)和Muscle數(shù)據(jù)集上,較果蠅早期發(fā)育數(shù)據(jù)集上的提升要大.我們推測可能的原因在于,果蠅早期發(fā)育數(shù)據(jù)集上較大的搜索空間中所存在的大量噪聲造成的,即結(jié)構(gòu)信息的引入,可能增強(qiáng)原本的弱信號(hào),引入了假陽性,抵消了一部分結(jié)構(gòu)信息所帶來的性能提升.檢查SegHMC在3個(gè)數(shù)據(jù)集上所找到的順式調(diào)控模塊,可以發(fā)現(xiàn),SegHMC能找出模擬數(shù)據(jù)集中所植入的70%以上共現(xiàn)模體對(duì);在Muscle數(shù)據(jù)集中,找出了Mef2-Myf[7]、Myf-Sp1[7]以及Tef-Mef2[47]等經(jīng)過生物實(shí)驗(yàn)驗(yàn)證的模體對(duì);在果蠅早期發(fā)育數(shù)據(jù)集上找出了共現(xiàn)的同型模體對(duì)[48],如Bcd-Bcd、Kr-Kr以及Hb-Hb等.

    圖7 引入結(jié)構(gòu)信息后SegHMC在所有數(shù)據(jù)集上性能的提升Fig.7 An effect of inclusion of structural information on the prediction performance of SegHMC for all datasets

    3 結(jié)論

    本文將順式調(diào)控模塊的模體聚集特征和內(nèi)部結(jié)構(gòu)保守性結(jié)合起來,抽象出一種復(fù)雜的順式調(diào)控模塊的調(diào)控結(jié)構(gòu)表示(或稱調(diào)控語法),并借助于具有更強(qiáng)表達(dá)能力的Segmental HMM來表達(dá),建立了一種識(shí)別順式調(diào)控模塊的概率模型.與其他基于概率類型的順式調(diào)控模塊識(shí)別方法相比,本文方法有如下的優(yōu)點(diǎn):

    1)我們不僅將順式調(diào)控模塊表示為模體的組合,還將模體共同出現(xiàn)的頻率、模體順序偏好以及順式調(diào)控模塊中的相鄰模體之間距離分布等特征引入到順式調(diào)控模塊的調(diào)控語法表示當(dāng)中,這些特征可以有效提高順式調(diào)控模塊的識(shí)別精度.

    2)Segmental HMM的狀態(tài)可以表示一個(gè)序列片段,這一特性使得我們可以將所抽象的順式調(diào)控模塊的結(jié)構(gòu)元素(模體、背景等)直接映射為Segmental HMM的一個(gè)狀態(tài).在模體而不是堿基的抽象層次上對(duì)順式調(diào)控模塊的調(diào)控結(jié)構(gòu)進(jìn)行建模,從而使得整個(gè)模型結(jié)構(gòu)的表示更清晰、更自然.此外,Segmental HMM并不限定片段長度的具體分布,使得我們可以根據(jù)模型假定選取特定的片段長度分布.

    3)通過預(yù)先搜索模體實(shí)例,提前標(biāo)定相應(yīng)片段的類型,顯式構(gòu)建Segmental HMM的狀態(tài)轉(zhuǎn)換圖,有效減少了待搜索空間,提高了算法效率,使得本文方法可以處理真核生物普遍的長調(diào)控區(qū)序列.此外, segmental HMM模型所固有的在線(Online)性質(zhì),使得本文模型可以經(jīng)過在一組具有相似調(diào)控結(jié)構(gòu)基因的調(diào)控區(qū)上進(jìn)行訓(xùn)練,然后用于搜索其他所要研究基因甚至整個(gè)基因組中的相似順式調(diào)控模塊.

    為進(jìn)一步提高識(shí)別順式調(diào)控模塊的性能,我們?cè)诤罄m(xù)的研究中打算從以下幾個(gè)方面著手,以進(jìn)一步完善我們的工作:一方面,我們將搜集更多標(biāo)注的順式調(diào)控模塊,并使用更系統(tǒng)的方法,比如交叉驗(yàn)證,來輔助模型參數(shù)的選擇.我們相信這些措施能使我們方法的性能得到進(jìn)一步的提升.另一方面,當(dāng)前新的生物測序技術(shù),如染色質(zhì)共沉淀后微陣列分析(Chip-chip)或染色質(zhì)共沉淀后測序(Chip-seq)等的快速發(fā)展,產(chǎn)生了DNA雙螺旋結(jié)構(gòu)Profile、染色質(zhì)結(jié)構(gòu)、組蛋白修飾、蛋白質(zhì)占位等的大量實(shí)驗(yàn)數(shù)據(jù),這些數(shù)據(jù)隱藏著基因表達(dá)調(diào)控規(guī)律,進(jìn)一步分析提取這些可用于順式調(diào)控模塊預(yù)測的生物信息,并將這些信息與現(xiàn)有的計(jì)算模型結(jié)合,構(gòu)建更有效的順式調(diào)控模塊識(shí)別方法也是我們今后進(jìn)一步努力的方向.

    1 Wasserman W W,Sandelin A.Applied bioinformatics for the identification of regulatory elements.Nature Reviews Genetics,2004,5(4):276?287

    2 Davidson E H.The Regulatory Genome:Gene Regulatory Networks in Development and Evolution.San Diego,California:Academic Press/Elsevier,2006.

    3 Wang Pei,Lv Jin-Hu.Control of genetic regulatory networks: opportunities and challenges.Acta Automatica Sinica,2013,39(12):1969?1979 (王沛,呂金虎.基因調(diào)控網(wǎng)絡(luò)的控制:機(jī)遇與挑戰(zhàn).自動(dòng)化學(xué)報(bào), 2013,39(12):1969?1979)

    4 Chen L N,Wang R S,Zhang X S.Biomolecular Networks: Methods and Applications in Systems Biology.Hoboken, New Jersey:Wiley,2009.

    5 Kleinjan D A,Seawright A,Mella S,Carr C B,Tyas D A,Simpson T I,Mason J O,Price D J,van Heyningen V. Long-range downstream enhancers are essential for Pax6 expression.Developmental Biology,2006,299(2):563?581

    6 Hardison R C,Taylor J.Genomic approaches towards finding cis-regulatory modules in animals.Nature Reviews Genetics,2012,13(7):469?483

    7 Matys V,Kel-Margoulis O V,Fricke E,Liebich I,Land S,Barre-Dirrie A,Reuter I,Chekmenev D,Krull M,Hornischer K,Voss N,Stegmaier P,Lewicki-Potapov B,Saxel H,Kel A E,Wingender E.TRANSFAC?and its module TRANSCompel?:transcriptional gene regulation in eukaryotes.Nucleic Acids Research,2006,34(Database issue): D108?D110

    8 Portales-Casamar E,Thongjuea S,Kwon A T,Arenillas D,Zhao X,Valen E,Yusuf D,Lenhard B,WassermanWW,Sandelin A.JASPAR 2010:the greatly expanded open-access database of transcription factor binding profiles.Nucleic Acids Research,2010,38(Database issue): D105?D110

    9 Klepper K,Sandve G K,Abul O,Johansen J,Drablos F.Assessment of composite motif discovery methods.BMC Bioinformatics,2008,9:123

    10 Su J,Teichmann S A,Down T A.Assessing computational methods of cis-regulatory module prediction.PLoS Computational Biology,2010,6(12):e1001020

    12 Suryamohan K,Halfon M S.Identifying transcriptional cisregulatory modules in animal genomes.Wiley Interdisciplinary Reviews:Developmental Biology,2015,4(2):59?84

    13 Thompson J A,Congdon C B.GAMI-CRM:using de novo motif inference to detect cis-regulatory modules.In:Proceedings of the 2014 IEEE Congress on Evolutionary Computation.Beijing,China:IEEE,2014.1022?1029

    14 Zheng Shu-Rui.Research of Cis-regulatory Module Discovery Method Based on HMM Model[Master dissertation], Xidian University,China,2012 (鄭樹銳.基于HMM模型的順式調(diào)控模塊識(shí)別方法的研究[碩士學(xué)位論文],西安電子科技大學(xué),中國,2012)

    15 Navarro C,Lopez F J,Cano C,Garcia-Alcalde F,Blanco A.CisMiner:genome-wide in-silico cis-regulatory module prediction by fuzzy itemset mining.PLoS One,2014,9(9): e108065

    16 Rouault H,Santolini M,Schweisguth F,Hakim V.Imogene:identification of motifs and cis-regulatory modules underlying gene co-regulation.Nucleic Acids Research,2014, 42(10):6128?6145

    17 Potier D,Seyres D,Guichard C,Iche-Torres M,Aerts S, Herrmann C,Perrin L.Identification of cis-regulatory modules encoding temporal dynamics during development.BMC Genomics,2014,15(1):534

    18 Thompson J A,Congdon C B.Initial results in using de novo motif inference to detect cis-regulatory modules.In: Proceedings of the 2013 International Conference on Bioinformatics,Computational Biology and Biomedical Informatics.Washington DC,USA:ACM,2013.687

    19 Lemnian I M,Eggeling R,Grosse I.Extended sunflower hidden Markov models for the recognition of homotypic cisregulatory modules.In:Proceedings of the 2013 German Conference on Bioinformatics.Gottingen,Germany,2013. 101?109

    20 Zhou Q,Wong W H.CisModule:de novo discovery of cisregulatory modules by hierarchical mixture modeling.Proceedings of the National Academy of Sciences of the United States of America,2004,101(33):12114?12119

    21 Gan Y L,Guan J H,Zhou S G,Zhang W X.Identifying cisregulatory elements and modules using conditional random fields.IEEE/ACM Transactions on Computational Biology and Bioinformatics,2014,11(1):73?82

    22 Alkema W B,Johansson O,Lagergren J,Wasserman W W. MSCAN:identification of functional clusters of transcription factor binding sites.Nucleic Acids Research,2004,32(Web Server issue):W195?W198

    23 Aerts S,Van Loo P,Thijs G,Moreau Y,De Moor B.Computational detection of cis-regulatory modules.Bioinformatics, 2003,19(Suppl 2):ii5?ii14

    24 Sharan R,Ovcharenko I,Ben-Hur A,Karp R M.CREME:a framework for identifying cis-regulatory modules in humanmouse conserved segments.Bioinformatics,2003,19(Suppl 1):i283?i291

    25 Arnold P,Erb I,Pachkov M,Molina N,van Nimwegen E. MotEvo:integrated Bayesian probabilistic methods for inferring regulatory sites and motifs on multiple alignments of DNA sequences.Bioinformatics,2012,28(4):487?494

    26 Sinha S,He X.MORPH:probabilistic alignment combined with hidden Markov models of cis-regulatory modules.PLoS Computational Biology,2007,3(11):e216

    28 Bailey T L,Noble W S.Searching for statistically significant regulatory modules.Bioinformatics,2003,19(Suppl 2):ii16?ii25

    29 Leoncini M,Montangero M,Pellegrini M,Tillan K P.CMStalker:a combinatorial tool for composite motif discovery. IEEE/ACM Transactions on Computational Biology and Bioinformatics,2015,12(5):1123?1136

    30 Chan B Y,Kibler D.Using hexamers to predict cisregulatory motifs in Drosophila.BMC Bioinformatics,2005, 6:262

    31 Kolbe D,Taylor J,Elnitski L,Eswara P,Li J,Miller W, Hardison R,Chiaromonte F.Regulatory potential scores from genome-wide three-way alignments of human,mouse, and rat.Genome Research,2004,14(4):700?707

    32 Sinha S,van Nimwegen E,Siggia E D.A probabilistic method to detect regulatory modules.Bioinformatics,2003, 19(Suppl 1):i292?i301

    33 Nikulova A A,Favorov A V,Sutormin R A,Makeev V J, Mironov A A.CORECLUST:identification of the conserved CRM grammar together with prediction of gene regulation. Nucleic Acids Research,2012,40(12):e93

    34 Durbin R,Eddy S R,Krogh A,Mitchison G.Biological Sequence Analysis:Probabilistic Models of Proteins and Nucleic Acids.Cambridge:Cambridge University Press,1998.

    35 Lin T H,Ray P,Sandve G K,Uguroglu S,Xing E P.Bay-Cis:a Bayesian hierarchical HMM for cis-regulatory module decoding in metazoan genomes.In:Proceedings of the 12th Annual International Conference on Research in Computational Molecular Biology.Singapore:Springer,2008.66?81

    36 Zhou Q,Wong W H.Coupling hidden Markov models for the discovery of Cis-regulatory modules in multiple species. Annals of Applied Statistics,2007,1(1):36?65

    37 Hu J F,Hu H Y,Li X M.MOPAT:a graph-based method to predict recurrent cis-regulatory modules from known motifs. Nucleic Acids Research,2008,36(13):4488?4497

    38 Russell M J.A segmental HMM for speech pattern modelling.In:Processing of the 1993 IEEE International Conference on Acoustics,Speech,and Signal Processing.Minneapolis,MN,USA:IEEE,1993.499?502

    39 Stormo G D.DNA binding sites:representation and discovery.Bioinformatics,2000,16(1):16?23

    40 Liu X,Brutlag D L,Liu J S.BioProspector:discovering conserved DNA motifs in upstream regulatory regions of co-expressed genes.In:Proceedings of the 6th Pacific Symposium on Biocomputing.Hawaii,USA,2001.127?138

    41 Wasserman W W,Fickett J W.Identification of regulatory regions which confer muscle-specific gene expression.Journal of Molecular Biology,1998,278(1):167?181

    42 Kulakovskiy I V,Makeev V J.Discovery of DNA motifs recognized by transcription factors through integration of different experimental sources.Biophysics,2009,54(6): 667?674

    43 Tweedie S,Ashburner M,Falls K,Leyland P,McQuilton P, Marygold S,Millburn G,Osumi-Sutherland D,Schroeder A,Seal R,Zhang H,The FlyBase Consortium.FlyBase: enhancing Drosophila Gene Ontology annotations.Nucleic Acids Research,2009,37(Database issue):D555?D559

    44 Gallo S M,Gerrard D T,Miner D,Simich M,Des Soye B, Bergman C M,Halfon M S.REDfly v3.0:toward a comprehensive database of transcriptional regulatory elements in Drosophila.Nucleic Acids Research,2011,39(Database issue):D118-D123

    45 Tompa M,Li N,Bailey T L,Church G M,De Moor B,Eskin E,Favorov A V,Frith M C,Fu Y T,Kent W J,Makeev V J,Mironov A A,Noble W S,Pavesi G,Pesole G,Rgnier M,Simonis N,Sinha S,Thijs G,van Helden J,Vandenbogaert M,Weng Z P,Workman C,Ye C,Zhu Z.Assessing computational tools for the discovery of transcription factor binding sites.Nature Biotechnology,2005,23:137?144

    46 Shaw W M Jr,Burgin R,Howell P.Performance standards and evaluations in IR test collections:cluster-based retrieval models.Information Processing&Management, 1997,33(1):1?14

    47 Maeda T,Gupta M P,Stewart A F R.TEF-1 and MEF2 transcription factors interact to regulate muscle-specific promoters.Biochemical and Biophysical Research Communications,2002,294(4):791?797

    48 Lifanov A P,Makeev V J,Nazina A G,Papatsenko D A. Homotypic regulatory clusters in Drosophila.Genome Research,2003,13(4):579?588

    郭海濤 西安電子科技大學(xué)計(jì)算機(jī)學(xué)院博士研究生.主要研究方向?yàn)樯镄畔W(xué)算法,并行算法.

    E-mail:ght_311@sina.com

    (GUO Hai-Tao Ph.D.candidate at the School of Computer Science and Technology,Xidian University.His research interest covers bioinformatics algorithms and pararllel algorithms.)

    霍紅衛(wèi) 博士,西安電子科技大學(xué)計(jì)算機(jī)學(xué)院教授.主要研究方向?yàn)榇髷?shù)據(jù)算法與壓縮數(shù)據(jù)結(jié)構(gòu),生物信息學(xué)算法,算法工程.本文通信作者.

    E-mail:hwhuo@mail.xidian.edu.cn

    (HUO Hong-Wei Ph.D.,professor at the School of Computer Science and Technology,Xidian University.Her research interest covers big data algorithms and compressed data structures,bioinformatics algorithms,algorithm engineering.Corresponding author of this paper.)

    于 強(qiáng) 博士,西安電子科技大學(xué)計(jì)算機(jī)學(xué)院講師.主要研究方向?yàn)樯镄畔W(xué)算法,并行算法.

    E-mail:qyu@mail.xidian.edu.cn

    (YU Qiang Ph.D.,lecturer at the School of Computer Science and Technology,Xidian University. His research interest covers bioinformatics algorithms and pararllel algorithms.)

    SegHMC:an Algorithm for Discovery of Cis-regulatory Module Based on Segmental HMM

    GUO Hai-Tao1HUO Hong-Wei1YU Qiang1

    Cis-regulatory module(CRM)plays a key role in metazoan gene transcriptional regulation,and the discovery of cis-regulatory module has been a crucial research topic recently.Many computational methods have been proposed to predict the cis-regulatory module,but it is still a main task to further improve the prediction accuracy for cis-regulatory modules.Combining multiple features of cis-regulatory module together can improve the prediction accuracy for cisregulatory module.Based on this,the paper presents an algorithm SegHMC(Segmental HMM model for discovery of cis-regulatory module)for the discovery of cis-regulatory module based on segmental HMM.The model further extends the representation of the structure of cis-regulatory module(or regulatory grammar),which not only describes a CRM as a combination of a group of motifs but also further introduces the frequency of the occurrence of motifs,the favour of the order of motifs,and the distance distribution between the adjacent motifs and other features.Experiments on the benchmark datasets demonstrate that the proposed algorithm outperforms the present main algorithms in the prediction accuracy.

    Gene transcriptional regulation,motif,segmental HMM,discovery of cis-regulatory module

    郭海濤,霍紅衛(wèi),于強(qiáng).SegHMC:一種基于Segmental HMM模型的順式調(diào)控模塊識(shí)別算法.自動(dòng)化學(xué)報(bào),2016, 42(11):1718?1731

    Guo Hai-Tao,Huo Hong-Wei,Yu Qiang.SegHMC:an algorithm for discovery of cis-regulatory module based on segmental HMM.Acta Automatica Sinica,2016,42(11):1718?1731

    2015-05-18 錄用日期2016-06-06

    Manuscript received May 18,2015;accepted June 6,2016

    國家自然科學(xué)基金(61173025,61373044,61502366),中國博士后科學(xué)基金(2015M582621)資助

    Supported by National Natural Science Foundation of China (61173025,61373044,61502366),the Chinese Postdoctoral Science Foundation(2015M582621)

    本文責(zé)任編委黃慶明

    Recommended by Associate Editor HUANG Qing-Ming

    1.西安電子科技大學(xué)計(jì)算機(jī)學(xué)院西安710071

    1.School of Computer Science and Technology,Xidian University,Xi'an 710071

    DOI 10.16383/j.aas.2016.c150309

    猜你喜歡
    模體集上調(diào)控
    基于Matrix Profile的時(shí)間序列變長模體挖掘
    Cookie-Cutter集上的Gibbs測度
    鏈完備偏序集上廣義向量均衡問題解映射的保序性
    如何調(diào)控困意
    經(jīng)濟(jì)穩(wěn)中有進(jìn) 調(diào)控托而不舉
    中國外匯(2019年15期)2019-10-14 01:00:34
    植入(l, d)模體發(fā)現(xiàn)若干算法的實(shí)現(xiàn)與比較
    復(fù)扇形指標(biāo)集上的分布混沌
    基于網(wǎng)絡(luò)模體特征攻擊的網(wǎng)絡(luò)抗毀性研究
    順勢而導(dǎo) 靈活調(diào)控
    基于模體演化的時(shí)序鏈路預(yù)測方法
    少妇猛男粗大的猛烈进出视频| 中文字幕精品免费在线观看视频 | 中国三级夫妇交换| 久久国产精品大桥未久av| 97精品久久久久久久久久精品| 男女边吃奶边做爰视频| 中国美白少妇内射xxxbb| 三级国产精品片| 久久ye,这里只有精品| 国产免费现黄频在线看| 亚洲精品自拍成人| 成年人免费黄色播放视频| 美女cb高潮喷水在线观看| 免费观看a级毛片全部| 久热这里只有精品99| 国产精品久久久久久精品古装| 国产亚洲欧美精品永久| av网站免费在线观看视频| 国产成人午夜福利电影在线观看| 人人澡人人妻人| 亚洲欧美精品自产自拍| 午夜免费男女啪啪视频观看| 麻豆乱淫一区二区| 夜夜骑夜夜射夜夜干| 哪个播放器可以免费观看大片| 欧美日韩一区二区视频在线观看视频在线| 人妻一区二区av| 国产色爽女视频免费观看| 曰老女人黄片| 激情五月婷婷亚洲| 免费大片黄手机在线观看| 在线观看人妻少妇| 午夜久久久在线观看| 国产精品久久久久久久久免| 又粗又硬又长又爽又黄的视频| 国产 一区精品| 97超碰精品成人国产| 久久影院123| 91久久精品电影网| 日韩电影二区| 精品午夜福利在线看| a级毛片在线看网站| 精品亚洲成a人片在线观看| 少妇 在线观看| 亚洲av.av天堂| 亚洲av日韩在线播放| 日韩欧美精品免费久久| 26uuu在线亚洲综合色| 欧美人与善性xxx| 一边亲一边摸免费视频| 久久国产精品大桥未久av| 欧美 亚洲 国产 日韩一| 日日啪夜夜爽| 十分钟在线观看高清视频www| 纵有疾风起免费观看全集完整版| 免费看光身美女| 狂野欧美白嫩少妇大欣赏| 色婷婷av一区二区三区视频| 欧美xxⅹ黑人| 热99久久久久精品小说推荐| 高清不卡的av网站| 人人妻人人澡人人看| 亚洲三级黄色毛片| 亚洲精品第二区| 国产黄片视频在线免费观看| 天天操日日干夜夜撸| 国产在线一区二区三区精| 欧美bdsm另类| 又黄又爽又刺激的免费视频.| 久久久精品94久久精品| 国产永久视频网站| 国产 精品1| 国产精品不卡视频一区二区| 天堂俺去俺来也www色官网| 人人妻人人添人人爽欧美一区卜| 日韩一本色道免费dvd| 久久国内精品自在自线图片| 母亲3免费完整高清在线观看 | 黄色毛片三级朝国网站| 如日韩欧美国产精品一区二区三区 | 亚洲精品美女久久av网站| av一本久久久久| 国产日韩欧美在线精品| 亚洲精品中文字幕在线视频| 制服人妻中文乱码| 亚洲综合色惰| 久久亚洲国产成人精品v| 精品人妻一区二区三区麻豆| 女人精品久久久久毛片| 人人妻人人澡人人看| 亚洲欧美成人综合另类久久久| 国产毛片在线视频| 国产精品久久久久久久久免| 女人精品久久久久毛片| 高清毛片免费看| 3wmmmm亚洲av在线观看| 青春草视频在线免费观看| 韩国高清视频一区二区三区| 一二三四中文在线观看免费高清| 国产精品人妻久久久影院| a级毛色黄片| 丝瓜视频免费看黄片| 永久免费av网站大全| 国产av国产精品国产| 成人漫画全彩无遮挡| 寂寞人妻少妇视频99o| 免费少妇av软件| 国产欧美日韩一区二区三区在线 | 三级国产精品片| 少妇的逼水好多| 精品国产国语对白av| 国产男人的电影天堂91| 国产精品久久久久久久电影| 精品亚洲成国产av| 免费观看无遮挡的男女| 蜜桃在线观看..| 成年av动漫网址| 黄色视频在线播放观看不卡| 91精品伊人久久大香线蕉| 亚洲综合色网址| 最近中文字幕2019免费版| 国产成人精品一,二区| 亚洲精品色激情综合| 夫妻性生交免费视频一级片| 欧美激情 高清一区二区三区| 久久99热6这里只有精品| 五月玫瑰六月丁香| 亚洲精品自拍成人| 亚洲精品成人av观看孕妇| 肉色欧美久久久久久久蜜桃| 国内精品宾馆在线| 久久ye,这里只有精品| 国产成人av激情在线播放 | 18+在线观看网站| 午夜福利视频在线观看免费| 亚洲人成网站在线播| 国产69精品久久久久777片| 少妇人妻精品综合一区二区| 18禁在线无遮挡免费观看视频| 黑人欧美特级aaaaaa片| 高清av免费在线| 国产成人91sexporn| 亚洲激情五月婷婷啪啪| 天堂俺去俺来也www色官网| 精品一区二区三区视频在线| 亚洲,一卡二卡三卡| 日韩成人av中文字幕在线观看| 欧美精品一区二区免费开放| 久久久久国产网址| 在线观看人妻少妇| 晚上一个人看的免费电影| 成年人午夜在线观看视频| 男男h啪啪无遮挡| 免费播放大片免费观看视频在线观看| 亚洲精品乱码久久久久久按摩| 久久久久久久久久久久大奶| 日韩,欧美,国产一区二区三区| 亚洲欧美日韩另类电影网站| 亚洲国产最新在线播放| 少妇猛男粗大的猛烈进出视频| 国产亚洲欧美精品永久| 国产在视频线精品| 亚洲国产精品一区三区| 一二三四中文在线观看免费高清| 亚洲国产最新在线播放| 国产视频内射| 乱码一卡2卡4卡精品| 亚洲成人一二三区av| 99热6这里只有精品| 欧美xxxx性猛交bbbb| 国产老妇伦熟女老妇高清| 免费观看无遮挡的男女| 日日摸夜夜添夜夜爱| 大片电影免费在线观看免费| 久久午夜综合久久蜜桃| 久久午夜综合久久蜜桃| 校园人妻丝袜中文字幕| 美女福利国产在线| 嫩草影院入口| 日本欧美国产在线视频| 香蕉精品网在线| 高清毛片免费看| 赤兔流量卡办理| 欧美+日韩+精品| 精品熟女少妇av免费看| 亚洲第一区二区三区不卡| 91aial.com中文字幕在线观看| 久久精品国产自在天天线| 亚洲av成人精品一二三区| 97在线视频观看| 日韩成人av中文字幕在线观看| 国产黄片视频在线免费观看| 国产老妇伦熟女老妇高清| 建设人人有责人人尽责人人享有的| 亚洲精品中文字幕在线视频| 亚洲美女搞黄在线观看| 夜夜爽夜夜爽视频| 亚洲性久久影院| 久久女婷五月综合色啪小说| 日本黄大片高清| 交换朋友夫妻互换小说| 亚洲av成人精品一二三区| 亚洲av.av天堂| 免费观看性生交大片5| 卡戴珊不雅视频在线播放| 一级毛片我不卡| 青春草国产在线视频| 女的被弄到高潮叫床怎么办| 精品亚洲成a人片在线观看| 秋霞伦理黄片| 色5月婷婷丁香| 国产在线一区二区三区精| 国产在线一区二区三区精| 日本91视频免费播放| 老女人水多毛片| kizo精华| 免费大片黄手机在线观看| 亚洲av日韩在线播放| 日本色播在线视频| 精品国产一区二区久久| 亚洲内射少妇av| 免费大片18禁| 91久久精品国产一区二区成人| 国产国拍精品亚洲av在线观看| 国产视频首页在线观看| 久久久久久久精品精品| 一级毛片我不卡| 国产69精品久久久久777片| 国产成人aa在线观看| 一区二区三区乱码不卡18| 成年人午夜在线观看视频| 亚洲精品乱码久久久v下载方式| 99热这里只有精品一区| 18禁观看日本| 精品一区二区三卡| 欧美老熟妇乱子伦牲交| 国产成人freesex在线| 国产免费现黄频在线看| 制服丝袜香蕉在线| 欧美亚洲 丝袜 人妻 在线| 日本wwww免费看| 777米奇影视久久| 国产精品三级大全| 免费av中文字幕在线| 夫妻午夜视频| 成人漫画全彩无遮挡| av卡一久久| 老司机影院成人| 啦啦啦啦在线视频资源| 免费黄网站久久成人精品| tube8黄色片| 午夜影院在线不卡| a级片在线免费高清观看视频| 成人免费观看视频高清| av网站免费在线观看视频| 精品久久久久久电影网| 一级毛片电影观看| 激情五月婷婷亚洲| 一二三四中文在线观看免费高清| 久久久久网色| 精品99又大又爽又粗少妇毛片| 久久这里有精品视频免费| 国产探花极品一区二区| 精品国产露脸久久av麻豆| 国产免费一级a男人的天堂| 欧美亚洲 丝袜 人妻 在线| 久久久久久久久大av| 男女国产视频网站| 日日摸夜夜添夜夜爱| 日本av手机在线免费观看| 一区二区三区免费毛片| 国产精品一二三区在线看| 丝袜喷水一区| 91精品三级在线观看| 国产免费福利视频在线观看| 婷婷色综合大香蕉| 一本色道久久久久久精品综合| 免费大片18禁| 中文字幕精品免费在线观看视频 | 肉色欧美久久久久久久蜜桃| 这个男人来自地球电影免费观看 | 一区二区三区免费毛片| 免费日韩欧美在线观看| 免费看av在线观看网站| 国产在线一区二区三区精| 免费人成在线观看视频色| 国产在线视频一区二区| 麻豆成人av视频| 日本黄大片高清| 午夜日本视频在线| 桃花免费在线播放| 中国国产av一级| 国语对白做爰xxxⅹ性视频网站| 久久久精品区二区三区| 涩涩av久久男人的天堂| 18+在线观看网站| 国产成人精品无人区| 如日韩欧美国产精品一区二区三区 | 热re99久久精品国产66热6| 亚洲图色成人| 免费少妇av软件| 色5月婷婷丁香| 伦精品一区二区三区| 在线天堂最新版资源| 欧美精品一区二区免费开放| 精品久久久久久电影网| √禁漫天堂资源中文www| 亚洲精品av麻豆狂野| 亚洲欧美中文字幕日韩二区| 国产一级毛片在线| 亚洲丝袜综合中文字幕| 涩涩av久久男人的天堂| 亚洲av成人精品一区久久| 在线播放无遮挡| 岛国毛片在线播放| 天天躁夜夜躁狠狠久久av| 18禁裸乳无遮挡动漫免费视频| 黄片播放在线免费| 曰老女人黄片| 久久精品国产a三级三级三级| 热99久久久久精品小说推荐| 男女免费视频国产| 日韩不卡一区二区三区视频在线| 国产视频首页在线观看| 亚洲av免费高清在线观看| 熟妇人妻不卡中文字幕| 国产成人一区二区在线| 22中文网久久字幕| 亚洲av成人精品一二三区| 大码成人一级视频| 亚洲成人手机| 亚洲精品乱码久久久v下载方式| 亚洲国产精品999| 精品少妇黑人巨大在线播放| 亚洲无线观看免费| 国产av精品麻豆| 欧美激情极品国产一区二区三区 | 亚洲欧美清纯卡通| 蜜桃国产av成人99| 一区二区日韩欧美中文字幕 | 国产成人免费无遮挡视频| 亚洲人与动物交配视频| 久久久久久久国产电影| 久热这里只有精品99| 九九在线视频观看精品| 九九久久精品国产亚洲av麻豆| 在线观看美女被高潮喷水网站| 少妇人妻精品综合一区二区| 免费少妇av软件| 精品国产露脸久久av麻豆| 哪个播放器可以免费观看大片| 久久久国产一区二区| 免费观看的影片在线观看| 18禁观看日本| 97超碰精品成人国产| 国产免费福利视频在线观看| 久久99蜜桃精品久久| 青春草亚洲视频在线观看| 99热6这里只有精品| 在线播放无遮挡| 欧美日韩成人在线一区二区| 国产精品三级大全| 久久精品夜色国产| 久热这里只有精品99| 亚洲av日韩在线播放| 啦啦啦在线观看免费高清www| 国产精品久久久久成人av| 亚洲精品一区蜜桃| 久久狼人影院| 一级毛片黄色毛片免费观看视频| 少妇丰满av| 久久青草综合色| 免费观看的影片在线观看| 乱码一卡2卡4卡精品| 免费观看在线日韩| 一本久久精品| 99热全是精品| 欧美激情 高清一区二区三区| 91精品三级在线观看| 国产午夜精品久久久久久一区二区三区| 伦理电影免费视频| 国产成人aa在线观看| 亚洲丝袜综合中文字幕| 一级毛片黄色毛片免费观看视频| 街头女战士在线观看网站| 黑人猛操日本美女一级片| av女优亚洲男人天堂| 精品一区二区三卡| 国产免费福利视频在线观看| 精品一区二区免费观看| 中文字幕av电影在线播放| 成人二区视频| 在线精品无人区一区二区三| 亚洲精品日韩在线中文字幕| 免费不卡的大黄色大毛片视频在线观看| 亚洲四区av| 伦精品一区二区三区| 中国三级夫妇交换| 欧美日韩av久久| 亚洲精品美女久久av网站| 777米奇影视久久| 999精品在线视频| 日本黄大片高清| 日本av手机在线免费观看| 亚洲欧美色中文字幕在线| 大香蕉久久成人网| 午夜福利影视在线免费观看| 蜜桃在线观看..| 大香蕉久久成人网| 在线观看国产h片| 久久久久国产精品人妻一区二区| 男人添女人高潮全过程视频| 亚洲精品aⅴ在线观看| 九色成人免费人妻av| 美女内射精品一级片tv| 99热这里只有精品一区| 国产在线视频一区二区| 一级,二级,三级黄色视频| 黑人猛操日本美女一级片| 国产精品久久久久久精品电影小说| 熟妇人妻不卡中文字幕| videosex国产| 菩萨蛮人人尽说江南好唐韦庄| 一二三四中文在线观看免费高清| 久久精品国产亚洲av天美| 人妻 亚洲 视频| 亚洲精品自拍成人| 伦理电影免费视频| 国产精品国产三级国产专区5o| 秋霞伦理黄片| 亚洲激情五月婷婷啪啪| 99九九线精品视频在线观看视频| 日韩一区二区三区影片| 另类亚洲欧美激情| 夜夜看夜夜爽夜夜摸| 五月玫瑰六月丁香| 大香蕉久久网| 七月丁香在线播放| 中文字幕精品免费在线观看视频 | 国产成人精品婷婷| 亚洲四区av| 又黄又爽又刺激的免费视频.| 99精国产麻豆久久婷婷| 国产精品人妻久久久影院| 黄色视频在线播放观看不卡| 日韩欧美一区视频在线观看| 大香蕉久久成人网| 亚洲欧美日韩卡通动漫| 国产日韩欧美亚洲二区| .国产精品久久| 少妇人妻久久综合中文| 王馨瑶露胸无遮挡在线观看| 亚洲国产毛片av蜜桃av| 国产女主播在线喷水免费视频网站| 亚洲欧洲精品一区二区精品久久久 | 一个人免费看片子| 国产av码专区亚洲av| 久久久久久伊人网av| 99久国产av精品国产电影| 新久久久久国产一级毛片| 久久精品国产自在天天线| 精品久久蜜臀av无| 欧美日韩国产mv在线观看视频| 一级爰片在线观看| 男女边摸边吃奶| 高清不卡的av网站| 精品一区二区三卡| 母亲3免费完整高清在线观看 | 国产伦精品一区二区三区视频9| 日韩免费高清中文字幕av| 国产淫语在线视频| 2022亚洲国产成人精品| 热re99久久国产66热| 97在线人人人人妻| 免费观看av网站的网址| 男女高潮啪啪啪动态图| 亚洲不卡免费看| 天天操日日干夜夜撸| 黄色怎么调成土黄色| 日日啪夜夜爽| 婷婷色麻豆天堂久久| 天美传媒精品一区二区| 午夜免费男女啪啪视频观看| 免费人妻精品一区二区三区视频| 久久国产亚洲av麻豆专区| 在线观看免费视频网站a站| 99久久中文字幕三级久久日本| 亚洲精品中文字幕在线视频| 久久国产精品大桥未久av| 日本色播在线视频| 在线免费观看不下载黄p国产| 久久久国产精品麻豆| 色视频在线一区二区三区| kizo精华| 黑人猛操日本美女一级片| 一级爰片在线观看| 日日撸夜夜添| 高清毛片免费看| 国产黄频视频在线观看| 美女国产高潮福利片在线看| 国产av码专区亚洲av| 一区二区三区四区激情视频| 伦精品一区二区三区| 国产亚洲午夜精品一区二区久久| 久久精品国产鲁丝片午夜精品| 国产午夜精品久久久久久一区二区三区| 国产黄色视频一区二区在线观看| 久久99一区二区三区| 国产免费一区二区三区四区乱码| 国产欧美亚洲国产| 精品久久国产蜜桃| 在线观看免费视频网站a站| 爱豆传媒免费全集在线观看| 飞空精品影院首页| 国产黄频视频在线观看| 亚洲精品成人av观看孕妇| 寂寞人妻少妇视频99o| 免费观看的影片在线观看| 国产av精品麻豆| 久久午夜福利片| 性高湖久久久久久久久免费观看| 成人影院久久| 在线观看免费视频网站a站| 九色成人免费人妻av| 亚洲av成人精品一区久久| 亚洲精品乱久久久久久| 又黄又爽又刺激的免费视频.| 桃花免费在线播放| 精品人妻一区二区三区麻豆| 视频区图区小说| 蜜桃久久精品国产亚洲av| 欧美日韩成人在线一区二区| 精品人妻熟女av久视频| 成人毛片a级毛片在线播放| 一本一本综合久久| 男女免费视频国产| 老司机影院成人| 国产成人精品无人区| 一本一本综合久久| 国产女主播在线喷水免费视频网站| 中文乱码字字幕精品一区二区三区| 夜夜看夜夜爽夜夜摸| 婷婷色综合www| 精品国产露脸久久av麻豆| 久久国内精品自在自线图片| 精品国产露脸久久av麻豆| 中文字幕精品免费在线观看视频 | 伦精品一区二区三区| 久久影院123| 国产精品不卡视频一区二区| 日韩强制内射视频| 国产av一区二区精品久久| √禁漫天堂资源中文www| 在线观看免费视频网站a站| 午夜久久久在线观看| 日韩一区二区视频免费看| 亚洲内射少妇av| av专区在线播放| 男女高潮啪啪啪动态图| 免费看不卡的av| 在线观看国产h片| 亚洲性久久影院| 亚洲精品久久午夜乱码| 看免费成人av毛片| 91精品国产国语对白视频| 黄色配什么色好看| 少妇猛男粗大的猛烈进出视频| videosex国产| 国产免费视频播放在线视频| 亚洲av二区三区四区| 大片免费播放器 马上看| 欧美日韩综合久久久久久| 免费观看的影片在线观看| 少妇熟女欧美另类| 纵有疾风起免费观看全集完整版| 成年人午夜在线观看视频| 精品人妻偷拍中文字幕| 亚洲人与动物交配视频| 这个男人来自地球电影免费观看 | 内地一区二区视频在线| 亚洲综合色惰| 中文字幕制服av| 亚洲国产av新网站| videos熟女内射| 大码成人一级视频| videos熟女内射| 大码成人一级视频| 80岁老熟妇乱子伦牲交| 青春草亚洲视频在线观看| 国产成人freesex在线| 99久久综合免费| 国产有黄有色有爽视频| av女优亚洲男人天堂| 亚洲精品乱码久久久久久按摩| 欧美丝袜亚洲另类| 日韩不卡一区二区三区视频在线| 亚洲精品一区蜜桃| 成人国产av品久久久| 2018国产大陆天天弄谢| 精品人妻在线不人妻| 国产日韩一区二区三区精品不卡 | 国产日韩欧美视频二区| 久久免费观看电影| 成人二区视频| 国产精品秋霞免费鲁丝片| 亚洲av成人精品一区久久| av一本久久久久| 国产成人一区二区在线| 美女福利国产在线| 肉色欧美久久久久久久蜜桃| 人妻少妇偷人精品九色| 中文欧美无线码| 丝瓜视频免费看黄片| av又黄又爽大尺度在线免费看| 少妇被粗大的猛进出69影院 |