La più lunga sottosequenza comune

Problema: Date due sequenze di elementi contigui $X=[x_0\dots x_m]$ e $Y=[y_0\dots y_n]$ trovare la più grande sequenza $Z=[z_0\dots z_p]$ di elementi comuni e ordinati tra $X$ e $Y$ ma non necessariamente contigui.

Esempio, con $X=ADEBACAEEB$ e $Y=CDBBACACEB$, $Z=DBACAEB$


Soluzione: Programmazione dinamica
Complessità: $O(mn)$ per lcs_lenght, $O(m + n)$ per print_lcs

Indicando con $X_i$ il prefisso di una sequenza con i primi $i$ caratteri della sequenza stessa (ad esempio, $X_3=ADE$ con $X$ definito come nell'esempio precedente) si può ricavare il numero di elementi (indicandolo con $c[i,j]$) della più lunga sottosequenza comune tra $X_i$ ed $Y_i$ (sottosequenza che indicheremo con $Z_{ij}$) attraverso la formula:

$c[i, j]:= \begin{cases} 0, & \text{if $i = 0$ or $j = 0$,}\\ c[i-1, j-1] + 1, & \text{if $i, j > 0$ and $x_i = y_j$,}\\ max(c[i, j-1], c[i-1, j]), & \text{if $i, j > 0$ and $x_i \ne y_j$.} \end{cases}$

In effetti, se una delle due sequenze è vuota ($i = 0$ o $j = 0$) una sottosequenza comune non esiste nemmeno ed il numero di elementi della sottosequenza comune sarà zero.
Se $x_i = y_j$ (l’i-esimo e il j-esimo elemento delle due sottosequenze sono uguali) allora basta prendere il numero di elementi della più grande sottosequenza comune tra $X_{i-1}$ e $Y_{j-1}$ ed aggiungere 1 (l’elemento corrente in comune).
Se, in ultimo, $x_i \ne y_j$, non basta prendere semplicemente il numero di elementi della più grande sottosequenza comune tra $X_{i-1}$ e $Y_{j-1}$ perché facendolo si potrebbe escludere un potenziale elemento della sottosequenza comune. Ad esempio con $X=ABCDE$ e $Y=ABEA$ si avrebbe $Z=ABE$ ma se si va a calcolare $c[5, 4]$ (avendo $x_5 \ne y_4$) e si usa semplicemente il valore ritornato da $c[i-1, j-1]$, allora avremo che l’ultimo elemento di $X_5$ viene escluso e si avrebbe $Z=AB$. E' necessario quindi selezionare il massimo tra i numeri di elementi delle più grandi sottosequenze comuni tra $X_{i}$ e $Y_{j-1}$ e tra $X_{i-1}$ e $Y_{j}$ (ovvero, e più semplicemente, il massimo tra $c[i, j-1]$ e $c[i-1, j])$.

Il codice usato per l’implementazione usa un approccio bottom-up piuttosto che uno ricorsivo. $c[i, j]$ e $b[i, j]$ sono due matrici (array di array in questo caso) che conservano, nel primo caso, il numero di elementi della più grande sottosequenza comune tra $X_{i}$ e $Y_{j}$ e, nel secondo caso, dettagli implementativi utili a ricostruire la più grande sottosequenza comune. In $b[i, j]$ indichiamo il movimento relativo alla posizione corrente con una ‘d’ (che sta per diagonal, movimento in diagonale $\nwarrow\,$) se $x_i = y_j,\,$ una ‘u’ (che sta per up, movimento in alto $\uparrow\,$) se $c[i-1, j]\ge c[i, j-1],\,$ altrimenti con una ‘l’ (che sta per left, movimento a sinistra $\longleftarrow\,$).  Per maggiori dettagli si veda [1].

#include <iostream>
 
using namespace std;
 
void lcs_lenght(const char*, const char*, const int, const int, int*, char*);
void print_lcs(const char*, const char*, int, const int);
 
void print_mat_b(const char*, const int, const int);
void print_mat_c(const int*, const int, const int);
 
int main(int argc, char **argv)
{
 char X[10] = { 'A', 'D', 'E', 'B', 'A', 'C', 'A', 'E', 'E', 'B' };
 char Y[10] = { 'C', 'D', 'B', 'B', 'A', 'C', 'A', 'C', 'E', 'B' };
 char mat_b[10][10];
 int mat_c[10][10];
 
 lcs_lenght(X, Y, 10, 10, &mat_c[0][0], &mat_b[0][0]);
 print_lcs(&mat_b[0][0], X, 10 * 10, 1);
 std::cout << "\n\n\n";
 print_mat_c(&mat_c[0][0], 10, 10);
 std::cout << "\n\n\n";
 print_mat_b(&mat_b[0][0], 10, 10);
 
 return 0;
}
 
void lcs_lenght(const char* x, const char* y, const int m, const int n, int* c, char* b)
{
 // imposta la prima riga della matrice c[i, j] a zero (j = 0 vuol dire che la sequenza Y è vuota)
 // imposta la prima riga della matrice b[i, j] ad 'n' (nessun movimento)
 for (int i = 0; i < m; ++i)
 {
  c[i * m + 0] = 0;
  b[i * m + 0] = 'n';
 }
 
 // imposta la prima colonna della matrice c[i, j] a zero (i = 0 vuol dire che la sequenza X è vuota)
 // imposta la prima colonna della matrice b[i, j] ad 'n' (nessun movimento)
 for (int j = 0; j < n; ++j)
 {
  c[0 * n + j] = 0;
  b[0 * n + j] = 'n';
 }
 
 for (int i = 1; i < m; ++i)
 {
  for (int j = 1; j < n; ++j)
  {
   if (x[i] == y[j])
   {
    c[i*n + j] = c[(i - 1)*n + (j - 1)] + 1;
    b[i*n + j] = 'd';
   }
   else if (c[(i - 1)*n + j] >= c[i*n + (j - 1)])
   {
    c[i*n + j] = c[(i - 1)*n + j];
    b[i*n + j] = 'u';
   }
   else
   {
    c[i*n + j] = c[i*n + (j - 1)];
    b[i*n + j] = 'l';
   }
  }
 }
}
 
// Usa ricorsione per stampare gli elementi in ordine inverso, cioè corretto.
void print_lcs(const char* b, const char* x, int index, const int offset)
{
 index -= offset;
 
 if (index <= 0)
  return;
 
 if (b[index] == 'd')
 {
  print_lcs(b, x, index, 11);
  std::cout << x[index / 11];
 }
 else if (b[index] == 'u')
  print_lcs(b, x, index, 10);
 else
  print_lcs(b, x, index, 1);
}
 
void print_mat_b(const char* b, const int m, const int n)
{
 for (int i = 0; i < m; ++i)
 {
  for (int j = 0; j < n; ++j)
  {
   std::cout << b[i*n + j] << ' ';
  }
 
  std::cout << '\n';
 }
}
 
void print_mat_c(const int* c, const int m, const int n)
{
 for (int i = 0; i < m; ++i)
 {
  for (int j = 0; j < n; ++j)
  {
   std::cout << c[i*n + j] << ' ';
  }
 
  std::cout << '\n';
 }
}

Riferimenti:
[1] Introduzione agli Algoritmi e Strutture Dati - Cormen, Leiserson, Rivest, Stein