Needleman–Wunsch算法

**尼德曼 - 翁施算法**（英語：Needleman-Wunsch Algorithm）是基於生物信息學的知識來匹配蛋白序列或者 DNA 序列的算法。這是將動態算法應用於生物序列的比較的最早期的幾個實例之一。該算法是由 Saul B. Needlman 和 Christian D. Wunsch 兩位科學家於 1970 年發明的。本算法高效地解決了如何將一個龐大的數學問題分解為一系列小問題，並且從一系列小問題的解決方法重建大問題的解決方法的過程。該算法也被稱為優化匹配算法和整體序列比較法。時至今日尼德曼 - 翁施算法仍然被廣泛應用於優化整體序列比較中。

## 新手指南

本算法可以用於比較任意兩個序列。本指南將會使用兩個短小的 DNA 序列作為例子來演示該算法

```
GCATGCU
GATTACA
```

### 初始化得分矩陣

首先建立如圖一的得分矩陣。從第一列第一行的位置起始。計算過程從 d 0 , 0 開始，可以是按行計算，每行從左到右，也可以是按列計算，每列從上到下。當然，任何計算過程，只要滿足在計算 d i , j 時 d i-1 , j、d i-1 , j-1 和 d i, j-1 都已經被計算這個條件即可。在計算 d i , j 後，需要保存 d i , j 是從 d i-1 , j、d i-1 , j-1 或 di, j-1 中的哪一個推進的，或保存計算的路徑，以便於後續處理。上述計算過程到 d m , n 結束。

### 選擇得分體系

接下來我們需要決定如何給每個獨立配對打分。從我們的序列中，你可能就能發現最好的序列配對之一：

```
GCATG-CU
G-ATTACA
```

我們可以看出每個位置配對都有三種可能情況：匹配， 不匹配與錯位（或插入）：

* 匹配： 兩個字母相同
* 不匹配：兩個字母不相同
* 錯位：一個字母與另一個序列中的間隔（空白）相匹配

這三種得分情況有很多打分標準，這些情況都總結在得分體系的部分。從現在開始，我們將使用 Needleman 和 Wunsch 創造的簡單體系來進行打分，即匹配得一分，不匹配得 - 1 分，錯位得 - 1 分.

### 填充得分矩陣

開始於第二行中的第二列，初始值為 0。通過矩陣一行一行移動，計算每個位置的分數。得分被計算為從現有的分數可能的最佳得分（即最高）的左側，頂部或左上方（對角線）。當一個得分從頂部計算，或從左邊這代表在我們的對準的插入缺失。當我們從對角線計算分數這表示兩個字母所得位置匹配的對準。定不存在「向上」或「左上」的位置對第二行，我們只能從現有單元向左添加。因此，我們添加 - 1 的權利，因為這代表了從以前的比分插入缺失。這導致在第一行是 0，-1，-2，-3，-4，-5，-6，-7。這同樣適用於第二列，因為我們只具有以上現有分數。因此，我們有：

第一種情況是存在向三個方向構築矩陣，周圍位置得分如圖

對接下來的位置，我們有不同的選擇。

在某些位置，可能有兩個或者三個方向得分相同，如圖所示：

完整的得分矩陣如下圖所示 

 

 

 

 

 

 
所以根據該矩陣就可以方便的出我們想要的匹配路徑了。 

```
结果：

序列         最佳匹配
---------    ----------------------
GCATGCU      GCATG-CU      GCA-TGCU      GCAT-GCU
GATTACA      G-ATTACA      G-ATTACA      G-ATTACA

解读：

你可以从最左位置开始解读匹配

匹配:  +GCATGCU
       +GATTACA
```
分數： 0 沒有任何分數

### 追溯本源

注意： 可能有多個路徑都能得到最高的分數，在這裡我們只展示其中一種

根據指向上一個位置的箭頭來回溯取得最佳匹配的路徑。然後根據從起始位置指向最終位置的路徑來完成匹配的建立。方法如下：

* 箭頭表示匹配或者不匹配，不管哪一種，都表明這兩個位置上的字母相關聯。如以下例子所示：

```
空--> G
空 --> G
```

* 如果有一個水平的箭頭表示會有一個縱向序列的間隔。如下圖例子所示：

