Articles

alineación por pares de secuencias de nucleótidos utilizando coincidencias exactas máximas


enfoque

en nuestro algoritmo propuesto, el primer paso hacia la alineación de secuencias es extraer MEMs entre secuencias comparándolas directamente. La figura 3a es un ejemplo que compara un objetivo y una secuencia de consulta donde CTC y AAA son dos MEMs identificados por la comparación. Cada grupo de símbolos idénticos continuos en la comparación, resulta en un MEM incluso si está compuesto de un solo símbolo coincidente., Para extraer todos los MEMs entre las secuencias, la secuencia de consulta debe ser desplazada hacia la derecha y hacia la izquierda un símbolo a la vez (ver Fig. 3b). Después de cada turno, el paso de comparación debe repetirse para identificar nuevos MEMs. Por ejemplo, la tercera línea de la Fig. 3b representa el caso en el que la secuencia de consulta se desplaza hacia el símbolo de la derecha y se compara con la secuencia de destino. El resultado de la comparación identifica a AAAAGC como un nuevo MEM. Todos los demás MEMs extraídos por las operaciones shift y compare también se resaltan en la Fig. 3b., Tres de los MEMs (Mx, My y Mz) están resaltados con diferentes colores que se utilizarán para una explicación posterior.

Fig. 3

extracción de MEM usando las operaciones shift y compare. a identificar MEMs por comparación directa de secuencias. b la consulta se desplaza a la izquierda hasta que el último símbolo de la secuencia de consulta se alinea con el primer símbolo de la secuencia de destino. A continuación, la secuencia de consulta se desplaza a la derecha hasta que el primer símbolo de la secuencia de consulta se alinea con el último símbolo de la secuencia de destino., Después de cada cambio, la parte superpuesta de la consulta y las secuencias de destino se comparan para identificar nuevos MEMs. Tres de los MEMs (Mx, My y Mz) se resaltan con diferentes colores para ser utilizados para una explicación posterior

en el modelo de puntuación afín-gap, la puntuación de alineación se calcula utilizando la EC., 1 donde Nm es el número de partidos cada uno recibiendo una puntuación de partido de Rm, Nx es el número de desajustes cada uno recibiendo una penalización de desajuste de Px,No es el número de aperturas de huecos cada uno recibiendo una penalización de apertura de hueco de Po y Ng es la longitud total de todos los huecos, cada hueco recibiendo una penalización de extensión de hueco de Pg. Se abriría una brecha para cada grupo de brecha continua. Por ejemplo, si hay dos huecos en la alineación, donde la longitud del primer hueco es tres y la longitud del segundo hueco es cuatro, entonces hay dos huecos (No=2) y la longitud total del hueco es siete (Ng=3+4=7).,

$$ {}COMO = (N_{m} \times R_{m}) – ((N_{x} \times P_{x}) + (N_{o} \times P_{o}) + (N_{g} \times P_{g})) $$
(1)

Dada la lista de todos los MEMs, la alineación puede ser calculada usando parcial de las alineaciones. Por ejemplo, considere MEMs Mx, My y Mz en la Fig. 3b. los alineamientos parciales realizados tomando diferentes combinaciones de Mx, My y Mz junto con el número de coincidencias, desajustes y huecos, así como los puntajes de alineación resultantes se muestran en la Fig. 4. La alineación que solo incluye Mx y Mz da como resultado la puntuación de alineación más alta., Tenga en cuenta que, My y Mz se superponen entre sí y cuando ambos se consideran en la misma alineación, la superposición se excluye de Mz. Considerando todos los MEMs en la Fig. 3b resulta en muchas más combinaciones donde ninguna de ellas logra una puntuación más alta.

Fig. 4

Todas las combinaciones posibles de los MEMs en la alineación

