El problema de volver a unir los reads

Lo que se obtiene del secuenciador son fragmentos de secuencia cortos (reads). Reads cortos de 100 pb y 150 pb, y reads largos de 10 a 100 kb. El genoma original es una única secuencia larga de millones a miles de millones de pb. El problema de volver a unir estos reads para reconstruir el genoma original se denomina ensamblaje (assembly).

Existen dos enfoques principales: OLC (Overlap-Layout-Consensus) y Grafo de De Bruijn. En esta sección se aborda OLC. En la siguiente sección (M23) se tratará el método de De Bruijn, y en la sección posterior (M24) se pasará al estándar actual para reads largos, hifiasm. Estas tres secciones constituyen los cimientos de la serie sobre ensamblaje.

Las tres etapas de OLC

Como su nombre indica, son tres etapas.

1. Overlap: Buscar los sufijos y prefijos comunes entre todos los pares de reads.
2. Layout: Construir un grafo donde los reads son nodos y las superposiciones son aristas, y encontrar la ruta óptima.
3. Consensus: Unir los reads a lo largo de la ruta y corregir errores mediante votación mayoritaria en las regiones de superposición.

Las tres etapas son difíciles, pero especialmente la etapa de Layout es el núcleo de este algoritmo. La ruta óptima es una ruta que visita cada read exactamente una vez, lo que se conoce como camino hamiltoniano. Este es un problema NP-completo. Por lo tanto, las herramientas reales utilizan varias heurísticas para obtener aproximaciones.

Cálculo de la superposición: ¿cómo se detecta la superposición?

Cuando tenemos el read A y el read B, buscamos una longitud de superposición (overlap length) o tal que el sufijo de A coincida con el prefijo de B en esa longitud. Por ejemplo:

A: ACGTACG
B:    TACGTAT

El sufijo de A TACG y el prefijo de B TACG coinciden en una longitud de 4. Los dos fragmentos se conectan con una superposición de 4.

El método estándar tiene una complejidad de O(N² × L), donde O(N²) corresponde a todas las parejas y O(L) al emparejamiento de sufijos y prefijos de cada pareja. Si N es del orden de millones, el cálculo resulta inviable. Por ello, se utiliza una indexación previa mediante un árbol de sufijos (M13) o un FM-index (M15), reduciendo la complejidad a O((N+L) log L). Aquí se evidencia por qué se estudiaron las estructuras de datos de indexación de cadenas en la fase previa al ensamblaje.

En la práctica, se establecen como umbrales una longitud mínima de superposición (por ejemplo, 50 pb) y una precisión mínima (por ejemplo, 95 %). Los fragmentos con errores pueden perder ligeramente la superposición, pero siguen conectándose con otros fragmentos.

Grafo de superposiciones

• Nodos: cada fragmento
• Aristas: superposición entre dos fragmentos (dirigida: A→B indica "el sufijo de A se solapa con el prefijo de B")
• Peso de la arista: longitud de la superposición (o puntuación de precisión)

Sobre este grafo, un camino que pasa exactamente una vez por cada nodo (camino hamiltoniano) reconstruye el orden del genoma original.

Ejemplo sencillo. 5 fragmentos.

R1: ACGTACG
R2:   GTACGT
R3:      CGTAT
R4: TTACGTA
R5: GTATGC

Tabla de superposiciones (para facilitar la lectura, solo se muestran las superposiciones de 3 o más):

Al representar esto en un gráfico, se obtiene una ruta dirigida: R4 → R1 → R2 → R3 → R5. Se visitan los cinco nodos. Al concatenar las lecturas siguiendo esta ruta:

R4:  TTACGTA
R1:      ACGTACG
R2:         GTACGT
R3:            CGTAT
R5:                GTATGC

Consensus: TTACGTACGTATGC

El genoma original se reconstruye. Este es el uso práctico de OLC.

Ruta de Hamilton: ¿por qué es difícil?

En el ejemplo anterior, con solo 5 nodos, es fácil resolverlo a mano. Pero con millones de lecturas, la situación cambia radicalmente. El problema de la ruta de Hamilton es NP-completo. No existen algoritmos conocidos de tiempo polinómico para resolverlo.

