LCS（最長共通部分列）問題を解く

$\newcommand{\dp}{\text{dp}}$

## 問題設定
与えられた2つの文字列の最長共通部分列(LCS)の長さを求めてください。
また、最長共通部分列(LCS)の例を一つ出力してください。
ただし各文字列の長さは1以上3000以下であるものとします。

## LCSの長さを計算する
最初に前半の問題を解きます。後半の問題は前半の結果を使って復元するという方針になります。
### アルゴリズム
具体例があったほうが分かりやすいので、文字列**atcoder**と**abcde**について考えてみます。
**atcoder**と**abcde**のLCSの長さを動的計画法で計算します。
次のようなテーブルを持ちます。
|  | - | a | b | c | d | e |  
|:-:|:-:|:-:|:-:|:-:|:-:|:-:|
| **-**  |0   | 0  | 0  | 0  | 0  | 0  |   
| **a**  | 0  |   |   |   |   |   |   
| **t**  | 0  |   |   |   |   |   |   
| **c**  |  0 |   |   |   |   |   |   
| **o**  |  0 |   |   |   |   |   |   
| **d**  | 0  |   |   |   |   |   |  
| **e**  | 0  |   |   |   |  ★ |   |   
| **r**  | 0  |   |   |   |   |   |     

例えば★には、**atcode**と**abcd**のLCSの長さが入るようにします。
以下、**atcoder**を文字列$S$、**abcde**を文字列$T$と呼称して一般化した状況を考えます。
$\dp[i][j]$を、「文字列$S$の$i$文字目までと、文字列$T$の$j$文字目までのLCSの長さ」とします。
上のテーブルで既に書いているように、$\dp[i][0]$と$\dp[0][j]$は全て$0$で初期化できます。なぜなら、どんな文字列も長さ$0$の文字列とLCSを取ると長さ$0$になるからです。
この状況では、次のような遷移により動的計画法が回ります。ただし、$i,j>0$であり、$S_i$は$S$の$i$文字目で、$T_j$は$T$の$j$文字目です。

```math
\dp[i][j] = 
\begin{cases}
\max(\dp[i-1][j],\dp[i][j-1]) && (S_i \neq T_j) \\
\dp[i-1][j-1] + 1 && (S_i = T_j)
\end{cases}
```

実際に遷移を実行してみると、上のテーブルは以下のように埋まります。

|  | - | a | b | c | d | e |  
|:-:|:-:|:-:|:-:|:-:|:-:|:-:|
| **-**  |0   | 0  | 0  | 0  | 0  | 0  |   
| **a**  | 0  | **1**  | 1  | 1  | 1  | 1  |   
| **t**  | 0  | 1  | 1  | 1  | 1  | 1  |   
| **c**  |  0 | 1  |  1 |  **2** |  2 |  2 |   
| **o**  |  0 | 1  | 1  | 2  | 2  | 2  |   
| **d**  | 0  | 1  | 1  | 2  |  **3** | 3  |  
| **e**  | 0  | 1  | 1  | 2  | 3  |  **4** |   
| **r**  | 0  |  1 | 1  |  2 |  3 | 4  |  

太字になっているマスは、$S_i = T_j$の方の処理が実行された部分です。
右下の数字$\dp[7][5]$を見ることで、LCSの長さ$4$が正しく求められました。


