12

問題: 2 つの文字列間の LCS の長さが必要です。文字列のサイズは最大 100 文字です。アルファベットはいつものDNAのもので、4文字の「ACGT」。動的なアプローチは十分に迅速ではありません。

私の問題は、たくさんのペア (私が見る限り数億のランク) を扱っていることです。LCS_length 関数の呼び出しを可能な限り最小限に減らしたと思うので、プログラムの実行を高速化する他の唯一の方法は、より効率的な LCS_Length 関数を使用することです。

通常の動的プログラミングのアプローチで実装することから始めました。それは正しい答えを与え、うまくいけば適切に実装されます。

#define arrayLengthMacro(a) strlen(a) + 1
#define MAX_STRING 101

static int MaxLength(int lengthA, int lengthB);

/* 
 * Then the two strings are compared following a dynamic computing
 * LCS table algorithm. Since we only require the length of the LCS 
 * we can get this rather easily.
 */
int LCS_Length(char *a, char *b)
{
    int lengthA = arrayLengthMacro(a),lengthB = arrayLengthMacro(b), 
        LCS = 0, i, j, maxLength, board[MAX_STRING][MAX_STRING];

        maxLength = MaxLength(lengthA, lengthB);

    //printf("%d %d\n", lengthA, lengthB);
    for (i = 0; i < maxLength - 1; i++)
    {
        board[i][0] = 0;
        board[0][i] = 0;
    }

    for (i = 1; i < lengthA; i++)
    {
        for (j = 1; j < lengthB; j++)
        {
/* If a match is found we allocate the number in (i-1, j-1) incremented  
 * by 1 to the (i, j) position
 */
            if (a[i - 1] == b[j - 1])
            {

                board[i][j] = board[i-1][j-1] + 1;
                if(LCS < board[i][j])
                {
                    LCS++;
                }
            }
            else
            {
                if (board[i-1][j] > board[i][j-1])
                {
                    board[i][j] = board[i-1][j];
                }
                else
                {
                    board[i][j] = board[i][j-1];
                }
            }
        }
    }

    return LCS;
}

それは O(mn) である必要があります (うまくいけば)。

次に、速度を求めて 、Myers によるO(ND) 論文 を提供するこの投稿Longest Common Subsequenceを見つけました。LCSを最短の編集スクリプト(SES)に関連付けるこれを試しました。彼らが与える関係は D = M + N - 2L で、D は SES の長さ、M と N は 2 つのストリングの長さ、L は LCS の長さです。しかし、これは私の実装では常に正しいとは限りません。私は私の実装を与えます(私は正しいと思います):

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define arrayLengthMacro(a) strlen(a)

int LCS_Length (char *a, char *b);
int MinLength (int A, int B);
int Max (int A, int B);
int snake(int k, int max, char *a, char *b, int lengthA, int lengthB);

int main(void)
{
    int L;
    char a[] = "tomato", b[] = "potato"; //must give LCS = 4
    L =  LCS_Length(a, b);
    printf("LCS: %d\n", L );    
    char c[] = "GTCGTTCGGAATGCCGTTGCTCTGTAAA", d[] = "ACCGGTCGAGTGCGCGGAAGCCGGCCGAA"; // must give LCS = 20
    L =  LCS_Length(c, d);
    printf("LCS: %d\n", L );
    char e[] = "human", f[] = "chimpanzee"; // must give LCS = 4
    L =  LCS_Length(e, f);
    printf("LCS: %d\n", L );
    char g[] = "howareyou", h[] = "whoareyou"; // LCS =8
    L =  LCS_Length(g, h);
    printf("LCS: %d\n", L );
    char i[] = "TTCTTTCGGTAACGCCTACTTTATGAAGAGGTTACATTGCAATCGGGTAAATTAACCAACAAGTAATGGTAGTTCCTAGTAGAGAAACCCTCCCGCTCAC", 
        j[] = "GCACGCGCCTGTTGCTACGCTCTGTCCATCCTTTGTGTGCCGGGTACTCAGACCGGTAACTCGAGTTGCTATCGCGGTTATCAGGATCATTTACGGACTC"; // 61
    L =  LCS_Length(i, j);
    printf("LCS: %d\n", L );


    return 0;
}

int LCS_Length(char *a, char *b)
{

    int D, lengthA = arrayLengthMacro(a), lengthB = arrayLengthMacro(b), 
        max, *V_, *V, i, k, x, y;

    max = lengthA + lengthB;
    V_ = malloc(sizeof(int) * (max+1));
    if(V_ == NULL)
    {
        fprintf(stderr, "Failed to allocate memory for LCS");
        exit(1);
    }
    V = V_ + lengthA;
    V[1] = 0;

    for (D = 0; D < max; D++)
    {
        for (k = -D; k <= D; k = k + 2)
        {
            if ((k == -D && V[k-1] < V[k+1]) || (k != D && V[k-1] < V[k+1]))
            {
                x = V[k+1];
            }
            else
            {
                x = V[k-1] + 1;
            }
            y = x - k;
            while ((x < lengthB) && (y < lengthA) && (a[x+1] == b[y+1]))
            {
                x++;
                y++;
            }
            V[k] = x;
            if ((x >= lengthB) && (y >= lengthA))
            {
                return (lengthA + lengthB - D)/2;
            }
        }
    }
    return  (lengthA + lengthB - D)/2;
}

