最长公共子序列

最大子序列、最长递增子序列、最长公共子串、最长公共子序列、字符串编辑距离

最长公共子序列

最长公共子序列与最长公共子串的区别在于最长公共子序列不要求在原字符串中是连续的，比如ADE和ABCDE的最长公共子序列是ADE。

我们用动态规划的方法来思考这个问题如是求解。首先要找到状态转移方程：

符号约定，C1是S1的最右侧字符，C2是S2的最右侧字符，S1‘是从S1中去除C1的部分，S2'是从S2中去除C2的部分。

LCS(S1,S2)等于下列3项的最大者：

（1）LCS（S1，S2’）

（2）LCS（S1’，S2）

（3）LCS（S1’，S2’）--如果C1不等于C2； LCS（S1'，S2'）+C1--如果C1等于C2；

边界终止条件：如果S1和S2都是空串，则结果也是空串。

下面我们同样要构建一个矩阵来存储动态规划过程中子问题的解。这个矩阵中的每个数字代表了该行和该列之前的LCS的长度。与上面刚刚分析出的状态转移议程相对应，矩阵中每个格子里的数字应该这么填，它等于以下3项的最大值：

（1）上面一个格子里的数字

（2）左边一个格子里的数字

（3）左上角那个格子里的数字（如果 C1不等于C2）； 左上角那个格子里的数字+1（ 如果C1等于C2）

举个例子：

 　　　　  G　　C　　T　　A

　　 0　　0　　0　　0　　0

G　 0　　1　　1　　1　　1

B　 0　　1　　1　　1　　1

T　 0　　1　　1　　2　　2

A    0　　1　　1　　2　　3

填写最后一个数字时，它应该是下面三个的最大者：

（1）上边的数字2

（2）左边的数字2

（3）左上角的数字2+1=3,因为此时C1==C2

所以最终结果是3。

在填写过程中我们还是记录下当前单元格的数字来自于哪个单元格，以方便最后我们回溯找出最长公共子串。有时候左上、左、上三者中有多个同时达到最大，那么任取其中之一，但是在整个过程中你必须遵循固定的优先标准。在我的代码中优先级别是左上>左>上。

下图给出了回溯法找出LCS的过程：

奉上代码：

 

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78

#include<iostream> #include<cstring> #include<stack> #include<utility> #define LEFTUP 0 #define LEFT 1 #define UP 2 using  namespace  std; int  Max( int  a, int  b, int  c, int  *max){            //找最大者时a的优先级别最高，c的最低.最大值保存在*max中      int  res=0;          //res记录来自于哪个单元格      *max=a;      if (b>*max){          *max=b;          res=1;      }      if (c>*max){          *max=c;          res=2;      }      return  res; } //调用此函数时请注意把较长的字符串赋给str1,这主要是为了在回溯最长子序列时节省时间。如果没有把较长的字符串赋给str1不影响程序的正确执行。 string LCS( const  string &str1, const  string &str2){      int  xlen=str1.size();               //横向长度      int  ylen=str2.size();               //纵向长度      if (xlen==0||ylen==0)                //str1和str2中只要有一个为空，则返回空          return  "" ;      pair< int , int > arr[ylen+1][xlen+1];    //构造pair二维数组，first记录数据，second记录来源      for ( int  i=0;i<=xlen;i++)         //首行清0          arr[0][i].first=0;      for ( int  j=0;j<=ylen;j++)         //首列清0          arr[j][0].first=0;      for ( int  i=1;i<=ylen;i++){          char  s=str2.at(i-1);          for ( int  j=1;j<=xlen;j++){              int  leftup=arr[i-1][j-1].first;              int  left=arr[i][j-1].first;              int  up=arr[i-1][j].first;              if (str1.at(j-1)==s)         //C1==C2                  leftup++;              int  max;              arr[i][j].second=Max(leftup,left,up,&arr[i][j].first); //          cout<<arr[i][j].first<<arr[i][j].second<<" ";          } //      cout<<endl;      }       /*矩阵构造完毕*/      //回溯找出最长公共子序列      stack< int > st;      int  i=ylen,j=xlen;      while (i>=0&&j>=0){          if (arr[i][j].second==LEFTUP){              if (arr[i][j].first==arr[i-1][j-1].first+1)                  st.push(i);              --i;              --j;          }          else  if (arr[i][j].second==LEFT){              --j;          }          else  if (arr[i][j].second==UP){              --i;          }      }      string res= "" ;      while (!st.empty()){          int  index=st.top()-1;          res.append(str2.substr(index,1));          st.pop();      }      return  res; } int  main(){      string str1= "GCCCTAGCG" ;      string str2= "GCGCAATG" ;      string lcs=LCS(str1,str2);      cout<<lcs<<endl;      return  0; }

