基因識(shí)別問題及其算法實(shí)現(xiàn)

1、基因識(shí)別問題及其算法實(shí)現(xiàn)一、背景介紹dna是生物遺傳信息的載體，只化學(xué)名稱為脫氧核糖核酸(deoxyribonucleic acid,縮 寫為dna)cdna分子是一種長(zhǎng)鏈聚合物,dna序列由腺喋吟adenine, 鳥瞟吟(guanine, g),胞i密噪(cytosine, c),胸腺n密唳thymine, t)這四種核昔酸(nucleotide)符號(hào)按一定 的順序連接1何成。其中帶有遺傳訊息的dna片段稱為基因(gm>)(見圖1第一行)。其他的 dna序列片段，冇些直接以自身構(gòu)造發(fā)揮作川，冇些則參與調(diào)控遺傳訊息的表現(xiàn)。在真核生物的dna序列屮，基因通常被劃分為許多間隔的片段(見圖1第

2、二行)，其屮編碼蛋白質(zhì)的部分，即編碼序列(coding sequence)片段，稱為外顯子(exon),不編碼的部 分稱為內(nèi)含子(buron)。外顯子在dva序列剪接(splicing)后仍然會(huì)被保存下來，并可在外顯子(exon) 內(nèi)含子(intron )圖1真核牛物dna用列(基因序列)結(jié)構(gòu)示意圖蛋白質(zhì)合成過程中被轉(zhuǎn)錄(transcription)復(fù)制(replication)而合成為蛋白質(zhì)(見圖2)。dna序列通過遺傳編碼來儲(chǔ)存信息，指導(dǎo)蛋白質(zhì)的合成，把遺傳信息準(zhǔn)確無誤地傳遞到蛋白質(zhì)(protein)上去并實(shí)現(xiàn)各種生命功能。dna序列蛋白質(zhì)序列圖2蛋白質(zhì)結(jié)構(gòu)示意圖對(duì)大最、復(fù)雜的基因序列的分

3、析，傳統(tǒng)生物學(xué)解決問題的方式是基于分子實(shí)驗(yàn)的方法, 其代價(jià)高昂。諾貝爾獎(jiǎng)獲得者w.吉爾伯特(waltergilbert, 1932；【美】,第一個(gè)制備 出混合脫氧核糖核酸的科學(xué)家)1991年曾經(jīng)指出：“現(xiàn)在，基于全部基因序列都將知曉，并 以電子可操作的方式駐超在數(shù)據(jù)庫中，新的生物學(xué)研究模式的出發(fā)點(diǎn)應(yīng)是理論的。一個(gè)科學(xué) 家將從理論推測(cè)出發(fā)，然后再回到實(shí)驗(yàn)中去，追蹤或驗(yàn)證這些理論假設(shè)?！彪S著世界人類基 因組工程計(jì)劃的順利完成，通過物理或數(shù)學(xué)的方法從人最的dna序列屮獲取豐富的生物信息，對(duì)生物學(xué)、醫(yī)學(xué)、藥學(xué)等諸多方而都具有重要的理論意義和實(shí)際價(jià)值，也是h前生物信 息學(xué)領(lǐng)域的-個(gè)研究熱點(diǎn)。二、數(shù)字序列

4、映射與頻譜3周期性：對(duì)給定的rwa序列，怎么去識(shí)別出其中的編碼序列（即外顯了），也稱為基因預(yù)測(cè)，是 一個(gè)尚未完全解決的問題，也是當(dāng)前生物信息學(xué)的一個(gè)最棊礎(chǔ)、最首婆的問題?；蝾A(yù)測(cè)問題的一類方法是基于統(tǒng)計(jì)學(xué)的。很多國(guó)際生物數(shù)據(jù)網(wǎng)站上也有“基因識(shí)別” 的算法。比如知名的數(shù)據(jù)網(wǎng)站http:/gcncs. mi t. odu/genscan. html提供的壟因識(shí)別軟件 genscan （由斯坦福人學(xué)研究人員研發(fā)的、可免費(fèi)使用的基因預(yù)測(cè)軟件），主要就是基于隱馬 爾科夫鏈（1imm）方法。但是，它預(yù)測(cè)人的基因組屮有45000個(gè)基因，相當(dāng)于現(xiàn)在普遍認(rèn)可 數(shù)目的兩倍。另外，統(tǒng)計(jì)預(yù)測(cè)方法通常需要將編碼序列信息

