最长公共子序列算法简介

最长公共子序列（Longest Common Subsequence，LCS）算法用于寻找两个序列中长度最长的共同子序列。本文对其基于动态规划（DP）的解法进行简单介绍。

我们首先说明什么是最长公共子序列。举例来讲，对于两个数组 [1, 2, 3, 4, 5]，和 [2, 4, 6]，其 LCS 就是 [2, 4]。

LCS 长度算法

我们首先解决一个更简单的问题：如何得到两个数组 $arr1$ 和 $arr2$ LCS 的长度。

这个问题可以采用 DP 解决，复杂度是 $O(m\times n)$，其中 $m$ 和 $n$ 分别是两个数组 $arr1$ 和 $arr2$ 的长度。

状态表定义

创建状态表 $dp[m+1][n+1]$，其中 $dp[i][j]$ ($0\leq i\leq m$，$0\leq j\leq n$) 代表 $arr1[0…i-1]$ 和 $arr2[0…j-1]$ 之间的 LCS 长度。换句话说，代表 $arr1$ 前 $i$ 个元素和 $arr2$ 前 $j$ 个元素的 LCS 长度。

基础情况

若 $i = 0$ 或者 $j = 0$，则至少一个用来计算 LCS 的子数组是空的，因此 $dp[i][j] = 0$。

递推公式

DP 的思路是利用已知的最优信息得到未知的最优信息。在一开始，我们：

• 已知的信息：$dp[0…m-1][0] = 0$，$dp[0][0…n-1] = 0$；
• 最后想要得到的信息：$dp[m][n]$。

容易想到 $i$ 和 $j$ 从 $0$ （已知）开始推导到 $m$ 和 $n$ （未知）。

那么根据状态表 $dp[i][j]$ 的定义，对于元素 $arr1[i]$ 和 $arr2[j]$，有两种情况：

$arr1[i] \neq arr2[j]$

这种情况下，我们不可以延长子序列。那么 $dp[i][j]$ 的候选集合为：

1. $dp[0…i-1][j]$
2. $dp[i][0…j-1]$
3. $dp[0…i-1][0…j-1]$

乍一看如果从这个集合选最大值会使算法出现平方复杂度。但是其实这个步骤可以用 $O(1)$ 完成，因为：

1. $MAX(dp[0…i-1][j]) = dp[i-1][j]$。缩短任一子数组的长度只会缩短 LCS 的长度。
2. 同理，$MAX(dp[i][0…j-1])=dp[i][j-1]$。
3. 同理，$MAX(dp[0…i - 1][0…j-1])=dp[i-1][j-1]$。
4. 同理，$MAX(dp[i-1][j-1], dp[i-1][j], dp[i][j-1])= MAX(dp[i-1][j], dp[i][j-1])$。

因此有 $$dp[i][j] = MAX(dp[i-1][j], dp[i][j-1]), arr1[i] \neq arr2[j]$$

$arr1[i] = arr2[j]$

这种情况下，我们可以延长子序列。延长子序列之前的 LCS 长度是 $dp[i-1][j-1]$，得到 $dp[i-1][j-1] + 1$。因此有 $$dp[i][j] = dp[i-1][j-1] + 1, arr1[i]=arr2[j]$$

代码

时间和空间复杂度均为 $O(m\times n)$。

template <typename T>
size_t getLongestCommonSequenceLength(const std::vector<T>& arr1,
                                      const std::vector<T>& arr2) {
  const size_t arr1Size = arr1.size();
  const size_t arr2Size = arr2.size();
  std::vector<std::vector<size_t>> dp(arr1Size + 1,
                                      std::vector<size_t>(arr2Size + 1, 0));
  for (size_t i = 1; i <= arr1Size; i++) {
    for (size_t j = 1; j <= arr2Size; j++) {
      if (arr1[i - 1] == arr2[j - 1]) {
        dp[i][j] = dp[i - 1][j - 1] + 1;
      } else {
        dp[i][j] = std::max(dp[i - 1][j], dp[i][j - 1]);
      }
    }
  }
  return dp[arr1Size][arr2Size];
}

LCS 内容算法

进一步地，我们可能需要直接得到 LCS 的内容。

暴力算法

