La contraparte hacia adelante del Coalescente

En M29 aprendimos a rastrear genealogías hacia atrás en el tiempo mediante el coalescente. La teoría en la que se asienta es el modelo de probabilidad hacia adelante conocido como el modelo Wright-Fisher, establecido por Sewall Wright y Ronald Fisher respectivamente en la década de 1930.

Wright-Fisher es la base de la genética de poblaciones. Representa cómo cambian las frecuencias alélicas entre generaciones en una población ideal de tamaño N mediante un proceso estocástico discreto. De este modelo se derivan la teoría neutral de Kimura, el coalescente e incluso los análisis más recientes de pangenomas.

Este capítulo es el último de la serie Micro Tier 30 y la culminación de la Fase A.

Definición del modelo

Tamaño poblacional N, dos alelos A y a. Frecuencia de A en el tiempo t: p_t.

Definición de Wright-Fisher: Cada individuo de la siguiente generación (t+1) selecciona aleatoriamente, de forma independiente y con reemplazo, un progenitor de la generación anterior (t). La probabilidad de que cada descendiente herede A es p_t.

Es decir, el número de individuos con alelo A en la siguiente generación, X, se extrae de una distribución binomial: X ~ Binomial(N, p_t). La nueva frecuencia es p_siguiente = X / N.

Este proceso estocástico discreto constituye todo el modelo. De una definición sorprendentemente simple se derivan resultados sumamente ricos.

Deriva genética

La predicción clave de Wright-Fisher es que las frecuencias alélicas fluctúan incluso sin presión selectiva. Se trata del error de muestreo aleatorio derivado de la distribución binomial. A esto se le llama deriva genética.

La varianza del cambio de frecuencia por generación es p_t * (1 - p_t) / N.

La deriva es mayor en poblaciones más pequeñas y menor cuando la frecuencia se acerca a los extremos (0 o 1). Finalmente, el alelo llegará a fijación (frecuencia 1) o a pérdida (frecuencia 0). Estos dos estados son estados absorbentes: una vez alcanzados, no se puede salir de ellos.

Probabilidad y tiempo de fijación

La probabilidad de que un alelo se fije eventualmente, partiendo de una frecuencia inicial p_0, es exactamente p_0.

Un resultado sorprendentemente simple: la frecuencia inicial es la propia probabilidad de fijación. Una nueva mutación aparece en un solo individuo durante una generación, por lo que su frecuencia inicial es 1 / (2N); por tanto, la probabilidad de fijación de una nueva mutación neutra es 1 / (2N).

El tiempo medio hasta la fijación (condicionado a que ocurra) es del orden de 4N generaciones. Si p_0 es pequeño, este tiempo se vuelve muy largo. Una nueva mutación que comienza en un solo individuo tarda una media de 4N generaciones en fijarse. Para humanos (con N aproximadamente 10,000), esto equivale a 40,000 generaciones, es decir, unos 1 millón de años.

Teoría neutral de la evolución de Kimura

Motoo Kimura (Nature 1968) amplió los resultados de este modelo a la teoría evolutiva en su conjunto con su teoría de la evolución neutral. La idea principal es:

La mayoría de las variaciones a nivel molecular se fijan por deriva neutral y no por selección natural.

Esta afirmación generó intensos debates cuando se publicó, al entrar en conflicto directo con la visión tradicional de que la selección natural darwiniana era el único motor de la evolución. Hoy en día, la perspectiva de Kimura se reconoce como la base de la evolución molecular.

Predicciones cuantitativas de la teoría neutral: La tasa de sustitución por sitio por generación es igual a la tasa de mutación por generación, $\mu$. En detalle: $(2N \times \mu) \times (1 / 2N) = \mu$. El tamaño de la población $N$ se cancela y solo queda $\mu$. El hecho de que las tasas de evolución molecular observadas entre especies coincidan con esta predicción fue una prueba sólida de la teoría neutral.

Simulación en R (20 líneas)

wright_fisher <- function(N, p0, generations, seed = NULL) {
  if (!is.null(seed)) set.seed(seed)
  p <- numeric(generations + 1)
  p[1] <- p0
  for (t in 2:(generations + 1)) {
    X <- rbinom(1, N, p[t-1])
    p[t] <- X / N
    if (p[t] == 0 || p[t] == 1) {
      p[t:(generations+1)] <- p[t]
      break
    }
  }
  return(p)
}

