Cómo funciona

Ligamiento: distinguir un haplotipo de una casualidad

Dos variantes en un mismo codón solo son una variante multinucleotídica si las mismas moléculas llevan ambas. get_MNV reporta con cuánta fuerza van juntas, como D' con su valor p, en las columnas Haplotype LD y Haplotype LD p (LD / LDP en la salida VCF). Requiere --bam.

Por qué un ratio de co-ocurrencia no basta#

La medida obvia es la fracción de moléculas que llevan una variante y llevan también la otra. Tiene un fallo que empeora cuanto más comunes son las variantes.

Coge dos sustituciones que están cada una en el 90% de las moléculas, sin nada que las relacione. Por pura aritmética se encuentran juntas en el 81% de ellas, y el ratio reporta 0.9000. Leído como ligamiento, eso dice "estas dos casi siempre van juntas". Lo que dice en realidad es "estas dos son comunes".

MNV Phasing Support es ese ratio y se sigue reportando, porque "con qué frecuencia co-ocurren" es una pregunta legítima con una respuesta útil. Pero no puede responder a la que hay debajo, que es si esa co-ocurrencia significa algo.

Qué mide D'#

D' es el exceso de co-ocurrencia sobre lo que predicen las dos frecuencias, normalizado por el mayor exceso que esas frecuencias podrían haber producido. Esa normalización es lo que lo hace comparable entre sitios, y deja el valor en [-1, 1].

Los mismos datos de arriba, puntuados como toca:

Las filas primera y quinta son la cuestión. Un ratio de 0.9 y uno de 0.5 parecen ligamiento fuerte y moderado; los dos son exactamente lo que predice la independencia.

El lado negativo#

Un D' cerca de -1 es un hallazgo, no un fallo. Los dos alelos están presentes en la muestra, a frecuencia apreciable, y ninguna molécula lleva ambos.

En un organismo haploide eso no es una variante con dos formas. Son dos linajes, y a nivel de codón es la firma a nivel de lectura de una infección mixta o de subpoblaciones distintas.

Un ratio de co-ocurrencia de 0.0000 sí lo reporta, y con honestidad: hubo moléculas que cruzaron el sitio y ninguna llevaba ambas. Lo que no puede hacer es decir cuánto se aleja eso del azar, así que se lee igual tanto si las dos se excluyen como si una es simplemente demasiado rara para esperar una molécula compartida. Un D' de -1 con un p pequeño distingue los dos casos.

El valor p#

Un D' de 1.0 con cuatro moléculas y un D' de 1.0 con cuatrocientas son el mismo número y no son la misma evidencia. Haplotype LD p es el test exacto de Fisher a dos colas sobre la tabla 2x2 de recuentos de moléculas, así que dice si el alejamiento de la independencia es más de lo que la profundidad podría producir por azar.

Léelos juntos. Un |D'| grande con un p cerca de 1 es una muestra pequeña hablando con seguridad de nada.

Cuándo está ausente#

- significa que la pregunta no se pudo responder, y eso es deliberadamente distinto de una respuesta de cero:

• Sin --bam, o una fila de una sola variante: no hay nada que correlacionar.
• Ninguna molécula observó las variantes juntas. Habitual cuando un codón cruza un intrón y los fragmentos son más cortos que él.
• Uno de los alelos está en todas las moléculas que observaron, o en ninguna. Sin variación en ese locus no hay correlación que medir, y el denominador de D' es cero. Reportar 0.0 ahí afirmaría que se midió independencia y se encontró.

Ese último caso es más común de lo que parece. Si doce moléculas llevan tres variantes y ocho llevan solo las dos primeras, esas dos primeras están en las veinte moléculas que observaron la ventana: la fila sale con sus ocho moléculas, y su ligamiento es -, porque dentro de esta ventana esas dos nunca varían.

Los haplotipos de indel también lo llevan#

El estadístico no es específico de codones. Una fila de haplotipo local de indel también lo lleva, calculado sobre las variantes que esa fila afirma que van juntas.

Esa fila antes solo tenía un recuento de lecturas, que no distingue una asociación real de dos variantes comunes que se encuentran. Dos ventanas, cada una con una inserción y una sustitución:

Las dos filas son reales: cinco y dieciocho moléculas llevan de verdad esas dos variantes juntas. En la primera la inserción está en el 40% de las moléculas y la sustitución en el 30%, y encontrarse en el 10% es lo que predice el azar. En la segunda las dos están en la mitad y se encuentran en el 45% donde el azar predice el 25%.

Leerlo junto al recuento de lecturas#

Responden a preguntas distintas y ninguno sustituye al otro.

• Event Reads / MNV Reads: cuántas moléculas son esa combinación.
• Haplotype LD: si sus variantes tienen algo que ver entre sí.

