1. 程式人生 > >動態規劃2-最長公共子序列

動態規劃2-最長公共子序列

參考

http://open.163.com/newview/movie/free?pid=M6UTT5U0I&mid=M6V2U1HL4

問題是給定字串x和y,求出兩個當中最長的公共子序列。比如x=abcdef y=acefg,那麼他們的最長公共子序列就是acef。就是求x的所有可能的子字串與y所有的子字串匹配,如果相同,那麼就是一個公共子序列,然後求最長的一個。

建議觀看上面的公開課,講的非常好。本文思路是根據上面的公開課總結實踐的。

我們先看看求一個字串的所有子序列,有什麼規律。比如給定一個字串abc,那麼有多少子序列呢?a b c ab ac bc abc還有一個空,就是2^3。按照上面公開課講的非常通俗易懂,就是字串一共有n個字元,那麼每一位的字元都有兩個狀態----0-不是子序列的一員 1-是子序列的一員。那麼總共的種類就是2^n種。這是最壞的情況,字串中沒有相同的字元,如果有,結果會比這個少。但是計算就是這麼計算的。

瞭解了上面,那麼我們進一步判斷,如果給定一個x的子序列x1,那麼怎麼判斷y中是否含有呢?簡單的方法就是從x1的起始位置和y的起始位置開始判斷,如果相同,每個索引加一,直到x1/y的字串判斷結束。如果不匹配,那麼y索引加一,再從x1的開頭開始計算。

這種是窮舉法,肯定能得到結果,但是時間和空間都是最差的。我們先實現看看(窮舉法排列組合可以參考排列組合)

struct MNODE
{
    string substr;
    MNODE* pnext;
    MNODE()
    {
        pnext = nullptr;
    }
};
MNODE* msubstr(string& srcstr, int index, int sublen)
{
    MNODE* tmpp = nullptr;
    MNODE* prep = nullptr;
    if (index >= sublen && sublen > 0)
    {
        if (sublen == 1)
        {
            for (int i = 0; i < index; i++)
            {
                MNODE* ftpn = nullptr;
                if (tmpp == nullptr)
                {
                    tmpp = new MNODE;
                    prep = tmpp;
                    ftpn = tmpp;
                }
                else
                {
                    ftpn = new MNODE;
                    prep->pnext = ftpn;
                    prep = ftpn;
                }
                ftpn->substr = srcstr[i];
            }
        }
        else if (sublen > 1)
        {
            MNODE* nsub = msubstr(srcstr, index - 1, sublen - 1);
            tmpp = nsub;
            while (nsub != nullptr)
            {
                nsub->substr = nsub->substr + srcstr[index - 1];
                prep = nsub;
                nsub = nsub->pnext;
            }
            nsub = msubstr(srcstr, index - 1, sublen);
            prep->pnext = nsub;
            while (nsub != nullptr)
            {
                nsub = nsub->pnext;
            }
        }
    }
    return tmpp;
}
bool comps(string& sstr, string& dstr)
{
    int sindex = 0;
    int dindex = 0;
    for (int i = sindex; i < sstr.size(); i++)
    {
        if (dindex == dstr.size())
        {
            break;
        }
        for (int j = dindex; j < dstr.size(); j++)
        {
            if (sstr[i] == dstr[j])
            {
                sindex++;
dindex = j + 1; break; } } } if (sindex == sstr.size() && dindex <= dstr.size()) { return true; } return false; } int main() { string a = "abcdef"; string b = "acefg"; for (int i = a.size(); i > 0; i--) { MNODE* psubstr = msubstr(a, a.size(), i); while (psubstr != nullptr) { if (comps(psubstr->substr, b)) { cout << psubstr->substr << endl; } MNODE* tmpp = psubstr; psubstr = tmpp->pnext; delete tmpp; } } char inchar; cin >> inchar; }

根據排列組合的方式獲取所有的子串,然後匹配,可以列出所有的公共子串,按照從長倒短列印。

 

