Автор: David Eppstein
Перевод: Кантор И.

Для начала определим, в чем заключается задача:

     Пусть у нас есть некоторый образец - короткая строка и текст - длинная.
Но, в отличие от точного поиска подстроки в строке, мы хотим знать, появляются ли буквы образца в том же порядке ( но, возможно, на разном расстоянии ) в тексте.
     Если да, то образец - подпоследовательность текста.

     Например: Является ли 'nano' подпоследовательностью 'nematode knowledge' ?
     Да, причем единственным образом.
Вот строка с выделенной подпоследовательностью: 'NemAtode kNOwledge'.

     Псевдокод этой задачи можно записать так:

    subseq(char * P, char * T)
    {
        while (*T != '\0')
            if (*P == *T++ && *++P == '\0')
                return TRUE;
        return FALSE;
    }
   

А теперь собственно наша задача.

     Что, если образец не появляется в тексте ? Тогда часто бывает нужно найти самую длинную подпоследовательность, которая появляется и там и там.

     Текст и образец теперь равноправны, поэтому обозначим их просто строками A и В. Длина A - m, B - n символов.

Рассмотрим рекурсивное решение.

     Попробуем решить нашу задачу, используя рекурсию, а далее улучшить его через динамическое программирование.

     Для начала - несколько простых наблюдений.

     Пусть у нас есть две строки, например 'nematode knowledge' и 'empty bottle'. Тогда мы можем представить подпоследовательность как способ записи этих строк так, что определенные буквы будут одна под другой:

    n e m a t o d e [пробел] k n o  w  l e d g e

      | |   |          |         |     | |

      e m p t y     [пробел]   b o t t l e

Рекурсивное решение:

    int lcs_length(char * A, char * B)

    {
	if (*A == '\0' || *B == '\0') return 0;
	else if (*A == *B) return 1 + lcs_length(A+1, B+1);
	else return max(lcs_length(A+1,B), lcs_length(A,B+1));
    }

Запоминание.

     Недостаток рекурсивного решения в том, что одинаковые подзадачи вызываются много раз в разное время.      Подзадача состоит в вызове lcs_length с суффиксами A и B в качестве аргументов, так что есть всего (m+1)(n+1) возможных подзадач. Это - весьма маленькое число. При ~2ⁿ рекурсивных вызовах некоторые из них будут решаться много раз заново.

     Решение, использующее динамическое программирование, состоит в проверке, в каком случае нам нужно решать подзадачу, а в каком - она уже решена. Тогда мы просто посмотрим вверх, на предыдущее решение, не вычисляя его еще раз.

     Делая это 'в лоб', мы просто добавим в наш рекурсивный алгоритм код, который будет смотреть вверх. Такая рекурсивная версия динамического программирования - 'сверху вниз' называется запоминанием.

     Подзадачи состоят из двух суффиксов наших строк. Чтобы хранить и заглядывать в решения подзадач было легче, будем представлять их начальными позициями в строках, а не указателями на символы (как было раньше).

    char * A;
    char * B;
    int lcs_length(char * AA, char * BB)
    {
	A = AA; B = BB;
	return subproblem(0, 0);
    }

    int subproblem(int i, int j)
    {
	if (A[i] == '\0' || B[j] == '\0') return 0;
	else if (A[i+1] == B[j+1]) return 1 + lcs_length(i+1, j+1);
	else return max(lcs_length(i+1, j), lcs_length(i, j+1));
    }
  

    char * A;
    char * B;
    array L;

    int lcs_length(char * AA, char * BB)
    {
	A = AA; B = BB;
	allocate storage for L;
	for (i = 0; i <= m; i++)
	    for (j = 0; j <= m; j++)
		L[i,j] = -1;
	return subproblem(0, 0);
    }

    int subproblem(int i, int j)
    {
	if (L[i,j] < 0) {
	    if (A[i] == '\0' || B[j] == '\0') L[i,j] = 0;
	    else if (A[i+1] == B[j+1]) L[i,j] = 1 + lcs_length(i+1, j+1);
	    else L[i,j] = max(lcs_length(i+1, j), lcs_length(i, j+1));
	}
	return L[i,j];
    }
                         

Динамическое программирование снизу вверх.

     То, что было описано выше - всего лишь более умный путь создания рекурсивного алгоритма. Но можно считать эту программу способом заполнения массива L. Последовательность заполнения контролирует рекурсия, однако те же результаты можно получить и при другом порядке.

     Мы можем придумать что-нибудь попроще, что будкт этим заниматься. Единственная проблема - для заполнения L[i,j] нам нужно знать то, от чего зависит его значение: в данном случае, это L[i+1,j], L[i,j+1] и L[i+1,j+1].

     А потому - будем идти по массиву назад, от последнего ряда к первому, и от последней колонке - вверх к первой. Наш алгоритм станет итерационным, используя циклы вместо рекурсии. Или алгоритмом снизу-вверх, так как мы будем заполнять массив от простых подзадач к более сложным.