Las herramientas prácticas utilizan varios trucos:

• Algoritmo voraz (greedy): En cada paso, se selecciona el borde de superposición más largo para continuar. Aunque converge solo hacia un óptimo local, funciona bien en la mayoría de los genomas reales.
• Reducción del grafo de cadenas: Se eliminan los bordes transitivos (si existe A→C y también A→B→C, se elimina A→C), lo que reduce drásticamente el tamaño del grafo. Esta es la idea central de Myers.
• Separación de repeticiones: Las secuencias repetitivas hacen que las lecturas provenientes de múltiples posiciones se agrupen en un solo nodo, enredando la topología del grafo. Esto se resuelve utilizando información de pares extremos y el alcance de lecturas largas.

Celera Assembler, utilizado por el equipo de Venter en el Proyecto Genoma Humano, fue una obra representativa de este enfoque. Posteriormente, herramientas para lecturas largas como miniasm, Canu y Falcon también pertenecen a la familia OLC.

¿Por qué OLC ha resurgido en la era de las lecturas largas?

Durante la era de las lecturas cortas (principios y mediados de la década de 2010), el grafo de De Bruijn (ver el siguiente capítulo) era el estándar. Dado que la longitud de las lecturas era de 100 pb, encontrar superposiciones era difícil en sí mismo, y los métodos basados en k-mer del grafo de De Bruijn eran mucho más eficientes.

Con la llegada de la era de las lecturas largas (década de 2020), la situación se invirtió. Una sola lectura puede tener entre 10 y 100 kb. Aunque la tasa de error es alta (5～15 %), proporciona abundante información sobre estructuras a larga distancia. La detección de superposiciones es estable y puede abarcar regiones repetitivas. Esta es la razón por la que OLC ha vuelto a convertirse en el estándar.

Herramientas modernas como hifiasm han incorporado las características de las lecturas HiFi (lecturas largas de alta precisión) sobre la estructura básica de OLC. Este tema se tratará en M24.

Implementación mínima de OLC en Python (~40 líneas)

Una implementación de juguete para verificar los principios. No sustituye a las herramientas prácticas.

def overlap(a, b, min_length=3):
    \"\"\"Devuelve la longitud máxima de superposición entre el sufijo de a y el prefijo de b.\"\"\"
    start = 0
    while True:
        start = a.find(b[:min_length], start)
        if start == -1:
            return 0
        if b.startswith(a[start:]):
            return len(a) - start
        start += 1

def build_overlap_graph(reads, min_length=3):
    \"\"\"Calcula las superposiciones para todos los pares de lecturas.\"\"\"
    graph = {}
    for a in reads:
        for b in reads:
            if a == b:
                continue
            olen = overlap(a, b, min_length)
            if olen > 0:
                graph.setdefault(a, []).append((b, olen))
    for a in graph:
        graph[a].sort(key=lambda x: -x[1])  # orden descendente por superposición
    return graph

def greedy_hamiltonian(reads, min_length=3):
    \"\"\"Camino hamiltoniano voraz (aproximación).\"\"\"
    graph = build_overlap_graph(reads, min_length)
    remaining = set(reads)
    # Comienza desde nodos sin aristas entrantes (si no hay, inicio arbitrario)
    in_deg = {r: 0 for r in reads}
    for a, edges in graph.items():
        for b, _ in edges:
            in_deg[b] += 1
    start = min(reads, key=lambda r: in_deg[r])
    path = [start]
    remaining.remove(start)

    current = start
    while remaining:
        candidates = graph.get(current, [])
        for nxt, _ in candidates:
            if nxt in remaining:
                path.append(nxt)
                remaining.remove(nxt)
                current = nxt
                break
        else:
            # Desconexión → reinicia desde nodos restantes
            current = remaining.pop()
            path.append(current)
    return path

def consensus(path, min_length=3):
    \"\"\"Une las lecturas recortando por la superposición a lo largo del camino.\"\"\"
    result = path[0]
    for i in range(1, len(path)):
        olen = overlap(path[i-1], path[i], min_length)
        result += path[i][olen:]
    return result