那麼我們有沒有更好的解法呢?按照上面的公開課,我們可以這樣計算。LCS(x,y)表示計算x和y的最長公共子串。那麼C(i,j)表示當前x的子串從0-i和y的子串從0-j的最長公共子串。那麼推導公式就是,如果x[i]=y[j],那麼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中0-i組成的字串與y中0-j組成的字串的最長公共子串,是有兩部分推導來的

如果x當前索引的字串等於y當前索引的字串,那麼最長公共子串就是i-1與j-1的最長公共子串加上當前的字元。這個很好理解,就是已知前面的最長公共子串,如果當前相等,那麼加上就好了。

如果不相等呢?我們需要分兩步判斷,就是先把i-1加一,算一下當前的最長公共子串,再把j-1加一算一下當前的最長公共子串,然後判斷哪一個長,就更新哪一個結果。

 

上面就是推導公式。

公開課中對其進行了證明

第一條證明,假設z(1,k)表示C(i,j)的最長公共子串,有k個元素,從1-k,當前x(i)==y(j),那麼z(k)==x(i)==y(j)。去掉當前的k,那麼LCS(i-1,j-1)就等於z(1,k-1)。對於這個,我們有一個問題,就是必須保證z(1,k-1)是LCS(i-1,j-1),不然我們後面的推論就是錯誤的,因為如果不是最長公共子串,那麼後續的加上1也不是。公開課中用了一個反證法,比如存在w是一個比z(1,k-1)更長的子串,那麼|w|(用|w|表示字串w的長度)一定大於k-1,那麼把w拼接上x(i)或y(j),那麼|w+x(i)|肯定大於k,這與我們一開始設定的z(1,k)是最長公共子串衝突。所以z(1,k-1)肯定是取出x(i)和y(j)之後剩下的字串的最長公共子串。同樣的方法可以證明第二條定理。

那麼動態規劃的特徵是什麼,對於什麼樣的問題可以使用動態規劃呢?

第一條就是存在最優子結構。也就是一個問題的最優解,包含了子問題的最優解。如果子問題不是最優解,那麼也不可能是總問題的最優解。

上面的公式最壞的情況是什麼?按照公開課說的,就是一直是第二條規則。因為第一條規則,不需要計算,直接加一,索引各加一;第二條規則是隻有一個索引加一,還需要分別計算兩次。

怎麼優化呢?按照公開課講的,我們可以看一下這個公式有沒有什麼規律,如果是第一條,那麼沒什麼可以省略的,因為已經很簡單了,直接更新資料就好了。如果是第二條呢?舉個例子,我們看一下第二部計算的分解

 

 

 我們沒有繼續往下分,這裡就可以看出,有重複,兩個[6,5],計算了兩遍,並且還是靠近樹的根的位置,就導致重複計算的內容更多,那麼如何優化呢?這裡大家估計都想到了,就是把一個[6,5]記錄下來,到下一個計算的時候,判斷一下,如果有,就不再重複計算了,這樣就用空間替換了時間,並且是非常可觀的,因為記錄結果花不了多少空間,但是確實能節省很多時間。

這裡就引出了動態規劃的第二個特徵,重疊子問題。就是我們在計算的過程中,因為分解問題導致會出現重複的子問題,我們可以使用備忘法,就是上面說的,記錄一下這個問題的答案,下次碰到直接獲取就好了。

但是這樣並沒有減少太多的運算,因為什麼呢?主要還是因為第二個規則,如果x(i)!=y(j),因為我們僅僅能利用i和j相同時的重複計算,但是並不能嚴格使用i-1和j-1的結果。比如x=abc123kmnp,y=kmnpabc456,那麼最長公共子串是kmnp,當我們計算到x和y的c的下一位時,所得的子串並不能通過前面的abc來判斷,因為隨著不斷的遍歷,最長公共子串是kmnp而不是abc,並且求kmnp與abc沒有任何關係。

這個方法就不實現了,每一小步都是上面的一個窮舉法,只不過如果碰到一樣的可以減少一次計算。