set.seed(42)
N <- 100
p0 <- 0.5
n_traj <- 20

trajectories <- replicate(n_traj, wright_fisher(N, p0, 500))
matplot(0:500, trajectories, type = 'l', lty = 1,
        xlab = "Generation", ylab = "Allele frequency",
        main = paste0("Wright-Fisher, N=", N, ", p0=", p0))
abline(h = c(0, 1), lty = 2, col = "gray")

Al graficar el trazado, las 20 trayectorias oscilan de manera distinta hasta que la mayoría se estabiliza en 0 o 1. La diversidad desaparece incluso sin presión selectiva: esta es la demostración visual de la deriva genética.

Simulación en Python y validación estadística

import random

def wright_fisher(N, p0, generations, seed=None):
    if seed is not None:
        random.seed(seed)
    p = [p0]
    for _ in range(generations):
        X = sum(1 for _ in range(N) if random.random() < p[-1])
        pt = X / N
        p.append(pt)
        if pt == 0 or pt == 1:
            p.extend([pt] * (generations - len(p) + 1))
            break
    return p

N = 100
p0 = 0.3
n_sims = 10000
fixed = 0
for i in range(n_sims):
    p = wright_fisher(N, p0, 500, seed=i)
    if p[-1] == 1:
        fixed += 1

empirical = fixed / n_sims
print("empirical fixation:", empirical)
print("theoretical (p0):", p0)

Al realizar 10,000 simulaciones para obtener la probabilidad fija empírica, se obtiene un valor muy cercano a $p_0 = 0.3$. La concordancia entre la teoría y la simulación se verifica estadísticamente.

Extensión — Selección, Migración, Mutación

Se pueden añadir tres fuerzas adicionales al modelo básico de Wright-Fisher.

Selección natural (Selection): Cada alelo tiene una aptitud $w_A, w_a$. El valor esperado de la frecuencia de $A$ en la siguiente generación se ajusta como $p_t$ multiplicado por $w_A$ dividido por ($p_t$ multiplicado por $w_A$ más ($1$ menos $p_t$) multiplicado por $w_a$), y este valor se utiliza como el parámetro $p$ de la distribución binomial.

Migración (Migration): Proporción $m$ de alelos que ingresan desde otra población en cada generación. La nueva frecuencia es $(1 \text{ menos } m)$ por $p_t$ más $m$ por $p_{\text{source}}$.

Mutación (Mutation): Tasa de mutación $\mu$ por generación de $A$ a $a$. La nueva frecuencia es $(1 \text{ menos } \mu)$ por $p_t$ más $\nu$ por $(1 \text{ menos } p_t)$.

El punto donde estas tres fuerzas y la deriva genética se equilibran constituye la distribución estacionaria, como el equilibrio mutación-deriva o el equilibrio selección-deriva. La mayor parte de la genética de poblaciones consiste en el análisis de estos estados de equilibrio.

Limitaciones de Wright-Fisher

• Supuesto de ausencia de solapamiento generacional (en poblaciones reales el solapamiento es alto).
• Supuesto de apareamiento aleatorio (en poblaciones reales existe una gran estructura geográfica y jerárquica).
• Supuesto de proporción sexual equilibrada.
• Supuesto de generaciones discretas (el modelo de generaciones continuas es el modelo de Moran).

A pesar de estas limitaciones, la razón por la que Wright-Fisher persiste es su simetría con el modelo coalescente. El hecho de que la dirección forward de Wright-Fisher y la dirección backward del modelo coalescente sean estadísticamente equivalentes es la columna vertebral de la teoría actual de la genética de poblaciones.

Modelo de Moran — Alternativa de generaciones continuas

Si Wright-Fisher utiliza generaciones discretas, el modelo de Moran utiliza generaciones continuas. En cada instante, un individuo muere y otro nace como descendiente de un progenitor elegido al azar. Debido a que posee características matemáticamente más manejables, se utiliza con frecuencia en investigación teórica. Para simulaciones prácticas, Wright-Fisher es el estándar.

SLiM — Herramienta práctica para simulaciones de Wright-Fisher

