El fin del GWAS a gran escala: "Asociación" e "Inferencia causal" no son lo mismo

Hasta S43 (REGENIE y SAIGE), hemos cubierto cómo encontrar rápidamente y con precisión SNPs estadísticamente significativos en biobancos de medio millón de individuos. Sin embargo, al abrir la tabla de resultados del GWAS, hay algo extraño: en lugar de un solo SNP significativo cerca de un locus (sitio genómico), aparecen decenas o incluso cientos agrupados.

Esto no significa que haya cientos de variantes causales en esa región. La mayoría son simplemente SNPs vecinos que están cerca de la verdadera variante causal y se heredan juntos (co-inherited). Este fenómeno se denomina desequilibrio de ligamiento (LD, Linkage Disequilibrium). Dos SNPs físicamente cercanos en el genoma no se separan fácilmente mediante recombinación meiótica, por lo que se heredan como un bloque. Como resultado, los SNPs "falsos compañeros" estadísticamente indistinguibles de la verdadera variante causal obtienen valores p igualmente significativos.

En esta sección, derivamos los principios de la estadística bayesiana del mapeo fino (Fine-mapping), que reduce las candidatas a variantes causales más probables dentro de una "región asociada", y examinamos dos herramientas implementadas: CAVIAR y SuSiE. Sin embargo, es importante aclarar desde el principio: la probabilidad proporcionada por el mapeo fino es la probabilidad posterior dentro del modelo estadístico dado y el conjunto de SNPs candidatos, no la relación causal biológica confirmada mediante experimentos funcionales.

La ilusión estadística creada por el LD y su resolución mediante factores bayesianos aproximados

La ilusión en situaciones de LD perfecto

Imaginemos un caso extremo: el SNP A es la verdadera variante causal, y el SNP B es un vecino con un coeficiente de correlación $r=1$ (LD perfecto) con A. Es decir, los genotipos de A y B son siempre idénticos en todos los individuos de esta cohorte. En este caso, si realizamos un análisis de regresión de un solo variante, A y B obtendrán z-scores completamente idénticos. Solo mirando las estadísticas, no hay forma de distinguir cuál es la verdadera causa.

Asignación de la "probabilidad de causalidad" para cada SNP mediante el Factor Bayesiano Aproximado (ABF)

El factor bayesiano aproximado (ABF) de Wakefield (2009) responde a este problema de manera probabilística. Si definimos la varianza de la estimación del efecto observado para cada SNP $j$ como $V_j$, y la distribución previa del tamaño del efecto verdadero como $\theta \sim \mathcal{N}(0, W)$, entonces para el factor de contracción (shrinkage factor) $r_j = \frac{W}{V_j + W}$:

$\text{ABF}_j = \sqrt{1 - r_j} \, \exp\!\left(\frac{r_j Z_j^2}{2}\right)$

y aquí $Z_j$ es la puntuación z del SNP $j$. Bajo el supuesto simplificador de que "hay exactamente una variante causal dentro de este locus", la probabilidad posterior de inclusión (PIP, Posterior Inclusion Probability) de que el SNP $j$ sea esa variante causal se normaliza como sigue:

$\text{PIP}_j = \frac{\text{ABF}_j}{\sum_{k=1}^{M} \text{ABF}_k}$

Ejemplo de cálculo manual: El PIP cambia drásticamente según la presencia de LD

Fijemos $r=0.5$ asumiendo que la varianza a priori y la varianza observada son iguales, y calculemos dos escenarios.

Escenario 1 (Sin LD, $r_{LD}=0$): El SNP1 es la verdadera variante causal con $Z_1 = 5$, y el SNP2 es un vecino no relacionado con $Z_2 = 0.1$.

$\text{ABF}_1 = \sqrt{0.5}\exp\left(\frac{0.5 \times 25}{2}\right) = 0.7071 \times e^{6.25} \approx 0.7071 \times 518.01 \approx 366.3$

$\text{ABF}_2 = \sqrt{0.5}\exp\left(\frac{0.5 \times 0.01}{2}\right) = 0.7071 \times e^{0.0025} \approx 0.7071 \times 1.0025 \approx 0.709$