5、已知的dna序列作為訓(xùn)練數(shù)據(jù) 集來確定模型中的參數(shù)，從而提高模型的預(yù)測(cè)水平。但在對(duì)基因信息了解不多的悄況下，基 因識(shí)別的準(zhǔn)確率會(huì)明顯下降。因此在目前基因預(yù)測(cè)研究屮，采用信號(hào)處理與分析方法來發(fā)現(xiàn)基因編碼序列也受到廣泛 重視。1. 數(shù)字序列映射在lwa序列研究中，首先需要把a(bǔ)、7 g、c四種核昔酸的符號(hào)序列，根據(jù)一定的規(guī) 則映射成相應(yīng)的數(shù)值序列，以便于對(duì)其作數(shù)字處理。令z = a,t,g,c,長(zhǎng)度（即核廿酸符號(hào)個(gè)數(shù)，又稱堿基對(duì)（base paid長(zhǎng)度，單位 記為切）為n的任意qm4序列，可表達(dá)為s = sn i s“， = 0丄2,n-1即a、t、g、c的符號(hào)序列s： s0,su,sw l?，F(xiàn)對(duì)于

6、任意確定的bwl ,令sn = bsn工 b稱之為vb刖映射，于是生成相應(yīng)的01序列（即二進(jìn)制序列）uhn：地0,旳j1,，,uhn -i (b e i )o例如，假設(shè)給定的一段dma序列片段為s = atcgtactg,則所牛成的四個(gè)0-1序列分 別為:uan: 1,0,0,0,0,1,0,0,0 ；%“： (),0,0丄0,0,0,0,1；ucn: 0,0,1,0,0,0,1,0,0 ；utn: 0,1,0,0,1,0,0,1,0。這樣產(chǎn)生的四個(gè)數(shù)7序列乂稱為dna序列的指示序列indicator sequence)。2. 頻譜3周期性為研究dn4編碼序列(外顯子)的特性，對(duì)指示序列分別做

7、離散fourier變換(dft)n-1_ 2兀nk5幻=工坷川幺7 n ,a=0,l,n-l(1)n=0以此可得到四個(gè)長(zhǎng)度均為n的復(fù)數(shù)序列uhk, bel.計(jì)算每個(gè)復(fù)序列uhk的平方功率譜，并相加則得到整個(gè)dna序列s的功率譜序列pr：pk = ua*|2 + utkf + ugkf 4-uckf, k = o,hn-(2)對(duì)于同一段dm4序列，其外顯子與內(nèi)含子序列片段的功率譜通常表現(xiàn)出不同的特性圖3編號(hào)為bk006948.2的酵母基因dva序列的功率譜(因?yàn)閷?duì)稱性，實(shí)際這里只給出了功率譜圖 的一半)。 上圖是基因上一段外顯了(區(qū)間為81787, 82920,長(zhǎng)1134切)對(duì)應(yīng)的指示序列映射的

8、功率 譜,它具有3周期性;(b)下圖是基因上一段內(nèi)含子(區(qū)間為96361,97551，長(zhǎng)1191切)的指示序列的功 率譜，它不具有3-周期性。n可以看到：外顯子序列的功率譜曲線在頻率k=處，具有較人的頻譜峰值3value),而內(nèi)含了則沒冇類似的峰值。這種統(tǒng)計(jì)現(xiàn)彖被稱為堿壟的3周期(3-base periodicity) 3o記dna序列s的總功率譜的平均值為n-1_工p幻e =旦(3)而將dna序列在特定位置，即k =n二丄處的功率譜值，與整個(gè)序列s的總功率譜的平均值的3e比率稱為dna序列的“信噪比” (signalnoise ratio, snr),即dna序列的信噪比值的大小，既表示頻譜

9、峰值(peak value)的相對(duì)高度，也反映編碼或 非編碼序列3周期性的強(qiáng)弱。信噪比/?大于某個(gè)適當(dāng)選定的閾值他(比如/?0=2), ina序列上編碼序列片段(外 顯了)通常滿足的特性，而內(nèi)含了則一般不具冇該性質(zhì)叫在dm4序列sn, n=0,l,2,n 1中，若n為3的倍數(shù)，將核廿酸符號(hào)b ei = a,t,g,c出現(xiàn)在該序列的 0,3,6,. n-3 與 l,4,7,.n-2 以及 2,5,&.n-1 等位置上n的頻數(shù)分別記為必和z方，則二處的總功率譜值即為c n 2/r/rn-1_ : *n只寧=工匕卻=z為訕“ =工£訕小£山-也=1be!be/be/w=0

