BIOINFORMÁTICA

El algoritmo de Needleman-Wunsch sirve para realizar alineamientos globales de dos secuencias. Se suele utilizar en el ámbito de la bioinformática para alinear secuencias de proteínas o de ácidos nucleicos. Fue propuesto por primera vez en 1970, por Saul Needleman y Christian Wunsch. Se trata de un ejemplo típico deprogramación dinámica. El algoritmo funciona del mismo modo independientemente de la complejidad o longitud de las secuencias y garantiza la obtención del mejor alineamiento.¹​

Las dos secuencias a alinear, llamadas A y B en los ejemplos, de longitud ${\displaystyle |A|=m}$ y ${\displaystyle |B|=n}$, están formadas por elementos de un alfabeto finito de símbolos. El algoritmo necesita saber qué símbolos son diferentes entre sí y cuáles son iguales. Podemos utilizar una matriz cuadrada (S) para este propósito, en la que cada elemento ${\displaystyle S_{ij}}$indique la similitud entre los elementos i y j del alfabeto usado. Si nuestro alfabeto de símbolos no fuese finito, en vez de una matriz podríamos usar una función ${\displaystyle R^{2}\rightarrow R}$ que tuviese como parámetros ambos símbolos a comparar y cuya salida fuese la similitud entre ambos. También se necesita otro parámetro (d) que nos indique cómo vamos a valorar que un símbolo no quede alineado con otro y que en su lugar se utilice un hueco.

Por ejemplo podemos definir la siguiente matriz:

${\displaystyle {\begin{array}{ccccc}-&A&G&C&T\\A&10&-1&-3&-4\\G&-1&7&-5&-3\\C&-3&-5&9&0\\T&-4&-3&0&8\end{array}}}$

Y entonces el siguiente alineamiento:

AGACTAGTTAC
CGA---GACGT

con una penalización por hueco de ${\displaystyle d=-5}$ nos devolvería como solución óptima:

${\displaystyle S(A,C)+S(G,G)+S(A,A)+3\times d+S(G,G)+S(T,A)+S(T,C)+S(A,G)+S(C,T)=-3+7+10-3\times 5+7+-4+0+-1+0=1}$

Para determinar la puntuación óptima y poder reconstruir el alineamiento que devolvería esa puntuación se necesita otra matriz, F, que almacena los resultados parciales de cada posible alineamiento. Las dimensiones de la matriz F son el número de elementos en la secuencia A y el de B (${\displaystyle |A|\times |B|}$).

En cada iteración del algoritmo recibe valor un elemento de la matriz F. El valor que recibe el elemento ${\displaystyle F_{ij}}$representa la puntuación obtenida al alinear de forma óptima los primeros i elementos de A y los primeros j de B. Cuando el algoritmo termine, el último elemento de F (${\displaystyle F_{mn}}$, con ${\displaystyle m=|A|}$ y ${\displaystyle n=|B|}$) contendrá la puntuación para el alineamiento óptimo de ambas secuencias.

  Inicio del algoritmo:
  ${\displaystyle F_{0j}=d*j}$
  ${\displaystyle F_{i0}=d*i}$
  Recursión para obtener el siguiente elemento de forma óptima:
  ${\displaystyle F_{ij}=\max(F_{i-1,j-1}+S(A_{i},B_{j}),F_{i,j-1}+d,F_{i-1,j}+d)}$

La matriz F se calcula con el siguiente algoritmo:

   for i=0 to length(A)-1
     F(i,0) <- span=""> d*i
   for j=0 to length(B)-1
     F(0,j) <- span=""> d*j
   for i=1 to length(A)
     for j = 1 to length(B)
     {
       Choice1 <- span=""> F(i-1,j-1) + S(A(i), B(j))
       Choice2 <- span=""> F(i-1, j) + d
       Choice3 <- span=""> F(i, j-1) + d
       F(i,j) <- span=""> max(Choice1, Choice2, Choice3)
     }