Examinar todas las combinaciones posibles de los MEMs sería exhaustiva., En la sección» algoritmo de alineación » describimos un nuevo algoritmo de programación dinámica DP-MEM que encuentra eficientemente la mejor combinación sin considerar todos los casos. DP-MEM necesita saber qué partes de las secuencias coinciden, pero no los símbolos reales en las secuencias. La entrada a DP-MEM es el posicionamiento de MEMs en el destino y en las secuencias de consulta que se obtienen durante el proceso de extracción de MEM descrito en la sección «Extracción de MEM»., En el resto de esta sección se explica cómo se representan los MEMs con sus posiciones y cómo se calcula el número de coincidencias, desajustes y huecos cuando los MEMs se combinan en una alineación. La figura 5 es otro ejemplo de alineación con seis MEMs (M1 A M6) que forman la alineación entre la secuencia de destino T y la secuencia de consulta Q. para simplificar, no hay superposición entre MEMs en este ejemplo. Cada MEM Mi se representa como un triplete de números enteros: las posiciones iniciales en T y en Q (STi y SQi respectivamente) y su longitud (Li)., Las posiciones finales en T y en Q se pueden calcular más tarde (Φ2E del algoritmo 2). La tabla 1 muestra la longitud y el posicionamiento de M1 A M6 en T y en Q.

Fig. 5

Un ejemplo de alineación con relieve MEM

Tabla 1 Inicio y final de la posición de MEMs en la Fig., 5
Tabla 2 Computación número de discrepancias y las diferencias entre los MEMs en la Fig. 5

en caso de que haya desajustes y brechas entre Mi y Mj, todas las brechas se consideran continuas para reducir la penalización de brecha abierta (solo se aplica una penalización de brecha abierta para una brecha continua). Por lo tanto, para todos los MEMs adyacentes que tienen espacios entre ellos, solo se aplica una penalización de apertura de espacio., La colocación de desajustes y la única brecha continua no es importante, ya que no afectaría la puntuación de alineación. Asumimos que la penalización de desajuste es constante (esto es habitual para las secuencias de ADN).

extracción de MEM

existen métodos para extraer coincidencias exactas máximas entre secuencias largas, como un genoma completo. Sin embargo, estos métodos se basan en el preprocesamiento y la indexación de una o ambas secuencias, lo que requiere mucho tiempo. Por ejemplo, en DNA read aligner, el genoma de referencia se indexa una vez, y el mismo índice se usa cada vez que se alinea una nueva lectura., Estamos buscando un algoritmo rápido para identificar MEMs entre secuencias relativamente cortas que cambian para cada alineación. Un método de fuerza bruta para este problema (archivo adicional 1: Sección II) es lento e ineficiente (con la complejidad de O(n3)). Proponemos un método paralelo rápido a nivel de broca para acelerar el proceso de extracción de MEM. Nuestro método de extracción MEM se basa en las operaciones de desplazamiento y comparación que se muestran en la Fig. 3b. el primer paso es representar secuencias con vectores de bits, donde A, C, T Y G están codificados como 00, 01, 10 y 11, respectivamente (archivo adicional 1: Sección III)., La figura 6 ilustra un par de secuencias de ejemplo, junto con las representaciones de vectores de bits correspondientes. En una computadora de productos básicos, la palabra máquina es generalmente de 64 bits que pueden acomodar 32 símbolos de nucleótidos. Dado que una secuencia es generalmente más grande que 32 Símbolos, el vector de bits correspondiente se almacena en varias palabras de máquina. Cada operación en vectores de bits de secuencias de símbolos de tamaño n actúa sobre palabras máquina \(\lceil \frac {n}{32} \rceil\).

Fig., 6

la Representación de secuencias de bits en vectores. Salida XOR (X) con MEMs resaltados. Edges bit-vector (e) identifica el inicio y el final de cada MEMs

con la representación de vectores de bits de secuencias, cambiar una secuencia por un símbolo es lo mismo que cambiar el vector de bits por dos bits, y la comparación de secuencias se puede hacer con la instrucción XOR (32 símbolos a la vez). En la salida XOR (X), 00 significa que los símbolos coinciden, y una secuencia de 00s muestra un MEM., Un conjunto de operaciones de desplazamiento y bits como se muestra en el algoritmo 1 calcula X y posteriormente el vector de bits de borde (E) en el que el inicio y el final de cada MEM se resaltan con bits establecidos (bits con un valor de uno). La figura 6 muestra los vectores de bits X Y E con MEMs resaltados. El posicionamiento de MEMs en secuencias se calcula a partir del vector de bits de borde (archivo adicional 1: Sección IV).,

algoritmo de alineación