10、w=03戻7易見，當(dāng)四種核昔酸符號(hào)b (/?/)在序列的上述第一、第二、第三個(gè)子序列上出現(xiàn)的頻n數(shù)風(fēng)小心越接近相等吋，一處的譜值也就越接近于零。所以，基因外顯子序列的功率譜n曲線，在二頻率處具有較人的頻譜峰值(peak value),反映了在基因外顯子片段上，四種核 3莒酸符號(hào)在序列的三個(gè)了序列上分布的“非均衡性”。通常認(rèn)為這種現(xiàn)象源于編碼基因序列 “密碼子”(code)使用的偏向性(bias)。雖然目前對(duì)此現(xiàn)象產(chǎn)生的“機(jī)理”還不是十分 地清楚，但是頻譜的3-周期性被普遍認(rèn)為是可用于識(shí)別基因編碼序列(外顯子)的一個(gè)重要 的特征信息。3. 基因識(shí)別頻譜峰值特征的發(fā)現(xiàn)，或者頻譜與信噪比概念的引入，

11、其最終冃的是要探測(cè)、預(yù)報(bào)一個(gè)尚未被注釋的完整的dma序列的所侑基因編碼序列(外顯了)片段。圖4基于序列頻譜3周期性的的基因預(yù)測(cè)方法流程圖已經(jīng)冇一些研究者提出了識(shí)別基因的算法（如參見及其后面的文獻(xiàn)）。冃前利川信噪 比的基因識(shí)別算法通常有兩種：一是固定長(zhǎng)度窗ii滑動(dòng)法叫另一是移動(dòng)信噪比曲線識(shí)別 法?；诠潭ㄩL(zhǎng)度滑動(dòng)窗口上頻譜曲線的基因識(shí)別方法:對(duì)一個(gè)dna序列s和它的指示序列他川 " w /皿二0丄2,n 1。収長(zhǎng)度m （通 常取為3的倍數(shù),例如m=99, 129,255,513等）作為固定窗口長(zhǎng)度。m -1 m -對(duì)任意n<o<n<n- ）,在以料為中心的長(zhǎng)度為m的序

12、列片段w+2 2上（當(dāng)n接近序列的兩端時(shí)，窗ii實(shí)際啟效長(zhǎng)度可能會(huì)小于m）,作四個(gè)指示序列的離散m-fourier 變換(dft)i=n+_ .imubk=工 uhie j m , k=o,l,m-l2并求出它在理處總頻譜/9（/7；）,即2r m22m以3+叫+噲)把這樣得到的頻譜值"理），二0,1,2,2-1,經(jīng)過標(biāo)準(zhǔn)化處理（即除以最大頻譜值m。酬"心；丁）并回出其頻譜曲線圖5固定長(zhǎng)度滑動(dòng)窗口的頻譜p = 1心;)曲線(人類線粒體棊因，nc_012920_l.fasta)圖中紅色水平細(xì)線條是dna序列實(shí)際的基因外顯子的區(qū)間?；瑒?dòng)窗i i頻譜p(/2；)曲線的 峰與皋因外

13、顯了區(qū)間具冇“對(duì)應(yīng)”關(guān)系?；赿na序列上“移動(dòng)序列”信噪比曲線的基因識(shí)別方法：設(shè)已知dna序列s和它的指示序列uhn, bwl ,=1。對(duì)任意n(0<n<7v-l ),通常取3的倍數(shù)并逐漸增大。在的左邊一個(gè)長(zhǎng)度為77的序列片段(), n-lj h,相應(yīng)的了序列s稱為dna序列s的“移動(dòng)子序列”，作該移動(dòng)了序列對(duì)應(yīng)的四 個(gè)指示序列的離散fourier變換(£>ft)i=n-_ j2jrik如幻=£唄弘昕，心0,1,1i=q并求出移動(dòng)子序列n = 0丄，n 1上的信噪比rnnr r r 1r r1 r r r 1r r r 1p匕才+ “丁匕+如匕+&quo

14、t;c§n-i其中可旳為移動(dòng)子序列sow_幻的功率譜的平均值£=o在坐標(biāo)系中呦出移動(dòng)序n列so*的信噪比曲線/?(稱為信噪比移動(dòng)曲線(snr walk curve),見圖6)1412c)ccccns7000030003500400045005000550060006500nuclootide position n圖6 dna移動(dòng)序列其指示序列的信噪比曲線。(人類線粒體基因，nc_012920_l.fasta)圖中紅色水平細(xì)線條是dna序列實(shí)際的基因外顯子的區(qū)間。dna序列的信噪比移動(dòng)曲線 的峰、谷與基因外顯子區(qū)間的端點(diǎn)也具有較“明顯的”的對(duì)應(yīng)關(guān)系。三、請(qǐng)研究的幾個(gè)問題：1.

