計算生物學講座

1、計算生物學講座李伍舉軍事醫(yī)學科學院基礎醫(yī)學研究所 2022/8/30基本概念主要類容 數(shù)據(jù)庫 序列比較 RNA二級結構預測 外源基因高效表達數(shù)學模型 蛋白質結構預測 進化樹構建 基因預測：coding region, noncoding RNA genes 基因表達譜分析 計算工具：數(shù)據(jù)挖掘、模式識別等 BioSun軟件介紹What is computational biology Computational biology is a new field of research which develops models and software implementations for com

2、putational problems in molecular biology, biotechnology, and genetics. Such means are needed both in the basic research and in the industrial applications of biotechnology. Computer-aided DNA sequencing, sequence comparison, prediction of protein structures, docking of molecules, and the interpretat

3、ion of electron microscopy or NMR data are examples of typical computationally intensive tasks. Bioinformatics refers to the development and use of the (molecular) biological databases 計算生物學定義 計算生物學是一門概念性學科，以生物信息為基礎，以計算為工具，解決生物學問題。 與生物信息學的定義類似，只是側重點有所不同。計算生物學側重于計算與問題，通過計算解決問題；生物信息學側重于數(shù)據(jù)的管理與數(shù)據(jù)庫的構建。組成

4、部分表現(xiàn)形式創(chuàng)新點相關背景數(shù)據(jù)數(shù)據(jù)庫數(shù)據(jù)數(shù)據(jù)庫管理系統(tǒng)算法程序高效算法 數(shù)學、物理、計算機解決問題理論問題提出新的問題生物學計算生物學/生物信息學定義解析數(shù)據(jù)庫 目前，各式各樣的生物數(shù)據(jù)庫很多，比較著名的有美國NCBI提供的Genbank，歐洲生物信息學研究所的EBI，日本國家遺傳學研究所的DDBJ等。 每年的Nucleic Acid Res.數(shù)據(jù)庫專刊了解目前的生物信息資源為計算生物學研究提供基礎序列比較計算分子生物學中基本技術、非常重要多方面應用：蛋白質結構預測、RNA二級結構預測、進化樹構建、進化譜構建、序列模式構 建、基因功能預測和基因預測等。比較形式多樣：兩個序列比較、多個序列比較；

5、全局比較、局部比較；最優(yōu)比較、次優(yōu)比 較；掃描數(shù)據(jù)庫等。流行的比較：blast，參見NCBI,EBI網(wǎng)頁目前的難點：多序列的最優(yōu)比較序列比較中的基本參數(shù)空格罰分：單獨空格、延伸空格 Wx=g+rx基本的相似性矩陣：殘基與殘基之間的 相似性，堿基與堿基之間相似性。最簡單的就是單位矩陣。Needleman-Wunsch算法ATTGCTTG0001000T0110011T0110011C0000100T0110011T0110011Seq1=ATTGCTT Seq2=GTTCTTATTGCTTG5434210T4543221T3443221C20222310T1221121T0110011ATTGC

6、TT : :GTT CTTRNA二級結構預測為什么要研究RNA二級結構預測tRNA的轉運氨基酸功能核酶的催化功能5和3端的非翻譯區(qū)(UTR) 的結構對基因表達的影響RNA二級結構中的基本結構螺旋區(qū)(helical regions):穩(wěn)定發(fā)卡環(huán)(hairpin loop):不穩(wěn)定內部環(huán)(interior loop):不穩(wěn)定膨脹圈(bulge loop):不穩(wěn)定多分支環(huán)(multibranch loop):不穩(wěn)定評價指標:自由能RNA二級結構自由能計算RNA二級結構自由能=各個基本結構的自由能之和基本結構的自由能由自由能計分系統(tǒng)確定具體計算方法見： 軍事醫(yī)學科學院院刊，1995，19：293自由能

7、計分系統(tǒng)Salser自由能系統(tǒng)Turner自由能系統(tǒng)http:/zukerm/rna/堿基對之間的堆積能+環(huán)區(qū)的不穩(wěn)定能自由能計算 5-AC-3 -1.8 5-CG-3 -3.4 3-UG-5 3-GC-5 5-GU-3 -1.8 H(4,10)=4.4 3-CA-5 TotalG=-1.8-3.4-1.8+4.4=-2.6 (Kcal/Mol)RNA二級結構預測方法分類Dynamic programming algorithmsKinetic folding algorithmsGenetic algorithmsComparative methodsHiggs. Quarterly Rev

8、iews of Biophysics 33:199-253,2000Dynamic programming algorithms 基于自由能計分系統(tǒng)，以最低自由能為目標函數(shù)，首先采用迭代方法求出所有可能子片斷的二級結構自由能，再用回歸方法求出RNA序列的最低自由能結構。Maximum matching model (Nussinov et al. 1980)Minimum free energy (Zuker, NAR 1981)Zukers Minimum free energy V(i,j)表示i 與j配對時的子片斷最低自由能;W(i,j)表示不論i 與j配對與否時的子片斷最低自由能；基于

9、W，尋找多分枝環(huán)，基于V，尋找分枝環(huán)的結構。Kinetic folding algorithms 鑒于RNA的二級結構最終是由若干個螺旋區(qū)組合而成，所以RNA二級結構的預測目標就是想辦法找出有關的螺旋區(qū)。動力學折疊算法的主要思想就是沿著自由能降低的折疊方向，以最低自由能為目標函數(shù)，去模擬RNA二級結構。 如果不考慮所有可能的折疊路徑，盡管動力學折疊算法以最低自由能為目標去模擬RNA二級結構，并不能保證最終獲得的結構為最低自由能結構。Genetic algorithms 遺傳算法，本質上是一種模仿生物進化的優(yōu)化算法，對一個復雜問題，如果沒有數(shù)學上嚴格的最優(yōu)解獲得辦法，可用遺傳算法來獲得比較優(yōu)化的

