生物信學(xué)第三章序列比對-課件

生物信息學(xué)第三章序列比對為什么要序列比對？尋找進化過程中的同源序列;基于同源物鑒定的功能預(yù)測;基本假設(shè)：序列的保守性功能的保守性注意：1.蛋白質(zhì)一般在三級結(jié)構(gòu)的層面上執(zhí)行功能；2.蛋白質(zhì)序列的保守性決定于其編碼DNA的保守性；通常本章內(nèi)容提要第一節(jié)：數(shù)學(xué)基礎(chǔ)：概率及概率模型第二節(jié)：雙序列比對算法的介紹Dotmatrix動態(tài)規(guī)劃算法(Needleman-Wunsch,Smith-Waterman算法)

FASTA和BLAST算法第三節(jié)：打分矩陣及其含義第四節(jié)：多序列比對第一節(jié)序列比對的數(shù)學(xué)基礎(chǔ)排列組合從N個物品中取出k個物品的排列數(shù)：從N個物品中取出k個物品的組合數(shù)：概率模型概率模型:一個能夠通過不同的概率產(chǎn)生不同結(jié)果的模型。概率模型可以模擬或者仿真某一類型的所有事件，并且對每個事件賦予一個概率。色子模型：一個色子存在6個概率值：p1,p2,…,p6，其中，擲出i的概率為pi(i=1,2,…,6)。因此：pi≥0，且考慮三次連續(xù)的擲色子，結(jié)果為[1，6，3]，則總概率為：p1p6p3概率分布考慮連續(xù)變量x，例如：物體的重量。則當(dāng)重量確切為1公斤時的概率，為0。變量的區(qū)間：P(x0≤x≤x1)

當(dāng)區(qū)間無限小->0時，上式：P(x-δx/2

≤x≤x+δx/2

)=f(x)δx f(x)稱為概率密度函數(shù)因此：且二項分布1.事件只有兩種可能出現(xiàn)的結(jié)果。例如擲硬幣，正面記為“1”，反面記為“0”。2.則擲硬幣N次，有k次是1的概率為：二項分布的期望值與標(biāo)準(zhǔn)方差期望值E(x)=μ方差VarX=σ2泊松分布

(Poissondistribution)1.稀有事件發(fā)生的概率：在一個連續(xù)的時間或空間中，稀有離散變量出現(xiàn)的概率2.N->∞,E(x)=μe=2.71828…泊松分布與二項分布的近似對于大的N及小的p值的二項分布，能夠相當(dāng)準(zhǔn)確地用一個參數(shù)為μ=Np的泊松分布近似。當(dāng)實驗次數(shù)很多而概率很小時：二項分布~泊松分布例1：鳥槍法的覆蓋率假設(shè)：需要測序的BAC長度200kbp;總共測序的序列數(shù)量：N;每次測序：500bp；每次測序的覆蓋率p：500/200kbp=0.0025因此：每個點平均覆蓋到的次數(shù):μ=N*pk:測序能夠覆蓋到點X的次數(shù)。鳥槍法：覆蓋率點X被覆蓋k次的概率：(二項分布~泊松分布)當(dāng)點X一次都不被覆蓋時，k=0;此時的概率為：覆蓋率vs.準(zhǔn)確性例2：泊松分布Prof.Gene發(fā)現(xiàn)一種序列上的調(diào)控信號，在人的基因組上平均每500kbp一個。那么，隨機給一條1mbp的序列，在上面發(fā)現(xiàn)5個這樣的信號，完全是隨機產(chǎn)生的概率是多少？本例中，E(x)=μ=2(1mbp/500kbp)統(tǒng)計顯著性：p-value<0.05超幾何分布與二項式分布的區(qū)別：不放回抽樣。例：有N個球，其中紅球M個，白球N-M個，每次拿出一個球再放回，總共n次，其中有m個球是紅球的概率為(二項式分布)：p=M/N超幾何分布(2)上例改為：有N個球，其中紅球M個，白球N-M個，每次拿出一個球不放回，總共n次，其中有m個球是紅球的概率為：并且，0≤m≤M<N超幾何分布右尾概率上例再改為：有N個球，其中紅球M個，白球N-M個，每次拿出一個球不放回，總共n次，其中至少有m個球是紅球的概率為：并且，0≤m≤M<N超幾何分布左尾概率上例再改為：有N個球，其中紅球M個，白球N-M個，每次拿出一個球不放回，總共n次，其中最多有m個球是紅球的概率為：并且，0≤m≤M<N超幾何分布雙尾概率所有出現(xiàn)概率<=觀察表概率的概率之和Fisher'sExactTest超幾何分布的精確概率計算。前提是固定邊際分布，即a+b、c+d、a+c與b+d的值不變。RAFisher,1935年文章示例：猜測先放的飲料牛奶茶合計實際先放的飲料牛奶a=3b=1a+b=4茶c=1d=3c+d=4合計a+c=4b+d=4n=8Fisher'sExactTest