reads = ["ACGTACG", "GTACGT", "CGTAT", "TTACGTA", "GTATGC"]
path = greedy_hamiltonian(reads, min_length=3)
print("Layout:", path)
print("Consensus:", consensus(path, min_length=3))

Al ejecutarlo, los resultados coinciden con el cálculo manual anterior. Estas 40 líneas constituyen el núcleo de OLC. Las herramientas prácticas son implementaciones complejas que añaden indexación, reducción transitiva y gestión de repeticiones a esta base, pero la estructura fundamental no se desvía de este código.

Complejidad

• Fase de superposición: Método directo O(N² × L) → Indexación para O((N+L) log L)
• Fase de diseño: Camino hamiltoniano = NP-completo. Aproximación voraz O(N)
• Fase de consenso: Longitud del camino O(N × L)

La superposición es el cuello de botella. Para millones de lecturas, sin indexación, el cálculo en sí mismo es imposible. Aquí es donde los árboles de sufijos y el índice FM, aprendidos previamente, dan sus frutos.

OLC vs. De Bruijn: ¿por qué coexisten estos dos modelos?

La diferencia clave entre los dos modelos es la definición de nodo. OLC opera a nivel de lectura, mientras que De Bruijn opera a nivel de k-mer. Esta única elección determina la complejidad algorítmica y el tipo de lecturas adecuadas. La comparación detallada continúa en el siguiente capítulo.

Mapeo CS: las dificultades de los algoritmos de grafos

El concepto de NP-completo tratado en DryBench aparece aquí como un caso práctico.

• Camino hamiltoniano: Visitar cada nodo exactamente una vez → NP-completo
• Camino euleriano: Recorrer cada arista exactamente una vez → Tiempo polinómico

Una simple diferencia (nodo vs. arista) invierte la complejidad computacional. La forma en que se define el problema cambia fundamentalmente la dificultad de su resolución. Aquí se revela la estética del diseño algorítmico. Al estudiar los grafos de De Bruijn en el siguiente capítulo, volveremos a observar este contraste.

Evaluación en Rosalind

Rosalind tiene varios problemas sobre OLC: OverlapGraphs(BA3C), LongOverlapMerge(BA3M), etc. Se resuelven directamente reutilizando la función overlap de la implementación en Python anterior. Escribir una función de superposición estable es la base práctica de esta serie.

Ramificaciones hacia el siguiente capítulo

• Siguiente capítulo (M23): Grafo de De Bruijn — Nodos k-mer + camino euleriano. El estándar de la era de las lecturas cortas. La belleza del tiempo polinómico.
• Dos capítulos después (M24): hifiasm: ensamblador diploide con PacBio HiFi. La implementación más avanzada de OLC. Ensamblaje específico para padres (trio binning).
• Ampliación: string graph (Myers 2005): la base teórica de la contracción de grafos en OLC. El fundamento teórico de hifiasm.

Para profundizar

• Myers (2005), The fragment assembly string graph, Bioinformatics: el artículo original sobre string graph. La esencia de la contracción de aristas transitivas.
• Chin et al. (2016), Phased diploid genome assembly with single-molecule real-time sequencing, Nature Methods: el artículo sobre Falcon. El flujo de trabajo práctico de OLC con lecturas largas.
• Stanford CS262: clase de ensamblaje del profesor Gill Bejerano (con subtítulos completos). Compare OLC y De Bruijn mediante anotaciones en la pizarra.
• Ben Langmead (JHU): hay un episodio sobre ensamblaje en su lista de reproducción. Visualice la topología de grafos.
• Compau & Pevzner Bioinformatics Algorithms: integrado con Rosalind. El Capítulo 3 trata sobre el ensamblaje.

Antes de pasar al siguiente capítulo, genere aleatoriamente entre 20 y 50 lecturas e insértelas en la implementación de Python anterior. El verdadero ejercicio de este capítulo consiste en verificar manualmente cuándo falla el algoritmo voraz y cómo se complica con las secuencias repetidas.