最长公共子序列(LCS, Longest Common Sequence)
前言
最长公共子序列是一类典型的动态规划问题，此类动态规划问题与之前提到的 钢条切割 和 矩阵链相乘的问题有不同的地方，不同的地方在于子问题的递归或者迭代需要根据条件进行选择。在钢条切割问题和矩阵链相乘问题中，我们对所有的条件下的子问题都进行相关的处理，没有任何的约束条件。而最长公共子序列问题中，需要根据不同的条件选择不同的子问题进行处理，我们可以称之为带条件的子问题处理。
与之前的处理动态规划问题的思路相类似，我们采用算法导论中的CRCC模板对这个问题进行处理。CRCC步骤包括：
- 表征最优解的结构(Characterize the structure of optimal solution)
- 递归定义解的值(Recursively define the value of optimal solution)
- 计算最优解的值(Compute the value of optimal solution)
- 构建最优解（Construct the optimal solution from the computed information)
利用这些步骤，我们可以把动态问题框架模块化，有助于帮助大家进行有效的思考和编程代码的撰写。我们就利用CRCC定义对最长公共子序列进行解析。
表征最优解的结构
最长公共子序列需要在两个序列中，寻找最长的公共子序列，子序列的索引在两个序列中单调递增。假定我们有两个序列Xm=<X1,X2,…Xm>和Y=<Y1,Y2…Yn>，假定Zk=<Z1,Z2…Zk>为这两个序列的子序列，并且长度为k。分三类情况进行讨论，
- 如果Xm=Yn, 我们可以得到Xm=Yn=Zk，并且Zk-1为Xm-1和Yn-1的子序列，Zk-1最长子序列的长度为k-1
- 如果Xm≠Zk,我们可以推断出，Zk为Xm-1和Yn两个序列的子序列，且最长子序列的长度为k
- 如果Yn≠Zk, 我们可以推断出,Zk为Xm和Yn-1两个序列的子序列，且最长子序列的长度为k
对于这三类情况，第一类情况，我们可以理解为两个子序列齐头并退（同步后退），第二/三种情况为单步后退，或者左脚后退或者右脚后退，后退后的子问题与当前的子问题形式相同，并且约束条件也相同，所以我们可以容易求出最优解的结构
递归定义解的值
假定Xm=<A,B,C>，并且m=3；假定Yn=<B,D,C>，并且n=3, 我们利用c[i,j]记录两个序列的最长子序列，i和j分别表示Xm和Yn元素中的下标，c[2,2]表示X2=<A,B>中和Y2=<B,D>两个序列的最长子序列，可以理解为子序列的最长子序列，有点拗口，上面表述中的两个子序列的含义不同，第一个子序列指的是原生子序列，第二子序列前面有定于最长，表示两个子序列的公共序列，也即公共最长子序列。如果下标为零，那么c[0,n]或c[m,0]的最长公共子序列的长度为0.
| C[i,j] | ∅(j=0) | B(j=1) | D(j=2) | C(j=3) |
|---|---|---|---|---|
| ∅(i=0) | 0 | 0 | 0 | 0 |
| A(i=1) | 0 | 0 | 0 | 0 |
| B(i=2) | 0 | c[1] [1]+1=1 | max{c[2] [1],c[1] [2]}=1 | max{c[1] [3],c[2] [2]}=1 |
| C(i=3) | 0 | max{c[2] [1], c[3] [0]}=1 | max{c[3] [1], c[2] [2]}=1 | c[2] [2]+1=2 |
鉴于上面的c[i] [j]矩阵的求值过程，很自然地递归定义解的值。
此递归定义是解决动态规划问题的核心和关键，如果此定义不明确或者无法准确描述问题，之后的递归或迭代程序便无从谈起。很多教程中称为状态转移方程，其实质是归纳表达式，假定子问题的解都已知，然后利用“假定”的已知子问题的解，求解父问题的解。在从子问题到父问题求解过程中，我们会经常看到求解最大或最小值的情况，而恰恰是这些最大或最小值的求解过程构成了递归或迭代核心元素，也体现了动态规划的最优解求解的特征，不断进行最大化或最小化过程，直至最终求得问题的最优解。
程序处理动态递归过程，可以选择利用递归（Top-down)，充分利用“大处着眼，小处着手”的原则，从n=upper bound开始，不断减少问题的复杂度，直至答案到递归的终止条件。如果此问题用递归的处理，那么递归的终止条件将是i=0 或j=0的时候，我们返回0值，此时表示最长公共子序列的长度为0，因为其中一个公共子序列不含有任何元素。
也可以选择利用bottom-up进行处理，解决问题从 最小子问题开始，不断利用迭代对后续问题进行处理，之前我们提到过DAG(有向无环图)的拓扑排序，利用bottom-up解决问题的过程，就是逆拓扑排序，不断向上求优的过程。
这是一张DAG的拓扑排序图，递归的过程一般④顶点开始，自上而下；迭代的过程一般从(O)顶点开始，自下而上。二者求解的复杂度一般情况下相同，求解效果很多时候并没有大的差异。
计算最优解的值
为了更好理解top-down和bottom-up的解题思路，我们利用两种方法解决本问题，同时会附加上解决问题过程的C语言代码，更好理解问题的本质。
a) Top-down计算最优解的值，本过程的核心是用记忆递归解决问题，自顶层之底层解决问题，充分利用递归的思想对问题进行充分解决。C语言的代码实现过程。
a-1) 头文件定义(longest_common_sequence.h)，数组c[i] [j]保存X1…Xi 和 Y1…Yj两个序列中的最长子序列的长度，数组b[i] [j]表示当前的c[i] [j]的值来自哪里，其可以来自三个不同的方向，↖对角线方向c[i-1] [j-1]还是←方向c[i] [j-1]亦或是↑c[i-1] [j]。
函数int lcs_length(char * X, char *Y, int m, int n, int (c)[N + 1], char ( b)[N + 1]) 返回X和Y的最长公共子序列的长度值，并进行相应的初始化工作。
/**
* @file longest_common_sequence.h
* @author your name (you@domain.com)
* @brief
* @version 0.1
* @date 2023-02-22
*
* @copyright Copyright (c) 2023
*
*/
#ifndef LONGEST_COMMON_SEQUENCE_H
#define LONGEST_COMMON_SEQUENCE_H
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <limits.h>
#define M 7 //length of string X
#define N 6 // length of string Y
/**
* @brief Get the longest common sequence of string A and string B
*
* @param A string X
* @param B string Y
* @param m length of string X
* @param n length of string Y
* @param c longest common sequence size
* @param b break symbol 'L'=left, 'D'=Diagonal, 'U'=Upper
*/
int lcs_length(char *X, char *Y, int m, int n, int (*c)[N + 1], char (*b)[N + 1]);
/**
* @brief Get the longest common sequence of string A and string B
*
* @param A string X
* @param B string Y
* @param m length of string X
* @param n length of string Y
* @param c longest common sequence size
* @param b breaking symbol,indication the path of longest commone sequence 'L'=left, 'D'=Diagonal, 'U'=Upper
*/
int lcs_length_aux(char *X, char *Y, int i, int j, int (*c)[N+ 1], char (*b)[N + 1]);
/**
* @brief print longest common sequence
*
* @param b char matrix(two dimension)
* @param X String A
* @param i index of A
* @param j index of B
*/
void print_lcs(char (*b)[N + 1], char *X, int i, int j);
#endif
a-2) 函数的实现,在递归实现的过程中，同样对递归函数进行了条件约束，当Xi=Yj我们调用子问题结果，当Xi≠Yj的时候，我们相当于调用两个递归函数，然后求取最大的值。
/**
* @file longest_common_sequence.c
* @author your name (you@domain.com)
* @brief
* @version 0.1
* @date 2023-02-22
*
* @copyright Copyright (c) 2023
*
*/
#ifndef LONGEST_COMMON_SEQUENCE_C
#define LONGEST_COMMON_SEQUENCE_C
#include "longest_common_sequence.h"
int lcs_length(char *X, char *Y, int m, int n, int (*c)[N+ 1], char (*b)[N + 1])
{
int i;
int j;
for(i=0;i<=m;i++)
{
for(j=0;j<=n;j++)
{
c[i][j]=INT_MIN;
}
}
return lcs_length_aux(X,Y,m,n,c,b);
}
int lcs_length_aux(char *X, char *Y, int i, int j, int (*c)[N + 1], char (*b)[N + 1])
{
if(i==0 || j==0)
{
c[i][j]=0;
}
if(c[i][j]>=0)
{
return c[i][j];
}
int temp1;
int temp2;
//from the top to down, that means from the bigger to smaller
//while
if(X[i]==Y[j])
{
c[i][j]=lcs_length_aux(X,Y,i-1,j-1,c,b)+1;
b[i][j]='D';
}
else
{
temp1 = lcs_length_aux(X, Y, i - 1, j, c, b);
temp2 = lcs_length_aux(X, Y, i, j - 1, c, b);
if(temp1>=temp2)
{
c[i][j]=temp1;
b[i][j] = 'U';
}
else
{
c[i][j]=temp2;
b[i][j] = 'L';
}
}
return c[i][j];
}
void print_lcs(char (*b)[N+1], char *X, int i, int j)
{
if(i==0 || j==0)
{
return;
}
if(b[i][j]=='D')
{
print_lcs(b, X, i - 1, j-1);
printf("%c\n",X[i]);
}
else
{
if(b[i][j]=='U')
{
print_lcs(b,X,i-1,j);
}
else if(b[i][j]=='L')
{
print_lcs(b,X,i,j-1);
}
}
}
#endif
a-3) 主函数，主函数提供X和Y的序列值，并对递归函数进行测试。
/**
* @file longest_common_sequence_main.c
* @author your name (you@domain.com)
* @brief
* @version 0.1
* @date 2023-02-22
*
* @copyright Copyright (c) 2023
*
*/
#ifndef LONGEST_COMMON_SEQUENCE_MAIN_C
#define LONGEST_COMMON_SEQUENCE_MAIN_C
#include "longest_common_sequence.c"
int main(void)
{
int m;
int n;
char X[M+1] = { '0','A','B','C','B','D','A','B'};
char Y[N+1]={'0','B','D','C','A','B','A'};
int c[M+1][N+1];
char b[M+1][N+1];
m=M;
n=N;
lcs_length(X,Y,m,n,c,b);
print_lcs(b,X,m,n);
getchar();
return EXIT_SUCCESS;
}
#endif
b) Bottom-up 求解最优解
由于篇幅的原因，我们就不对bottom-up的C语言实现进行详细解析，只贴出其实现的核心函数
void lcs_length(char *X, char *Y, int m, int n, int (*c)[N+1], char (*b)[N+1])
{
int i;
int j;
for(i=0;i<=m;i++)
{
c[i][0]=0;
}
for(j=0;j<=n;j++)
{
c[0][j]=0;
}
for(i=1;i<=m;i++)
{
for(j=1;j<=n;j++)
{
if (X[i] == Y[j]) // X =<A,B,C,B,D,A,B> and Y =<B,D,C,A,B,A>
{
c[i][j]=c[i-1][j-1]+1;
b[i][j]='D';
}
else
{
if(c[i-1][j]>=c[i][j-1])
{
c[i][j] = c[i - 1][j];
b[i][j]='U';
}
else
{
c[i][j] = c[i][j - 1];
b[i][j] = 'L';
}
}
}
}
}
总结
由于动态规划问题的求解很烧脑，状态方程的归纳是实现算法的核心，状态方程需要使用归纳法，归纳法属于高中的范畴，大学工科高等数学中很少提到归纳法，现在拾起来感觉非常吃力，很多时候感觉智力不足，无法突破问题分析，但是相信多加练习，假以时日，常规的动态规划解题就没有大问题。
参考书籍:
-《Introduction to algorithm 4ed edition》
1004
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