10、解，但并沒有嚴格的數(shù)學證明，保證獲得的解是最優(yōu)解，從這個意義上講，是Monte-Carlo模擬方法的特殊情況。 對RNA二級結構預測來說，由于已有動態(tài)編程算法來求RNA的最低自由能結構，所以，可視遺傳算法為動力學折疊算法的一種。Comparative methods 基本假定：來自不同物種的RNA序列如Phe-tRNA序列，如果它們擁有相同的功能，那么這些序列就應該擁有相同或類似的結構。因此，必須有多個相關的序列，方可使用此方法來研究RNA二級結構預測。收集來自不同物種的具有相同功能的RNA序列多序列比較識別所有可能的螺旋區(qū)，并按一定的計分方法排序RNA二級結構的組裝Juan(JMB,1999

11、,289:935只要有5個序列即可。排序的指標：自由能和保守性；Likelihood-ratio test.將排序后的螺旋區(qū)依次加入到RNA二級結構中去，直至形成一個穩(wěn)定結構比較方法預測RNA二級結構流程兩種最低自由能預測方法比較 Zuker的動態(tài)規(guī)劃算法與Pipas的螺旋區(qū)組合算法均可以用來求最低自由能結構，在多數(shù)情況下，結果是一致的，導致不一致的原因有：在Pipas算法中，通常只考慮全長的螺旋區(qū)，而對其子螺旋區(qū)不予考慮，在有些情況下，可能是其子螺旋區(qū)在最終的最低自由能結構中。在Zuker算法中，是以單點配對為基礎來求最低自由能結構，而在最終的最低自由能結構中，要將之過濾掉。哪一個結構最好

12、無論是改進的Zuker算法（MFOLD），還是螺旋區(qū)組合類算法，對一個特定的RNA序列來說，均可以求出多個RNA二級結構，如果不考慮實驗數(shù)據(jù)，究竟哪一個結構比較合理？目前，一般認為由比較方法求出的結構比較合理。但如果沒有相關的序列數(shù)據(jù)，如何尋找合理的結構？我們的工作螺旋區(qū)每個螺旋區(qū)H(S,E,L)用三個參數(shù)表示：S:螺旋區(qū)起點；E:螺旋區(qū)終點; L:螺旋區(qū)長度一級螺旋區(qū)概念設有一個RNA序列，長度為N,Hk（Sk,Ek,Lk）(k=1,2,n)為區(qū)間i,j（1ijN）上的n個螺旋區(qū)，并滿足下列關系：iS1E1S2E2SkEkSnEnj則稱這n個螺旋區(qū)為區(qū)間i,j上的一級螺旋區(qū)。螺旋區(qū)之間的關系

13、RNA二級結構預測 根據(jù)上述概念，提出了下列算法： 基于螺旋區(qū)隨機堆積的RNA二級結構預測 生物物理學報，12:213-218；1996 Prediction of RNA secondary tructure based on helical regions distribution Bioinformatics 14(8):700-706, 1998） RNA二級結構預測系統(tǒng)構建 生物化學與生物物 理進展，23:449-453；1996。中心思想 RNA分子在溶液中可以有多種結構與之對應，但每種結構出現(xiàn)的頻率不同，有的結構出現(xiàn)頻率較大，起主導作用。為了求出主導結構，以最低自由能為目標函數(shù)，

14、采用非決定性Monte-Carlo模擬方法獲得一定數(shù)目的二級結構，然后進行統(tǒng)計分析。 RNA二級結構預測的主要步驟給定RNA序列求出所有可能的螺旋區(qū)列表進行隨機堆積，獲得一定數(shù)目的RNA二級結構統(tǒng)計每個結構的出現(xiàn)頻率 出現(xiàn)頻率0.4，獲得主導結構出現(xiàn)頻率0.4，求每個螺旋區(qū)的出現(xiàn)頻率，然后反復迭代，最終獲得主導結構。Phe-tRNA的三葉草結構 Phe-tRNA的最低自由能結構 自由能：-20.8Kcal/Mol基于1161tRNA序列的RS, HD和MFold比較主要結論RS和HD:不容許螺旋區(qū)端點的GU配對RS方法優(yōu)于HD如果螺旋區(qū)端點的GU配對容許,T3.0;如果螺旋區(qū)端點的GU配對不容

15、許,T25螺旋區(qū)最佳長度為3 bpRS的預測精度為54.65% (523/957)，HD 為52.14 (499/957)，Zuker方法為32.92 (315/957) 。BJRNAFold程序 鑒于RS與HD方法較慢，為此，基于上述比較獲得的優(yōu)化參數(shù)，對Zuker方法進行了改進，新的程序命名為BJRNAFold對1139tRNA序列來說，BJRNAFold優(yōu)于RS方法。對其它家簇來說，BJRNAFold與MFold預測結果相當。外源基因高效表達數(shù)學模型構建 實現(xiàn)外源基因在原核系統(tǒng)或真核系統(tǒng)的高效表達具有重要的理論和實際意義。如一些重要的細胞因子，在體內含量甚微，單靠提取方法，難以獲得足量的

16、細胞因子，從而阻礙了結構與功能關系研究，更談不上臨床的大規(guī)模應用。 目前，有一些指導性原則可用于外源基因的高效表達設計，如利用表達系統(tǒng)的優(yōu)勢密碼子替換外源基因中的稀有密碼子或改變RNA二級結構以提高表達水平，但是，這些原則都是定性的，其次，這些因素是分別考慮的，沒有綜合考慮這些因素的共同作用。 當時，國內的情況是：國內科學家自己組建的原核高效表達載體pBV220，自構建以來，得到了廣泛應用，利用該載體已成功實現(xiàn)了多種細胞因子等外源基因的高效表達，但也有表達水平比較低的情況，給我們提供了豐富的數(shù)據(jù)資源。所以，我們以之為基礎，來定量研究外源基因表達水平與其他因素之間的關系。數(shù)據(jù)收集 Table 1

17、 Related data of 22 foreign genes carried by pBV220 vectorNo. Name Level Class D Ref. No. Name Level Class D Ref. 1 PCG12 1.96 1 9 3 12 PCG14 23.41 2 11 3 2 PCG11 5.66 1 15 3 13 HIFN 24 2 5 1 3 HIL4_CD 5-10 1 6 4 14 MIL4 25-30 2 5 9 4 HIV1 8 1 5 5 15 PCG18 27.20 2 8 3 5 RB 10-15 1 6 6 16 HIL6_T 28 2

18、 5 10 6 PCG16 11.32 1 11 3 17 PCG13 28.93 2 7 3 7 PCG17 11.92 1 6 3 18 HBV 30.41 2 8 11 8 HCV_NS3 14 1 8 7 19 PCG15 30.78 2 10 3 9 NAPIL8 18.5 1 5 8 20 HIL4 30-40 2 6 1210 HIL2 20 2 5 1 21 GMCSF_W 40 2 6 1311 GMCSF_Z 20 2 5 2 22 HIL6 71 2 5 14方法 RNA二級結構預測:采用基于螺旋區(qū)隨機堆積的RNA二級結構預測方法。密碼子偏性:采用CAI指標。判別分析:采

19、用Bayes判別分析方法。RNA二級結構與高(低)表達關系 原核啟動子含有兩個保守區(qū)域：-10與-35區(qū),當然也包含了Gold所證實的富含核糖體結合位點信息的-2521這個區(qū)域,那究竟以哪個區(qū)域為標準進行計算呢 ?另外, 3端二級結構又用哪個區(qū)域進行計算呢? 為了確定與高(低)表達具有顯著性統(tǒng)計學意義的區(qū)域，我們圍繞起始密碼子AUG和終止密碼子TAA附近，隨機選取了多個片段，并運用基于螺旋區(qū)隨機堆積的RNA二級結構預測方法，分析了每個片段的二級結構，然后運用判別分析方法考察每個片段二級結構自由能與高 (低)表達之關系，結果發(fā)現(xiàn)三個區(qū)間 21,89,23,95和28,93的自由能與高(低)表達具

20、有顯著的統(tǒng)計學意義,結果見下表。 表: 區(qū)間選取與高(低)表達之間關系 區(qū)間 F值 概率 判別符合率 21,89 9.7195 0.0012 90.5% 23,95 10.4115 0.0009 86.4% 28,93 7.2998 0.0044 73.5% 5和3端聯(lián)合判別函數(shù)其中以區(qū)間21,89判別效果最好,由此得到5和3端的聯(lián)合判別函數(shù)為:LES=-10.8036-0.4732*G5-1.8649*G3 (1)HES=-17.1970+0.1559*G5-2.6214*G3 (2)E=HES-LES上式中G5表示5端-3039區(qū)域(即區(qū)間21,89)的二級結構自由能,G3表示3端30 -

21、39區(qū)域的二級結構自由能,判別準則是:當E0時,判斷樣品為低表達;當E0時,判斷樣品為高表達。5端判別函數(shù)5端判別函數(shù):LES5=-2.4489-0.9339*G5 (3)HES5=-0.6790-0.4918*G5 (4)統(tǒng)計量F=5.8383,P0.02540.05,判別符合率為:59.1%3端判別函數(shù):LES3=-0.0535*G32-1.0971*G3- 6.7433 (5)HES3=-0.2130*G32-5.6572*G3-37.9814 (6)統(tǒng)計量F=9.4318, P0.006,判別符合率為68.2%。3端判別函數(shù)5和3端條件要實現(xiàn)外源基因的高效表達,根據(jù)判別準則, 下式必須

22、成立;LES5HES5 LES3HES3由上面的兩個不等式可得到G5,G3的取值范圍為:G5-4.0 (Kcal/mol) (7)-17.21G3-11.38 (kcal/mol) (8)SD序列、ATG與TAA在二級結構中位置與高(低)表達關系 表：SD序列、ATG與TAA在二級結構中參與配對的堿基數(shù)目 Class AUG SD TAA Class AUG SD TAA 1 3 0 0 2 1 1 0 1 0 3 0 2 0 0 0 1 1 3 3 2 1 0 0 1 0 3 0 2 3 0 0 1 0 0 3 2 0 0 0 1 0 3 0 2 3 0 0 1 2 0 0 2 0 0 0

23、1 3 0 0 2 1 0 0 1 0 0 0 2 1 3 0 2 0 2 3 2 1 0 0 2 0 4 2 2 3 3 0 通過SAS軟件分析,AUG,SD與TAA 在二級結構中參與配對的堿基數(shù)目與高(低)表達無顯著關系。不過，由于所選數(shù)據(jù)均是表達數(shù)據(jù)，通過分析上表可以看出:AUG, SD與TAA最好為0。局部密碼子偏性與高(低)表達關系 為了探討5與3端局部密碼子偏性與高(低)表達關系,利用Goldkey軟件計算下列值:X1,X2,.，X30；Y1，Y2，.，Y30。并且運用判別分析方法考察Xi與Yj的所有可能900種組合與高(低) 表達的關系,結果發(fā)現(xiàn):有8種組合,P值較小,詳見下表。

24、 表：5和3端局部密碼子偏性聯(lián)合作用與高(低)表達關系 5偏性 3偏性 F值 概率 判別符合率 X1 Y3 5.4572 0.0362 73.3% X3 Y3 3.0462 0.0851 73.3% X4 Y3 3.0686 0.0839 73.3% X9 Y3 3.2001 0.0769 73.3% X10 Y2 3.1687 0.0785 80.0% X10 Y3 3.4554 0.0653 66.7% X11 Y3 3.0782 0.0834 66.7% X12 Y3 2.9887 0.0885 66.7% 由上表可以看出: P值均在0.05左右擺動,以Y3的P值最小(因X1=1), 并

25、得到下列判別方程:LEC3=-0.44192+30.73696*Y3 (9)HEC3=-2.35737+70.99093*Y3 (10)P0.0362,由判別準則可得到外源基因高效表達條件是: Y 3 0.04758,即要求外源基因3端的3個密碼子(包括TAA)必須是大腸桿菌的優(yōu)勢密碼。SD序列與起始密碼子ATG之間堿基數(shù)(D)與高(低)表達關系由于外源基因與載體連接后,距離D就是定值,因此,只能根據(jù)表 1中的數(shù)值作判別分析,得到的判別函數(shù)為: LED=-4.4464+1.1273*D (11) HED=-3.1267+0.9453*D (12)統(tǒng)計量F=1.2325,P0.2801,判別符合

26、率為59.1%,因此,外源基因的高(低)表達與D值無顯著關系。由表1可以看出,D值在83范圍內較為合適。一般判別函數(shù)構建與回顧性分析通過綜合分析，最終得出X6,G5與G3聯(lián)合作用的判別函數(shù):LESC=-13.6401+12.8546*X6-0.3612*G5-2.0377*G3 (13)HESC=-21.8205+16.4293*X6+0.2990*G5-2.8423*G3 (14)統(tǒng)計量F=6.8386，P0.0029，下表為22個外源基因利用方程(13)和(14)的判別情況。從表中可以看出:22個外源基因中只有第1個被誤判。 通過分析有關數(shù)據(jù)可知：將第1個樣品誤判的原因是5端自由能太高了。

27、 表：22個外源基因回顧性分析 No. LESC HESC DC OC No. LESC HESC DC OC 1 18.4134 19.8864 2 1 12 19.1358 19.2885 2 2 2 19.5331 18.9596 1 1 13 28.6516 36.1672 2 2 3 6.6124 2.3394 1 1 14 14.4364 15.8918 2 2 4 12.1912 10.8165 1 1 15 18.9552 19.4380 2 2 5 4.4814 0.8667 1 1 16 13.8560 14.7747 2 2 6 19.6776 18.8400 1 1 1

28、7 18.5578 19.7668 2 2 7 18.1501 15.0622 1 1 18 22.2221 26.1855 2 2 8 18.4046 14.9546 1 1 19 18.9552 19.4380 2 2 9 5.2970 1.7473 1 1 20 19.4498 20.2457 2 2 10 18.3370 22.4037 2 2 21 18.6168 22.9754 2 2 11 26.5447 27.6698 2 2 22 18.0869 19.4216 2 2一般判別函數(shù)驗證表：4個外源基因判別結果 No. Name G5 G3 X6 Level OC DC 1 H

29、CV_L -4.8 -11.4 0.22918 10 1 1 2 HCV_W -4.8 -12.8 0.22918 11 1 1 3 JSNA 1.0 -14.0 0.58962 60 2 2 4 TNFA -6.5 -14.9 0.26419 30 2 2 有關論文 李伍舉，吳加金：pBV220載體中外源基因表達水平 定量分析 病毒學報，13:126-133；1997。李伍舉，吳加金：pBV220載體中外源基因二級結構 與表達水平關系 生物技術通訊,7:149-151,1996。Li Wuju, Wu Jiajin: GeneDn: for high-level expression des

30、ign of heterologous genes in a prokaryotic system Bioinformatics 14:884-885,1998。 模型的實驗驗證 進行Ricin-A鏈在E.coli 中的高效表達設 計，結果獲得了高效表達（表達水平20%） 裴吳紅 沈倍奮 李伍舉等：細胞與分子免 疫學雜志，1998，14（1）：33人FKBP12在E.coli中的高效表達設計，結 果獲得了高效表達（表達水平20%） 裴武紅 胡美茹 李伍舉等：中國生物化學 與分子生物學報， 2000，16（3）：322人SCF在E.coli 中的高效表達設計，結果獲得了高效表達（表達水平30%）

31、 洪海燕等，待發(fā)表。宋曉國等構建的高效原核融合表達載體pBVIL1， 在保持5和3端結構的情況下，改變編碼區(qū)內部 結構以實現(xiàn)不同基因的融合，按照我們的數(shù)學模 型，這些基因的表達水平均應在20%以上，結果 得到了實驗驗證。細胞與分子免疫學雜志， 2001,17(3):231與大連醫(yī)科大學合作，進行人NMDA受體靶片斷在E.coli中的高效表達設計，結果獲得了高效表達，表達水平從未改造前的6%上升到29%。軍事醫(yī)學科學院院刊，2002，第三期。其它情況：有兩例設計與預期不符，一例是關于人Insulin，一例是關于NATO基因，它們的表達水平均在15%左右，而沒有達到預期的20%以上，推測其原因，有

32、一些因素我們可能沒有考慮到：如表達序列的長短，被表達蛋白質的穩(wěn)定性，編碼區(qū)內部的稀有密碼子串聯(lián)作用等。其它多例實驗數(shù)據(jù)的驗證：表達水平均較低，通過我們分析，均找出了原因。 思考：表達模型對原核中的其它載體？對真核(Yeast system)中的載體？我們的目標：對某一個特定的外源基因， 能夠在實驗之前,就能從理論上知道該基因的表達水平。分子生物學中的其它實驗能否達到計算 機模擬?蛋白質結構預測二級結構預測： 神經(jīng)網(wǎng)絡方法、Chou-Fasman、Garnier等，目前最好的方法是Garnier方法、預測精 度在78%左右，不到80%。三級結構預測： 從頭預測，同源模建等應用：抗原表位分析、小分

33、子藥物設計等方面http:/進化樹構建收集一組相關的序列多序列比較，計算序列之間相似性基于相似性矩陣，可分別采用下列方法構建Maximum parsimony method (Minimum evolution method)Distance methods (Neighbor-joining methods)Maximum likelihood approachhttp:/phylip.html基于全基因組的物種進化樹構建從統(tǒng)計學上講：對物種進化樹的真實情況是無知的，因此，只能通過樣本情況推測總體情況：基于不同基因家簇給出的進化樹進行疊加基于全基因含量物種進化樹構建基于若干保守的蛋白家簇聯(lián)合

34、基因預測編碼區(qū)預測: GeneScan，Geneie，nc RNA 基因預測 ncRNA基因預測的專用方法 ncRNA基因預測的通用方法ncRNA基因預測的專用方法 主要是為識別某個特定的ncRNA基因家簇的新成員而設計的一類方法，一般采用啟發(fā)式算法，即根據(jù)特定的ncRNA基因家簇的一級結構和二級結構特征，發(fā)現(xiàn)一些規(guī)則，然后根據(jù)這些規(guī)則掃描基因組序列，并從基因組序列中發(fā)現(xiàn)符合這些規(guī)則的基因片段，這些基因片段即為該ncRNA基因家簇的可能新成員，如能通過比較基因組學在相近物種的基因組中發(fā)現(xiàn)類似的保守片段，便可進一步確證這些新成員的身份，當然，最終要通過實驗來驗證。 一級結構和二級結構特征分析 為

35、了尋找某一特定的ncRNA基因家簇的一級結構和二級結構特征，通常對該家簇的所有已知成員（序列）或部分已知成員（序列）進行多序列比較，然后，基于多序列比較形式提取該ncRNA基因家簇的一級結構和二級結構信息，通常有三種方法：加權矩陣即通常的頻數(shù)矩陣、模式分析和Eddy提出的用于RNA序列分析的共變化模型。ncRNA基因預測的通用方法 基本的出發(fā)點是希望能夠找出一種通用的預測方法，不依賴于某個特定的ncRNA基因家簇信息，將ncRNA基因從基因組序列中識別出來。目前已發(fā)展了下列四種方法：堿基組成方法、神經(jīng)網(wǎng)絡方法、比較基因組學方法和轉錄起始位點與終止位點預測法。堿基組成方法 基本設想是將一個基因組

36、中ncRNA基因序列當作信號，基因組序列當作背景，然后利用ncRNA基因的堿基組成與基因組中的堿基組成的差別將ncRNA基因識別出來，這種差別越大，越有利于ncRNA基因的識別。 神經(jīng)網(wǎng)絡方法 基本思想是將整個基因組序列分為三個部分：編碼蛋白質的基因部分、ncRNA的基因部分（第一類）和其它的非編碼基因間區(qū)（第二類），并進一步假設第二類中只有小部分含有ncRNA基因，然后以大腸桿菌為例，運用神經(jīng)網(wǎng)絡方法探討了第一類與第二類的區(qū)分問題。鑒于第一類的長度要遠遠小于第二類，于是從第二類中隨機選出一個與第一類大小相當?shù)牟糠謥磉M行訓練，然后，運用獲得的神經(jīng)網(wǎng)絡預測第二類中的其余部分，從而識別出第二類中的