計算公式：=統(tǒng)計顯著性

假設(shè)檢驗中的P值(Pvalue)Pvalue：一種在原假設(shè)為真的前提下出現(xiàn)觀察樣本以及更極端情況的概率。顯著性水平A：認(rèn)為預(yù)先設(shè)定的顯著性水平閾值，P<A為顯著。一般以P<0.05為顯著，P<0.01為非常顯著，其含義是樣本間的差異由抽樣誤差所致的概率小于0.05或0.01。假設(shè)檢驗本例中，零假設(shè)H0：該女同事只是隨便亂猜答案；備擇假設(shè)H1：該女同事所言不虛；

Pvalue計算：

P(a=3|a+b=c+d=a+c=b+d=4)=0.229P(a=4|a+b=c+d=a+c=b+d=4)=0.014例3：超幾何分布Prof.Gene從人的26873個蛋白質(zhì)中預(yù)測了2264個能結(jié)合某類金屬離子X?，F(xiàn)已知，人的26873個蛋白質(zhì)中有421個蛋白質(zhì)具有某種功能結(jié)構(gòu)域D，而在預(yù)測的2264個X金屬蛋白中，有94個具有結(jié)構(gòu)域D。問：結(jié)構(gòu)域D在2264個X金屬蛋白中是顯著出現(xiàn)，顯著不出現(xiàn)，還是隨機出現(xiàn)？例3：超幾何分布問題轉(zhuǎn)化：在26873個蛋白質(zhì)的體系中，取出2264個蛋白質(zhì)，其中至少有94個蛋白質(zhì)具有功能結(jié)構(gòu)域D的概率是多少？N=26873;n=2264;M=421;m=94;例3：超幾何分布非X金屬蛋白X金屬蛋白合計不含結(jié)構(gòu)域DN-M+m-nM-mN-n含結(jié)構(gòu)域Dn-mmn合計N-MMN例3：超幾何分布a+b+c+d=26873c+d=2264b+d=421d=94langsrud/fisher.htm第二節(jié)，雙序列比對算法1.DotMatrix，點陣法2.動態(tài)規(guī)劃算法：Global:Needleman-WunschLocal:Smith-Waterman3.Wordork-tuple算法：FASTA,BLAST1.點陣法1970年，Gibbs&McIntyre；尋找兩條序列間所有可能的比對；發(fā)現(xiàn)蛋白質(zhì)或者DNA序列上正向或者反向的重復(fù)；發(fā)現(xiàn)RNA上可能存在的互補區(qū)域。工具：myhits.isb-sib.ch/cgi-bin/dotlet/molkit/dnadot/點陣法：自身的比對AKGFKCADEA100000100K10010000G1000000F100000K10000C1000A100D10E1點陣法：重復(fù)序列AKGFDKGFEA100000000K10001000G1000100F100010D10000K11000G1100F110E1點陣法：反向重復(fù)/回文AUGCACGUCA100010000U10000010G1000100C101000A10000C11001G1100U110C1點陣法：不同序列的比對PKDFCKALVP100000000K10001000F0100000T00000K11000A100I00V011:PKDFCKALV2:PK-FTKAIVSeq1Seq2點陣法的序列比對Sequence1#1nSequence2#1m“-”Insertion“-”Insertion計算效率用CPU的計算時間和內(nèi)存占用量來衡量；對于需要解決的問題，其單位數(shù)量n在某算法下運算的基本操作重復(fù)執(zhí)行次數(shù)表示為f(n)；時間復(fù)雜度:T(n)=O(f(n))；如果需要解決的問題的大小與單位數(shù)量n的平方成正比，則O(n2)對于算法來說：O(1)>O(log(n))>O(n)>O(n2)>O(an)>O(n!)NP問題1.一般的，O(nk),當(dāng)k≤3時，為多項式時間，較為容易處理。2.當(dāng)O(an)，則難以處理。3.NP完全問題（NPC）：無法找到能夠在多項式時間復(fù)雜度內(nèi)解決方法的問題；4.近似算法/優(yōu)化算法，求近似解。P/NP問題-千禧年大獎難題之一1900年，德國數(shù)學(xué)家DavidHilbert提出的23個歷史性數(shù)學(xué)難題。千禧年大獎難題美國克雷數(shù)學(xué)研究所（ClayMathematicsInstitute,CMI）于2000年5月公布七個世界數(shù)學(xué)難題。千禧年大獎難題P/NP問題：P=NP？霍奇猜想龐加萊猜想（已證明）黎曼猜想楊-米爾斯存在性與質(zhì)量間隙納維-斯托克斯存在性與光滑性貝赫和斯維訥通-戴爾猜想P/NP/NPC問題P問題：PolynomialProblems可以在多項式(polynomial)時間內(nèi)解決的問題;NP:“Non-deterministicPolynomial”,并非