en la sección» Approach», mostramos que al considerar diferentes combinaciones de MEMs y calcular la puntuación de alineación para la alineación correspondiente, se puede identificar la combinación de MEMs que resulta en la puntuación de alineación máxima. Sin embargo, examinar todas las combinaciones posibles de MEMs es una solución ingenua. Una forma más sistemática de encontrar la alineación de manera eficiente es utilizar la programación dinámica.

la Programación Dinámica es el método de abordar la solución a un problema definiendo y resolviendo subproblemas más pequeños., Las soluciones para subproblemas se utilizan para resolver un problema más grande en cada paso. El proceso se repite hasta que se resuelvan todos los subproblemas. Eventualmente, la solución a uno de los subproblemas sería la solución al problema inicial. Cuando se resuelven todos los subproblemas, un proceso de retroceso identifica una serie de soluciones que contribuyen a la solución final. En la programación dinámica, debe haber un orden de los datos de entrada a lo largo de los cuales procede el procedimiento de recursión.

ordenamos todos los MEMs de acuerdo con la posición de su final en la secuencia de consulta (EQ)., Los MEMs que terminan en la misma posición se ordenan de manera arbitraria. El subproblema J-ésimo es encontrar la alineación de las subsecuencias de T Y Q que terminan en J-ésimo MEM Mj (T Y Q respectivamente). Demostraremos que este orden de MEM es suficiente para soportar la recursión correcta.

en la lista ordenada de MEMs, EQI = EQj indica que uno de Mi o Mj se solapa completamente con el otro MEM en la secuencia de consulta. Dado que en Φ2B del algoritmo 2 se excluye la región de superposición, Mi y Mj no pueden estar en la misma alineación., Por lo tanto, los subproblemas ith y jth se resuelven independientemente entre sí y el orden de I y j en la lista ordenada podría ser arbitrario. Si EQk>EQj (k> j en la lista ordenada), Mk no podría ser parte de la alineación que termina en Mj. Por lo tanto, los subproblemas jth se pueden resolver independientemente de la solución al subproblema kth. Tenga en cuenta que también es posible ordenar MEMs en función de su posición final en la secuencia de destino (ET) utilizando una justificación similar.

nuestro algoritmo de programación dinámica propuesto (DP-MEM) se elabora en el algoritmo 2., Para el ejemplo MEMs extraído en la Fig. 3b, la tabla de programación dinámica y el valor intermedio calculado en el algoritmo se muestran en las Figs. 7 y 8 respectivamente. La entrada a DP-MEM es la lista de MEMs donde cada MEM (Mj) es un triplete de enteros . La segunda entrada n es el número de MEMs en la lista. La salida S es la puntuación de alineación de las secuencias. El algoritmo imprime los índices de todos los MEMs que forman la alineación, donde el primero y el último número impreso son los índices de los MEMs más a la derecha y los más a la izquierda en la alineación respectivamente., Todos los pasos del algoritmo 2 se comentan en lo siguiente:

Fig. 7

tabla de programación dinámica utilizada en el algoritmo 2 para procesar MEMs extraídos en la Fig. 3b. La celda i y j representan el valor de \(S_{i}^{j}\). Las celdas vacías no se evalúan en Φ2. La evaluación de las células con la marca cruzada se salta en Φ2A. el valor inicial de Sj se calcula en Φ1. El valor Final de Sj y su fuente (lo que maximiza Sj) se resaltan para cada fila. El Sj más alto (S13) es la puntuación de alineación., M13 es el último MEM en la alineación y el MEM antes es MW = M3. Desde W = -1, M3 es el primer MEM en la alineación. El sistema de puntuación para esta alineación es Rm=2,Px=3,Po=4 y Pe=1

Fig. 8