$\text{PIP}_1 = \frac{366.3}{366.3 + 0.709} \approx 0.998, \quad \text{PIP}_2 \approx 0.002$

La señal se concentra casi al 100% en el SNP1: al no haber LD, el modelo estadístico reduce la lista de candidatos a variante causal únicamente al SNP1 (esta conclusión se basa en el supuesto del modelo de que "exactamente uno entre los SNPs candidatos es causal", y en la premisa de que las puntuaciones z observadas reflejan bien los valores verdaderos).

Escenario 2 (LD perfecto, $r_{LD}=1$): Como se explicó anteriormente, el SNP1 y el SNP2 tienen genotipos completamente idénticos, por lo que $Z_1 = Z_2 = 5$.

$\text{ABF}_1 = \text{ABF}_2 \approx 366.3 \implies \text{PIP}_1 = \text{PIP}_2 = \frac{366.3}{366.3 + 366.3} = 0.5$

Ambos SNPs tienen un PIP de 0.5, empatando. Si construimos un conjunto creíble del 95% (un conjunto que incluye SNPs ordenados por PIP decreciente hasta que la PIP acumulada supere 0.95), en este caso obtenemos un conjunto de tamaño 2 que incluye tanto al SNP1 como al SNP2. Esto refleja fielmente el hecho de que "uno de estos dos SNPs es la variante causal, pero los datos por sí solos no permiten estar seguros de cuál", lo cual no es un fallo del método.

CAVIAR: Exploración combinatoria de configuraciones causales

En un locus real, hay cientos o miles de SNPs y puede haber más de una variante causal. CAVIAR (Hormozdiari et al., 2014) modela el vector completo de puntuaciones z como una distribución normal multivariante que sigue una matriz de correlación de LD $\mathbf{R}$ como covarianza, y asigna probabilidades posteriores a todas las combinaciones posibles de variantes causales (configuraciones) bajo la restricción de que "como máximo $k$ de ellas son variantes causales". Cuando hay $M$ SNPs, el número de combinaciones es $\sum_{i=0}^{k}\binom{M}{i}$, lo que provoca una explosión combinatoria; por ello, en la práctica se limita $k$ a valores como 2 o 3.

SuSiE: suma iterativa de efectos individuales

SuSiE (Sum of Single Effects, Wang et al., 2020) no explora todas las combinaciones posibles, sino que descompone el vector completo de efectos $\boldsymbol{\beta}$ en la suma de $L$ vectores de efecto individual.

$\boldsymbol{\beta} = \sum_{l=1}^{L} \boldsymbol{\beta}^{(l)}, \quad \boldsymbol{\beta}^{(l)} = \gamma^{(l)} b^{(l)}$

donde $\gamma^{(l)}$ es un vector one-hot con exactamente un elemento igual a 1 y el resto iguales a 0 (es decir, se asume que la señal $l$ proviene de un solo SNP), y $b^{(l)}$ es la magnitud del efecto de ese SNP. Se entrena todo el modelo mediante el algoritmo IBSS (Selección Bayesiana Paso a Paso Iterativa), que ajusta iterativamente cada $\boldsymbol{\beta}^{(l)}$ mediante regresión bayesiana de efecto individual. Una iteración de IBSS tiene una complejidad aproximada de $O(nML)$ para un tamaño de muestra $n$, un número de SNPs $M$ y un número de señales $L$. Dado que se obtiene un conjunto creíble independiente para cada una de las $L$ señales, SuSiE puede manejar loci con múltiples señales sin la explosión combinatoria típica del enfoque de CAVIAR.

Práctica en R: cálculo de PIP y conjuntos creíbles

# Cálculo de PIP basado en ABF de Wakefield — Reproducción del ejemplo manual con código
compute_pip <- function(z, r_shrink = 0.5) {
  abf <- sqrt(1 - r_shrink) * exp(r_shrink * z^2 / 2)
  abf / sum(abf)
}

# Escenario 1: Sin LD
z_no_ld <- c(5, 0.1)
cat("PIP sin LD:", round(compute_pip(z_no_ld), 4), "\n")

# Escenario 2: LD perfecto (dos SNP con mismo puntaje z)
z_perfect_ld <- c(5, 5)
cat("PIP con LD perfecto:", round(compute_pip(z_perfect_ld), 4), "\n")