Cuando el algoritmo acaba tenemos calculada la matriz F; el resultado es la puntuación devuelta por el mejor alineamiento posible, de acuerdo a los parámetros que hemos definido. Para obtener la secuencia se necesita ejecutar el siguiente algoritmo, que hace uso de la matriz F. Este algoritmo comienza por el último elemento, ${\displaystyle F_{mn}}$, y va retrocediendo hasta llegar a un elemento de la primera fila o la primera columna de F. En cada paso se comparan 3 elementos de F para ver cuál de ellos es el que se ha seguido en la solución óptima. Para cada ${\displaystyle F_{ij}}$debemos comparar ${\displaystyle F_{i-1,j},F_{i,j-1}}$ y ${\displaystyle F_{i-1,j-1}}$. Si el elemento usado es ${\displaystyle F_{i-1,j}}$, entonces ${\displaystyle A_{i}}$ se ha alineado con un hueco; si es ${\displaystyle F_{i,j-1}}$, entonces ${\displaystyle B_{i}}$ se ha alineado con un hueco; y si no, si el elemento elegido es ${\displaystyle F_{i-1,j-1}}$, los elementos ${\displaystyle A_{i}}$ y ${\displaystyle B_{i}}$ han sido alineados. Es importante destacar que el que dos elementos sean alineados no implica necesariamente que sean iguales; significa que entre esa posibilidad, alinear con huecos o alinear símbolos diferentes, esa era la mejor opción. El pseudo-algoritmo que permite obtener el alineamiento correcto es el siguiente:

   AlignmentA <- span=""> ""
   AlignmentB <- span=""> ""
   i <- span=""> length(A) 
   j <- span=""> length(B)
   while (i > 0 AND j > 0)
   {
     Score <- span=""> F(i,j)
     ScoreDiag <- span=""> F(i - 1, j - 1)
     ScoreUp <- span=""> F(i, j - 1)
     ScoreLeft <- span=""> F(i - 1, j)
     if (Score == ScoreDiag + S(A(i), B(j)))
     {
       AlignmentA <- span=""> A(i-1) + AlignmentA
       AlignmentB <- span=""> B(j-1) + AlignmentB
       i <- span=""> i - 1
       j <- span=""> j - 1
     }
     else if (Score == ScoreLeft + d)
     {
       AlignmentA <- span=""> A(i-1) + AlignmentA
       AlignmentB <- span=""> "-" + AlignmentB
       i <- span=""> i - 1
     }
     otherwise (Score == ScoreUp + d)
     {
       AlignmentA <- span=""> "-" + AlignmentA
       AlignmentB <- span=""> B(j-1) + AlignmentB
       j <- span=""> j - 1
     }
   }
   while (i > 0)
   {
     AlignmentA <- span=""> A(i-1) + AlignmentA
     AlignmentB <- span=""> "-" + AlignmentB
     i <- span=""> i - 1
   }
   while (j > 0)
   {
     AlignmentA <- span=""> "-" + AlignmentA
     AlignmentB <- span=""> B(j-1) + AlignmentB
     j <- span=""> j - 1
   }

Se puede demostrar formalmente que tanto el tiempo de ejecución como el espacio necesario para ejecutar el algoritmo son de orden O(nm). Para alguna aplicaciones, sobre todo en bioinformática, el requerimiento de espacio es prohibitivo, puesto que se alinean secuencias muy largas. Existe una optimización de este algoritmo, denominada algoritmo de Hirschberg, que sólo necesita espacio del orden O(m+n), pero a costa de incrementar el tiempo de computación.

El algoritmo Needleman-Wunsch realiza un alineamiento global de dos secuencias (aquí llamadas A y B).

Es comúnmente usado en bioinformática para alinear secuencias de nucleótidos o proteínas. Fue propuesto en 1970 por Saul Needleman y Christian Wunsch en su paper A general method applicable to the search for similarities in the amino acid sequence of two proteins, J Mol Biol. 48(3):443-53.

El algoritmo Needleman–Wunsch es un ejemplo de programación dinámica, y está garantizado que encuentre el alineamiento con el puntaje máximo. Needleman–Wunsch fue la primera aplicación de programación dinámica para la comparación de secuencias biológicas.

Los puntajes para caracteres alineados son especificados por una matriz de similitud. Aquí, S(i,j) es la similitud de los caracteres i y j. Esta usa una penalidad por hueco (gap) lineal, aquí llamada d.

Por ejemplo, si la matriz de similitud era:

Entonces el alineamiento es:

AGACTAGTTAC

CGA—GACGT

Con una penalidad por hueco de -5, tendríamos el siguiente puntaje…

Para encontrar el alineamiento con el puntaje más alto, una matriz bidimensional es asignada. Esta matriz a menudo es llamada matriz F y su (i,j)ésima entrada frecuentemente es denotada Fij . Allí hay una columna para cada carácter de la secuencia A, y una fila para cada caracter de la secuencia B. Así si estamos alineando secuencias de tamaños n y m, el tiempo de ejecución del algoritmo es O(nm) y la cantidad de memoria utilizada es O(nm).

Sin embargo hay una versión modificada del algoritmo que usa solo O(n+m) espacio, al costo de un tiempo de ejecución más grande. Esta modificación es de hecho una técnica general que aplicamos a muchos algoritmos de programación dinámica; este método fue introducido en el algoritmo de Hirschberg para resolver el problema de la subcadena más larga en común.

Cuando el algoritmo progresa, Fij puede ser asignada para ser el puntaje óptimo para el alineamiento de los primeros i caracteres en A y los primeros j caracteres en B. El principio de optimización es entonces aplicado como sigue:

Base:

F0j = d * j

Fi0 = d * i

Recursión, basada en el principio de optimización:

Fij = max(Fi − 1,j − 1 + S(Ai − 1,Bj − 1),Fi,j − 1 + d,Fi − 1,j + d)

El pseudo-código para el algoritmo que computa la matriz A por lo tanto luce algo así (los índices de arreglos y secuencias empiezan en 0):

for i=0 a largo(A)-1

F(i,0) <- d="" i="" p="">

for j=0 a largo(B)-1

F(0,j) <- d="" j="" p="">

for i=1 a largo(A)

for j = 1 a largo(B)

{

Elección1 <- b="" f="" i-1="" j-1="" p="" s="">

Elección2 <- d="" f="" i-1="" j="" p="">

Elección3 <- d="" f="" i="" j-1="" p="">

F(i,j) <- elecci="" lecci="" max="" n1="" n2="" n3="" p="">

}

Una vez que la matriz F es computada, note que la esquina inferior derecha de la matriz es el máximo puntaje para cualquier alineamiento. Para computar que alineamiento da actualmente ese puntaje, puedes empezar desde la celda que se encuentra al fondo a la derecha, y comparar el valor con las tres posibles fuentes (Elección1, Elección2, Elección3) para ver de donde proviene.

Si era Elección1, entonces A(i) y B(i) están alineadas, si era Elección2 entonces A(i) está alineado con un gap, y si era Elección3, entonces B(i) está alineada con un hueco (gap).

AlineamientoA <- p="">

AlineamientoB <- p="">

i <- largo="" p="">

j <- largo="" p="">

while (i > 0 AND j > 0)

{

Score <- f="" i="" j="" p="">

ScoreDiag <- 1="" f="" i="" j="" p="">

ScoreUp <- 1="" f="" i="" j="" p="">

ScoreLeft <- 1="" f="" i="" j="" p="">

if (Score == ScoreDiag + S(A(i-1), B(j-1)))

{

AlineamientoA <- a="" alineamientoa="" i-1="" p="">

AlineamientoB <- alineamientob="" b="" j-1="" p="">

i <- 1="" i="" p="">

j <- 1="" j="" p="">

}

else if (Score == ScoreLeft + d)

{

AlineamientoA <- a="" alineamientoa="" i-1="" p="">

AlineamientoB <- alineamientob="" p="">

i <- 1="" i="" p="">

}

otherwise (Score == ScoreUp + d)

{

AlineamientoA <- alineamientoa="" p="">

AlineamientoB <- alineamientob="" b="" j-1="" p="">

j <- 1="" j="" p="">

}

}

while (i > 0)

{

AlineamientoA <- a="" alineamientoa="" i-1="" p="">

AlineamientoB <- alineamientob="" p="">

i <- 1="" i="" p="">

}

while (j > 0)

{

AlineamientoA <- alineamientoa="" p="">

AlineamientoB <- alineamientob="" b="" j-1="" p="">

j <- 1="" j="" p="">