37、含有ncRNA基因的序列片段，為了預測結果可靠，將上述過程多次重復進行。 比較基因組學方法 基本假設是ncRNA基因在相近物種的基因組中，不僅一級結構有一定的保守性，更重要的是其二級結構也非常保守，根據(jù)這個設想，Rivas和Eddy提出了一個識別ncRNA基因的自動方法，其相應程序為QRNA。該方法的核心部分是提出了三個概率模型，它們分別是編碼區(qū)模型COD、ncRNA基因模型RNA和零假設模型OTH，然后基于序列的比較形式，采用這三種模型分別計算 Bayes后驗概率，并根據(jù)數(shù)值的大小判定被比較序列片段是編碼區(qū)或ncRNA基因或其它的序列。 轉錄起始位點與終止位點預測法 在已經(jīng)注釋的基因組序列基

38、礎上，通過預測轉錄起始位點與終止位點來識別ncRNA基因的一種方法，主要包含四個步驟。首先，從待分析的基因組中找出“空白”的基因間區(qū)，在這空白的基因間區(qū)中，不包含任何已注釋的基因（ORF，tRNA，rRNA等）；其次，進行轉錄起始位點預測；再次，進行轉錄終止位點預測；最后，采用序列比較方法進行保守性分析,找出可能的ncRNA基因。但對真核生物來說,轉錄起始位點的預測并非易事。小結 上述四種方法均可用于ncRNA基因的識別問題,但是預測結果的可靠性均不及蛋白質編碼區(qū)預測算法。通過對大腸桿菌基因組序列測試表明，從精度與被證實的ncRNA基因的數(shù)目來說，以比較基因組學方法最好。 基于基因表達譜的生物

