這篇筆記是我蹭課《生物信息學(xué)》的課堂筆記，本人為非生信專業(yè)，主要做實(shí)驗(yàn)，所以此文僅供個(gè)人學(xué)習(xí)罐监。

序列搜索算法的進(jìn)化

雖然現(xiàn)代高通量生物分子技術(shù)產(chǎn)生了各種類型的數(shù)據(jù)，但生物序列數(shù)據(jù)（biosequence data）仍然是生物信息學(xué)分析的核心瞒爬。

NCBI在2007年底推出了SRA（Sequence Read Archive）數(shù)據(jù)庫弓柱，用于存儲(chǔ)二代測序數(shù)據(jù)，目前包含大約15PB堿基侧但，其中幾乎一半在開放訪問中矢空。SRA數(shù)據(jù)庫地址https://www.ncbi.nlm.nih.gov/sra/。SRA是國際核苷酸測序數(shù)據(jù)庫合作（INSDC）項(xiàng)目的一部分禀横，INSDC包括NCBI的SRA屁药，歐洲生物信息研究所（EBI）和日本DNA數(shù)據(jù)庫（DDBJ），SRA數(shù)據(jù)庫是美國國立衛(wèi)生研究院（NIH）存儲(chǔ)高通量測序數(shù)據(jù)的主要數(shù)據(jù)庫柏锄。數(shù)據(jù)上傳到以上任何一個(gè)組織者祖，其他兩個(gè)組織都是可以共享的立莉。

序列比對(duì)的意義：相似的序列可能會(huì)有相似的結(jié)構(gòu)绢彤，從而有相似的功能七问。序列之間的相似性可以幫助我們推斷一個(gè)未知新序列可能的功能。不同物種中相似的序列往往意味著其具有共同的祖先茫舶，即同源械巡。事實(shí)上，序列間的相似性是在演化分析中用來構(gòu)建演化樹的重要依據(jù)之一饶氏。（序列同源是指兩條序列有共同的祖先讥耗，然后經(jīng)過突變變得不同。同源無法考證疹启，只能通過相似性來推斷同源古程。同源不一定相似，相似也不一定同源喊崖。）

相似性矩陣（similarity matrix）或打分矩陣：兩條序列做全局比對(duì)挣磨，然后計(jì)算全局比對(duì)中一致字符的個(gè)數(shù)和相似字符的個(gè)數(shù)，再除以全局比對(duì)的長度荤懂，就可以得到兩序列的一致度（identity）和相似度（similarity）茁裙。替換記分矩陣是反映殘基之間相互替換率的矩陣，也就是說它描述了殘基兩兩相似的量化關(guān)系节仿。DNA的替換計(jì)分矩陣有等價(jià)矩陣（相同1分晤锥，不同0分）、轉(zhuǎn)換-顛換矩陣（轉(zhuǎn)換-1分廊宪，顛換-5分）矾瘾、BLAST矩陣（相同5分，不同-4分箭启，是實(shí)踐經(jīng)驗(yàn)所得）壕翩。蛋白質(zhì)的常用替換計(jì)分矩陣PAM 矩陣基于進(jìn)化原理，如果兩種氨基酸替換頻繁册烈，說明自然界容易接受這種替換戈泼，那么這一對(duì)氨基酸替換的得分就高∩蜕基礎(chǔ)的PAM-1矩陣反應(yīng)的是進(jìn)化產(chǎn)生的每一百個(gè)氨基酸平均發(fā)生一個(gè)突變的量值大猛，由統(tǒng)計(jì)方法得到。PAM-1自乘n次淀零，可以得到PAM-n挽绩，表示兩序列的親緣關(guān)系更遠(yuǎn)。目前最常用的蛋白序列比對(duì)打分矩陣是BLOSUM 62矩陣驾中，這個(gè)矩陣是由一致度≥62%的序列計(jì)算而來的唉堪∧Ａ≥0 的得分代表對(duì)應(yīng)的一對(duì)氨基酸為相似氨基酸，<0 的是不相似的氨基酸唠亚。

早期算法：動(dòng)態(tài)規(guī)劃（dynamic programing）

檢測兩條短序列（如基因）的相似性（similarities）链方，通常只考慮點(diǎn)突變（point mutations），即字符替換（character substitutions）灶搜、插入或刪除（insertions or deletions祟蚀，indel）。最簡單的兩序列比對(duì)方法有打點(diǎn)法割卖，用對(duì)角線表示相同區(qū)域前酿。序列比對(duì)就是運(yùn)用特定的算法找出兩個(gè)或者多個(gè)序列之間產(chǎn)生最大相似度。根據(jù)比對(duì)序列的個(gè)數(shù)可以把序列比對(duì)分為雙序列比對(duì)和多序列比對(duì)鹏溯。根據(jù)序列比對(duì)的算法不同罢维，雙序列比對(duì)又分為全局比對(duì)和局部比對(duì)。

