分享

【串和序列处理 6】LCS 最长公共子序列

 静听沙漏 2012-01-08

LCS:又称最长公共子序列。其中子序列(subsequence)的概念不同于串的子串。它是一个不一定连续但按顺序取自字符串X中的字符序列。例如:串"AAAG"就是串“CGATAATTGAGA”的一个子序列。

字符串的相似性问题可以通过求解两个串间的最长公共子序列(LCS)来得到。当然如果使用穷举算法列出串的所有子序列,一共有2^n种,而每个子序列是否是另外一个串的子序列又需要O(m)的时间复杂度,因此这个穷举的方法时间复杂度是O(m*(2^n))指数级别,效率相当的可怕。我们采用动态规划算法来解决这个问题。

动态规划算法解决最长公共子序列

假如我们有两个字符串:X=[0,1,2....n] Y=[0,1,2...m]。我们定义L(i, j)为X[0...i]与Y[0...j]之间的最长公共子序列的长度。通过动态规划思想(复杂问题的最优解是子问题的最优解和子问题的重叠性质决定的)。我们考虑这样两种情况:

(1) 当X[i]=Y[j]时, L(i, j)=L(i-1, j-1)+1。证明很简单。

(2) 当X[i]!=Y[j]时,说明此事X[0...i]和Y[0...j]的最长公共子序列中绝对不可能同时含有X[i]和Y[j]。那么公共子序列可能以X[i]结尾,可能以Y[j]结尾,可以末尾都不含有X[i]或Y[j]。因此

L(i, j)= MAX{L(i-1 , j), L(i, j-1)}

LCS动态规划Java源代码

Java代码
  1. package net.hr.algorithm.string;
  2. /**
  3. * 最长公共子串LCS
  4. * @author heartraid
  5. */
  6. public class LCS {
  7. /**字符串X的字符数组*/
  8. private char[] charArrayX=null;
  9. /**字符串Y的字符数组*/
  10. private char[] charArrayY=null;
  11. public LCS(String sa,String sb){
  12. charArrayX=new char[sa.length()+1];
  13. System.arraycopy(sa.toCharArray(),0,charArrayX,1,sa.length());
  14. charArrayY=new char[sb.length()+1];
  15. System.arraycopy(sb.toCharArray(),0,charArrayY,1,sb.length());
  16. }
  17. /**
  18. * 得到最长公共子序列的长度
  19. */
  20. public void getLCS(){
  21. int[][] length=new int[charArrayX.length+1][charArrayY.length+1];
  22. for(int m=1;m<charArrayX.length;m++)
  23. for(int n=1;n<charArrayY.length;n++){
  24. if(charArrayX[m]==charArrayY[n]){
  25. length[m][n]=length[m-1][n-1]+1;
  26. }
  27. else
  28. length[m][n]=max(length[m-1][n],length[m][n-1]);
  29. }
  30. //打印最长公共子序列
  31. String lcstr="";
  32. int x=charArrayX.length-1;
  33. int y=charArrayY.length-1;
  34. while(x>=1&&y>=1){
  35. if(charArrayX[x]==charArrayY[y]){
  36. lcstr=charArrayX[x]+lcstr;
  37. x--;
  38. y--;
  39. }else{
  40. if(length[x-1][y]<=length[x][y-1])
  41. y--;
  42. else
  43. x--;
  44. }
  45. }
  46. System.out.println("最长公共子序列为:"+lcstr+" [length="+lcstr.length()+"]");
  47. }
  48. /**
  49. * 取最大值
  50. */
  51. private int max(int m,int n){
  52. return m>n?m:n;
  53. }
  54. /**
  55. * 测试
  56. */
  57. public static void main(String[] args) {
  58. LCS lcs=new LCS("GTTCCTAATA","CGATAATTGAGA");
  59. lcs.getLCS();
  60. }
  61. }

这里解释一下上面的代码,其中getLength()的作用是递归获取最长公共子串的长度,并得到递归过程中每个子串之间最长公共子串长度的状态表lcs[][],这张表运行的结果如下:

C G A T A A T T G A G A

0 0 0 0 0 0 0 0 0 0 0 0 0
G 0 0 1 1 1 1 1 1 1 1 1 1 1
T 0 0 1 1 2 2 2 2 2 2 2 2 2
T 0 0 1 1 2 2 2 3 3 3 3 3 3
C 0 1 1 1 2 2 2 3 3 3 3 3 3
C 0 1 1 1 2 2 2 3 3 3 3 3 3
T 0 1 1 1 2 2 2 3 4 4 4 4 4
A 0 1 1 2 2 3 3 3 4 4 5 5 5
A 0 1 1 2 2 3 4 4 4 4 5 5 6
T 0 1 1 2 3 3 4 5 5 5 5 5 6
A 0 1 1 2 3 4 4 5 5 5 6 6 6

红色数字的位置就揭示了最长公共子序列: CTATTA。也就是说分析这张表就可以得到最长公共子序列字符串。

为什么呢?因为通过表的回溯过程,从后向前重构了一个最长公共子序列。对于任何位置lcs[i][j],确定是否X[i]=Y[j]。如果是,那么X[i]必是最长公共子序列的一个字符。如果否,那么移动到lcs[i,j-1]和lcs[i-1, j]之间的较大者。

动态规划方法LCS效率:

动态规划方法构造最长公共子序列需要O(m*n)的代价,另外,如果想要得到最长公共子序列,又需要O(m+n)的时间来读取csl[][]数组。尽管如此,其时间复杂度仍然比蛮力穷举的指数级别要强的多。