39、信息學計算機輔助寡核苷酸微陣列探針設計基于基因表達譜的分類系統(tǒng)構建基于基因表達譜的分型系統(tǒng)構建基于基因表達譜的調控網(wǎng)絡構建綜合性數(shù)據(jù)庫的構建等李伍舉：基因表達譜的生物信息學軍事醫(yī)學科學院院刊,26:73, 2002為什么要研究基于基因表達譜的樣本分型與分類問題人類基因組測序模式生物基因組測序。Sequence Database 這些序列的功能是什么？基因芯片技術蛋白芯片技術其它高通量技術？基因芯片技術Gene 1Gene 2Gene 3Gene mTissues基因表達譜 基于基因表達譜，可以同時考察在特 定生理或病理過程中細胞內基因群的動態(tài)表達水平，從而將基因的活動狀態(tài)比較完整地展現(xiàn)出來，使

40、研究人員能夠在基因組水平上以系統(tǒng)的、全局的觀念去研究生命現(xiàn)象及其本質。但是，從數(shù)據(jù)處理的角度來看，通過基因芯片實驗直接獲得的是一個基因表達譜，相關的數(shù)據(jù)分析和挖掘已經(jīng)無法僅僅通過簡單的計算來進行，輔助數(shù)據(jù)分析系統(tǒng)的建立已經(jīng)成為基因芯片技術的必要工具。 基于基因表達譜的樣本分型研究聚類分析：譜系聚類、K-平均值變量選擇：標準差、比值法樣本分型研究樣本分型方案的評價聚類分析中的變量選擇 Lukashin（2001）采用下列規(guī)則來對酵母細胞的周期數(shù)據(jù)進行聚類：要求所有17個時間點上的表達數(shù)據(jù)的絕對值大于或等于100，并且表達水平至少有2.5倍的變化 Welsh（2001）采用的基因變量選擇方法是要求