編輯距離（edit distance）是指兩個(gè)字符串之間丙挽，由一個(gè)轉(zhuǎn)成另一個(gè)所需的最少編輯操作次數(shù)肺孵。可進(jìn)行替換取试、插入或刪除操作悬槽。一般來說，編輯距離越小瞬浓，兩個(gè)字符串的相似度越大初婆。簡單的編輯距離不足以反映生物的距離，被加權(quán)編輯距離（weighted edit distance）所取代猿棉，通過相似性評(píng)分（similarity score）來定義磅叛，其中不同的字符替換可能有不同的罰分（penalties）。然后將兩個(gè)序列的整體比較形式化為最優(yōu)全局比對(duì)（optimal global alignment）問題萨赁，并用生物信息學(xué)中的經(jīng)典動(dòng)態(tài)規(guī)劃算法Needleman-Wunsch算法求解弊琴。

動(dòng)態(tài)規(guī)劃是求解最優(yōu)化問題的一種途徑，其基本思想是將過程分成若干個(gè)互相聯(lián)系的階段杖爽，在它的每一階段都作出最優(yōu)決策敲董，從而使整個(gè)過程達(dá)到最好的活動(dòng)效果。動(dòng)態(tài)規(guī)劃算法對(duì)于序列長度m慰安、n具有時(shí)間復(fù)雜度（time complexity）O(nm)腋寨，因此對(duì)于大規(guī)模序列搜索的計(jì)算要求太高。雖然全局比對(duì)的空間復(fù)雜度（space complexity）可以通過一個(gè)簡潔的分而治之的技巧（divide-and-conquer trick）使之線性化化焕。

最優(yōu)全局比對(duì)算法——Needleman-Wunsch算法：使用先前對(duì)于較小子序列的最優(yōu)比對(duì)結(jié)果萄窜，來進(jìn)行下一步的最優(yōu)比對(duì)。按照線性空位罰分（linear gap penalty）（1個(gè)gap -2分，k個(gè)gap -2k分）計(jì)算序列比對(duì)的得分查刻。

舉個(gè)栗子：右下圖比對(duì)键兜，繪制多一行多一列的表格，從左上角0開始穗泵，向下移動(dòng)（在橫向序列引入一個(gè)gap）得-2分普气，向右移動(dòng)（在縱向序列引入一個(gè)gap）得-2分，對(duì)角線移動(dòng)得到match得1分火欧，mismatch得-1分棋电。每個(gè)位置得分由它的上一個(gè)位置得分加上移動(dòng)得分計(jì)算得到。最終右下角最高得分為S(i,j)=-1苇侵，有3種最優(yōu)全局比對(duì)方式（從右下角到左上角回溯的路徑）。

最優(yōu)局部序列比對(duì)（Local sequence alignment）算法——Smith-Waterman算法：尋找輸入的兩個(gè)序列V和M的子序列（subsequences）的最優(yōu)全局比對(duì)企锌。適用于尋找相似的小片段榆浓，如蛋白質(zhì)的共有結(jié)構(gòu)域、不同序列的保守區(qū)域撕攒。計(jì)算方法與全局比對(duì)相似陡鹃，只是整個(gè)表格中沒有負(fù)分，最小為0分抖坪。在整個(gè)計(jì)分表中找到最大的分?jǐn)?shù)萍鲸，并向前追溯到0停止。

Gotoh’s算法進(jìn)一步將Needleman-Wunsch算法和Smith-Waterman算法擴(kuò)展到仿射空位罰分（affine gap penalties）的情況擦俐，對(duì)k個(gè)連續(xù)（consecutive）插入或刪除的罰分是空位開始罰分g（gap opening penalty）與空位延伸罰分e（k gap extension penalties）之和脊阴，penalty = g + (k ? 1)e，|e|≤|g|蚯瞧，往往第一個(gè)空位罰分較多嘿期。gap開頭罰分大，gap延長罰分小的時(shí)候埋合，gap集中連成長串出現(xiàn)备徐。

長度分別為m和n的兩條序列的所有可能的比對(duì)個(gè)數(shù)為

Filtration-based heuristics and database search 基于過濾的試探法和數(shù)據(jù)庫搜索

基本局部比對(duì)搜索工具BLAST（Basic Local Alignment Search Tool）是一種啟發(fā)式的（heuristic）序列比對(duì)方法甚颂，它遵循兩種一般的思想：過濾和索引（filtration and indexing）蜜猾。與給出最優(yōu)解的動(dòng)態(tài)規(guī)劃算法不同，BLAST是一種近似算法振诬。