Итерационное LCS ( Longest Common Subsequence ):

    int lcs_length(char * A, char * B)

    {
	allocate storage for array L;
	for (i = m; i >= 0; i--)
	    for (j = n; j >= 0; j--)
	    {
		if (A[i] == '\0' || B[j] == '\0') L[i,j] = 0;
		else if (A[i+1] == B[j+1]) L[i,j] = 1 + L[i+1, j+1];
		else L[i,j] = max(L[i+1, j], L[i, j+1]);
	    }
	return L[0,0];
    }                        

Сама подпоследовательность.

     Что, если нам нужна сама подпоследовательность, а не ее длина, как раньше ?
Заполнив массив L, как показано выше, мы можем ее найти, просто просмативая его.

    sequence S = empty
    i = 0
    j = 0

    while min(i,j) > 0
    {
	if (A[i]==B[j])
	{
	    add A[i] to end of S;
	    i--; j--;
	}
	else if (L[i+1,j] >= L[i,j+1]) i++;
	else j++;
    }
  

        n  e  m  a  t  o  d  e  _  k  n  o  w  l  e  d  g  e

    e   7  7  6  5  5  5  5  5  4  3  3  3  2  2  2  1  1  1  0

    m   6  6  6  5  5  4  4  4  4  3  3  3  2  2  1  1  1  1  0

    p   5  5  5  5  5  4  4  4  4  3  3  3  2  2  1  1  1  1  0

    t   5  5  5  5  5  4  4  4  4  3  3  3  2  2  1  1  1  1  0

    y   4  4  4  4  4  4  4  4  4  3  3  3  2  2  1  1  1  1  0

    _   4  4  4  4  4  4  4  4  4  3  3  3  2  2  1  1  1  1  0

    b   3  3  3  3  3  3  3  3  3  3  3  3  2  2  1  1  1  1  0

    o   3  3  3  3  3  3  3  3  3  3  3  3  2  2  1  1  1  1  0

    t   3  3  3  3  3  2  2  2  2  2  2  2  2  2  1  1  1  1  0

    t   3  3  3  3  3  2  2  2  2  2  2  2  2  2  1  1  1  1  0

    l   2  2  2  2  2  2  2  2  2  2  2  2  2  2  1  1  1  1  0

    e   1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  0

        0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0

     Чтобы найти наибольшую общую подпоследовательность, посмотрим на первую ячейку L[0,0]. Там 7, значит в последовательности 7 символов. L[0,0] была вычислена как max( L[0,1] , L[1,0] ), которые отвечают подзадачам при удалении 'n' из первой строки или 'e' из второй.
     Удаление 'n' дает нам подпоследовательность длины L[0,1]=7, в то время как удаление 'e' - только L[1,0]=6, так что мы можем удалить лишь 'n'. Посмотрим на ячейку L[0,1], оcтающуюся после этого удаления. A[0]=B[1]='e', так что мы можем без проблем включить это 'e' в подпоследовательность и перейти к L[1,2]=6. Аналогично, это дает 'm'... Продолжая в том же духе ( и делая заметки, как в алгоритме выше, двигаясь вниз, а не наискосок ) мы получаем общую подпоследовательность: 'emt ole'.

     Итак, мы нашли LCS за время O(mn). Вообще, взглянув на матрицу, можно заметить, что она имеет весьма странную структуру: большие блоки с одинаковыми значениями - и маленькое количесво 'углов', где значения меняются. Если это использовать, то можно получить асимптотику O(n log s + c log log min(c,mn/c)), где с - число таких углов, а s - число символов, появляющихся в двух строках.

Уменьшаем требоватия к памяти.

     Один из недостатков описанного метода динамического программирования, по сравнению с исходной рекурсией - в том, что он требует O(mn) памяти на массив L ( рекурсия - O(n+m) ).
Но итерационная версия вполне может быть подредактирована, чтобы использовать меньше пространства: после вычисления i-го ряда нам уже не нужны значения в ряду i+1.

Не требовательное к памяти LCS:

    int lcs_length(char * A, char * B)
    {
	allocate storage for one-dimensional arrays X and Y
	for (i = m; i >= 0; i--)
	{
	    for (j = n; j >= 0; j--)
	    {
		if (A[i] == '\0' || B[j] == '\0') X[j] = 0;
		else if (A[i+1] == B[j+1]) X[j] = 1 + Y[j+1];
		else X[j] = max(Y[j], X[j+1]);
	    }
	    Y = X;
	}
	return X[0];
    }