问题拓展:设A,B,C是三个长为n的字符串,它们取自同一常数大小的字母表。设计一个找出三个串的最长公共子串的O(n^3)的时间算法。(来自《Algorithm Design》(中文版:算法分析与设计) - Chapter9 - 文本处理 - 创新题C-9.18)

分析:可以通过《LCS最长公共子序列》的动态规划算法,设计Java源代码如下:

Java代码 复制代码 收藏代码
  1. public class TriLCS{
  2. char[] charA=null;
  3. char[] charB=null;
  4. char[] charC=null;
  5. public TriLCS(String sa,String sb,String sc){
  6. charA=new char[sa.length()+1];
  7. System.arraycopy(sa.toCharArray(),0,charA,1,sa.length());
  8. charB=new char[sb.length()+1];
  9. System.arraycopy(sb.toCharArray(),0,charB,1,sb.length());
  10. charC=new char[sc.length()+1];
  11. System.arraycopy(sc.toCharArray(),0,charC,1,sc.length());
  12. }
  13. public void getTriLCS(){
  14. int[][][] length=new int[charA.length][charB.length][charC.length];
  15. for(int a=1;a<charA.length;a++)
  16. for(int b=1;b<charB.length;b++)
  17. for(int c=1;c<charC.length;c++){
  18. if(charA[a]==charB[b]&&charA[a]==charC[c]){
  19. length[a][b][c]=length[a-1][b-1][c-1]+1;
  20. }
  21. else if(charA[a]==charB[b]&&charA[a]!=charC[c]){
  22. length[a][b][c]=max(length[a-1][b-1][c],length[a][b][c-1]);
  23. }
  24. else if(charA[a]==charC[c]&&charA[a]!=charB[b]){
  25. length[a][b][c]=max(length[a-1][b][c-1],length[a][b-1][c]);
  26. }
  27. else if(charB[b]==charC[c]&&charA[a]!=charB[b]){
  28. length[a][b][c]=max(length[a][b-1][c-1],length[a-1][b][c]);
  29. }
  30. else{
  31. length[a][b][c]=max(length[a-1][b][c],length[a][b-1][c],length[a][b][c-1],length[a-1][b-1][c],length[a-1][b][c-1],length[a][b-1][c-1]);
  32. }
  33. }
  34. //打印最长公共子序列
  35. String lcstr="";
  36. int a=charA.length-1;
  37. int b=charB.length-1;
  38. int c=charC.length-1;
  39. while(a>=1&&b>=1&&c>=1){
  40. if(charA[a]==charB[b]&&charA[a]==charC[c]){
  41. lcstr=charA[a]+lcstr;
  42. a--;
  43. b--;
  44. c--;
  45. }
  46. else if(charA[a]==charB[b]&&charA[a]!=charC[c]){
  47. if(length[a-1][b-1][c]<=length[a][b][c-1])
  48. c--;
  49. else{
  50. a--;
  51. b--;
  52. }
  53. }
  54. else if(charA[a]==charC[c]&&charA[a]!=charB[b]){
  55. if(length[a-1][b][c-1]<=length[a][b-1][c])
  56. b--;
  57. else{
  58. a--;
  59. c--;
  60. }
  61. }
  62. else if(charB[b]==charC[c]&&charA[a]!=charB[b]){
  63. if(length[a][b-1][c-1]<=length[a-1][b][c])
  64. a--;
  65. else{
  66. b--;
  67. c--;
  68. }
  69. }
  70. else{
  71. int maxize=max(length[a-1][b][c],length[a][b-1][c],length[a][b][c-1],length[a-1][b-1][c],length[a-1][b][c-1],length[a][b-1][c-1]);
  72. if(maxize==length[a-1][b][c])
  73. a--;
  74. if(maxize==length[a][b-1][c])
  75. b--;
  76. if(maxize==length[a][b][c-1])
  77. c--;
  78. if(maxize==length[a-1][b-1][c]){
  79. a--;
  80. b--;
  81. }
  82. if(maxize==length[a-1][b][c-1]){
  83. a--;
  84. c--;
  85. }
  86. if(maxize==length[a][b-1][c-1]){
  87. b--;
  88. c--;
  89. }
  90. }
  91. }
  92. System.out.println("最长子串为:"+lcstr+"(length="+length[charA.length-1][charB.length-1][charC.length-1]+")");
  93. }
  94. /**
  95. * 取最大值
  96. */
  97. private int max(int m,int n){
  98. return m>n?m:n;
  99. }
  100. /**
  101. * 取最大值
  102. */
  103. private int max(int x,int y,int z,int k,int m,int n){
  104. int maxizen=0;
  105. if(maxizen<x) maxizen=x;
  106. if(maxizen<y) maxizen=y;
  107. if(maxizen<z) maxizen=z;
  108. if(maxizen<k) maxizen=k;
  109. if(maxizen<m) maxizen=m;
  110. if(maxizen<n) maxizen=n;
  111. return maxizen;
  112. }
  113. public static void main(String[] args){
  114. TriLCS tri=new TriLCS("aadsbbbcs","adsabcs","adbsbsdcs");
  115. tri.getTriLCS();
  116. }
  117. }

    本站是提供个人知识管理的网络存储空间,所有内容均由用户发布,不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息,谨防诈骗。如发现有害或侵权内容,请点击一键举报。
    转藏 分享 献花(0

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多