前言

學(xué)過生信的肯定知道Needleman–Wunsch算法和Smith–Waterman 算法，一個用來進(jìn)行全局比對，一個用來進(jìn)行局部比對。單純看算法抽象后的算法公式的話也不復(fù)雜，用兩個短序列來套公式計算的話，畫過箭頭的同學(xué)都知道不難。本文簡單講述Needleman–Wunsch算法是如何利用動態(tài)規(guī)劃來尋找兩序列最大相似度。

Python代碼實(shí)現(xiàn)：Needleman–Wunsch 算法代碼

動態(tài)規(guī)劃

動態(tài)規(guī)劃是數(shù)學(xué)上的優(yōu)化方法，也是計算機(jī)編程的一種方法。對于一名生物狗，從數(shù)學(xué)層面去看動態(tài)規(guī)劃，就太過硬核了。我們就從計算機(jī)編程角度理解動態(tài)規(guī)劃就好了。
動態(tài)規(guī)劃是一種解決問題的方法，將一個復(fù)雜大問題拆分為幾個小的問題，從幾個小問題的最優(yōu)解推出大問題的最優(yōu)解。當(dāng)然這不是說將一個大問題拆分了就能用動態(tài)規(guī)劃了。同時問題也需要滿足兩個性質(zhì)，這兩個性質(zhì)后面再講。先回到我的主要探討問題全局比對上：對于兩條蛋白質(zhì)或者核酸序列來說，尋找它們的最大相似度，如何拆分成幾個小問題呢？（提示：矩陣、算法公式）

拆分問題

就這樣憑空想，如何拆分，我覺得我是想破腦袋也想不了T_T. 所以我們按照，學(xué)全局比對的套路手動來一遍序列比對操作，從中看看如何拆分問題。

為方便計算，得分矩陣就弄個簡單點(diǎn)的：

也就說完全匹配得10分，不全得5分，匹配到空位-5分，空位對空位一開始用來初始化時，才有意義。在之后比對中空位對空位沒有意義。
先給出公式（d為gap罰分，即堿基比對空格的得分）：

來自高歌老師PPT

現(xiàn)假設(shè)我們有兩條序列Seq1: CATCCCCCAAA 和 Seq2: CAAAA。長度分別為11，5。根據(jù)它們的長度畫一個合適的二維矩陣，填充好初始數(shù)據(jù)：

image.png

加粗的分別是 Seq1,Seq2，其中 -代表空格。一開始空格對空格，我們初始化，記，同時記空格分別對應(yīng)的行和列是第0行，第0列。然后第1行，第1列分別是兩條首字母對應(yīng)的行列。用代表（行，列），能有多大，取決于序列對的長度。的分?jǐn)?shù)是一條序列中第個堿基和另一條序列中第個堿基比對的分?jǐn)?shù)。放在矩陣中就是了第行的堿基和第列的堿基比對結(jié)果。

按照函數(shù)公式，的值由左上的值(),上方值()和左方值()確定。
（注:這里的左上，上，左是一個相對位置，在當(dāng)前設(shè)定矩陣坐標(biāo)系中才有意義，這里是將位置放在了左上方，向右和向下延伸。但其實(shí)，你也可以將出放在其它位置，比如右下方）

極端比對

我們看下第0行，它們只能由左方的值得到。我們算幾個值：



第0列，只能由上方的值得到，計算同理，不再演示了。
因?yàn)榈?行和0列，很容易就得出來了，所以一般畫矩陣時，就順便一起畫出來了。

你應(yīng)該知道，往右或往下延申，則比對到了空格。那這是為啥？
就拿第0行的比對情況舉個例子，第0行中，已經(jīng)固定了，它代表一條序列（后面稱為序列）的第0元素。,分別代表另一條序列（后面稱為序列）第0，1，2，3...元素。序列第一個元素是序列的第一個堿基。第0個是啥玩意，第0個啥也不是，那我們就把它當(dāng)做空格吧。

1. 序列比對從兩條序列的第元素開始，兩個空格比對，得分為0。
2. 固定，一直表示序列第0個元素，但是是一直增加的，所以時，就是用序列的的第個元素和序列的第個元素比對，但是序列第元素已經(jīng)和序列的第元素比對過了，你不能重復(fù)比對啊，你也不能私自和序列的元素比對，所以沒辦法，就只能和一個空格比對，產(chǎn)生空位罰分，在之前得分基礎(chǔ)上上加上空位罰分。
3. 序列的第元素都是如此，只能比對空格，產(chǎn)生空位罰分。

比對往右延申的過程就是行數(shù)固定了，沒有增加，但是列數(shù)增加，又兩條序列的中的元素不能重復(fù)比對，也不能比對序列其它的元素，所以就只能比對空格了。（向下延申同理。）
我們看兩個極端的例子：

image.png

##黑色路線
CATCCCCCAAA-----
-----------CAAAA
##紅色路線
-----CATCCCCCAAA
CAAAA-----------

你們說這第0行，0列的有啥用處，放這么多負(fù)分而產(chǎn)生極端比對。要排除這些不合適比對結(jié)果也要花費(fèi)許多資源。

比對開始

上面進(jìn)行了一次小探討。比對直接向右或向下會產(chǎn)生gap。因?yàn)樾辛杏幸环奖还潭?，沒法提供可比對元素。如果要正確比對的話，要使得和同時增加，即向右下方延申。除了第0行和0列的元素，只有單一來源，矩陣其它的元素都有三個來源，即上方，左方，左上。如何確定的值，按照公式，選擇其中最大的值。(為兩條序列的長度)為最終序列比對的相似得分。
既然我們想使最大，肯定不能過多的產(chǎn)生空位罰分吧。那么我們從一開始就一直向右下方延申，直至不能繼續(xù)，能到得到最大的相似得分呢？我們試一下吧：