```
G --> GC
G --> G-
```

* 如果有一個垂直的箭頭表示會有一個橫向序列的間隔。如下圖例子所示：

```
GCA --> GCA-
G-A --> G-AT
```

* 注意可能會有多個箭頭指向同一個位置，這代表了匹配的路徑有很多不同的分支，如果這些路徑的得分相同，代表這些匹配方法實質上是等效的

根據以上規則我們可以得出其中一種可能的匹配路徑如下:

```
G →  GC → GCA → GCA- → GCA-T → GCA-TG → GCA-TGC →  **GCA-TGCU**
G →  G- → G-A → G-AT → G-ATT → G-ATTA → G-ATTAC →  **G-ATTACA**
```

## 得分體系

### 基礎得分標準

最簡單的打分體系只考慮匹配，不匹配與間隔三種情況。在新手指南中使用的打分體系是匹配 = 1，不匹配 =-1，間隔 =-1。在這種情況下序列越長，得分越低，在本體系中我們期望好的匹配取得高分。另一種打分體系如下：

* 匹配 = 0
* 間隔 = 1
* 不匹配 = 1

在本體系中得分將會代表兩序列中匹配的距離。不同的得分體系適用於不同的情況。例如，如果在你的匹配中間隔被認為有很大的負面影響，你可以使用間隔扣分很多的打分體系：

* 匹配 = 0
* 間隔 = 10
* 不匹配 = 1

### 相似度矩陣

更複雜的打分體系不僅僅考慮選擇的不同類型，也考慮每個不同的字母的影響。例如，字母 A 與 A 匹配將得到 4 分。為了表示所有字母匹配得分的情況，我們將會使用相似度矩陣來幫助計算

每個分數代表從一個字母轉換到下一個位置的所有可能的情況的得分。

不同的得分矩陣適用於不同的情況。該得分體系對於蛋白序列的匹配比對尤為重要。以下為兩種廣為引用的蛋白序列相似度矩陣

* PAM
* BLOSUM

### 間隔補償

當比對序列往往存在差距（即插入缺失），有時路數。生物學上，有很大的差距更容易，而不是多個單一缺失的發生為一體的大型刪除。因此，我們應該進球兩個小的插入缺失是有過之而無不及一個大的。簡單和常見的方式做，這是通過一個大的差距，比分開始一個新的插入缺失和一個較小的間隙延伸得分為每個擴展了插入缺失信。例如，新插入缺失可能會花費 - 5 和擴展，插入缺失可能會花費 - 1。以這種方式取向，例如：

```
GAAAAAAT
G--A-A-T
```

該序列有多種匹配方式，類似的匹配方式有：

```
GAAAAAAT
GAA----T
```

或者任意有四個間隔的匹配方式

## 進階矩陣

例如，如果相似度矩陣如下圖所示：

那麼匹配結果如下：

```
AGACTAGTTAC
CGA---GACGT
```

如果間隔懲罰因子是 - 5，那麼該匹配的得分為：

　
: _S_(A,C) + _S_(G,G) + _S_(A,A) + (3 × _d_) + _S_(G,G) + _S_(T,A) + _S_(T,C) + _S_(A,G) + _S_(C,T)

　
: = -3 + 7 + 10 - (3 × 5) + 7 + (-4) + 0 + (-1) + 0 = 1

計算序列相似程度時還應該考慮序列長度的影響。令 S (s, t) 表示兩條長度分別為 m 和 n 的序列的相似性得分，假設 S (s, t) = 99，兩條序列有 99 個字符一致。如果 m=n=100 ，則可以說這兩條序列非常相似，幾乎一樣。但是，如果 m=n=1000，則僅有 10% 的字符相同。所以，在實際序列比較時，使用相對長度的得分就更有意義，定義：
（4）
以 sim（s，t）作為衡量序列相似性的指標。
從所使用的資料結構 d i , j 本身及其計算過程來看，序列兩兩比對基本算法的空間複雜度和時間複雜度都是 O (mn)。如果兩條序列的長度相近，空間和時間複雜度就為 O (n2)。

