本發(fā)明涉及基因編碼序列(外顯子)的預測，提出了一種周期三臨近特征外顯子預測方法。
背景技術(shù)：
人類基因組計劃的順利完成，使得剖析大量的DNA序列中豐富的生物信息成為可能，并對醫(yī)學、藥學、生物學等諸多領(lǐng)域具有重要的理論和實用價值。其中DNA序列信息的分析已逐漸成為分析化學的研究熱點?；蚴且唤M包含遺傳信息的DNA序列，其能夠通過遺傳編碼指導蛋白質(zhì)的合成，把遺傳信息準確無誤地傳遞到蛋白質(zhì)并實現(xiàn)生命的各種活動。存在于真核基因中的蛋白質(zhì)編碼序列由多不連續(xù)的片斷組成，其中外顯子編碼蛋白質(zhì)，內(nèi)含子不編碼蛋白質(zhì)，通過剪切后去掉內(nèi)含子后多個外顯子連在一起，最后翻譯成具有生命現(xiàn)象的蛋白質(zhì)。
杜竹青的文章《一種提高外顯子預測的改進周期3消噪策略》側(cè)重于采用濾波算法實現(xiàn)周期三的去噪，邵建峰的《DNA序列信號3周期特性》則分析了周期三法閾值的確定方法，但上述文獻均未考慮周期三行為的存在性。
田元新《外顯子周期三行為特征的研究》對21個基因序列的103條外顯子進行研究，發(fā)現(xiàn)許多外顯子不具有周期三行為，但將它們連在一起編碼蛋白質(zhì)時大都表現(xiàn)出周期三行為。其研究結(jié)果表明：單個外顯子的周期三行為與蛋白質(zhì)編碼序列長度、堿基在密碼子偏好性以及氨基酸密碼子的使用均有密切關(guān)系。文中給出的人類8個基因共43條外顯子功率譜最大值對應的頻率并非都現(xiàn)出周期三行為，但其功率譜的最大值均出現(xiàn)在1/3頻率處或該頻率臨近。田元新的方法是基于三聯(lián)碼、蛋白質(zhì)的翻譯次序以及密碼子使用偏好性對外顯子序列進行調(diào)整和預測，但是該方法實際操作的難度大且密碼子使用偏好性并不是有效識別外顯子的指標。
常規(guī)FASTA和BLAST等基于序列相似性對比的預測方法對于短長度的外顯子容易漏檢。基于機器學習的外顯子預測方法諸如神經(jīng)網(wǎng)絡、動態(tài)規(guī)劃等方法對迭代階數(shù)要求較高，階數(shù)過低不能反映完整信息，階數(shù)過高增加訓練時間。
技術(shù)實現(xiàn)要素：
針對現(xiàn)有的DNA序列周期三外顯子預測方法存在的易漏檢、檢測效率低的問題，本專利旨在提出一種基于周期三臨近特征的外顯子預測方法，以達到更好的預測進度和 預測效率。
本發(fā)明旨在提供一種基于周期三臨近特征的外顯子預測方法，包括如下步驟：
步驟1、訓練階段：
選取已知的DNA序列，記錄長度區(qū)間，建立DNA堿基序列的功率譜，并分別選取3倍數(shù)的外顯子和內(nèi)因子序列計算模板功率譜特征值；
步驟2、識別階段：
取待識別DNA序列，滑動與訓練階段同樣的長度區(qū)間，計算得到每組識別功率譜特征值；
步驟3、計算訓練階段模板的功率譜特征值與識別階段的每組識別功率譜特征值的歐式距離，對1/3頻率的臨近區(qū)間進行距離判斷；
步驟4、根據(jù)模板功率譜特征值與識別功率譜特征值的歐式距離，判斷外顯子位置。
上述的基于周期三臨近特征的外顯子預測方法，所述的步驟1具體為：
對DNA序列中每種堿基，分別將其位置記錄為序列，得到4個堿基A、T、G、C的位置的序列和n為時間的離散值，式中，l表示DNA序列中外顯子、內(nèi)因子和外顯子內(nèi)因子混合的標記，記l＝1表示外顯子序列，記l＝2表示內(nèi)因子序列，舍棄內(nèi)因子和外顯子混合序列，將每一個序列按照如下公式進行快速傅里葉變換得到每個序列傅里葉變換值和
式中，j為虛數(shù)單位，n為時間的離散值，k為頻域的離散值，N為周期長度，
采用如下公式，選取3倍數(shù)的外顯子和內(nèi)因子序列計算模板功率譜特征值，記為Sl[q]，其中q為3倍頻域的離散值
3倍頻域的離散值q＝…N/3-2,N/3-1,N/3,N/3+1,N/3+2…。
上述的基于周期三臨近特征的外顯子預測方法，所述的步驟2具體為：
取待識別DNA序列，滑動與訓練階段同樣的長度區(qū)間，對待識別DNA序列中每種堿基，分別將其位置記錄為序列，得到4個堿基A、T、G、C的位置的待識別序列ya(n)，yt(n)，yg(n)和yc(n)，將每一個待識別序列按照如下公式進行快速傅里葉變換得到每個待識別序列傅里葉變換值Ya[k]、Yt[k]、Yg[k]和Yc[k]：
式中，j為虛數(shù)單位，n為時間的離散值，k為頻域的離散值，N為周期長度，
采用如下公式，計算得到每組識別功率譜特征值，記為Zp[k]；
Zp[k]＝|Ya[k]|2+|Yt[k]|2+|Yg[k]|2+|Yc[k]|2
其中p為滑動數(shù)，p＝1,2,…,N-M+1，M為滑動窗大小，N＞＞M。
上述的基于周期三臨近特征的外顯子預測方法，步驟3具體為：根據(jù)如下公式計算歐式距離dl：
其中ωb為權(quán)重，b為特征值的數(shù)量，3≤b≤9且b為奇數(shù)，ω(b+1)/2為周期三位置。
上述的基于周期三臨近特征的外顯子預測方法，步驟3中所述的b＝7，權(quán)重ω1＝ω2＝ω3＝ω5＝ω6＝ω7＝0.1，ω4＝1-6×0.1＝0.4。
上述的基于周期三臨近特征的外顯子預測方法，所述的步驟4，判定的具體方法為，若d1＜d2，則該位置為外顯子位置，若d1＞d2，則該位置為內(nèi)因子位置。
本發(fā)明具有如下有益效果。
真核基因蛋白質(zhì)編碼序列中的外顯子通過剪切翻譯成具有生命現(xiàn)象的蛋白質(zhì)，但是由于周期三預測法具有特征的不確定性，因此本發(fā)明首次提出了一種具有通用性的周期三臨近特征外顯子預測方法，提取了周期三臨近功率譜特征值作為有效指標。
本發(fā)明則是提出了一種具有通用性的周期三臨近特征外顯子預測方法，該方法針對 周期三行為出現(xiàn)位置不確定性，提取了周期三臨近功率譜特征值作為有效指標，因此該算法更具有通用性，且預測精度更高。本發(fā)明時首次提出在DNA外顯子的識別中，采用對1/3頻率的臨近區(qū)間進行距離判斷的方法，將“臨近區(qū)間”用于判斷DNA外顯子位置。
具體實施方式
為進一步了解本專利的內(nèi)容，下面結(jié)合實施例對本專利作詳細描述。
實施例1
本實施例為采用具體方法的基于周期三臨近特征的外顯子預測方法，具體步驟如下。
步驟1、訓練階段：
選取已知的DNA序列，記錄長度區(qū)間，建立DNA堿基序列的功率譜：
對DNA序列中每種堿基，分別將其位置記錄為序列，得到4個堿基A、T、G、C的位置的序列和n為時間的離散值，式中，l表示DNA序列中外顯子、內(nèi)因子和外顯子內(nèi)因子混合的標記，記l＝1表示外顯子序列，記l＝2表示內(nèi)因子序列，舍棄內(nèi)因子和外顯子混合序列，將每一個序列按照如下公式進行快速傅里葉變換得到每個序列傅里葉變換值和
式中，j為虛數(shù)單位，n為時間的離散值，k為頻域的離散值，N為周期長度。
分別選取3倍數(shù)的外顯子和內(nèi)因子序列計算模板功率譜特征值：
采用如下公式，選取3倍數(shù)的外顯子和內(nèi)因子序列計算模板功率譜特征值，記為Sl[q]，其中q為3倍頻域的離散值
3倍頻域的離散值q＝…N/3-2,N/3-1,N/3,N/3+1,N/3+2…。
步驟2、識別階段：
取待識別DNA序列，滑動與訓練階段同樣的長度區(qū)間，計算得到每組識別功率譜特征值：
取待識別DNA序列，滑動與訓練階段同樣的長度區(qū)間，對待識別DNA序列中每種堿基，分別將其位置記錄為序列，得到4個堿基A、T、G、C的位置的待識別序列ya(n)，yt(n)，yg(n)和yc(n)，將每一個待識別序列按照如下公式進行快速傅里葉變換得到每個待識別序列傅里葉變換值Ya[k]、Yt[k]、Yg[k]和Yc[k]：
式中，j為虛數(shù)單位，n為時間的離散值，k為頻域的離散值，N為周期長度，
采用如下公式，計算得到每組識別功率譜特征值，記為Zp[k]；
Zp[k]＝|Ya[k]|2+|Yt[k]|2+|Yg[k]|2+|Yc[k]|2
其中p為滑動數(shù)，p＝1,2,…,N-M+1，M為滑動窗大小，N＞＞M。
步驟3、計算訓練階段模板的功率譜特征值與識別階段的每組識別功率譜特征值的歐式距離，對1/3頻率的臨近區(qū)間進行距離判斷，具體為根據(jù)如下公式計算歐式距離dl：
其中ωb為權(quán)重，b為特征值的數(shù)量，3≤b≤9且b為奇數(shù)，ω(b+1)/2為周期三位置。
可優(yōu)選b＝7，權(quán)重ω1＝ω2＝ω3＝ω5＝ω6＝ω7＝0.1，ω4＝1-6×0.1＝0.4。
步驟4、根據(jù)模板功率譜特征值與識別功率譜特征值的歐式距離，判斷外顯子位置，判定的具體方法為，若d1＜d2，則該位置為外顯子位置，若d1＞d2，則該位置為內(nèi)因子位置，記判斷結(jié)果為Rp，且
本實施例的方法中，q＝N/3，即為周期三預測法。采用臨近特征法采用1/3頻率處臨近的多特征值，以多維特征向量結(jié)合權(quán)重進行歐式距離的判斷和預測，能提取更多的 有效信息，因此該方法更具有通用性，周期三法僅為此通用算法的特例。
采用本實施例方法實際使用，并進行評估，以敏感度和特異度的平均值作為預測精度的衡量指標：
敏感度：Sn＝TP/(TP+FN)，特異度：Sp＝TN/(TN+FP)，精度：Ac＝(Sn+Sp)/2。
式中，TP表示被正確判為外顯子的個數(shù)，TN表示被正確判為內(nèi)含子的個數(shù)，F(xiàn)N表示被錯誤地判為內(nèi)含子的個數(shù)，F(xiàn)P表示被錯誤地判為外顯子的個數(shù)。
對DNA序列F56F11.4(編號：AF099922)五個外顯子的預測效果。長度為8000bp的序列F56F11.4中五個外顯子的位置分布分別為：928-1039、2528-2857、4114-4377、5465-5644和7255-7605。
為提取有效的多維基因序列特征，分別采用平均特征值法和本實施例的方法進行實驗。
平均特征值法，在幀長360時，分別等間隔選取不同長度的外顯子和內(nèi)因子序列功率譜中特征值作為模板特征，將實驗序列通過滑動窗口選取有限序列并提取特征值與模板特征作歐式距離比對，并計算序列位置的預測精度。
采用本發(fā)明的方法，幀長360時，分別等間隔選取不同長度的外顯子和內(nèi)因子序列功率譜中周期三臨近的特征值作為模板特征，將實驗序列通過滑動窗口選取有限序列并提取特征值與模板特征作權(quán)值平均的歐式距離比對，并計算序列位置的預測精度。
實驗結(jié)果表明，特征值數(shù)增加時，“平均特征值”法精度較低且均在0.6以下，本實施例方法的精度均達到0.84以上。周期三臨近特征值法充分利用了外顯子和內(nèi)含子周期三臨近特征之間的差異，在相同幀長時其比周期三法的0.8261的精度高，遠高于平均特征值法。
本實施例共選取了NCBI GenBank數(shù)據(jù)庫3組基因序列作為測試集，其中還包括長度為1667bp面包麥的DNA序列(編號：AB166873)和長度為2275bp山羊的DNA序列(編號：K01714)。表格中間列為傳統(tǒng)的周期三預測法預測的精度，右側(cè)列為本實施例方法的周期三臨近特征預測法，取特征值b＝7，ω1＝ω2＝ω3＝ω5＝ω6＝ω7＝0.1，ω4＝0.4。
表中可看出本實施例的算法敏感度較周期三法有大幅度提升，而敏感度主要表征外顯子正確判斷的個數(shù)，特異度稍稍下降，總體預測精度較高。
表1兩種方法外顯子區(qū)間預測精度比較