LCS兩個字串最長公共子串
LCS問題就是求兩個字串最長公共子串的問題。解法就是用一個矩陣來記錄兩個字串中所有位置的兩個字元之間的匹配情況，若是匹配則為1，否則為0。然後求出對角線最長的1序列，其對應的位置就是最長匹配子串的位置。
下面是字串21232523311324和字串312123223445的匹配矩陣，前者為X方向的，後者為Y方向的。不難找到，紅色部分是最長的匹配子串。通過查詢位置我們得到最長的匹配子串為：21232
0 0 0 1 0 0 0 1 1 0 0 1 0 0 0
0 1 0 0 0 0 0 0 0 1 1 0 0 0 0
1 0 1 0 1 0 1 0 0 0 0 0 1 0 0
0 1 0 0 0 0 0 0 0 1 1 0 0 0 0
1 0 1 0 1 0 1 0 0 0 0 0 1 0 0
0 0 0 1 0 0 0 1 1 0 0 1 0 0 0
1 0 1 0 1 0 1 0 0 0 0 0 1 0 0
1 0 1 0 1 0 1 0 0 0 0 0 1 0 0
0 0 0 1 0 0 0 1 1 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
但是在0和1的矩陣中找最長的1對角線序列又要花去一定的時間。通過改進矩陣的生成方式和設定標記變數，可以省去這部分時間。下面是新的矩陣生成方式：
0 0 0 1 0 0 0 1 1 0 0 1 0 0 0
0 1 0 0 0 0 0 0 0 2 1 0 0 0 0
1 0 2 0 1 0 1 0 0 0 0 0 1 0 0
0 2 0 0 0 0 0 0 0 1 1 0 0 0 0
1 0 3 0 1 0 1 0 0 0 0 0 1 0 0
0 0 0 4 0 0 0 2 1 0 0 1 0 0 0
1 0 1 0 5 0 1 0 0 0 0 0 2 0 0
1 0 1 0 1 0 1 0 0 0 0 0 1 0 0
0 0 0 2 0 0 0 2 1 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 1 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
不用多說，你大概已經看出來了。當字元匹配的時候，我們並不是簡單的給相應元素賦上1，而是賦上其左上角元素的值加一。我們用兩個標記變數來標記矩陣中值最大的元素的位置，在矩陣生成的過程中來判斷當前生成的元素的值是不是最大的，據此來改變標記變數的值，那麼到矩陣完成的時候，最長匹配子串的位置和長度就已經出來了。
這樣做速度比較快，但是花的空間太多。我們注意到在改進的矩陣生成方式當中，每生成一行，前面的那一行就已經沒有用了。因此我們只需使用一維陣列即可。最終的程式碼如下：
Private Function LCS(ByVal str_1 As String, ByVal str_2 As String) As String
If str_1 = "" Or str_2 = "" Then Return ""
Dim c(str_1.Length) As Integer
Dim max, maxj, i, j As Integer
maxj = 0 : max = 0 '這兩個是標誌變數
For i = 0 To str_2.Length - 1
For j = str_1.Length - 1 To 0 Step -1
If str_2.Chars(i) = str_1.Chars(j) Then
If i = 0 Or j = 0 Then
c(j) = 1
Else
c(j) = c(j - 1) + 1
End If
Else
c(j) = 0
End If
If c(j) > max Then '把>改成>=則返回最後一個最長匹配子串
max = c(j) : maxj = j '更新標誌變數
End If
Next
Next
If max = 0 Then Return ""
Return str_1.Substring(maxj - max + 1, max) '直接從標誌變數得出結果
End Function首先將要看到如何運用動態程式設計查詢兩個 DNA 序列的最長公共子序列（longest common subsequence，LCS）。發現了新的基因序列的生物學家通常想知道該基因序列與其他哪個序列最相似。查詢 LCS 是計算兩個序列相似程度的一種方法：LCS 越長，兩個序列越相似。
子序列中的字元與子字串中的字元不同，它們不需要是連續的。例如，ACE 是 ABCDE 的子序列，但不是它的子字串。請看下面兩個 DNA 序列：
- S1 = DE>GCCCTAGCGDE>
- S2 = DE>GCGCAATGDE>
這兩個序列的 LCS 是 GCCAG。（請注意，這僅是一個 LCS，而不是唯一的
首先，考慮如何遞迴地計算 LCS。令：
- C1 是 S1 最右側的字元
- C2 是 S2 最右側的字元
- S1' 是 S1 中 “切掉” C1 的部分
- S2' 是 S2 中 “切掉” C2 的部分
有三個遞迴子問題：
- L1 = LCS(S1', S2)
- L2 = LCS(S1, S2')
- L3 = LCS(S1', S2')
結果表明（而且很容易使人相信）原始問題的解就是下面三個子序列中最長的一個：
- L1
- L2
- 如果 C1 等於 C2，則為 L3 後端加上 C1 ，如果 C1 不等於 C2，則為 L3。
（基線條件（base case）是 S1 或 S2 為長度為零的字串的情形。在這種情況下，S1 和 S2 的 LCS 顯然是長度為零的字串。）
但是，就像計算斐波納契數的遞迴過程一樣，這個遞迴解需要多次計算相同的子問題。可以證明，這種遞迴解法需要耗費指數級的時間。相比之下，這一問題的動態程式設計解法的執行時間是 Θ(mn)，其中 m 和 n 分別是兩個序列的長度。
為了用動態程式設計有效地計算 LCS，首先需要構建一個表格，用它儲存部分結果。沿著頂部列出一個序列，再沿著左側從上到下列出另一個序列，如圖 2 所示：
圖 2. 初始 LCS 表格
這種方法的思路是：將從上向下、從左到右填充表格，每個單元格包含一個數字，代表該行和該列之前的兩個字串的 LCS 的長度。也就是說，每個單元格包含原始問題的一個字問題的解。例如，請看第 6 行第 7 列的單元格：它在 GCGCAATG 序列第二個 C 的右側，在 GCCCTAGCG 的 T 的下面。這個單元格最終包含的數字就是 GCGC 和 GCCCT 的 LCS 的長度。
首先看一下表格的第二行中應該是什麼條目。這一行的單元格儲存的 LCS 長度對應的是序列 GCGCAATA 的零長前端和序列 GCCCTAGCG 的 LCS。顯然，這些 LCS 的值都是 0。類似的，沿著第二列向下的值也都是 0，這與遞迴解的基線條件對應。現在表格如圖 3 所示：
圖 3.填充了基線條件的 LCS 表格
接下來，要實現與遞迴演算法中遞迴子問題對應的情景，但這時使用的是表格中已經填充的值。在圖 4 中，我已經填充了一半左右的單元格：
圖 4.填充了一半的 LCS 表格
在填充單元格時，需要考慮以下條件：
- 它左側的單元格
- 它上面的單元格
- 它左上側的單元格
下面三個值分別對應著我在前面列出的三個遞迴子問題返回的值。
- V1 = 左側單元格的值
- V2 = 上面單元格的值
- V3 = 左上側單元格的值
在空單元格中填充下面 3 個數字中的最大值：
- V1
- V2
- 如果 C1 等於 C2 則為 V3 + 1，如果 C1 不等於 C2，則為 V3 ，其中 C1 是當前單元格上面的字元，C2 是當前單元格左側的字元
請注意，我在圖中還添加了箭頭，指向當前單元格值的來源。後面的 “回溯” 一節將用這些箭頭建立實際的 LCS（與僅僅發現 LCS 長度相反）。
現在填充圖 4 中接下來的空單元格 — 在 GCCCTAGCG 中第三個 C 下面和 GCGCAATG 第二個 C 的右側的單元格。它上面的值是 2，左側的值是 3，左上側的值是 2。這個單元格上面的字元和左側的字元相等（都是 C），所以必須選擇 2、3 和 3（左上側單元格中的 2 + 1）的最大值。所以，這個單元格的值為 3。繪製一個箭頭，從該單元格指向其中的值的源單元格。在這個示例中，新的數值可能來自不止一個單元格，所以可以任選一個：例如左上側單元格。
作為練習，您可以嘗試填充表格的餘下部分。如果在關聯過程中，一直按照左上側-上側-左側的順序選擇單元格，那麼會得到如圖 5 所示的表格。（當然，如果在關聯過程中做了不同的選擇，那麼箭頭就會不同，但是數字是相同的。）
圖 5.填充好的 LCS 表格
回想一下，任何單元格中的數字都是該單元格所在行之上和列之前的字串的 LCS 長度。所以，表格右下角的數字就是字串 S1 和 S2 （在本例中是 GCCCTAGCG 和 GCGCAATG）的 LCS 的長度。所以，這兩個序列的 LCS 長度是 5。
這是在所有動態程式設計演算法中都要牢記的關鍵點。表格中的每個單元格都包含單元格所在行上面和所在列左側序列前端問題的解。
接下來要做的就是尋找實際的 LCS。使用單元格箭頭進行回溯可以完成。在構建表格的時候，請記住，如果箭頭指向左上側的單元格，那麼當前單元格中的值要比左上側單元格的值大 1，這意味著左側單元格和上面單元格中的字元相等。構建 LCS 時，這會將相應的字元新增到 LCS 中。所以，構建 LCS 的途徑就是從右下角的單元格開始，沿著箭頭一路返回。每當沿著對角箭頭回到左上角的單元格而且 該單元格的值比當前單元格的值小 1 時，就要將對應的公共字元新增到 正在構建的 LCS 的前端。請注意，之所以將字元放在 LCS 前端，是因為我們是從 LCS 末端開始的。（在 圖 5 的示例中，右下角的 5 與要新增的第 5 個字元對應。）
依此類推，繼續構建 LCS。從右下側的單元格開始，看到單元格指標指向左上側的單元格，而且當前單元格的值（5）比其左上側單元格的值（4）大 1。所以將字元 G 新增到最初的零長度的字串之前。下一個箭頭，從包含 4 的單元格開始，也指向左上側，但是值沒有變。接著這個箭頭也是如此。下一個單元格的箭頭還是指向左上側，但是這次值從 3 變為 4。這意味著需要將這一行和這一列中的公共字元 A 新增到 LCS 中。所以現在的 LCS 是 AG。接下來，沿著指標向左（對應著跳過上面的 T）到達另一個 3。然後有一個對角指標指向 2。因此，又添加了在當前行和當前列中的公共字元 C，生成的 LCS 為 CAG。繼續使用這種方式，直到最後到達 0。圖 6 顯示了整個回溯過程：
圖 6.在填滿的 LCS 表格上進行回溯
通過這種回溯方法，得到的 LCS 為 GCCAG
python的程式碼實現：
def lcs_similar(x, y): x = unicode(x, 'utf-8', 'ignore') y = unicode(y, 'utf-8', 'ignore')
m = len(x) + 1 n = len(y) + 1 if m == 1 and n == 1: return 100
elif m==1 or n == 1: return 0
lcs_len = [[0 for a in range(n)] for b in range(m)]
for i in range(1, m): for j in range(1, n): if(x[i-1] == y[j-1]): lcs_len[i][j] = lcs_len[i-1][j-1] + 1 else: max_len = lcs_len[i-1][j] if lcs_len[i][j-1] > lcs_len[i-1][j]: max_len = lcs_len[i][j-1] lcs_len[i][j] = max_len i = m - 2 j = n - 2 common_substr = u'' while True: if i == -1 or j == -1: break if x[i] == y[j]: common_substr = u"%s%s" % (x[i], common_substr) i = i - 1 j = j -1 else: if lcs_len[i][j+1] > lcs_len[i+1][j]: i = i-1 else: j = j-1 print common_substr //////////////////////////////////////////////////////////////////
新增一個C++版本的：
string LCS_2(string s1,string s2) { if(s1==""||s2=="") return ""; int m=s1.size()+1; int n=s2.size()+1; int lcs[m][n],i,j; for(i=0;i<m;i++) for(j=0;j<n;j++) lcs[i][j]=0; for(i=1;i<m;i++) for(j=1;j<n;j++) { if(s1[i-1]==s2[j-1]) lcs[i][j]=lsc[i-1][j-1]+1; else lcs[i][j]=lcs[i-1][j]>=lcs[i][j-1]?lcs[i-1][j]:lcs[i][j-1];//取上面或左側最大值 } i=m-2; j=n-2; string ss="" while(i!=-1&&j!=-1) { if(x[i]==y[j]) { ss+=s[i]; i--; j--; } else { if(lcs[i+1][j+1]==lcs[i][j]) { i--; j--; } else { if(lcs[i][j+1]>=lcs[i+1][j]) i--; else j--; } } } reverse(ss);//逆轉ss return ss; }