valores Intermedios para calcular \(S_{i}^{j}\) en la Fig. 7. Tenga en cuenta que Sij en esta figura se refiere a \(s_{i}^{j}\)

  • Φ1: puntuación de cada MEM para todos sus símbolos coincidentes., Tenga en cuenta que hay símbolos LJ coincidentes en Mj. Sj representa la puntuación de alineación más alta para la alineación que termina en Mj. Inicializar Sj en este paso es similar a calcular la puntuación de alineación parcial cuando solo se incluye Mj en la alineación. W se utiliza para retroceder. El valor de -1 indica que el Sj actual se obtiene considerando Mj solo en la alineación.

  • Φ2: cálculo de Sj para cada MEM (Mj)., Para calcular Sj, para cada MEM Mi Donde Mi aparece antes de Mj en la lista, el algoritmo agrega Mj a la alineación que termina en Mi (extendiendo alineaciones previamente encontradas) y busca la extensión que maximiza Sj (resolviendo un subproblema más grande usando subproblemas previamente resueltos).

  • Φ2A: Omitir la extensión cuando no es posible. Si ETi> ETj entonces Mi contiene parte de la secuencia de destino que está más allá de la alineación que termina en Mj y la extensión no es posible. Si EQi=Eqj o ETi = ETj o SQI≥SQJ o STI≥STj, entonces uno de los MEMs se solapa completamente con el otro MEM., En este caso, Mi y Mj no pueden estar alineados juntos.

  • Φ2B: calculando la longitud de la superposición entre Mi y Mj. Si \({MO}_{i}^{j}\) es menor o igual a cero, entonces no existe superposición.

  • Φ2C: mantener una copia de Mj antes de excluir la superposición.

  • Φ2D: si existe superposición, excluyendo la superposición de Mj

  • Φ2E: calcular la posición final de Mj en T y Q.

  • Φ2F: calcular la distancia (número de símbolos) entre Mi y Mj en T y Q.

  • Φ2G: calcular el número de desajustes y huecos entre Mi y Mj.,

  • Φ2H: calculando la penalización por los desajustes y huecos entre Mi y Mj (\(P_{i}^{j}\)). Si existe brecha, solo se resta una penalización de apertura de brecha.

  • Φ2I: cálculo de la puntuación de alineación \(\left (s_{i}^{j}\right)\) cuando se agrega Mj a la alineación que termina en Mi. La puntuación de todos los símbolos coincidentes en Mj (Lj×Rm) se añade a la puntuación de alineación para la alineación que termina en Mi (Si). Luego se resta la penalización por los huecos y desajustes entre Mi y Mj\(\left (P_{i}^{j}\right)\).,

  • Φ2J: si extender Mj a la alineación que termina en Mi resulta en una puntuación \(\left (s_{i}^{j}\right)\) mayor que la puntuación actual para Mj (Sj), entonces la nueva puntuación se almacena en Sj. También W se establece en i para realizar un seguimiento de la Mi que maximizar la puntuación de Mj.

  • Φ2K: restaurar el valor de Mj antes de la exclusión para que Mj pueda usarse en otras extensiones de alineación.

  • Φ3: buscando el MEM con el Sj más alto. Este MEM es el último MEM en la alineación (Me)., La puntuación más alta (Se) se devuelve como S, que es la puntuación de alineación más alta para las secuencias dadas. El índice del MEM que maximiza Sj se almacena en e Para comenzar a retroceder de mí.

  • Φ4: en la alineación, el MEM anterior inmediato para mí es el que maximiza la puntuación de alineación para mí. El índice de dichos MEM se almacena en W. como resultado, la iteración de F←W visita el índice de todos los MEMs en la alineación. Cuando W es igual a -1, Mf es el primer MEM en la alineación y la iteración se detiene.,

en nuestro algoritmo, no penalizamos desajustes y huecos antes del primer MEM y después del último MEM en la alineación. Esto resulta en un algoritmo de alineación local. Al considerar estas penalizaciones el algoritmo genera una alineación global (archivo adicional 1: Sección V).

la ecuación para calcular \(p_{i}^{j}\) en Φ2H del algoritmo 2 asume que no hay un símbolo coincidente entre T Y Q en el área entre Mi y Mj (todos los símbolos se cuentan como desajustes o huecos)., Aunque esta suposición no es cierta para todos los Mi, siempre es cierta para el Mi que conduce al máximo \(S_{i}^{j}\) que anula el efecto de que la suposición sea incorrecta para otros Mi. Como prueba, supongamos que hay un símbolo coincidente en el área entre Mi y Mj. El símbolo correspondiente sería un MEM (Mk). Mk ya está extendido a la alineación que termina en Mi. Por lo tanto, al extender Mj A Mk se logra una puntuación más alta cuando se compara con extender Mj a Mi.

