More than 5 years have passed since last update.

最長部分列問題

Posted at 2019-12-06

初めまして。とりという名前でTwitterをしています。
研究で共通文字列問題(LCS)を実装する必要があり、理解度を深める目的で記事を書いてみました。
初投稿ですが分からないなりに頑張ってみようと思います。

最長共通部分列問題(LCS)とは

Wikipediaによると以下のように記述されています。応用先としてDNA配列の解析でLCSを求めるアルゴリズムが使われているみたいです。

最長共通部分列問題（さいちょうきょうつうぶぶんれつもんだい、英: Longest-common subsequence problem, LCS）とは、与えられた列の集合（しばしば、2つの列からなる集合）の最長共通部分列を見つけ出す問題である。（ここで部分列(subsequence)は、部分文字列(substring)とは異なることに注意する。前者は元の列の連続した項からなる必要はない。）

ここで、部分列と部分文字列は連続した文字列である必要がないという点が異なっています。例えばABCDEFとABBCCDのLCSはABCDになります。
競プロだとLCSの文字列長を求める問題が多いですが、今回はLCSまで求めるプログラムを実装していきたいと思います。

LCSのアルゴリズム

LCSを求めることを考えていきます。まず、入力される2つの文字列の組をそれぞれX=(x₁,x₂...x_m)、Y=(y₁,y₂...y_n)とします。説明のため具体的にX=ABCDEF,Y=ABBCCDとおいてみます。ここでx₂はBになります。そして大きさ(m+1)×(n+1)のDPテーブルdp[m+1][n+1]を用意します。ここで配列の大きさを文字列より1つ多めに取っています。以上で準備は終わりです。ここからLCSを求めていきます。

1. dpの1行目と1列目の値を0に初期化する。
2. dp[i][j](i=1,2,...,m ,j=1,2,...,n)の値を以下の式から求める。