41、每個基因表達譜的標準差大于或等于250，最終選出1243個基因對上皮細胞卵巢癌表達數(shù)據(jù)進行聚類分析 Perou(1999)則采用R/G比值方法，要求在26個上皮或乳癌樣本中，至少有三個樣本的表達水平在3以上，最終選出1247個基因來對樣本進行聚類。然而，上述基因變量的選擇方法并不具有普遍適用的特點 為了較好地解決聚類分析中的變量選擇問題，Xing（2001）提出了一個通用的用于聚類分析的基因變量選擇方法CLIFF，其主要思想為：在假定樣本表型未知情況下，采用某種聚類方法獲得初步的樣本分型方案，然后以之為基礎，采用監(jiān)督學習方法，選擇一定數(shù)目的基因來對樣本進行進一步的分型，該過程反復迭代，直到最后

42、樣本的分型方案沒有變化為止。在Xing方案的監(jiān)督學習中，基因變量的個數(shù)是靠經(jīng)驗來確定的，最終的樣本分型方案與選中的基因變量個數(shù)有關，其次，并沒有從理論上給出最佳的樣本分型個數(shù)。因此，到目前為止，基于基因表達譜的樣本分型問題仍沒有得到很好解決。 我們在詳細研究他人算法的基礎上，嘗試了多種變量選擇方法后發(fā)現(xiàn)，以變異系數(shù)和t-檢驗為基礎進行基因變量選擇是一個行之有效的方法，通過這種基因變量選擇方法與譜系聚類、K-平均值方法和自組圖方法的整合，構建了基于基因表達譜的樣本分型的整合系統(tǒng)SamCluster，較好地解決了聚類分析中的變量選擇問題?；诨虮磉_譜的樣本分類研究分類方法：Fisher、距離判別