下面给一个Java版本

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31

public  static  <E> List<E> longestCommonSubsequence(E[] s1,E[] s2){          int [][] num= new  int [s1.length+ 1 ][s2.length+ 1 ];          for ( int  i= 1 ;i<s1.length+ 1 ;i++){              for ( int  j= 1 ;j<s2.length+ 1 ;j++){                  if (s1[i- 1 ].equals(s2[j- 1 ])){                      num[i][j]= 1 +num[i- 1 ][j- 1 ];                  }                  else {                      num[i][j]=Math.max(num[i- 1 ][j],num[i][j- 1 ]);                  }              }          }          System.out.println( "lenght of LCS= " +num[s1.length][s2.length]);          int  s1position=s1.length,s2position=s2.length;          List<E> result= new  LinkedList<E>();          while (s1position> 0  && s2position> 0 ){              if (s1[s1position- 1 ].equals(s2[s2position- 1 ])){                  result.add(s1[s1position- 1 ]);                  s1position--;                  s2position--;              }              else  if (num[s1position][s2position- 1 ]>=num[s1position- 1 ][s2position]){                  s2position--;              }              else {                  s1position--;              }          }          Collections.reverse(result);          return  result;      }

std::endl是一个特殊值，称为操纵符（manipulator），将它写入输出流时具有输出换行的效果，并刷新与设备相关联的缓冲区(buffer)。通过刷新缓冲区用户可立即看到写入到流中的输出。

我在调试以上代码时在45行（cout<<endl;）处设置断点，结果发现“43行（cout<<arr[i][j].first<<arr[i][j].second<<" ";） 没有执行”，这就是因为43行末尾没有加endl，所以用户没有立即看到输出流中的数据。

字符串编辑距离

 要想把字符串S1变成S2，可以经过若干次下列原子操作：

1.删除一个字符

2.增加一个字符

3.更改一个字符

字符串S1和S2的编辑距离定义为从S1变成S2所需要原子操作的最少次数。

解法跟上面的最长公共子序列十分相似，都是动态规划，把一个问题转换为若干个规模更小的子问题，并且都借助于一个二维矩阵来实现计算。

约定：字符串S去掉最后一个字符T后为S'，T1和T2分别是S1和S2的最后一个字符。

则dist(S1,S2)是下列4个值的最小者：

1.dist(S1',S2')--当T1==T2

2.1+dist(S1',S2)--当T1!=T2，并且删除S1的最后一个字符T1

3.1+dist(S1,S2')--当T1!=T2，并且在S1后面增加一个字符T2

4.1+dist(S1',S2')--当T1!=T2，并且把S1的最的一个字符T1改成T2

 

把问题转换为二维矩阵：

arr[i][j]表示S1.sub(0,i)和S2.sub(0,j)的编辑距离，则

arr[i][j]=min{1+arr[i][j-1], 1+arr[i-1][j], 1+arr[i-1][j-1](当S1[i]!=S2[j]), arr[i-1][j-1](当S1[i]==S2[j])}

边界情况：arr[0][j]=j, arr[i][0]=i

代码就不写了，跟最长公共子序列很相似。

 

计算两个字条串的相似度除了Edit Distance，还有一种方法是计算Jaro Distance。具体怎么算读者可以搜一下。