LCS最长公共子串

问题介绍

LCS问题(longest common subsequence problem)指的是求解两个字符串最长公共子序列问题。这里的子序列是可以不连续的。LCS问题广泛地出现在计算生物学中（DNA序列、系统生成树等等）。这里介绍如何解决LCS问题，以及算法的正确性证明和性能分析。

解决方案

假设需要求解串X,Y的LCS，其中|X|=n,|Y|=m,c[i][j]表示X[1…i]和Y[1…j]的LCS长度,z[1…k]表示X[1…i]和Y[1…j]的LCS,k=c[i][j]，则问题就是求解$c[n][m]和z[c[n][m]]$

朴素的想法是，按照问题的要求，我们可以得到串X的所有子串，共有$2^n$个，然后判断该子串是否出现在串Y中，每次判断都需要遍历串Y，因此时间复杂度为$O(2^n\ast m)$，这显然是我们不能接受的复杂度。

为了解决这个问题，我们需要得到该问题的一些性质：

定理:

$$ifX[i]==Y[j]:c[i][j]=c[i-1][j-1]+1$$

$$otherwise：c[i][j]=max(c[i-1][j],c[i][j-1])$$

引理1:如果X[i]==Y[j]，则$z[c[i][j]]=X[i]$

证明：如果$z[c[i][j]]\ne X[i]$，且X[i]==Y[j]，那么不妨将X[i]加入到LCS中，$c[i][j]$加一，因此$z[1..c[i][j]]$不是LCS，与条件矛盾，证毕。

引理2:如果X[i]==Y[j]，则X[i-1]和Y[j-1]的LCS是$z[1..c[i][j]-1]$

证明：X[i-1]和Y[j-1]的LCS不是$z[1..c[i][j]-1]$，则使用X[i-1]和Y[j-1]的LCS替换$z[1..c[i][j]-1]$后再加上X[i]会得到X[i]和Y[j]的一个更长的一个LCS，与条件矛盾，证毕。引理2这里展示了问题的最优子结构

引理3:如果两个串的LCS包含两个串的末尾元素X[i]和Y[j]，则这两个元素相等

证明：如果X[i]不等于Y[j]，则在LCS中，X[i]对应Y[j’],j’<j，Y[j]对应X[i’],i’<i，这不符合公共子串保留原串顺序的性质，矛盾。证毕。

当X[i]==Y[j]时，由引理1保证了此时的LCS是串X[i]对应串Y[j’],$j^{\prime }⩽j$，对于$j^{\prime }<j$的情况，我们不妨用$j$来替换$j^{\prime }$，这样也不会对LCS的长度有什么影响，然后由引理2，$c[i][j]=c[i-1][j-1]+1$

对于第二个条件，因为X[i]$\ne$Y[j]，那么此时的LCS要么属于$c[i-1][j]$，要么属于$c[i][j-1]$。

假设都不属于，那么此时的LCS一定包含了$c[i-1][j]$中没有的元素X[i]和$c[i][j-1]$中没有的元素Y[j]，由引理3，矛盾。

在证明了上述定理以后我们可以根据该式子计算LCS:

1. 递归法

int LCS(string& x, string &y,int n,int m)
{
    if(-1==n || -1==m) return 0;
    if(x[n]==y[m]) return LCS(x, y, n-1, m-1)+1;
    else return max(LCS(x, y, n-1, m), LCS(x, y, n, m-1));
}

性能分析

分析递归树，最坏的情况是每次x[n]!=y[m]，那么会得到一个高度为n+m的二叉树，时间复杂度为$O(2^{n+m})$，空间复杂度为$O(n+m)$

2. 备忘录方法（记忆化搜索）

分析上面时间复杂度我们发现，在搜索过程中很多的子问题都是一模一样的，也就是具有重叠子问题性质，因此我们不妨每计算出一个子问题的结果就进行一次记录，后面再次需要求解结果的时候就不需要再计算，而是直接返回结果。

int LCS(string& x, string &y,int n,int m,int *c)
{
    if(-1==n || -1==m) return 0;
    if(c[n*y.size()+m]) return c[n*y.size()+m];
    int ret;
    if(x[n]==y[m]) ret = LCS(x, y, n-1, m-1, c)+1;
    else ret = max(LCS(x, y, n-1, m, c), LCS(x, y, n, m-1, c));
    return c[n*y.size()+m]=ret;
}

int main()
{
    string X,Y;
    cout<<"请输入字符串X:"; cin>>X;
    cout<<"请输入字符串Y:"; cin>>Y;
    int* c = new int[X.size()*Y.size()+10]();
    cout<<"字符串X和Y的LCS的长度为"<<LCS(X, Y, X.size(), Y.size(), c)<<endl;
    delete c;
    return 0;
}

性能分析

我们可以把对结果是否已经计算出的判断和返回答案的耗费记录在调用该状态答案的耗费上，把实际结果的计算记录在该状态中，则最坏情况下每种状态都要计算出来，因此时间复杂度为$O(nm)$，空间复杂度为$O(nm)$

3. 自底向上计算（动态规划法）

我们观察计算的过程，如果我们对状态空间按照从左向右从上向下进行求解，就可以计算出所有的答案

int LCS(string& x, string &y,int n,int m,int *c)
{
    for(int i=1; i<=x.size(); ++i)
    {
        for(int j=1;j<=y.size(); ++j)
        {
            if(x[i-1] == y[j-1])
                c[i*y.size()+j]=c[(i-1)*y.size()+j-1]+1;
            else
                c[i*y.size()+j]=max(c[(i-1)*y.size()+j], c[(i)*y.size()+j-1]);
        }
    }
    return c[n*y.size()+m];
}

性能分析

时间复杂度$O(nm)$，因为访问的是连续的内存空间，因此这里的$O(nm)$应该比上面小。空间复杂度$O(nm)$，如果使用滚动数组还能够将空间复杂度降低为$O(min(n,m))$，如果不使用滚动数组，想要得到完整的LCS串需要在计算的时候设置指针，最后进行回溯。如果使用滚动数组，需要使用分治法得到LCS串。

回溯如图：