公開課最後給出了這個問題的最終方法。這個方法我在看一些演算法題題解的時候有人用到,當時大家都驚為天人,原來一切都還是出自於國外,可惜沒有一個作者給其他人指出這種思路出自哪裡,如何使用,都自己默默的接受了別人的讚許。

 上面是給定的一個例子求ABCBDAB與BDCABA的最長公共子序列。有多個解,比如BCBA和BDAB。看上面的圖,我們把兩個字串用一個二維陣列表示,每個元素是x(i)與y(j)的對應關係。陣列中灰色的是預設填充0,表示的是,x和y字串對應空字串的時候,那麼肯定都是沒有相同的元素。然後我們一行一行的比較(每行顏色我用了淡橙色和淡藍色區分),第一行對於x(0)和y(0),我們知道A與B不相等,所以他們的相同的子串的個數是0,也就是最長公共字序列是0,也就是空字串。x(1)與y(0)比較發現B與B是相等的,所以就在上一次計算得出的公共子序列基礎上加1(這樣我們前面的每次計算都是有用的,沒有額外的開銷),到了x(2)和y(0),因為C與B不相等,所以當前就等於前面已經求得的公共子序列,也就是1.這裡也用到了上面的推導公式,如果相同就加1,不相同就等於LCS(x(i-1),y(j-1))分別加上i和j,求得的兩個公共子序列的最大值。那好,到了x(3)和y(0),這時有一個B和B相等,為什麼還是1呢?根據公式推導如果x(i)==y(j),那麼當前的最長公共子串等於LCS(x(i-1),y(j-1))+1,因為要退兩個,也就是x(2)和y(-1),是0,所以LCS(x(3),y(0))還是1.從圖示中可以看出,如果相等,那麼就是左上角頂格的數字加一,如果不等,就是前一格或是上一格的最大值。前一格就是LCS(x(i-1), y(j)),上一格就是LCS(x(i), y(j-1));同樣相等的話就是左上角的加一,左上角就是LCS(x(i-1),y(j-1))。真是完美的解決方案,非常形象,也很好理解,也很好記。把數學和幾何完美的結合在一起。那麼最長公共子串怎麼求呢?很簡單,一行一行的掃描或是一排一排的掃描,在上一次結果所在的i,j位置往後,第一次數字增加的地方就是公共子串的位置,為了求最長的公共子串,我們可以得到最大的數值,然後逆向操作,第一次數字減少的地方就是公共子串的元素。

void lcs()
{
    string a = "abc123kmnp";
    string b = "kmnpabc456";
    int* lcsmap = new int[a.size() * b.size()]();
    int maxnum = 0;
    int maxi = 0;
    int maxj = 0;
    for (int j = 0; j < b.size(); j++)
    {
        for (int i = 0; i < a.size(); i++)
        {
            if (a[i] == b[j])
            {
                if (i > 0 && j > 0)
                {
                    *(lcsmap + a.size() * i + j) = *(lcsmap + a.size() * (i - 1) + j - 1) + 1;
                }
                else
                {
                    *(lcsmap + a.size() * i + j) = 1;
                }
            }
            else
            {
                if (i > 0 && j > 0)
                {
                    *(lcsmap + a.size() * i + j) = max((*(lcsmap + a.size() * (i - 1) + j)), (*(lcsmap + a.size() * i + j - 1)));
                }
                else if (i > 0)
                {
                    *(lcsmap + a.size() * i + j) = *(lcsmap + a.size() * (i - 1) + j);
                }
                else if (j > 0)
                {
                    *(lcsmap + a.size() * i + j) = *(lcsmap + a.size() * i + j - 1);
                }
            }
            if (*(lcsmap + a.size() * i + j) > maxnum)
            {
                maxnum = *(lcsmap + a.size() * i + j);
                maxi = i;
                maxj = j;
            }
        }
    }
    string slcs;
    for (int j = maxj; j >= 0; j--)
    {
        for (int i = maxi; i >= 0; i--)
        {
            if (a[i] == b[j] && *(lcsmap + a.size() * i + j) == maxnum)
            {
                maxi = i;
                maxnum--;
                slcs = a[i] + slcs;
                break;
            }
        }
    }
    cout << slcs;
}