算法导论-----最长公共子序列LCS(动态规划)

目录


一.概念梳理

  1. 子序列(subsequence): 一个特定序列的子序列就是将给定序列中零个或多个元素去掉后得到的结果(不改变元素间相对次序)。例如序列<A,B,C,B,D,A,B><A,B,C,B,D,A,B>的子序列有:<A,B><A,B><B,C,A><B,C,A><A,B,C,D,A><A,B,C,D,A>等。 
  2.公共子序列(common subsequence): 给定序列X和Y,序列Z是X的子序列,也是Y的子序列,则Z是X和Y的公共子序列。例如X=<A,B,C,B,D,A,B>X=<A,B,C,B,D,A,B>Y=<B,D,C,A,B,A>Y=<B,D,C,A,B,A>,那么序列Z=<B,C,A>Z=<B,C,A>为X和Y的公共子序列,其长度为3。但ZZ不是XXYY的最长公共子序列,而序列<B,C,B,A><B,C,B,A><B,D,A,B><B,D,A,B>也均为XXYY的最长公共子序列,长度为4,而XXYY不存在长度大于等于5的公共子序列。 
  3.最长公共子序列问题(LCS:longest-common-subsequence problem):In the longest-common-subsequence problem, we are given two sequences X=<x1,x2,...,xm>X=<x1,x2,...,xm> and Y=<y1,y2,...,yn>Y=<y1,y2,...,yn>and wish to find a (not “the”) maximum-length common subsequence of XX and YY . This section shows how to efficiently solve the LCS problem using dynamic programming.


最长公共子序列


二.最长公共子序列解决方案

方案1:蛮力搜索策略

蛮力搜索策略的步骤如下:

  1. 枚举序列XX里的每一个子序列xixi
  2. 检查子序列xixi是否也是YY序列里的子序列;
  3. 在每一步记录当前找到的子序列里面的最长的子序列。

蛮力策略也叫做暴力穷举法,是所有算法中最直观的方法,但效率往往也是最差的。在第1步枚举XX中所有的子序列有2m2m个,每个子序列在YY中检查时间复杂度为O(n)O(n)。因此蛮力策略的最坏时间复杂度为O(n2m)O(n2m),这是指数级算法,对较长的序列求LCS是不适用的。


方案2:动态规划策略

  • LCS问题具有最优子结构 
    令 X=<x1,x2,...,xm>X=<x1,x2,...,xm> 和 Y=<y1,y2,...,yn>Y=<y1,y2,...,yn> 为两个序列,Z=<z1,z2,z3,...,zk>Z=<z1,z2,z3,...,zk>XXYY的任意LCS。则

如果xm=ynxm=yn,则zk=xm=ynzk=xm=ynZk1Zk−1Xm1Xm−1Yn1Yn−1的一个LCS。 
如果xmynxm≠yn,那么zkxmzk≠xm,意味着ZZXm1Xm−1YY的一个LCS。 
如果xmynxm≠yn,那么zkynzk≠yn,意味着ZZXXYn1Yn−1的一个LCS。

  从上述的结论可以看出,两个序列的LCS问题包含两个序列的前缀的LCS,因此,LCS问题具有最优子结构性质。在设计递归算法时,不难看出递归算法具有子问题重叠的性质。 
  设C[i,j]C[i,j]表示XiXiYjYj的最长公共子序列LCS的长度。如果i=0i=0j=0j=0,即一个序列长度为00时,那么LCS的长度为0。根据LCS问题的最优子结构性质,可得如下公式:

C[i,j]=0C[i1,j1]+1MAX(C[i,j1],C[i1,j])i=0j=0i,j>0xi=yji,j>0xiyjC[i,j]={0,当i=0或j=0C[i−1,j−1]+1,当i,j>0且xi=yjMAX(C[i,j−1],C[i−1,j])当i,j>0且xi≠yj

根据上述的递归公式和初值,有如下伪代码和实现。

  • 伪代码1(递归):计算LCS的长度

伪代码1

  • 伪代码2(非递归):计算LCS的长度

伪代码2

  • 伪代码:构造一个LCS

LCS

三、C代码实现

C代码实现1

  实现1是完全按照《算法导论》中的伪代码编写而成。整个过程中表c和b的内容如下图所示: 

这里写图片描述

#include <stdio.h>
#include <string.h>
#define MAXLEN 50

void LCSLength(char *x, char *y, int m, int n, int c[][MAXLEN], int b[][MAXLEN])
{
    int i, j;

    for(i = 0; i <= m; i++)
        c[i][0] = 0;
    for(j = 1; j <= n; j++)
        c[0][j] = 0;
    for(i = 1; i<= m; i++)
    {
        for(j = 1; j <= n; j++)
        {
            if(x[i-1] == y[j-1])
            {
                c[i][j] = c[i-1][j-1] + 1;
                b[i][j] = 1;                    //如果使用'↖'、'↑'、'←'字符,会有警告,也能正确执行。
            }                                   //本算法采用1,3,2三个整形作为标记
            else if(c[i-1][j] >= c[i][j-1])
            {
                c[i][j] = c[i-1][j];
                b[i][j] = 3;
            }
            else
            {
                c[i][j] = c[i][j-1];
                b[i][j] = 2;
            }
        }
    }
}

void PrintLCS(int b[][MAXLEN], char *x, int i, int j)
{
    if(i == 0 || j == 0)
        return;
    if(b[i][j] == 1)
    {
        PrintLCS(b, x, i-1, j-1);
        printf("%c ", x[i-1]);
    }
    else if(b[i][j] == 3)
        PrintLCS(b, x, i-1, j);
    else
        PrintLCS(b, x, i, j-1);
}

int main()
{
    char x[MAXLEN] = {"ABCBDAB"};
    char y[MAXLEN] = {"BDCABA"};

    int  b[MAXLEN][MAXLEN];                        //传递二维数组必须知道列数,所以使用MAXLEN这个确定的数
    int  c[MAXLEN][MAXLEN];

    int m, n;

    m = strlen(x);
    n = strlen(y);

    LCSLength(x, y, m, n, c, b);
    PrintLCS(b, x, m, n);

    return 0;
}
  •   最终屏幕输出的结果为:B C B A 。完全正确。

实现2(空间优化)

  实现1中用了两个二维的表b和c,在时空开销上有改进的余地。我们完全可以去掉表b。因为每个c[i][j]c[i][j]只依赖于c[i1][j]c[i−1][j]c[i][j1]c[i][j−1]c[i1][j1]c[i−1][j−1]三项,当给定c[i][j]c[i][j]时,我们可以在O(1)O(1)的时间内判定出c[i][j]c[i][j]是使用了三项中的哪一项。从而节约了Θ(mn)Θ(mn)的空间。

#include <stdio.h>
#include <string.h>
#define MAXLEN 50

void LCSLength(char *x, char *y, int m, int n, int c[][MAXLEN]) {
    int i, j;

    for(i = 0; i <= m; i++)
        c[i][0] = 0;
    for(j = 1; j <= n; j++)
        c[0][j] = 0;
    for(i = 1; i<= m; i++) {
        for(j = 1; j <= n; j++) {
            if(x[i-1] == y[j-1]) {                          //仅仅去掉了对b数组的使用,其它都没变
                c[i][j] = c[i-1][j-1] + 1;
            } else if(c[i-1][j] >= c[i][j-1]) {
                c[i][j] = c[i-1][j];
            } else {
                c[i][j] = c[i][j-1];

            }
        }
    }
}
/*
void PrintLCS(int c[][MAXLEN], char *x, int i, int j) {         //非递归版PrintLCS
    static char s[MAXLEN];
    int k=c[i][j];
    s[k]='\0';
    while(k>0){
        if(c[i][j]==c[i-1][j]) i--;
        else if(c[i][j]==c[i][j-1]) j--;
        else{
            s[--k]=x[i-1];
            i--;j--;
        }
    }
    printf("%s",s);
}
*/
void PrintLCS(int c[][MAXLEN], char *x, int i, int j) {
    if(i == 0 || j == 0)
        return;
    if(c[i][j] == c[i-1][j]) {
        PrintLCS(c, x, i-1, j);
    } else if(c[i][j] == c[i][j-1])
        PrintLCS(c, x, i, j-1);
    else {
        PrintLCS(c, x, i-1, j-1);
        printf("%c ",x[i-1]);
    }
}

int main() {
    char x[MAXLEN] = {"ABCBDAB"};
    char y[MAXLEN] = {"BDCABA"};
    //char x[MAXLEN] = {"ACCGGTCGAGTGCGCGGAAGCCGGCCGAA"}; //算法导论上222页的DNA的碱基序列匹配
    //char y[MAXLEN] = {"GTCGTTCGGAATGCCGTTGCTCTGTAAA"};

    int  c[MAXLEN][MAXLEN];     //仅仅使用一个c表

    int m, n;

    m = strlen(x);
    n = strlen(y);

    LCSLength(x, y, m, n, c);
    PrintLCS(c, x, m, n);

    return 0;
}

  第一组测试序列,最终屏幕输出的结果为:B C B A 。完全正确。 

  第二组测试序列,最终屏幕输出的结果为:G T C G T C G G A A G C C G G 。和教材提供的结果相同,完全正确。

转自: https://blog.youkuaiyun.com/so_geili/article/details/53737001

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值