43、變量選擇：逐步優(yōu)化、t-檢驗樣本分類研究樣本分類方案的評價主要應用：基于基因表達譜的腫瘤分類研究。 目前，樣本類型預測研究主要應用于腫瘤分類，利用基因表達譜在腫瘤與正常組織中的差異對腫瘤進行分類與診斷已形成共識。當前的腫瘤分類技術高度依賴于病理學工作者對腫瘤組織的主觀判斷，而基于基因芯片技術，即使一些組織沒有顯著變化，利用基因表達譜也可以對之做出早期診斷；另外，特別重要的一點是可以根據(jù)基因表達譜的變化來區(qū)分形態(tài)學上相似的腫瘤，這樣對腫瘤類型的精確識別有助于制定配套的最佳方案，從而達到增加療效、降低毒性的目的；另外，基于基因表達譜的腫瘤分類對腫瘤發(fā)生機制的理解以及征服這些腫瘤提供了重要思路。 究

44、竟采用多少個基因以及如何選擇這些基因，特別是對多類腫瘤情況如何選擇這些基因，到目前為止，還沒有定論；不過，采用最少的基因達到最高的預測精度將是追求的目標。為此，以我們提出的分類穩(wěn)定性概念為基礎，探討了Fisher線性判別分析方法與Monte-carlo模擬和逐步優(yōu)化等基因變量選擇方法的整合，以便對某個特定的基因表達譜，尋找最好的基因集合，這將對腫瘤發(fā)生機制與藥物作用機制的理解具有一定意義，并且，所找出的最佳基因或基因集合可用于腫瘤診斷與藥物靶基因識別等方面，為進一步的分子生物學實驗提供思路。第一部分基于基因表達譜的樣本分型系統(tǒng)構建 數(shù)據(jù)與方法四個數(shù)據(jù)集說明 為了說明樣本分型系統(tǒng)Samclust

45、er的性能，并與其它分型方法進行比較 ，主要運用四個數(shù)據(jù)集：COLON、 OVARIAN 、LEUKEMIA72、LEUKEMIA38 。1、COLON：該基因表達譜包含2000個基 因和62例樣本，其中有22例正常組織， 其余40例樣本為結腸癌組織。2、LEUKEMIA72: 該基因表達譜包含6817 個基因和72例樣本，其中有47例急性淋 巴細胞白血?。ˋLL），其余25例為急 性粒細胞白血?。ˋML）。3、LEUKEMIA38: 該基因表達譜包含6817個基因和38例樣本，在38例樣本中，有 27例急性淋巴細胞白血?。ˋLL）和11 例急性粒細胞白血病（AML）, 此數(shù)據(jù)由LEUKEMI

46、A72提取而 來。4、OVARIAN: 該基因表達譜包含7129個基因和 36例樣本，在36例樣本中，有27例卵巢癌組 織、5例正常組織和4例惡性的卵巢癌細胞系。變異系數(shù)計算（CV） 設G=gij(mn) 表示基因表達矩陣，其中m和n分別表示基因和樣本的個數(shù)，gij表示第i個基因在第j個樣本中的表達水平，對第i個基因來說，用gi=gi1, gi2, ,gin表示該基因在n個樣本中的表達向量。 兩點特別考慮1、如果在一個基因表達矩陣中，由于標準 化或對數(shù)化等原因導致某個基因的最小 表達水平(gmin)為負值，則要對該基因 的表達向量進行調整，調整的方式為該 基因的表達向量加上gmin的絕對值，調

