算法实践--最长公共子序列(Longest Common Subsquence)

什么是最长公共子序列

X=ACCG

Y=CCAGCA

长度为1的公共子序列: {A} {C} {G}

长度为2的公共子序列：{AC} {CC} {CG} {AG}

长度为3的公共子序列：{ACG}

长度为4的公共子序列

最长公共子序列即为 {ACG}

问题：长度为N和M的两个序列如何求他们的最长公共子序列？

X = ACCGGGTTACCGTTTAAAACCCGGGTAACCT

 

Y = CCAGGACCAGGGACCGTTTACCAGCCTTAAACCA

 

简单算法

for (int i=N; i>0; --i) {

　　找到X所有长度为i的子序列;

　　找到Y所有长度为i的子序列;

　　如果存在公共子序列则终止

}

复杂度O(2^N)

动态规划算法

定义C[i][j]为序列X[i..i]和Y[1..j]的最长子序列的长度

C[i][0]和C[0][j]==0

递推公式

C[i][j]=C[i-1][j-1]+1,　　X[i]==Y[j]

C[i][j]=MAX(C[i-1][j], C[i][j-1]),　　X[i]!=Y[j]

C++实现

template<typename T>
void prtM(vector< vector<T> >& arr, string msg = "") {
    cout << msg << " " << endl;
    for  (auto i: arr) {
        for (auto j: i) {
            cout << j << " ";
        }
        cout << endl;
    }
    cout << endl;
}

template <typename T>
void prt(vector<T>& arr, string msg = "") {
    cout << msg << " ";
    for  (auto i: arr) {
        cout << i << " ";
    }
    cout << endl;
}


void prt_LCS(vector< vector<char> >& S, string& X, int i, int j) {　// 递归调用
    //cout << "i=" << i << " j=" << j << "  " << S[i][j] << endl;
    if (i==0 || j == 0)
        return;
    if ('s' == S[i][j]) {
        prt_LCS(S, X, i-1, j-1);
        cout << X[i];
    } else if ('j' == S[i][j])  {
        prt_LCS(S, X, i, j-1);
    } else {
        prt_LCS(S, X, i-1, j);
    }
}


void calc_LCS(string& X,  string& Y) {
    cout << "X: " << X << endl;
    cout << "Y: " << Y << endl;

    vector< vector<int> > C;  // 序列X[0..i] and Y[0..j]的公共子序列的长度
    vector< vector<char> > S;　　//最长公共子串的位置
    for (int i=0; i<X.size(); i++) {
        C.push_back( vector<int>(Y.size()) );　　
        S.push_back( vector<char>(Y.size()) );
    }

    for (int i=0; i<X.size(); i++) 
        C[i][0] = 0;
    for (int j=0; j<Y.size(); j++) 
        C[0][j] = 0;

    for (int i=1; i<X.size(); i++)
        for (int j=1; j<Y.size(); j++) {
            if (X[i] == Y[j]) {
               C[i][j] = C[i-1][j-1] + 1;
               S[i][j] = 's';  // 相同，可以打印
            } else if ( C[i][j-1] > C[i-1][j] ) {
                C[i][j] = C[i][j-1];
                S[i][j] = 'j';　　//Y的最后一个可以去掉，不影响LCS
            } else {
                C[i][j] = C[i-1][j] ;
                S[i][j] = 'i';　　//X的最后一个可以去掉，不影响LCS
            }
        }

    prtM(C);
    prtM(S);
    prt_LCS(S, X, X.size()-1, Y.size()-1);
    cout<< endl;
}

int main() {

    string S1 = " ACCGGGTTAC";  // 为了方便表示，忽略第一个字符
    string S2 = " AGGACCA";

    calc_LCS(S1, S2);

    return 0;
}