“Non-Polynomial”

可以在多項式的時間里驗證一個解的問題;NPC:NP-complete2.動態(tài)規(guī)劃算法1.打分模型、替代矩陣以及空位罰分。2.比對算法：遞歸及動態(tài)規(guī)劃算法；3.全局優(yōu)化比對：Needleman-Wunsch4.局部優(yōu)化比對：Smith-Waterman5.工具資源：.au/course/lectures2019/Likic.pdfsnippets.dzone/posts/show/2199/NW-align/jaligner.sourceforge/打分模型1.字符相同：identity2.字符替代：similarity，相似性，氨基酸/堿基之間的替代和突變3.插入和缺失4.空位罰分BLOSUM62替代矩陣空位罰分1.線性罰分：d,每次罰分的分?jǐn)?shù)；g，空位數(shù)2.修正的罰分：d,第一次罰分的分?jǐn)?shù)；g，空位數(shù)；e,修正后的參數(shù)遞歸和動態(tài)規(guī)劃算法兩條序列的比較，無空位：時間復(fù)雜度為O(n2);

兩條序列比對，允許空位，時間復(fù)雜度為：因此，有空位的雙序列比對，時間復(fù)雜度為：O(22n)，指數(shù)增加，NPC問題！遞歸和動態(tài)規(guī)劃算法(2)

數(shù)學(xué)上保證提供最優(yōu)解。動態(tài)規(guī)劃算法：比較所有可能的字符對，考慮匹配、錯配以及空位罰分，并且將比對次數(shù)控制在多項式時間內(nèi)。替代矩陣：BLOSUM62，空位罰分：11延伸的空位罰分：1(BLAST工具)例：全局比對序列1：VDS–CY序列2：VESLCY替代矩陣中的分?jǐn)?shù)：424-1197兩序列比對的總分：Score=Σ(AApairscores)–gappenalty=15動態(tài)規(guī)劃算法：全局比對GapVDSCYGap01gap2gap…V1gapE2gapS…LCY本例：線性罰分全局比對(2)GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66要求解Sij的分?jǐn)?shù)，我們必須先知道Si-1,j-1,Si-1,j,以及Si,j-1的分?jǐn)?shù)，這種方法叫做遞歸算法；采用這種方法，可以把大的問題分割成小的問題逐一解決，即動態(tài)規(guī)劃算法；需要存儲如何得到Sij分?jǐn)?shù)的過程。全局比對(3)GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66ijNeedleman-Wunsch算法；Sij=maxofSi-1,j-1+σ(xi,yj)

Si-1,j-d(從左到右)Si,j-1-d(從上到下)全局比對(4)GapVDSCYGap0-11-22-33-44-55V-11SijE-22S-33L-44C-55Y-66ij4-11-11Needleman-Wunsch算法；Sij=maxofSi-1,j-1+σ(xi,yj)

Si-1,j-d(從左到右)Si,j-1-d(從上到下)全局比對(5)GapVDSCYGap0-11-22-33-44-55V-114E-22S-33L-44C-55Y-664-11-11全局比對(6)GapVDSCYGap0-11-22-33-44-55V-114SijE-22S-33L-44C-55Y-66-3-11-11Needleman-Wunsch算法；Sij=maxofSi-1,j-1+σ(xi,yj)

Si-1,j-d(從左到右)Si,j-1-d(從上到下)全局比對(7)GapVDSCYGap0-11-22-33-44-55V-114-7E-22S-33L-44C-55Y-66-3-11-11全局比對(8)GapVDSCYGap0-11-22-33-44-55V-114-7-18-29-40E-22-76-5-16-27S-33-18-510-1-12L-44-29-16-19-3C-55-40-27-1287Y-66-51-38-23-31542回溯：比對結(jié)果GapVDSCYGap0-11-22-33-44-55V-114-7-18-29-40E-22-76-5-16-27S-33-18-510-1-12L-44-29-16-19-3C-55-40-27-1287Y-66-51-38-23-31542VDS–CYVESLCY局部優(yōu)化比對下例：局部優(yōu)化打分兩條序列如下：LDS–CHGESLCK目標(biāo)：使用局部優(yōu)化算法尋找比對的結(jié)果局部優(yōu)化比對(2)1.Smith-Waterman算法；2.時間復(fù)雜度O(n2)；3.Sij=maxof0 Si-1,j-1+σ(xi,yj) Si-1,j-d(從左到右) Si,j-1-d(從上到下)本例中：gap:12，線性罰分模型。局部優(yōu)化比對(3)GapLDSCHGap00