主に例があります。例えば。"tomato" と "potato" (コメントしないでください) の LCS の長さは 4 です。実装では、SES (コード内の D) が 4 ではなく 2 として出力されるのは 5 サイスであることがわかります ( 「t」を削除、「p」を挿入、「m」を削除、「t」を挿入)。おそらく O(ND) アルゴリズムも置換をカウントするのではないかと思う傾向がありますが、これについてはよくわかりません。

実装可能なアプローチ (私は膨大なプログラミング スキルを持っていません) であれば、どんなアプローチでも大歓迎です。(たとえば、小さなアルファベットを利用する方法を誰かが知っている場合)。

編集: いつでも任意の 2 つの文字列間で LCS 関数を使用することを何よりも強調すると便利だと思います。したがって、数百万の他の文字列と比較して、s1 などの文字列だけではありません。s1000 で s200、s10000 で s0、s100000 で 250 の可能性があります。ほとんどの使用文字列も追跡できない可能性があります。近似アルゴリズムを実装しており、余分なエラーを追加したくないため、LCSの長さは近似ではないことが必要です。

編集: callgrind を実行しました。10000 個の文字列の入力の場合、文字列のさまざまなペアに対して、lcs 関数を約 50,000,000 回呼び出すようです。(10000 文字列は実行したい最小値であり、最大値は 100 万です (実行可能な場合))。

4

3 に答える 3

1

私は、LCS を計算するための動的プログラミングよりも手の込んだアルゴリズムに精通していませんが、いくつかのことを指摘したいと思います。

まず、O(ND) アプローチは、非常に大きく、非常に類似した文字列を比較する場合にのみ意味があります。これはあなたには当てはまらないようです。

第二に、LCD_Length 関数の漸近的なパフォーマンスを高速化することは、文字列がかなり短いため、おそらく注目すべきことではありません。類似または非類似のペア (およびすべてのペアの正確な LCS ではない) を見つけることのみに関心がある場合は、Yannick が言及した BK ツリーが有望な方法のように見えます。

最後に、DP の実装について気になることがいくつかありました。コード内の「board[i][j]」の正しい解釈は、「文字列 a[1..i]、b[1..j] の最長のサブシーケンス」です (この表記では 1 インデックスを使用しています)。 )。したがって、for ループには i = lengthA と j = lengthB を含める必要があります。arrayLengthMacro(a) を導入することでコード内のこのバグをハッキングしたようですが、そのハックはアルゴリズムのコンテキストでは意味がありません。「ボード」が満たされると、ボード[長さA][長さB]でソリューションを検索できます。つまり、不要な「if (LCS < ボード[i][j])」ブロックを取り除き、ボード[を返すことができます。長さA][長さB]。また、ループ境界は初期化で間違っているように見えます ((i = 0; i <= maxLength; i++) ここで maxLength = max(lengthA, lengthB)) である必要があると確信しています)。

于 2011-07-02T09:03:05.917 に答える
1

計算を高速化するには、いくつかの方法があります。

  • 単純な動的プログラミングの代わりに、A* 検索を使用できます ((i,j) までの部分的なアラインメントには必ず |ij| の削除または挿入が必要であるというヒューリスティックを使用して)。
  • 1 つのシーケンスを他のシーケンス全体と比較する場合、そのプレフィックスに至る部分の動的計画法テーブル (または A* 検索状態) を計算し、計算の一部を再利用することで時間を節約できます。 . 動的プログラミング テーブルに固執する場合は、文字列のライブラリを辞書順に並べ替え、変更された部分のみを破棄できます (たとえば、「バナナ」があり、それを「パナマ」や「パナメリカニズム」と比較したい場合、 「panam」をカバーするテーブルの部分を再利用できます)。
  • ほとんどの文字列が非常に似ている場合は、共通の接頭辞を探して計算から共通の接頭辞を除外することで時間を節約できます (たとえば、「パナマ」と「パナメリカニズム」の LCS は、共通の接頭辞「パナム」に「パナム」の LCS を加えたものです)。 「a」と「エリカニズム」)
  • 文字列が非常に似ていない場合は、記号カウントを使用して編集回数の下限を取得できます (たとえば、"AAAB" から "ABBB" への変換には少なくとも 2 回の編集が必要です。もう一方)。これは、A* 検索のヒューリスティックでも使用できます。

編集:同じセットの文字列と比較する場合、 サンプルサイズが大きい場合に文字列の類似度スコアを計算する効率的な方法で BK-Tree データ構造を提案した人がいますか?

于 2011-07-02T08:29:35.733 に答える
0

文字列、ツリー、シーケンスに関するガスフィールドのアルゴリズムのコピーを入手することをお勧めします。これは、計算生物学の文字列操作に関するものです。

于 2011-07-02T08:37:23.917 に答える