### 証明
以下このことの正当性を証明します。漸化式が正しいことを確かめれば十分です。
$S_i \neq T_j$のときを考えます。$k = \max(\dp[i-1][j],\dp[i][j-1])$とおきます。$\dp[i][j]$が$k$以上であることは明らかです。高々$k+1$なので、$\dp[i][j] \neq k+1$であることを背理法で示します。
$k = \dp[i-1][j]$であるとしましょう。文字列$S$の$i-1$文字目までと、文字列$T$の$j$文字目までのLCSの末尾が、$T_q ~ (1 \le q \le j)$に対応していたとします。仮定$\dp[i][j] = k+1$より、$S_i$は文字列$T$の$q+1$文字目以降のどれかに対応していなければなりません。$S_i \neq T_j$だったので、$S_i$に対応する文字は$T_r ~ (q+1 \le r < j)$と表せるはずです。すると、文字列$S$の$i$文字目までと、文字列$T$の$j$文字目までのLCSを考えるためには、文字列$S$の$i$文字目までと、文字列$T$の$r$文字目までのLCSを考えれば十分だったということになります。特に、$k+1 = \dp[i][j]$は$\dp[i][j-1]$以下でなくてはなりません。ところが、$k$の取り方より$\dp[i][j-1] \le k$であるので、矛盾します。以上の議論により、$S_i \neq T_j$のときは$\dp[i][j] = k$であることが分かります。
$S_i = T_j$のときを考えます。$k = \dp[i-1][j-1]$とおきます。明らかに$k+1 \le \dp[i][j] \le k+2$です。もし$\dp[i][j] = k+2$なら、$\dp[i-1][j] = \dp[i][j-1] = k+1$です。すなわち、文字列$S$の$i-1$文字目までと、文字列$T$の$j-1$文字目までのLCSの末尾を$S_p = T_q ~ (1 \le p \le i-1, ~ 1 \le q \le j-1)$とおくと、ある$p \le p' \le i-1, ~ q \le q' \le j-1$が存在して $U = S_i = T_j$が$S_{p'} = T_{q'} = U$を満たしていなければならないということになります。この時点で、文字列$S$の$i-1$文字目までと、文字列$T$の$j-1$文字目までのLCSに$S_{p'} = T_{q'} = U$を追加して長さ$k+1$の共通部分列が作れてしまうので矛盾します。以上の議論により、$S_i = T_j$のときは$\dp[i][j] = \dp[i-1][j-1] + 1$であることが示されました。

### ソースコード（C++）
```C++
#include <iostream>
#include <string>
#include <vector>
using namespace std;

int main(){
    
    string S,T;
    cin >> S >> T;
    int s = S.size(), t = T.size();
    vector<vector<int>> dp(s+1,vector<int>(t+1));
    for (int i=1; i<=s; i++){
        for(int j=1; j<=t; j++){
            dp[i][j] = (S[i-1] == T[j-1]) ? dp[i][j] = dp[i-1][j-1] + 1 : max(dp[i-1][j],dp[i][j-1]);
        }
    }
    cout << dp[s][t] << endl;

}
```

## LCSを出力する
後半の問題にとりかかります。LCSを具体的に一つ出力せよという問題です。

### アルゴリズム
動的計画法で作ったテーブルを右下から辿るという方針で解きます。

|  | - | a | b | c | d | e |  
|:-:|:-:|:-:|:-:|:-:|:-:|:-:|
| **-**  |0   | 0  | 0  | 0  | 0  | 0  |   
| **a**  | 0  | **1**  | 1  | 1  | 1  | 1  |   
| **t**  | 0  | 1  | 1  | 1  | 1  | 1  |   
| **c**  |  0 | 1  |  1 |  **2** |  2 |  2 |   
| **o**  |  0 | 1  | 1  | 2  | 2  | 2  |   
| **d**  | 0  | 1  | 1  | 2  |  **3** | 3  |  
| **e**  | 0  | 1  | 1  | 2  | 3  |  **4** |   
| **r**  | 0  |  1 | 1  |  2 |  3 | 4  | 