15、功率譜與信噪比的快速算法對(duì)于很長(zhǎng)的dna序列，在計(jì)算其功率譜或信噪比時(shí)，離散fourier變換（dft）的總體計(jì) 算量仍然很人，會(huì)影響到所設(shè)計(jì)的基因識(shí)別算法的效率。大家能否對(duì)映射，探求功率 譜與信噪比的某種快速計(jì)算方法？在基因識(shí)別研究中，為了通過引入更好的數(shù)值映射而獲取dna序列更多的信息，除了 上面介紹的wss映射外，實(shí)際上人們還研究過許多不同的數(shù)值映射方法。例如，著名的 z-curve映射（參見或者附件1）。試探討z-curve映射的頻譜為信噪比和voss映射下的頻 譜與信噪比之間的關(guān)系；此外，能否對(duì)實(shí)數(shù)映射，如：4to,ct1,gt2,tt3,也給出功率譜與信噪比 的快速計(jì)算公式？2對(duì)不

16、同物種類型基因的閾值確定對(duì)特定的基因類型的序列，將其信噪比的判別閾值取為r°=2,帶有一定的主 觀性、經(jīng)驗(yàn)性。對(duì)不同的基因類型，所選取的判別閾值也許應(yīng)該是不同的。附件小給出了來 口于苦名的牛:物數(shù)據(jù)網(wǎng)站：hup://guidc/的兒個(gè)基因序列數(shù)據(jù)，另外也 給出了帶有編碼外顯子信息的10（）個(gè)人和鼠類的，以及200個(gè)哺乳動(dòng)物類的基因序列的樣本 數(shù)據(jù)集合。大家還町以從生物數(shù)據(jù)庫下載更多的數(shù)據(jù)，找你們認(rèn)為具冇代表性的棊因序列， 并對(duì)每類基因研究其閾值確定方法和閾值結(jié)果。此外，對(duì)按照頻譜或信噪比特征將編碼與非 編碼區(qū)間分類的有效性，以及分類識(shí)別時(shí)所產(chǎn)住

17、的分類錯(cuò)誤作適當(dāng)分析。3.基因識(shí)別算法的實(shí)現(xiàn)我們的目的是要探測(cè)、預(yù)報(bào)尚未被注釋的、完整的dna序列的所有基因編碼序列（外 顯了）。冃前基因識(shí)別方面的多數(shù)篦法結(jié)果還不是很充分。例如前面所列舉的某些基因識(shí)別 算法，由于dna序列隨機(jī)噪聲的影響等原因，還很難“精確地”確定基因外顯子區(qū)間的兩 個(gè)端點(diǎn)。對(duì)此，你的建模i才i隊(duì)有沒有更好的解決方法？請(qǐng)對(duì)你們所設(shè)計(jì)的基因識(shí)別算法的準(zhǔn)確率 做出適當(dāng)評(píng)估，并將算法用于對(duì)附件中給出的6個(gè)未被注釋的dna序列（gene6）的編碼 區(qū)域的預(yù)測(cè)。4. 延展性研究在基因識(shí)別研究中，還冇很多問題冇待深入探討。比如（1）采用頻譜或信噪比這樣單一的判別特征，也許是影響、限制基因

18、識(shí)別正確率的一 個(gè)重要原因。人們發(fā)現(xiàn)，對(duì)某些dn4序列而言，其部分編碼序列（外顯子），尤其是短的（長(zhǎng) 度小于1 00/7/2 ）的編碼序列，就可能不具有頻譜或者信噪比顯著性。你們團(tuán)隊(duì)能否總結(jié)，甚 至獨(dú)口提出一些識(shí)別基因編碼序列的其它特征指數(shù)，并對(duì)此做相關(guān)的分析？（2）“基因突變”是生物醫(yī)學(xué)等方面的一個(gè)關(guān)注熱點(diǎn)?；蛲蛔儼╠na序列中單個(gè)核 廿酸的替換，刪除或者插入等。那么，能否利用頻譜或信噪比方法去發(fā)現(xiàn)基因編碼序列可能 存在的突變呢？上而捉出的基于頻譜3-周期性的基因預(yù)測(cè)四個(gè)方而問題中，“快速算法”與“閾值確定” 是為設(shè)計(jì)基因預(yù)測(cè)算法做準(zhǔn)備的。此外，在最后的延展性研究中，各隊(duì)也可以對(duì)你們口己認(rèn) 為冇價(jià)值的其它相關(guān)問題展開探討。參考文獻(xiàn)：【1 】burge, c., karlin, s., 1997. prediction of complete gene structures in human genomic dna. j. mol. biol. 268, 78-94.【2】 anastassiou, d., 2000. frequency-domain analysis of biomolecular sequences. bioinformatics 16, 1073-1081.3 kotlar, d., lavner, y., 20