encadenar semillas colineales como se discute en han sido ampliamente utilizados en la alineación de grandes secuencias, es decir, la alineación genoma-a-genoma., También se ha utilizado para identificar regiones candidatas para una lectura dada una serie de MEMs en BWA. Los Algoritmos de Encadenamiento con puntuación son similares al algoritmo de programación dinámica que propusimos (DP-MEM). Sin embargo, hay diferencias que hacen que DP-MEM sea adecuado para la alineación por pares de secuencias cortas. DP-MEM tiene en cuenta que todos los MEMs dentro de un cierto tamaño de espacio se proporcionan en la entrada y optimiza el número de iteración en el algoritmo. DP-MEM también implementa un enfoque heurístico para compensar el efecto de MEMs cortos eliminados de la lista de entrada resultantes de las brechas entre MEMs.,

Optimization

dadas las secuencias de longitud n, El algoritmo para extraer MEMs (proporcionado en la sección «Extracción de MEM») requiere 2(n−1) shift y 2n−1 compare operaciones en vectores de bits (cada acto en palabras máquina \(\lceil \frac {n}{32} \rceil\)) que resultan en un algoritmo con complejidad de O(n2) para producir vectores de bits de borde para el par de secuencias dado. La complejidad de la función que calcula el posicionamiento de MEMs desde el vector de bits de borde y los ordena en función de EQ aún no se ha agregado. Además, si se extraen M MEMs, Φ2 del algoritmo 2(DP-MEM) tiene la complejidad de O (m2)., Teniendo en cuenta la complejidad de la extracción de MEM y DP-MEM, MEM-Align es mucho más lento que los Algoritmos de alineación disponibles. Para acelerar el proceso, presentamos varias optimizaciones para MEM-Align que logran un tiempo de ejecución más rápido al sacrificar la precisión. Por otro lado, introducimos métodos para mejorar la precisión con una pérdida de rendimiento mínima.

alineación de bandas

La alineación de bandas es un método heurístico conocido para acelerar el proceso de alineación. Esta técnica limita el patrón de los huecos en la alineación (archivo adicional 1: Sección VI)., En consecuencia, si la alineación entre dos secuencias no sigue este patrón, el algoritmo no identificará la alineación. En la programación dinámica tradicional, La alineación se obtiene después de calcular el valor de todas las celdas de la tabla. Sin embargo, con la optimización de alineación de bandas, solo se evalúan las celdas en el diámetro y cerca de la diagonal. gl es el ancho de la banda en la alineación de bandas donde las células más lejos que gl al diámetro no se evalúan. Células más oscuras en la Fig. 1 mostrar el caso donde gl = 1.,

a diferencia del enfoque de programación dinámica tradicional, MEM-Align no tiene una tabla similar para aplicar la alineación de bandas. Sin embargo, encontramos que podemos simular el mismo efecto limitando el número de operaciones de cambio en el proceso de extracción de MEM. Por ejemplo, si desplazamos la secuencia de consulta hasta gl a la derecha y a la izquierda, alcanzamos la alineación de bandas con la banda de gl. La alineación de bandas reduce la complejidad de la extracción de MEM de O(n2) A O(n.(2gl+1)) donde gl es un valor pequeño y fijo. Por lo tanto, la complejidad de la extracción MEM es O(n) cuando se aplica alineación de bandas., Además, con dicha alineación de bandas, es probable que se extraigan menos MEMs, lo que acelera los pasos algorítmicos posteriores.

eliminación de MEM cortos

observamos que la mayoría de los Mem extraídos son cortos y son el resultado de símbolos coincidentes al azar. Para acelerar la alineación MEM, los MEMs más cortos que sl se filtran durante el proceso de extracción MEM. Esto reduce el número de MEMs a extraer y procesar (posteriormente acelerando el algoritmo). El filtrado de MEM corto se realiza extendiendo el algoritmo 1 con un conjunto de operaciones de desplazamiento y bits (archivo adicional 1: Sección VII).,

por otro lado, hay casos raros en los que MEMs cortos son parte de la alineación; por ejemplo, un símbolo coincidente rodeado de desajustes. Sin tener todos los MEMs en la lista de entrada, DP-MEM no es capaz de encontrar la misma alineación que encuentra cuando todos los MEMs existen en la lista de entrada. Para compensar la pérdida de MEMs cortos en la entrada, modificamos Φ2H de DP-MEM para considerar la posibilidad de tener coincidencias cortas entre MEMs (archivo adicional 1: Sección VIII).