首先將要看到如何運用動態編程查找兩個 DNA 序列的 最長公共子序列（longest common subsequence，LCS）。發現了新的基因序列的生物學家通常想知道該基因序列與其他哪個序列最相似。查找 LCS 是計算兩個序列相似程度的一種方法：LCS 越長，兩個序列越相似。
子序列中的字符與子字符串中的字符不同，它們不需要是連續的。例如，ACE 是 ABCDE 的子序列，但不是它的子字符串。請看下面兩個 DNA 序列：

　
: •S1= GCCCTAGCG

　
: •S2= GCGCAATG

* 

　
: $F_{ij}=\max(F_{i-1,j-1}+S(A_{i},B_{j}),\;F_{i,j-1}+d,\;F_{i-1,j}+d)$

計算程序如下

```
**for** i=0 **to** **length**(A)
  F(i,0) ← d*i
**for** j=0 **to** **length**(B)
  F(0,j) ← d*j
**for** i=1 **to** **length**(A)
  **for** j=1 **to** **length**(B)
  {
    Match ← F(i-1,j-1) + S(Ai, Bj)
    Delete ← F(i-1, j) + d
    Insert ← F(i, j-1) + d
    F(i,j) ← **max**(Match, Insert, Delete)
  }
```

這兩個序列的 LCS 是 GCCAG。（請注意，這僅是 一個 LCS，而不是 唯一的 LCS，因為可能存在其他長度相同的公共子序列。這種最優化問題和其他最優化問題的解可能不止一個。）

```
AlignmentA ← ""
AlignmentB ← ""
i ← **length**(A)
j ← **length**(B)
**while** (i > 0 **or** j > 0)
{
  **if** (i > 0 **and** j > 0 **and** F(i,j) == F(i-1,j-1) + S(Ai, Bj))
  {
    AlignmentA ← Ai + AlignmentA
    AlignmentB ← Bj + AlignmentB
    i ← i - 1
    j ← j - 1
  }
  **else** **if** (i > 0 **and** F(i,j) == F(i-1,j) + d)
  {
    AlignmentA ← Ai + AlignmentA
    AlignmentB ← "-" + AlignmentB
    i ← i - 1
  }
  **else**
  {
    AlignmentA ← "-" + AlignmentA
    AlignmentB ← Bj + AlignmentB
    j ← j - 1
  }
}
```

## 歷史與算法發展

通過 Needleman 和 Wunsch 的描述的算法的最初目的是找到在兩種蛋白質的胺基酸序列的相似性。

Needleman 和 Wunsch 的描述他們的算法為明確的情況下，當對準被匹配和不匹配僅僅懲罰，並差距沒有懲罰（D =0）。從 1970 年最初發布提示遞歸。 $F_{ij}=\max _{h<i,k<j}\{F_{h,j-1}+S(A_{i},B_{j}),F_{i-1,k}+S(A_{i},B_{j})\}$ .

對應的動態規劃算法需要立方時間。該報還指出，遞歸可容納任意差距處罰公式：
懲罰因子，減去對於由每個間隙，可被評估為阻擋在允許的間隙。懲罰因子可以是間隙的大小和 / 或方向的函數。 [頁 444]

一個更好的動態規劃算法二次運行時間為同樣的問題（無間隙罰款）最早是由大衛 Sankoff1972 年類似的二次時間算法推出了由 TK Vintsyuk 於 1968 年獨立發現了一種語音處理（「時間扭曲」），由羅伯特・華格納和麥可・J・菲舍爾在 1974 年的字符串匹配。

Needleman 和 Wunsch 的配製最大化相似性的問題。另一種可能性是最小化序列之間的編輯距離，弗拉基米爾的 Levenshtein 引入。彼得・塞勒斯在 1974 年發現的兩個問題是等價的。

Needleman-Wunsch 算法仍是廣泛使用的最佳全局比對，特別是當全局比對的質量即結果是最重要的。然而，該算法是相對於時間和空間的昂貴，正比於兩個序列，因而不適合長度特別長的序列比對。

全局性配對和局部性配對的比較
全局性配對和局部性配對的基本區別是在局部性配對中，你將嘗試用你訴求的序列的片段或局部來進行配對。而在全局性配對中，你將嘗試使用整個序列（從起始端到結束端）來進行配對，也即是說，在全局性配對中，當你的比對序列長度不一致時，在配對過程中可能會產生間隔或者不匹配的情況。注意，在局部性配對中也有可能產生間隔。