SLiM(Bioinformatics 2019) es la herramienta estándar para simulaciones complejas de Wright-Fisher. Admite selección, migración, estructura poblacional y recombinación cromosómica. Los scripts se escriben en el lenguaje Eidos.

// Ejemplo de SLiM: Expansión y selección natural tras un cuello de botella
initialize() {
    initializeMutationRate(1e-7);
    initializeMutationType("m1", 0.5, "f", 0.0);
    initializeMutationType("m2", 0.5, "f", 0.05);
    initializeGenomicElementType("g1", c(m1, m2), c(0.95, 0.05));
    initializeGenomicElement(g1, 0, 999999);
    initializeRecombinationRate(1e-8);
}

1 { sim.addSubpop("p1", 10000); }
1000 { p1.setSubpopulationSize(1000); }
2000 { p1.setSubpopulationSize(10000); }
3000 late() { sim.outputFull(); }

Al ejecutar este script durante varias horas, se pueden observar los patrones de diversidad genética de un escenario específico bajo las predicciones de Wright-Fisher.

Correspondencia CS — La esencia de las cadenas de Markov

Wright-Fisher es una cadena de Markov con un espacio de estados que va de (0/N) a (N/N). El siguiente estado depende únicamente del estado actual.

• Matriz de transición: (N + 1) x (N + 1). Cada elemento es una probabilidad de la distribución binomial.
• Estados absorbentes: estado 0 (pérdida) y N (fijación).
• Distribución estacionaria: no existe debido a los escenarios de absorción; en su lugar, se definen las probabilidades y los tiempos de absorción.

El HMM de M18 también es una cadena de Markov. La diferencia radica en los estados ocultos frente a los estados observados. Una sola cadena de Markov puede adoptar múltiples formas. Modelos de secuencias (HMM), genética de poblaciones (Wright-Fisher), aprendizaje por refuerzo (MDP), PageRank. Al dominar este concepto, se pueden analizar instantáneamente los procesos de Markov en nuevos dominios.

Evaluación en Rosalind

Problema WFMD de Rosalind. Cálculo de la probabilidad generacional de Wright-Fisher. Implementación directa de la versión probabilística teórica de la función Python anterior.

Finalización de la Fase A y guía para los próximos episodios

M30 es el último de los 30 episodios del Micro Tier (de M01 a M30). Se completa la Fase A (matemáticas de los algoritmos centrales). Repasemos lo aprendido hasta ahora.

• Algoritmos de secuencias y alineamiento (Needleman-Wunsch, Smith-Waterman, BLAST)
• Índices de cadenas (árbol de sufijos, BWT, FM-index, minimap2, Diamond)
• HMM (Viterbi, Forward-Backward, Baum-Welch, Profile HMM)
• Ensamblaje (OLC, De Bruijn, hifiasm)
• MSA, filogenia y genética de poblaciones (MAFFT, NJ, Fitch, IQ-TREE, coalescencia, Wright-Fisher)

Cada uno de los 30 episodios es una obra original en coreano de entre 8 y 14 KB. A partir de ahora, nos preparamos para la Fase B — System S1 Tier (25 episodios de pipelines de NGS). La teoría del Micro Tier dará sus frutos en los pipelines prácticos que comienzan desde FASTQ.

Si desea profundizar más

• Kimura (1968), Evolutionary rate at the molecular level, Nature — Artículo original de la teoría neutral de la evolución.
• Ewens (2004), Mathematical Population Genetics — Texto estándar. Cubre Wright-Fisher, coalescencia y selección.
• Haller y Messer (2019), SLiM 3: forward genetic simulations beyond the Wright-Fisher model, MBE — artículo original de SLiM.
• Harvard HST.508 — Clase sobre Wright-Fisher (con subtítulos completos). Notas de Jonathan Pritchard.
• Stanford BIOMEDIN 214 — Clase de genética de poblaciones. Análisis de datos prácticos.

Ejecute el código de simulación en R y observe las trayectorias de la deriva genética variando el tamaño de la población (N) a 10, 100, 1000 y 10000. El último ejercicio práctico de esta sección consiste en confirmar visualmente la predicción de que la deriva es más lenta cuanto mayor es la población.