# Función para extraer el conjunto creíble al 95%
credible_set <- function(pip, threshold = 0.95) {
  ord <- order(pip, decreasing = TRUE)
  cum <- cumsum(pip[ord])
  ord[1:which(cum >= threshold)[1]]
}

cat("Índices del conjunto creíble al 95% con LD perfecto:", credible_set(compute_pip(z_perfect_ld)), "\n")

Mapeo de CS

• Factor de Bayes (Bayes Factor): Es la razón de las verosimilitudes de que se observen los datos bajo cada una de dos hipótesis (nula frente a alternativa). El fine-mapping se reduce al problema de calcular el factor de Bayes por hipótesis para la afirmación "este SNP es causal" y normalizarlo como una probabilidad posterior.
• Regresión de efecto único (Single-Effect Regression) + suma iterativa: Aproximar un problema complejo de señales múltiples como la suma de varios modelos simples donde "solo hay un componente activo", es una forma de pensamiento aditivo iterativo similar a cómo el boosting por gradiente construye un modelo fuerte sumando secuariamente aprendices débiles (weak learners).
• Búsqueda combinatoria y poda: Que CAVIAR limite el número de variantes causales a $k$ es una estrategia típica de optimización combinatoria que realiza una poda para mantener el rango computable, limitando la profundidad y el ancho de ramificación en un árbol de búsqueda que crece exponencialmente.

Defectos comunes

• No uso del LD in-sample: Si se utiliza el LD del panel de referencia (por ejemplo, 1000 Genomes) para realizar fine-mapping con estadísticas resumidas, la estructura de LD de la cohorte real del estudio puede no coincidir, lo que podría distorsionar el credible set. Siempre que sea posible, se debe utilizar el LD calculado directamente desde la cohorte original.
• Subestimación del número de variantes causales $L$·$k$: Si en un locus con dos señales reales se fija $L=1$, las dos señales se fusionarán en una sola, lo que hará que el credible set sea excesivamente grande o incluya SNPs incorrectos.
• No verificación de la pureza (Purity): Cada credible set de SuSiE incluye un indicador de pureza que combina valores como el mínimo y el promedio de la correlación (LD) entre los SNPs internos. En particular, una correlación absoluta mínima baja es una señal de que la confianza en ese credible set disminuye. No se debe juzgar solo por el PIP; también hay que verificar la pureza.
• Confusión del credible set con garantía de muestreo repetido: Un "credible set del 95%" es la cobertura teórica al aplicar el mismo modelo de manera repetida, no una garantía incondicional de que "el verdadero variante está contenido con un 95% de probabilidad en este conjunto específico de datos".

Para profundizar más

El texto ha sido reconstruido directamente por el equipo de investigación de BPD. Adentremosnos más en el artículo original y la documentación oficial.

• Artículo original de CAVIAR: Hormozdiari et al. (2014), Identifying Causal Variants at Loci with Multiple Signals of Association, Genetics, 198(2):497–508.
• Artículo original de SuSiE: Wang et al. (2020), A simple new approach to variable selection in regression, with application to genetic fine-mapping, Journal of the Royal Statistical Society: Series B, 82(5):1273–1300.
• Artículo original de FINEMAP: Benner et al. (2016), FINEMAP: efficient variable selection using summary data from genome-wide association studies, Bioinformatics, 32(10):1493–1501.
• Artículo original del factor de Bayes de Wakefield: Wakefield (2009), Bayes factors for genome-wide association studies: comparison with P-values, Genetic Epidemiology, 33(1):79–86.

Hemos aprendido a desentrañar las ilusiones estadísticas creadas por el LD mediante factores de Bayes. Sin embargo, identificar la variante causal no es el final del camino; surge otra pregunta: ¿en qué tejido y a través de qué gen modifica la expresión génica para influir en la enfermedad?

En el siguiente capítulo S45, exploraremos el mapeo de eQTL (locus de rasgo cuantitativo de expresión) que conecta directamente las variantes genéticas con los niveles de expresión génica, y examinaremos la especificidad de la regulación génica que varía entre tejidos utilizando datos de GTEx.