Un haplotipo con pocas moléculas puede estar perfectamente ligado, y uno con muchas puede ser una casualidad. El primero es una especie minoritaria real que merece reportarse; el segundo son dos variantes que resultan ser comunes en la misma muestra.

Detalles del cálculo#

• La tabla se construye sobre moléculas, no lecturas: los dos mates de un fragmento paired-end son una molécula. Mira Indels y haplotipos locales para qué cuenta como una.
• Solo se incluyen las moléculas que observaron todas las posiciones de la fila, en los dos lados de la tabla. Una molécula que paró entre dos posiciones no vio evidencia sobre el par.
• Con tres o más variantes se puntúa cada par y decide el más débil, porque la fila afirma que una molécula las lleva todas, y esa afirmación vale lo que su par menos ligado.
• r² se calcula a la vez y está disponible en la API. Se reporta D' porque la pregunta aquí es "están en las mismas moléculas", que es a lo que responde D'; r² refleja además cuánto se separan las dos frecuencias, que es otra cosa que saber.

Cómo se comprueba esto#

tests/phasing/ construye alineamientos molécula a molécula, eligiendo qué lecturas llevan qué variantes, y compara contra D' y un test exacto derivados por segunda vez en la propia suite. Barre tablas 2x2 que cubren la independencia, las dos direcciones de asociación, y las tablas degeneradas donde D' no tiene valor, junto a los barridos de soporte de fase.

Cuando las lecturas desmienten una lectura#

Un codón con dos sustituciones tiene dos lecturas, y una coincidencia del catálogo que hace solo una de ellas afirma algo de las moléculas. Una coincidencia a nivel de codón en la columna Catalogue afirma que una misma molécula lleva todas las sustituciones del codón; una coincidencia por SNV, que una lleva una sustitución sin las demás. Con --bam las moléculas pueden desmentir cualquiera de las dos. Una entrada que coincide solo por una lectura desmentida nombra un aminoácido que nada en la muestra lleva, y get_MNV la retiene: la entrada pasa a la columna Catalogue Withheld (la clave INFO CATW en el VCF), con el motivo.

Catalogue           -
Catalogue Withheld  RIF: Assoc w R (rpoB p.Ser450Tyr, codon reading only: no molecule carries the whole codon)

La entrada se sigue listando, así que no se oculta nada, pero deja de afirmarse como encontrada. Una coincidencia que hacen las dos lecturas no se toca, ni tampoco una coincidencia de nucleótido (c.1349C>A, g.761155C>A): esa base está en las moléculas, acompañada o no.

Deciden las moléculas, D' lo confirma. La lectura de codón queda desmentida cuando ninguna de las moléculas que cubren el codón lleva todas las sustituciones (MNV Reads es 0) y D' está en -0,5 o por debajo; la lectura propia de una sustitución, cuando toda molécula que la lleva lleva el codón entero (su SNP Reads es 0 y MNV Reads no) y D' está en 0,5 o por encima. D' solo no basta: con las dos sustituciones en más de la mitad de las moléculas, alguna tiene que compartirlas, así que D' llega a -1 con un 30% de las moléculas llevando el codón entero, y ese aminoácido está en la muestra. Requiere --bam, porque sin moléculas no hay nada que contar y nada cambia. Dos umbrales gobiernan el lado de D', y hay que cumplir los dos:

Exigir ambos es deliberado. El p-valor solo dejaría que una muestra profunda desmintiera una lectura con D' de -0,05, una exclusión medible y sin sentido; la magnitud sola dejaría que lo hicieran tres moléculas. --keep-unphased-codon-matches lista todas esas coincidencias en Catalogue como cualquier otra.

Por qué se retiene en vez de borrarse#

Dos muestras pueden dar VCF idénticos byte a byte y significar lo contrario: un 30% de moléculas con S450Y, frente a una parada prematura y un cambio sinónimo al 30% cada uno en moléculas distintas. La primera tiene una subpoblación resistente; la segunda no, y llamarla S450Y inventa una variante que nadie lleva. El VCF no las distingue y las lecturas sí, así que la salida honesta nombra la entrada, dice que se preguntó a las moléculas y dice qué contestaron.

Límites#

• El estadístico es dentro de una ventana: un codón, o una ventana de haplotipo local. get_MNV no infiere fase de largo alcance, y una variante fuera de la ventana no forma parte de la tabla.
• Describe las moléculas que se secuenciaron. El ligamiento a nivel de lectura en un esquema de amplicones está acotado por el amplicón, y un par separado por más de un fragmento simplemente no tiene moléculas en común que contar.
• No dice cuál de dos linajes en competencia es cuál. Un D' cerca de -1 reporta que existen dos; identificarlos es trabajo de una herramienta de fasing o de deconvolución de cepas.

Editar esta página en GitHub