局部性配對

 

 

 

 

全局性配對

 

 

 

 
以下是些基於兩種算法不同點的比較。
在 Needleman-Wunsch（全局性）算法中，分數跟蹤系統是基於得分矩陣的右下角開始的。然而在 Smith-Waterman（局部性）算法中，得分矩陣的起始點是基於分數最高的坐標位置。
全局性配對主要應用於比較同源性的基因的相似性而局部性配對主要應用於在非同源性的基因序列中尋找具有相似性的同源的基因域。

全局性配對是當你將全部的序列納入考慮範圍之內後進行的序列配對，而局部性配對則是將序列的一小部分進行優先比對。以下是為了幫助理解該概念而舉得一個例子：

 
假設你有一個較長的序列作為參考，假設該序列長為 2000 鹼基。你所要配對的序列長度較短，約為 100 鹼基左右。假設在你的參考序列中有與配對序列極為接近的鹼基序列片段出現，當你進行局部性配對時，你將得到一個非常好的匹配結果和相當高的匹配分數。然而，當你進行全局性配對時，匹配程度將會很低。因為該配對方法要求將兩個序列的所有鹼基都進行配對，所以當兩序列長度相差很大時，將會造成配對結果含有很長的間隔。儘管在某處你的配對序列和參考序列的匹配程度相當高，在全局性配對時，將不會把這兩個片段直接進行配對，而是儘量嘗試讓參考序列內的所有鹼基都參與匹配，造成結果與局部性配對方法有極大差別。

 
因此，在選擇配對算法的時候，必須參考兩序列的具體長度和配對的實際意義和需要，選擇適當的方法，才能得到有生物學意義的結果。

## 程序代碼

本算法可以在多個平台上實現運算，常見的程式語言均有公開的代碼可以參考使用。
以下提供在 C 語言中算法的實現
usage: ./bin/needleman\_wunsch [OPTIONS] [seq1 seq2]

```
     Needleman-Wunsch optimal global alignment (maximises score).  
     Takes a pair of sequences on the command line, or can read from a
     file and from sequence piped in.  Can read gzip files, FASTA and FASTQ.
```

```
     OPTIONS:
       --file <file>        同时读取比对的两序列文件
       --files <f1> <f2>    
       --stdin              读取标准文件头
```

```
       --case_sensitive     使用大小写字母区分
       --match <score>      [默认: 1]
       --mismatch <score>   [默认: -2]
       --gapopen <score>    [默认: -4]
       --gapextend <score>  [默认: -1]
```

```
       --scoring <PAM30|PAM70|BLOSUM80|BLOSUM62>
       --substitution_matrix <file>  
       --substitution_pairs <file>   
```

 

```
       freestartgap       间隔不扣分
       freeendgap         序列结束不扣分
```

```
       printscores        打印优化后的比对分数
       zam               
       printmatrices      打印动态矩阵
       printfasta         打印fasta格式的文件头
       pretty             打印描述行
       colour             打印颜色
```

```
     可选优化项：
       nogapsin1          被比对序列不允许有间隔
       nogapsin2          比对序列不允许有间隔
       nogaps             所有序列均不允许有间隔
       nomismatches       不允许有误差配对
```

## 使用 Needleman-Wunsch 算法的工具

* EMBOSS Needle and EMBOSS Stretcher Global Alignment Tools （頁面存檔備份，存於網際網路檔案館）
* Needleman-Wunsch alignment for two nucleotide sequences （頁面存檔備份，存於網際網路檔案館）
* MathWorks - Globally align two sequences using Needleman-Wunsch algorithm （頁面存檔備份，存於網際網路檔案館）

應用 —— 如何使用該算法
生物信息學國家中心網站提供了包括各種生物信息學計算工具的優秀資源，其中就有進行全局性配對計算的網頁連結。

NCBI 全局性配對

1 進入 NCBI 主頁 (ncbi.nlm.nih.gov)

2 尋找 Popular Resources, 點擊 Blast

3. 選擇 Needleman-Wunsch 算法
你可以選擇自行輸入任意序列對，下滑頁面至底部並點擊 Align，比對的結果將提供鹼基的配對結果，相似的百分比和最終得分等信息。