47、 整后的基因向量的各個分量均為非負 值，從而便于均值與標準誤的計算。2、由于實驗誤差等原因，常常導致某個基因在 某個樣本中的表達水平呈現(xiàn)異常（偏大或偏 ?。瑸榱丝朔@種異常對變異系數(shù)的影 響，在計算某個基因的變異系數(shù)時，不考慮 該基因向量的最大與最小分量。 CV=cv1, cv2, , cvm 用MCV和SCV表示變異系數(shù)的均值與標準誤。 t-檢驗譜系聚類 聚類，就是物以類聚之意，為了刻畫樣本之間的相似性并進行聚類，必須對樣本之間的相似性或距離進行定量，然后，按特定的方式進行聚類。1、標準化變換 標準化變換，首先對基因變量進行標準化，然后對樣本變量進行標準化，標準化之后其均值為0，方差為1。

48、 2、構建相關系數(shù)矩陣 構建相關系數(shù)矩陣，就是采用標準的相關系數(shù)計算方法，計算n個樣本之間的兩兩相關系數(shù)，從而獲得相關系數(shù)矩陣。 3、譜系聚類 首先將n個樣本視作n個類，從相關系數(shù)矩陣中尋找最大值元素即最相似的兩個類，并將它們合并成一類，由此原來的n個類便轉化為n-1類，對新的n-1類繼續(xù)構建相關系數(shù)矩陣，此過程反復進行，直到最后，所有的n個樣本在1類為止。 類與類之間的距離（或相關系數(shù)）定義有多種方式，比較常用的有如下6種方式：最短距離法、最長距離法、中間距離法、重心法、類平均法和變差平方和法，這在一般的統(tǒng)計學課本上均有介紹。在構建Samcluster過程中，我們采用了類平均法。基于基因表達

49、譜的樣本類型發(fā)現(xiàn)的整合方案 兩個假設1、要求基因在不同樣本中的表達水平必須有波動。 變異系數(shù)是一個很好的指標，利用它可以刻畫不 同數(shù)據(jù)集中的各個基因的表達水平的波動情況。為了 選出用于樣本類型發(fā)現(xiàn)的基因變量，必須要求有關基 因在各個樣本中的表達水平有波動，且其變異系數(shù)大 于某個指定的閾值,如閾值太低，選出的基因變量集合 中將包含一定數(shù)目的噪聲基因變量（即不利于樣本類 型發(fā)現(xiàn)的基因變量），如閾值太高，選出的基因變量 集合中將缺少一定數(shù)目的信號基因變量（即有利于樣 本類型發(fā)現(xiàn)的基因變量）。為此，我們采用下列公式 來選擇變異系數(shù)的閾值。 CVth = MCV+CiSCV 2、要求基因的表達水平在不同

50、的樣本分型中（指推定的 樣本分型），其表達水平的差異具有一定的顯著性。 由于在我們的研究中，總是假定兩種可能的樣本 類型，因此，可用t-檢驗來刻畫。通過第一個假設， 我們可以獲得某個特定的基因變量集合，以之為基礎 進行聚類分析，可得到兩個推定的樣本類型,對此進行 t-檢驗分析，可找出表達水平呈現(xiàn)一定顯著性差異的 一些基因（可假定P=0.01、0.05或0.1等），從而將 表達水平不顯著的一些噪聲基因變量去除，在保留的 基因變量集合基礎上，進行進一步的聚類分析，此過 程反復迭代，直到最后沒有可剔除的基因變量為止， 最終獲得一個穩(wěn)定的樣本分型方案。一致樣本類型構建 在固定P值的情況下，多次改變標準

51、誤的系數(shù)因子Ci的值，對每一次改動，我們將獲得一個推定的樣本分型方案，并由此構建樣本之間的關系矩陣Snn。我們稱以此為基礎進行聚類分析所獲得的樣本分型方案為樣本的一致分型方案。為了獲得最佳的樣本分型方案和對應的基因集合，考慮了樣本的一致分型方案與各個推定的樣本分型方案之間的距離，并稱距離最小者為最佳的樣本分型方案。 變異系數(shù)計算基因表達矩陣，給定P值選定CVth聚類分析兩個推定的樣本類型t-檢驗獲得推定的樣本分型方案改變CVth，獲得多個樣本分型方案構建樣本關系矩陣一致樣本分型獲得最佳的樣本分型方案與對應的基因變量集合計算一致樣本分型與推定的樣本分型之間的距離程序設計 運用Matlab程序設計

52、語言，構建了基于基因表達譜的樣本類型發(fā)現(xiàn)系統(tǒng)Samcluster， 結果基于原始基因表達數(shù)據(jù)的樣本分型沒有進行基因變量選擇基于整合系統(tǒng)的樣本分型研究進行基因變量選擇表1：基于Samcluster系統(tǒng)的結腸癌基因表達譜COLON的樣本分型情況圖1：基于Samcluster系統(tǒng)的結腸癌基因表達譜COLON的樣本分型情況作者 方法 變量選擇方法 基于臨床分型 誤判數(shù) Alon 雙向聚類 t-檢驗 是 8Xiong 譜系聚類 Fisher方法 是 5Samcluster 譜系聚類 CV計算與t-檢驗 否 6Samcluster 譜系聚類 CV計算與t-檢驗 是 5表2：基于Samcluster系統(tǒng)的白

53、血病基因表達譜LEUKEMIA72的樣本分型情況圖2：基于Samcluster系統(tǒng)的白血病基因表達譜LEUKEMIA72的樣本分型情況系統(tǒng) 方法 變量選擇方法 基于臨床分型 誤判數(shù) CLIFF 標準化分割 Bayer誤差 是 3 聚類算法 Information gain ranking Markov blanket filteringCLIFF 標準化分割 Bayer誤差 否 ？ 聚類算法 Information gain ranking Markov blanket filteringSamcluster 譜系聚類 CV計算與t-檢驗 否 2Samcluster 譜系聚類 CV計算與t-檢