我们可以继续用 LCS 长度算法，只不过把记录长度改成纪录 $arr1[0…i-1]$ 和 $arr2[0…j-1]$ 的 LCS，写成代码：

template <typename T>
std::vector<T> getLongestCommonSequence(const std::vector<T>& arr1,
                                        const std::vector<T>& arr2) {
  const size_t arr1Size = arr1.size();
  const size_t arr2Size = arr2.size();
  std::vector<std::vector<std::vector<T>>> dp(
      arr1Size + 1, std::vector<std::vector<T>>(arr2Size + 1));
  for (size_t i = 1; i <= arr1Size; i++) {
    for (size_t j = 1; j <= arr2Size; j++) {
      if (arr1[i - 1] == arr2[j - 1]) {
        dp[i][j] = dp[i - 1][j - 1];
        dp[i][j].push_back(arr1[i - 1]);
      } else {
        dp[i][j] = dp[i - 1][j].size() > dp[i][j - 1].size() ? dp[i - 1][j]
                                                             : dp[i][j - 1];
      }
    }
  }

  return dp[arr1Size][arr2Size];
}

这个算法是不可用的。因为这样做会复制和存储大量数组，时间和空间复杂度都是三次方级别。

改进算法

我们可以直接得到 LCS 而不是复制和存储大量的无用数据。简而言之，我们可以从 LCS 长度算法的状态表 $dp$ 倒推出 LCS。

我们用数组 [1, 2, 3, 4, 5] 和 [2, 4, 6] 举例。下图是 LCS 长度算法在这两个数组上的状态表，并且标出了每个状态值的来源。如果 $arr1[i] \neq arr2[j]$ 时 $dp[i-1][j] = dp[i][j-1]$，我们始终取 $dp[i-1][j]$。

我们需要倒推的是红色路线。容易发现：

1. 在值来源是左上的时候（即 $arr1[i]=arr2[j]$），LCS 添加元素。
2. 值的来源方向取决于 $arr1[i]$ 和 $arr2[j]$ 是否相等。如果相等，值来源左上。如果不相等，值来源取决于上值和左值哪个更大。如果上值和左值一样大，我们优先取上值。

那么算法思路就很明确了：

1. 算法开始时，$i=m$ 且 $j = n$。
2. 若 $arr1[i-1]=arr2[j-1]$，LCS 需要在前面插入 $arr1[i-1]$，并且向左上方移动（即 $i$ 和 $j$ 均减一）。
3. 若 $arr1[i-1] \neq arr2[j-1]$，我们需要判断是上移还是左移：
4. 若 $dp[i-1][j] \geq dp[i][j-1]$，我们应该上移，即 $i$ 减一。
5. 反之，我们应该左移，即 $j$ 减一。
6. 若 $i = 0$ 或者 $j = 0$，算法结束。

写成代码：

template <typename T>
std::vector<T> getLongestCommonSequence(const std::vector<T>& arr1,
                                        const std::vector<T>& arr2) {
  const size_t arr1Size = arr1.size();
  const size_t arr2Size = arr2.size();
  std::vector<std::vector<size_t>> dp(arr1Size + 1,
                                      std::vector<size_t>(arr2Size + 1, 0));
  for (size_t i = 1; i <= arr1Size; i++) {
    for (size_t j = 1; j <= arr2Size; j++) {
      if (arr1[i - 1] == arr2[j - 1]) {
        dp[i][j] = dp[i - 1][j - 1] + 1;
      } else {
        dp[i][j] = std::max(dp[i - 1][j], dp[i][j - 1]);
      }
    }
  }

  std::vector<T> lcs;
  lcs.reserve(dp[arr1Size][arr2Size]);
  size_t i = arr1Size;
  size_t j = arr2Size;
  while (i > 0 && j > 0) {
    if (arr1[i - 1] == arr2[j - 1]) {
      // Push back and reverse lcs later
      lcs.push_back(arr1[i - 1]);
      i--;
      j--;
    } else if (dp[i - 1][j] >= dp[i][j - 1]) {
      i--;
    } else {
      j--;
    }
  }

  std::reverse(lcs.begin(), lcs.end());
  return lcs;
}

容易得出时间和空间复杂度仍然是 $O(m\times n)$。