BLASTp是用蛋白質(zhì)序列搜索蛋白質(zhì)序列數(shù)據(jù)庫蹭睡，BLASTn是用核酸序列搜索核酸序列數(shù)據(jù)庫，BLASTx是將核酸序列按6條鏈翻譯成蛋白質(zhì)序列后搜索蛋白質(zhì)序列數(shù)據(jù)庫贷揽。tBLASTn 是用蛋白質(zhì)序列搜核酸序列數(shù)據(jù)庫棠笑，核酸數(shù)據(jù)庫中的核酸序列要按6條鏈翻譯成蛋白質(zhì)序列后再被搜索。tBLASTx是將核酸序列按6條鏈翻譯成蛋白質(zhì)序列后搜索核酸序列數(shù)據(jù)庫禽绪，核酸數(shù)據(jù)庫中的所有核酸序列也要按6條鏈翻譯成的蛋白質(zhì)序列后再被搜索。除了按照搜索內(nèi)容分類萨蚕，BLAST還可以根據(jù)搜索算法不同分為標(biāo)準(zhǔn)BLAST敷搪，PSI-BLAST（位點(diǎn)特異性迭代，Position-Specific Iterated）和PHI-BLAST（模式識(shí)別斩例，Pattern-Hit Initiated）等。

BLAST的思想是真正匹配的比對(duì)很有可能包含得分非常高的配對(duì)片段从橘，即seeds念赶。片段對(duì)是指兩個(gè)給定序列中的一對(duì)子序列，它們的長度相等恰力，且可以形成無空位的完全匹配叉谜。首先在搜索的序列中選取固定長度的小片段（如3個(gè)氨基酸或11個(gè)核苷酸），稱為seeds/anchors踩萎，也就是在序列中從頭到尾一個(gè)殘基一個(gè)殘基地移動(dòng)停局，每3個(gè)作為一個(gè)seeds。然后從打分矩陣計(jì)算各種3氨基酸組合的得分香府，設(shè)置一個(gè)閾值董栽，大于閾值的3氨基酸作為比對(duì)上的，數(shù)據(jù)庫中不能匹配的部分都從考慮中丟棄（過濾filtered out）企孩。然后從seeds開始向兩端擴(kuò)展锭碳，得到高分值片段對(duì)（high scoring segment pairs，HSPs）勿璃。注意擒抛，低復(fù)雜度的序列（重復(fù)序列）應(yīng)該首先過濾掉。這樣計(jì)算復(fù)雜度是序列長度n的一次方蝗柔，如果做雙序列比對(duì)話需要構(gòu)建一個(gè)n乘以n的表格闻葵，計(jì)算復(fù)雜度是n的二次方。

搜索使用數(shù)據(jù)庫的預(yù)構(gòu)建索引（pre-built index）癣丧，通常是哈希表（hash table）槽畔，加快定位seed patterns的位置。這種方法的性能被定義為靈敏性（Sensitivity）與特異性（specificity）之間的權(quán)衡（trade-off）胁编，并由用于過濾的種子類型決定厢钧。由于分?jǐn)?shù)小于閾值的seeds被直接過濾掉了，所以其靈敏性不高嬉橙。靈敏性是指真實(shí)為陽性的樣本中早直，用這種方法能夠檢測出來的比率。

由于連續(xù)11個(gè)核苷酸完全相同（k consecutive matching nucleotides）要求太高市框，又提出了一種改進(jìn)的算法霞扬，選擇根據(jù)給定模式的k個(gè)間隔的匹配核苷酸作為種子，即間隔的種子（spaced seed），可以改善靈敏性與特異性之間的權(quán)衡喻圃。

索引是為了加速對(duì)表中數(shù)據(jù)行的檢索而創(chuàng)建的一種分散的存儲(chǔ)結(jié)構(gòu)萤彩。索引是針對(duì)表而建立的，它是由數(shù)據(jù)頁面以外的索引頁面組成的斧拍，每個(gè)索引頁面中的行都會(huì)含有邏輯指針雀扶，以便加速檢索物理數(shù)據(jù)。哈希表（Hash table肆汹，也叫散列表）愚墓，是根據(jù)關(guān)鍵碼值(Key value)而直接進(jìn)行訪問的數(shù)據(jù)結(jié)構(gòu)，通過把關(guān)鍵碼值映射到表中一個(gè)位置來訪問記錄昂勉，以加快查找的速度浪册，這個(gè)映射函數(shù)叫做哈希函數(shù)或散列函數(shù)。Hash（散列硼啤、哈希）议经，就是把任意長度的輸入通過算法變成固定長度的輸出。直觀解釋起來谴返，就是對(duì)一串?dāng)?shù)據(jù)m進(jìn)行雜糅，輸出另一段固定長度的數(shù)據(jù)h咧织，作為這段數(shù)據(jù)的特征（指紋）嗓袱。MD5和SHA都是歷史悠久的Hash算法。

（參考文獻(xiàn)：Kucherov,G. (2019) Evolution of biosequence search algorithms: a brief survey. Bioinformatics, 35, 3547–3552.）