54、驗 是 1表3：基于Samcluster系統(tǒng)的白血病基因表達譜LEUKEMIA38的樣本分型情況圖3：基于Samcluster系統(tǒng)的白血病基因表達譜LEUKEMIA38的樣本分型情況 對LEUKEMIA38數(shù)據(jù)集來說，最好的樣本分型方案給出的分型精度94.7% (36/38)。如果結合樣本的臨床表型來選擇樣本的分型方案，則所有樣本都被正確分型。在Golub的結果中，他們使用自組圖方法SOM對LEUKEMIA38進行分型研究，在38個樣本分為兩個類型時，有4個樣本被誤判（1 AML，3 ALLs）。因此，對LEUKEMIA38數(shù)據(jù)集來說，Samcluster系統(tǒng)要優(yōu)于SOM。 從上表的D值分布來

55、看，當CVi位于區(qū)間0.5，1.0時，樣本分型往往達到最佳效果，因此，我們可以認為，當CVi值太高時，有許多有助于樣本分型的基因變量將被舍棄；當CVi值太低時，有許多不利于樣本分型的基因變量將被引進。而且，三個P值（0.01，0.05和0.1）均可以用來作為基因表達水平是否具有顯著性差異的檢驗指標，其中以0.05和0.1的效果較好。圖4：基于Samcluster系統(tǒng)的卵巢癌基因表達譜OVARIAN的樣本分型情況討論小結 通過這一部分工作，我們構建了基于基因表達譜的樣本類型自動識別系統(tǒng)Samcluster，其中心思想是將基于變異系數(shù)和t-檢驗的變量選擇方法整合到常用的譜系聚類算法中，并提出了一致

56、樣本類型的概念，在此基礎上獲得了基于基因表達譜的樣本的最佳分型方案。通過對4個基因表達數(shù)據(jù)的測試，結果表明，獲得的基因水平上的樣本分型與樣本的臨床分型具有較高的一致性。因此，這一工作對以基因芯片技術為基礎的相關研究（如中藥現(xiàn)代化等）具有較大意義。與CLIFF的基于基因表達譜的樣本分型系統(tǒng)比較比較項目 Samcluster CLIFF 結果聚類方法： 譜系聚類 標準化的分割聚類 非監(jiān)督學習 CV計算 Bayer error 不能確定最佳值監(jiān)督學習 t-檢驗 Information gain ranking 不能確定最佳值 Markov blanket filtering 不能確定最佳值基因表達分

57、布 無 混合Gaussian分布一致分型 有 無 參數(shù)敏感性 不敏感 敏感測試數(shù)據(jù)集 4個 1個LEUKEMIA72 1個誤判 3個誤判 第二部分基于基因表達譜的樣本分類系統(tǒng)構建 數(shù)據(jù)與方法1、基因表達數(shù)據(jù)集：COLON2、Fisher線性判別分析方法類內離差矩陣W和類間離差矩陣B M類樣本的基因表達矩陣計算第 l 類樣本和總體樣本的平均值 基因變量選擇方法 1、主成分分析2、t-檢驗3、鄰近關系分析中P值4、Monte Carlo模擬方法5、逐步優(yōu)化方法6、全局優(yōu)化方法分類精度的穩(wěn)定性分析 在基于基因表達譜的樣本分類研究中，用不同的基因變量選擇方法所選出的基因集合常常不一致，那究竟選擇哪一個

58、基因集合？即使用同一種方法，常常選出具有相同或相近分類精度的基因集合，那如何選擇較好的基因集合？另外，在分類研究中，常常將整個樣本按一定比例分為訓練組與試驗組，為此，我們不得不考慮樣本分配情況對分類精度的影響?；谶@些考慮，我們提出了分類穩(wěn)定性概念。 給定樣本的分配比例待評價的基因集合訓練組測試組訓練組分類精度分類規(guī)則測試組分類精度交叉有效性分類精度 分類精度均值A1、A2、An 考慮（A1+A2+An）/n的波動情況程序設計結果、基因變量個數(shù)與分類精度之間的關系主成分個數(shù)與分類精度關系基于t-檢驗的基因個數(shù)與分類精度之間關系基于P值的基因集合的分類精度基于Monte-Carlo模擬的基因集合

59、的分類精度基于逐步優(yōu)化的基因集合的分類精度多種基因變量選擇方法的比較分析基于主成分的樣本分類，盡管利用了基因表達譜的全部信息，但實際分類精度并不高；另外，為了對新樣本進行判別，必須要知道2000個基因的表達水平。因此，無論是實際應用，還是理論分析，此種方法價值不大。另外,通過上述結果可以看出：無論是采用T-檢驗，還是P值方法，隨著基因數(shù)目的增加，分類精度呈下降趨勢，而用較少的基因，其分類精度卻不是很高。因此，在基于基因表達譜的樣本(腫瘤)分類研究中，這些方法不是很好。Monte-Carlo模擬方法雖然運用較少的基因而達到較高的分類精度，但是，仍存在下列缺點：計算機運行時間較長，不能保證所得的基

60、因集合是最佳的，特別嚴重的是，每次計算所得的基因集合并不固定，因此，此方法也不宜使用。然而，逐步優(yōu)化方法與全組合方法卻能克服上述缺點，既保證運用較少的基因獲得較高的分類精度，又能保證每次運行獲得相同的基因集合，因此，值得推廣使用。在實際應用中，可以采用下列策略：如果基因表達譜中基因數(shù)目較少，采用全組合方法，尋求最佳基因集合；如基因數(shù)目很多，應用逐步優(yōu)化方法尋找最優(yōu)基因集合。小結通過這一部分工作，我們采用Fisher線性判別分析方法和多種變量選擇方法構建了基于基因表達譜的樣本分類系統(tǒng)Tclass, 并以結腸癌基因表達譜數(shù)據(jù)為例說明了Tclass系統(tǒng)的功能, 結果表明：運用較少的基因就可以達到較高