Puede haber casos más difíciles donde en la alineación, existen múltiples MEMs cortos entre dos MEMs(ver Fig. 9)., La única manera de identificar correctamente la puntuación para el área entre Mi y Mj en Φ2H es aplicar una alineación global a esta región. Sin embargo, Φ2H es una operación frecuente y debe permanecer rápido. En consecuencia, decidimos superar parcialmente el problema limitando los posibles casos (un método heurístico).

Fig., 9

Un ejemplo que muestra múltiples MEM en el área pequeña entre Mi y Mj en una alineación

Si hay lagunas en la zona entre Mi y Mj, suponemos que sólo hay una continua brecha, ya sea a la izquierda o a la derecha de la zona. Entonces, solo dos alineaciones son posibles para el área., El número de símbolos coincidentes se cuenta para ambos casos de manera secuencial y el que resulta en coincidencias máximas se toma como el número de coincidencias entre Mi y Mj (archivo adicional 1: Sección IX). La comparación secuencial es una operación costosa y diseñamos un método para evitar la comparación secuencial cuando sea posible (archivo adicional 1: Sección X).

cualquier otro caso que no se ajuste a la suposición anterior resulta en una alineación con una puntuación más baja., Sin embargo, teniendo en cuenta la baja tasa de brechas y desajustes, la posibilidad de tener múltiples brechas y desajustes en un área pequeña es baja.

extensión de alineación eficiente

en Φ2 de DP-MEM, Mj extiende todas las alineaciones que terminan en {M1} Mj−1} (si es posible). Sin embargo, para cada Mj hay un subconjunto más pequeño Ωj {{M1 M Mj−1} tal que al extender Mj a todas las alineaciones que terminan en un Mi∈Ωj se encuentra la alineación que termina en Mj (EC. 2). En otras palabras, habría menos \(s_{i}^{j}\) que evaluar. Definición del conjunto Ωj y la prueba de EC. 2 se proporcionan en el archivo adicional 1: Sección XI., La definición de Ωj se ve afectada cuando se aplica la optimización de eliminación de MEM corto (archivo adicional 1: Sección XII).

$$ \max\limits_{M_{i} \in \Omega_{j}}{S_{i}^{j}} = \max\limits_{1 \leq i \leq j-1}{S_{i}^{j}} $$
(2)

Híbrido alineación

Para mantener la precisión del algoritmo, se decidió utilizar un método híbrido que es una combinación de MEM-Alinear y Smith-Waterman algoritmo. Definimos tres casos en los que MEM-Align puede ser inexacto., Si la alineación de un par de secuencias cae en uno de estos casos, usamos el algoritmo Smith-Waterman para alinear secuencias. Estos casos son:

  • cuando las secuencias son repetitivas, y el número de MEMs extraídos supera el umbral TM. Encontramos que MEM-Align es probable que produzca una alineación inexacta al alinear secuencias repetitivas. Un valor de TM apropiado disminuye la posibilidad de reportar una alineación inexacta con un aumento insignificante en el tiempo promedio de procesamiento.,

  • Cuando la puntuación de alineación calculada para la alineación generada por MEM-Align es inferior a un umbral TS. Este caso ocurre principalmente cuando hay una brecha en la alineación que no se puede identificar debido a la alineación de bandas.

  • Cuando no existe MEM más largo que sl para ser extraído (un caso raro)., Si sl se establece en un valor alto y la similitud entre secuencias es baja,

aunque enviar pares de secuencias a un algoritmo externo resulta en un cálculo adicional, el número de secuencias enviadas al algoritmo externo permanece pequeño si se eligen valores apropiados para TM y TS.

omitir MEMs distantes

Cuando la distancia entre Mi y Mj es grande, no es probable que tenga Mi y Mj como MEMs adyacentes en la alineación., Por lo tanto, el algoritmo omite la extensión si la distancia entre Mi y Mj es más larga que un umbral TD (reduciendo aún más el número de \(s_{i}^{j}\) a evaluar). Esta optimización mejora ligeramente el rendimiento con un efecto secundario insignificante en la precisión.