或者， 你可以選擇搜索基因名稱來選擇要比對的序列。

1 進入 NCBI 主頁
2 尋找到 Popular Resources， 點擊 Gene
3 輸入你要搜索的基因名稱和種族名稱
4 點擊結果中顯示的基因名稱
5 選擇 Needleman-Wunsch 算法並運行

## 生物信息學以外的應用

### 三維立體畫

三維立體畫是利用人眼立體視覺現象製作的繪畫作品。普通繪畫和攝影作品，包括電腦製作的三維動畫，只是運用了人眼對光影、明暗、虛實的感覺得到立體的感覺，而沒有利用雙眼的立體視覺，一隻眼看和兩隻眼看都是一樣的。充分利用雙眼立體視覺的立體畫，將使你看到一個精彩的世界。人有兩隻眼，兩隻眼有一定距離，這就造成物體的影象在兩眼中有一些差異，見右圖，由圖可見，由於物體與眼的距離不同，兩眼的視角會有所不同，由於視角的不同所看到是影象也會有一些差異，大腦會根據這種差異感覺到立體的景象。三維立體畫就是利用這個原理，在水平方向生成一系列重複的圖案，當這些圖案在兩隻眼中重合時，就看到了立體的影象。參見下圖，這是一幅不能再簡單的立體畫了。圖中最上一行圓最遠，最下一行圓最近，請注意：最上一行圓之間距離最大，最下一行圓之間距離最小。重複圖案的距離決定了立體影象的遠近，生成三維立體畫的程序就是根據這個原理，依據三維影象的遠近，生成不同距離的重複圖案。
立體畫有兩種形式：第一種是由相同的圖案在水平方向以不同間隔排列而成，看起來是遠近不同的物體。　另一種立體畫較複雜，在這種立體畫上你不能直接看到物體的形象，畫面上只有雜亂的圖案，製作這樣的立體畫只有使用程序了。兩種作品看法是一樣的，原理都是使左眼看到左眼的影象，讓右眼看到右眼的影象。

 

### 電腦立體視覺效果

立體匹配是從一對立體圖像的三維重建進程的一個重要步驟。當圖像已被糾正，打個比方可以對準核苷酸 / 蛋白序列和匹配像素屬於掃描線，因為這兩個任務，著眼於人物的兩個字符串之間建立最佳的對應關係繪製。實際上，一對立體聲的「正確」的圖像可以被看作是「左」圖像的突變版本：噪聲和單個攝像機的靈敏度改變的像素值（即字符取代）; 和不同的視角揭示以前閉塞的數據，並引入新的閉塞（即插入和字符刪除）。作為結果，中的 Needleman-Wunsch 算法的小的修改，使其適合於立體匹配。 儘管該算法在準確性方面性能都不算是最先進的，該但算法的相對簡單性允許其在嵌入式系統上執行，使其有廣泛的應用.

.

## 其他連接

* Smith-Waterman 算法
* 序列比較

## 參考文獻

## 外部連結

* NW-align: A protein sequence-to-sequence alignment program by Needleman-Wunsch algorithm (online server & source code) （頁面存檔備份，存於網際網路檔案館）
* Needleman-Wunsch Algorithm as Haskell Code （頁面存檔備份，存於網際網路檔案館）
* A live Javascript-based demo of Needleman-Wunsch （頁面存檔備份，存於網際網路檔案館）
* An interactive Javascript-based visual explanation of Needleman-Wunsch Algorithm Archive.is 的存檔，存檔日期 2016-04-27
* B.A.B.A. （頁面存檔備份，存於網際網路檔案館） — an applet (with source) which visually explains the algorithm.
* A clear explanation of NW and its applications to sequence alignment
* Sequence Alignment Techniques at Technology Blog （頁面存檔備份，存於網際網路檔案館）
* Biostrings （頁面存檔備份，存於網際網路檔案館） R package implementing Needleman–Wunsch algorithm among others
* Parallel Needleman-Wunsch Algorithm for Grid Implementation by Tahir Naveed, Imitaz Saeed Siddiqui and Shaftab Ahmed - Bahria University
* Needleman-Wunsch Algorithm as Haxe Code （頁面存檔備份，存於網際網路檔案館）