动态规划——最长公共子序列LCS

【算法导论 ch15-4】

首先确定什么是最长公共子序列：

最长公共子串（longest common substring） ：其中的子串是原串的一个连续的部分；

最长公共子序列（longest common subsequence）：是不改变序列的顺序，从序列中得到任意的元素从而获得的新序列

二者区别=>子串中字符位置必须连续，子序列可以不必连续

（由于序列元素顺序的不同，得到的LCS可能不止一个）

用动态规划方法求解，需要先寻找最优子结构

定义c[i,j]为序列Xi和Yj的LCS长度，递归求解，将其保存到二维数组C中。可得递归式：

这样得到的表中最后一个值就是Xm和Yn的LCS；

还需要一个二维数组b，用于保存c[i,j]中的值是从c[i-1,j-1],c[i-1,j]还是c[i,j-1]中计算得到的，便于反向查找LCS元素时找到元素位置

伪代码：

便于反向查找LCS元素时找到元素位置

算导上的数组示意图：

构造LCS序列时，用b数组，从b[m,n]开始，延箭头反向追踪下去，每当遇到斜向尽头都证明该处对应的Xi和Yj是相同的，且是属于LCS的。一直追踪直到i或j到达最小。

为了保存下来每一个LCS，增加了一个数组用于保存结果，完整代码如下：

#include <stdio.h>#include <string>#include <math.h>using namespace std;/*可以考虑不适用b数组，直接从c[i-1,j-1],c[i-1,j],c[i,j-1]中计算得到c[i,j]运算方向*/const string str_x = "BDCABA";//"ACCGGTCGAGTGCGCGGAAGCCGGCCGAA";//const string str_y = "ABCBDAB";//"GTCGTTCGGAATGCCGTTGCTCTGTAAA";//int **c;int **b;char *z; /*建立的二维数组c和b的维数应当分别比序列x和y多1，用于存放初始的0值；所以在c和b数组中存放的值，对应序列x和y要大1，下标容易弄错代码输出一种LCS，但LCS不一定只有1种*/int count;int lenLCS;int getsublength(string x,string y,int m,int n){int i,j;//初始化for(i = 0;i<m;i++){c[i][0] = 0;b[i][0] = 0;}for(i = 0;i<n;i++){c[0][i] = 0;b[0][i] = 0;}//填写表格 ==> 表格c[i,j]是序列Xi和Yj的一个LCS的长度，用于得到LCS最终长度//b[i,j] ==>用来记录数组中箭头的方向，便于反向查找LCS元素时找到元素位置for (i=1;i<m ;i++ ){for (j = 1;j<n ;j++ ){if(x[i-1]==y[j-1])//斜角==3{c[i][j] = c[i-1][j-1] + 1;//LCS长度加上1b[i][j] = 3;}else if(c[i-1][j]==c[i][j-1])////左和上相同==4 ==> 需要当做两种情况处理{c[i][j] = c[i-1][j];b[i][j] = 4;}else if(c[i-1][j]>c[i][j-1])//上==2{c[i][j] = c[i-1][j];b[i][j] = 2;}else if(c[i-1][j]<c[i][j-1])//左==1{c[i][j] = c[i][j-1];b[i][j] = 1;} }}return c[m-1][n-1];}void getLCS(int x,int y,int k,char **ans)//增加一个数组保存结果,存下每一次运行的结果{if( (x==0||y==0)&& k == 1 ) {int t;count++;for (t = 1;t<lenLCS ;t++ ){//printf("%c",z[t]);ans[count][t] = z[t];}//printf("%d",k);//printf("\n");return;}if(b[x][y] == 4){//getLCS(x-1,y,k,ans);getLCS(x,y-1,k,ans);}else if (b[x][y] == 3){k--;z[k] = str_x[x-1];getLCS(x-1,y-1,k,ans);}else if(b[x][y] == 2) getLCS(x-1,y,k,ans);else if(b[x][y] == 1) getLCS(x,y-1,k,ans); //return;}void printarray(int **A,int x,int y){int i,j;for (i = 0;i<x ;i++ ){for (j = 0;j<y ;j++ ){printf ("%3d",A[i][j]);}printf("\n");}}int main(){int i,j;count = 0;int len_x = str_x.length();int len_y = str_y.length();c = new int* [len_x+1];b = new int *[len_x+1];for (i = 0;i<len_x+1 ;i++ ){c[i] = new int [len_y+1];b[i] = new int [len_y+1];}int sublen = getsublength(str_x,str_y,len_x+1,len_y+1);lenLCS = sublen + 1;printf("length of LCS = %d\n",sublen);printarray(c,len_x+1,len_y+1);printf("\n");printarray(b,len_x+1,len_y+1);z = new char [sublen+1];char** ans = new char *[len_x+1];for (i = 0;i<len_x ;i++ ){ans[i] = new char [sublen+1];}getLCS(len_x,len_y,sublen+1,ans);/*for (i = 1;i<sublen+1 ;i++ ){printf("%c",z[i]);}printf("\n");*/printf("\n%d kinds of LCS:\n",count);for (i = 1;i<=count ; i++){for (j = 1;j<=sublen ;j++ ){printf("%c",ans[i][j]);}printf("\n");}/* printf("%s\n",z);*/delete [] z;for (i=0;i<len_x+1 ;i++){delete [] c[i];delete [] b[i];delete [] ans[i];}delete [] c;delete [] b;delete [] ans;}