image.png

image.png

貪心算法：紅色路線選擇的策略，從一開始，三個方向的選擇，每一步，它都選擇三個方向中得分最大的那個，但最后結(jié)果卻不是最好的。這種在每一步都做出當(dāng)前看來是最好選擇的方法，而不考慮整體最優(yōu)的方法，叫做貪心算法。這種方法在這里得不到最優(yōu)結(jié)果。
窮舉搜索：這是一種很直觀的方法，就是把所有可能的情況都列舉出來，從中找最優(yōu)的。但是用于序列比對是是非常差的一種方法。從矩陣左上方到達(dá)矩陣右下方，有多少種路線？矩陣每一個小格有三種可能路線，矩陣大小,即有種路線。其中稍微取大一點(diǎn)的值，數(shù)值就爆炸了。完全不可取。

emm，寫到這里，還沒講到全局比對是如何拆分問題的，現(xiàn)在我們就開始討論這個問題吧。。
首先，在矩陣中當(dāng)比對到第格時，代表著比對終止。因?yàn)閮蓷l序列的每個元素都參加了比對過程，無論序列中的元素是比對到了空格，還是另一條序列的某個元素。因?yàn)樾蛄兄械脑囟急葘ν炅耍诒葘ο氯ゾ椭挥锌崭駥崭窳?，沒什么意義。所以是終點(diǎn)。
我們就從終點(diǎn)開始探討如何分割大問題為小問題吧。

其中和是知道的,那問題則可以通過求出,和的分?jǐn)?shù)得到。即一個大問題，分成了3個子問題了。三個子問題也不知道分?jǐn)?shù)，就繼續(xù)分子問題唄。如圖：

image.png

分治算法：分治的思想也是將一個問題分為幾個問題來解決。這與動態(tài)規(guī)劃的分離子問題有啥區(qū)別？從上面的圖我們可以看出來由分割的子問題們是存在重復(fù)問題的。分治分離的子問題一般獨(dú)立的不會相互重疊，而動態(tài)規(guī)劃的子問題則一般發(fā)生重疊現(xiàn)象。

兩個性質(zhì)

之前空著沒講的兩個性質(zhì)，現(xiàn)在來講一下。一個是最優(yōu)子結(jié)構(gòu)，一個是無后效性
最優(yōu)子結(jié)構(gòu)：如果一個問題的最優(yōu)解，可以將其分為幾個子問題，然后能從子問題的最優(yōu)解推出大問題的最優(yōu)解。就稱其有最優(yōu)子結(jié)構(gòu)性質(zhì)。這一性質(zhì)很容易從函數(shù)公式中看出來。
無后效性：給定某一階段的狀態(tài)，這一階段以后的發(fā)展過程不受這階段以前各段狀態(tài)的影響（來自阮行止的回答）這話是什么意思呢，以我們序列比對為例。比如：當(dāng)你確定了前面某個階段的比對，比對情況如下面的比對方案的前面一部分，這時比對分?jǐn)?shù)已經(jīng)確定了。無后效性就是前面的比對情況對后面的比對結(jié)果的分?jǐn)?shù)不會有任何影響。前面部分的比對方案的任意變動也不會影響后面部分的結(jié)果。

##1
CAT  CCCCC-----AAA
CA-  ----------AAA
##2
CAT  CCCCCAAA-----
CA-  AA----------A
##3
CAT  CCCCCAAA-----
CA-  ----------AAA

無后效性，可以讓我們放心大膽的在矩陣每個格子存放最優(yōu)的結(jié)果，舍棄許多不是最優(yōu)的分?jǐn)?shù)。因?yàn)椴淮嬖谇懊骐A段的最差比對方案，還能與后面階段的最差方案組合成為全局最優(yōu)比對方案。

復(fù)雜度

我們來看一下Needleman–Wunsch算法的復(fù)雜度。從二維矩陣中，可以看到我們需要計算每一格的所占最優(yōu)分值，同時也為了之后的計算方便也得儲存每一個的分值。故時間復(fù)雜度和空間復(fù)雜度都是。

其它

為了理解動態(tài)規(guī)劃，我在知乎上這個問題什么是動態(tài)規(guī)劃（Dynamic Programming）？動態(tài)規(guī)劃的意義是什么？的答案下看了許多回答。許多答主都給了自己的理解。他們的回答有相似之處，也有不同。不能說誰對誰錯（小聲bb:也不敢說誰對誰錯）
就全局比對我也談?wù)勎业目捶ǎ?br> 動態(tài)規(guī)劃和其它方法有許多相似的地方。要是將其它方法的定義定廣一點(diǎn)，也能把動態(tài)規(guī)劃說成其它方法。比如：
窮舉搜索： 之前上文曾提到過，不過那時我們說要窮舉的對象是每一個比對方案?，F(xiàn)在我們換一個觀點(diǎn)想一想：窮舉每一個矩陣中的格子，即窮舉每一個最優(yōu)的。這里我們的對象就從比對路徑而換成了（序列比對到第的最大相似度）
貪心算法: 之前提到的貪心算法，我們狹隘的將最優(yōu)，認(rèn)為是向下、向右或向右下而能得到的分?jǐn)?shù)?，F(xiàn)在我們將最優(yōu)看成相鄰中最大的那個。即按照這個貪心選取，它會可能不斷的對比相鄰的而選擇最優(yōu)。直至，到最后剩下三條路徑到達(dá),選中最優(yōu)的那個。
總之，不管黑貓白貓，能幫助解決問題就是好貓！選擇你認(rèn)為最好的一種理解。在以后多次遇到動態(tài)規(guī)劃后，在慢慢完善你的認(rèn)知。