まず`res = ""`として空な文字列を持っておきます。
このテーブルの右下$\dp[7][5] = 4$からスタートして、「どのマス目から遷移してきたのか？」を考えていきます。
まず$S_7 \neq T_5$です。左と上を見てみると、上の$\dp[6][5] = 4$と値が一致しています。すなわち、ここでは$\dp[7][5] = \max(\dp[6][5],\dp[7][4]) = \dp[6][5]$と遷移してきたのだと考えられます。
もし$\dp[6][5] = \dp[7][4]$であった場合どちらに遡れば良いのかという問題が生じますが、どちらでも良いです。実際、漸化式を見ればわかるように「マス目の左も上も現在いるマスと異なる数字である」という状況なら必ず$S_i = T_j$で遷移して来ているので、左にも上にも左上にも遡れないといった行き止まりの状況に遭遇することは絶対にありません。従って、矛盾のない範囲で好きにマス目を移動して構いません。
$\dp[6][5] = 4$は、$S_6 = T_5$より$\dp[6][5] = \dp[5][4] + 1$で遷移してきているはずです。すなわち$S_6 = T_5$がLCSの末尾に追加されながら遷移してきています。そこで、`res`に$S_6 = T_5$を**左から**追加して`res = "e"`とします。
以降同様に処理します。$S_5 = T_4$より$\dp[5][4] = \dp[4][3] + 1 = 3$という遷移が特定されるので`res`に$S_5 = T_4$を**左から**追加して`res = "de"`となります。$S_4 \neq T_3$なので左と上を見ると、$dp[4][3] = 2$は上の数字$\dp[3][3] = 2$に一致します。$S_3 = T_3$より$\dp[3][3] = \dp[2][2] + 1 = 2$という遷移が特定され、`res`に$S_3 = T_3$を**左から**追加して`res = "cde"`となります。$S_2 \neq T_2$です。左と上を見るとどちらも$\dp[2][2] = 1$に一致しているので、この場合はどちらに移っても良いです。$\dp[2][1]$に移りましょう。$S_2 \neq T_1$なのでまた左と上を見ます。一致している上$dp[1][1] = 1$に移ります。$S_1 = T_1$なので$\dp[1][1] = \dp[0][0] + 1 = 1$という遷移が特定され、`res`に$S_1 = T_1$を**左から**追加して`res = "acde"`となります。マス目の添字のどちらかが$0$に到達した時点で、空な文字列からLCSの要素を見出すことはできないので終了し、`res = "acde"`を出力します。

### ソースコード（C++）
```C++
#include <iostream>
#include <string>
#include <vector>
using namespace std;

int main(){
    
    string S,T;
    cin >> S >> T;
    int s = S.size(), t = T.size();
    vector<vector<int>> dp(s+1,vector<int>(t+1));
    for (int i=1; i<=s; i++){
        for(int j=1; j<=t; j++){
            dp[i][j] = (S[i-1] == T[j-1]) ? dp[i][j] = dp[i-1][j-1] + 1 : max(dp[i-1][j],dp[i][j-1]);
        }
    }
    cout << dp[s][t] << endl;
    int i = s, j = t;
    string res = "";
    while (i > 0 && j > 0){
        if (S[i-1] == T[j-1]){
            res = S[i-1] + res;
            i--; j--;
        }
        else{
            if (dp[i][j] == dp[i][j-1]) j--;
            else i--;
        }
    }
    cout << res << endl;

}
```
## 入出力例
### 入力
```
abracadabra
avadakedavra
```
### 出力
```
7
aaadara
```


## 計算量について
このアルゴリズムの時間計算量は$O(|S|\cdot|T|)$です（$\dp$テーブルの広さに一致）。よって、$|S|,|T|\le3000$の下では最大でも約$10^7$の計算量となり十分高速に動作します。
しかし驚くべきことに、$O(|S|+|T|+\delta(S,T)^2)$（ただし$S,T$のLCSを$L$としたとき$\delta(S,T) = (|S|-|L|) + (|T|-|L|)$）の解法が存在するようです[3]。

## 参考文献
[1] [Educational DP Contest / DP まとめコンテスト F-LCS](https://atcoder.jp/contests/dp/tasks/dp_f)
[2] [個人的な競プロメモ - F-LCS](https://scrapbox.io/pocala-kyopro/F_-_LCS)
[3] [えびちゃんのメモ - Output-sensitive algorithm の話 + diff や LCS など](https://rsk0315.github.io/jekyll-blog/2022/03/26/diff-algo.html)