dp_{i,j} = max\left\{
\begin{array}{ll}
dp_{i-1,j-1}+ s(x_i,y_j)\hspace{1mm},\\
dp_{i-1,j}\hspace{1mm}, \\
dp_{i,j-1}\hspace{1mm}.
\end{array}
\right.

\\
s(x_i,y_j) = \left\{
\begin{array}{ll}
1 & (x_i = y_j)\hspace{1mm},\\
0 & (x_i \neq y_j)\hspace{1mm}.
\end{array}
\right.

3. dp[i][j]がどの値を使って更新されたかがわかるように更新元のインデックスを覚えておく。値が等しい場合、優先度は左上>上>左とする。
4. i=m,j=nと初期化する。またLCSを求めるためにスタックを用意する。
5. dp[i][j]と更新元のdpテーブルの要素の値を比較し異なっていた場合x_i-1をスタックに pushする。また、iとjをそれぞれ更新元のインデックスに対応させる。
6. dpテーブルの要素の値が0ではなくなるまで5を繰り返す。
7. スタックからpopしていきLCSを求める。

この処理を簡単に説明していきます。ここで、入力される文字列をX=ABBCCD,Y=ABCDEFとしてm,nをそれぞれX,Yの文字列長とします。dpテーブルは以下のように初期化されます。

	A	B	C	D	E	F
-	0	0	0	0	0	0
A
B
B
C
C
D

次にi=1,j=1としてi,jが指すdpテーブルdp[i][j]の左上、左、上から最大値を求めます。
x_iとy_jが等しければ1を足してます。以下はi=2,j=3の時のdpテーブルの中身です。この時、x_i=B、y_j=Cなので参照元はdp[2][2]になります。

	A	B	C	D	E	F
-	0	0	0	0	0	0
A	1	1	1	1	1	1
B	1	2	2
B
C
C
D

これを繰り返していき、最終的にdpテーブルの中身は以下のようになります。

	A	B	C	D	E	F
-	0	0	0	0	0	0
A	1	1	1	1	1	1
B	1	2	2	2	2	2
B	1	2	2	2	2	2
C	1	2	3	3	3	3
C	1	2	3	3	3	3
D	1	2	3	4	4	4

最後にLCSを求めていきます。dpテーブルの右下から要素をたどっていきます。ここで注目する点は値が変わる(+1される)のはx_iとy_jが等しいときなのでdpテーブルの更新元と比較して値が変わっている時が文字列が一致した時とわかります。簡単な説明でしたがこのような手順でLCSを求めることができます。これをコードに落とし込んでいきます。

lcs.cpp

#include <bits/stdc++.h>
using namespace std;

struct LcsCell{
    int score;
    int row;
    int col;
    LcsCell* prev_cell;
};

// セルのメンバに値を代入する
void FillInCell(LcsCell *cell, int score, int i, int j, LcsCell *prev_cell){
    cell->score = score;
    cell->row = i;
    cell->col = j;
    cell->prev_cell = prev_cell;
}

// あるセルの左、上、左上にあるセルのスコアのうち最大値を返す
// 3つの値が同値の場合、返り値の優先度は左上>上>左とする
// また、最大値を持つセルのインデックスも求める
int CalcMaxAndPrevCell(int above_score, int left_score, int above_left_score,
                         int current_row, int current_col, int* prev_cell_row, int* prev_cell_col){
    int tmp_max;
    // 左上と上のセルを比較
    if(above_left_score >= above_score){
        *prev_cell_row = current_row - 1;
        *prev_cell_col = current_col - 1;
        tmp_max = above_left_score;
    }
    else{
        *prev_cell_row = current_row - 1;
        *prev_cell_col = current_col;
        tmp_max = above_score;
    }

    // 左上と上のセルを比較したものと左のセルを比較
    if(tmp_max < left_score){
        *prev_cell_row = current_row;
        *prev_cell_col = current_col - 1;
        tmp_max = left_score;
    }
    return tmp_max;
}

// LCSを求める関数
string FindLcs(LcsCell *cell, string x){
    LcsCell tmp = *cell;
    string lcs = "";

    printf("---trace route---\n");
    while(tmp.score){
        if(tmp.score != tmp.prev_cell->score){
            lcs.insert(lcs.begin(), x[(tmp.row)-1]);
            printf("i: %d, j: %d\n", tmp.row,tmp.col);
        }
        tmp = *(tmp.prev_cell);
    }
    putchar('\n');
    
    return lcs;
}


int main(void){
    string x = "ABBCCD";
    string y = "ABCDEF";
    int m = x.length();
    int n = y.length();

    cout << "String1 : " << x << endl;
    cout << "string2 : " << y << endl;
    putchar('\n');

    int i,j;

    //文字列が短い場合は LcsCell dp[m+1][n+1]; でも可能
    LcsCell **dp;
    dp = (LcsCell **)malloc(sizeof (LcsCell*) * (m+1));
    for(i = 0;  i <= m; i++){
        dp[i] = (LcsCell *)malloc(sizeof (LcsCell) * (n+1));
    }

    // dpテーブルの初期化
    for(i=0; i<=m; i++){
        FillInCell(&dp[i][0], 0, i, 0, NULL);
    }
    for(j=0; j<=n; j++){
        FillInCell(&dp[0][j], 0, 0, j, NULL);
    }

    // dpループ
    int tmp_prev_row, tmp_prev_col;
    tmp_prev_row = tmp_prev_col = 0;
    int max_score;
    for(i=1; i<=m; i++){
        for(j=1; j<=n; j++){
            if(x[i-1] == y[j-1]){
                max_score = CalcMaxAndPrevCell(dp[i-1][j].score, dp[i][j-1].score, dp[i-1][j-1].score + 1, i, j, &tmp_prev_row, &tmp_prev_col);
            }
            else{
                max_score = CalcMaxAndPrevCell(dp[i-1][j].score, dp[i][j-1].score, dp[i-1][j-1].score, i, j, &tmp_prev_row, &tmp_prev_col);
            }
            FillInCell(&dp[i][j], max_score, i, j, &dp[tmp_prev_row][tmp_prev_col]);
        }
    }

    // dpテーブルの表示
    printf("---dp table---\n");
    for(i=0; i<=m; i++){
        for(j=0; j<=n; j++){
            printf("%d ",dp[i][j].score);
        }
        putchar('\n');
    }
    putchar('\n');

    // LCSの表示
    string str = FindLcs(&dp[m][n], x);
    cout << "LCS : " << str << endl;
    return 0;
}

実行結果

String1 : ABBCCD
string2 : ABCDEF

---dp table---
0 0 0 0 0 0 0 
0 1 1 1 1 1 1 
0 1 2 2 2 2 2 
0 1 2 2 2 2 2 
0 1 2 3 3 3 3 
0 1 2 3 3 3 3 
0 1 2 3 4 4 4 

---trace route---
i: 6, j: 4
i: 5, j: 3
i: 3, j: 2
i: 1, j: 1

LCS : ABCD

まとめ

今回はC++で実装したのですが競プロのいい練習になって楽しかったです。
1ヶ月ぶりにC++を弄ったので自己参照構造体すら忘れていて...となりました。ポインタの扱いもぐだってしまってダメだなと思いました。練習していきたいと思います。

You get articles that match your needs
You can efficiently read back useful information
You can use dark theme

What you can do with signing up