Empirical Bayes data integreation for multi-response regression

1. Planteamiento del Problema

El artículo aborda el desafío de integrar datos provenientes de múltiples fuentes (por ejemplo, diferentes tejidos biológicos en estudios de asociación a nivel de transcriptoma o TWAS) para estimar coeficientes de regresión.

• Contexto: En estudios como TWAS, se busca predecir la expresión génica en múltiples tejidos utilizando datos de genotipos (SNPs). Tradicionalmente, cada tejido se analiza por separado, perdiendo información compartida, o se utilizan métodos de regresión multivariada que asumen estructuras específicas en los parámetros.
• Limitaciones de los métodos existentes:
  • Los métodos de regresión de rango reducido y los basados en penalización (sparsity) asumen que la matriz de coeficientes es de bajo rango o dispersa (sparse).
  • Estas suposiciones estructurales son difíciles de verificar en la práctica y a menudo no se cumplen en datos biológicos reales (donde los efectos pueden ser densos o de rango mixto).
  • Los enfoques de Bayes completo son computacionalmente intensivos y poco escalables para grandes conjuntos de datos.
• Objetivo: Desarrollar un método que integre datos de múltiples fuentes, sea escalable computacionalmente y funcione bien independientemente de la estructura subyacente de los parámetros (ya sea densa, dispersa, de bajo rango o no).

2. Metodología Propuesta

Los autores proponen un enfoque de Bayes Empírico basado en la contracción (shrinkage) de valores singulares.

A. Marco Teórico y Modelo

• Se asume un modelo de regresión lineal para cada fuente de datos (tejido $t$): $y_i^{(t)} = x_i^T \beta^{(t)} + \epsilon_i^{(t)}$.
• Los estimadores de Mínimos Cuadrados Ordinarios (OLS), $\hat{\beta}^{(t)}$, se tratan como datos observados bajo un modelo aditivo con ruido gaussiano: $\hat{\beta}^{(t)} = \beta^{(t)} + u^{(t)}$, donde $u^{(t)} \sim N(0, Q)$.
• El objetivo es estimar la matriz de coeficientes $B$ (donde las filas son $\beta^{(t)}$) minimizando la pérdida de Frobenius.

B. Estimador de Contracción Lineal Global

El método conecta la estimación de la matriz de coeficientes con la estimación de la matriz de covarianza marginal de los datos estandarizados.

1. Estimación de Covarianza: Se desarrolla un estimador óptimo asintótico para la matriz de covarianza $\Sigma$ bajo una pérdida de ahorro relativo (relative savings loss).
2. Regla de Contracción: Se propone un estimador de contracción lineal de la forma:
   $$\bar{B} = \hat{B}(I - Q^{1/2}\hat{\Sigma}_n Q^{1/2})$$
   Donde $\hat{\Sigma}_n$ es un estimador de la covarianza que contrae los autovalores de la matriz de covarianza muestral.
3. Optimalidad: Se demuestra que esta regla es asintóticamente óptima dentro de la clase de estimadores de contracción lineal, sin asumir que la matriz de coeficientes es de bajo rango o dispersa.

C. Estimación de la Matriz de Covarianza

• Se utilizan técnicas de análisis espectral y teoría de matrices aleatorias de alta dimensión (Ley de Marcenko-Pastur).
• Se deriva una regla de contracción óptima para los autovalores de la matriz de covarianza muestral $S_n$.
• Selección de Parámetros: Se utiliza un estimador de riesgo casi no sesgado (SURE - Stein's Unbiased Risk Estimate) para seleccionar automáticamente el parámetro de suavizado (bandwidth) necesario para la contracción, haciendo el método dependiente de los datos y robusto.
• Casos $p > n$: El método se extiende para manejar situaciones donde el número de predictores ($p$) es mayor que el número de muestras ($n$), utilizando una regla de contracción específica para los autovalores nulos.

D. Estimador de Contracción Lineal Local (Mezcla)

Para mayor flexibilidad, se extiende el método a un prior de mezcla de Gaussianas.

• En lugar de una contracción global única, se asume que los coeficientes provienen de una mezcla de distribuciones Gaussianas.
• Esto resulta en una regla de contracción localmente adaptativa (Local Linear Shrinkage - LLS), donde el grado de contracción depende de la magnitud de la señal observada.
• Se implementa mediante un algoritmo de MCMC (Data Augmentation) para estimar los componentes de la mezcla y las covarianzas asociadas.

3. Contribuciones Clave

1. Independencia de Estructura: A diferencia de los métodos de rango reducido o sparse, el enfoque propuesto no requiere que la matriz de coeficientes tenga una estructura específica (bajo rango o dispersa), lo que lo hace aplicable a escenarios donde estas suposiciones fallan.
2. Optimalidad Asintótica: Se demuestra teóricamente que el estimador de contracción lineal es óptimo asintóticamente bajo la pérdida de Frobenius, basándose en la estimación óptima de la matriz de covarianza bajo pérdida de ahorro relativo.
3. Escalabilidad Computacional: Al ser un enfoque de Bayes Empírico con soluciones analíticas y algoritmos eficientes, es mucho más escalable que los procedimientos de Bayes completo (MCMC intensivo) o métodos de optimización convexa complejos.
4. Selección Automática de Parámetros: El uso de SURE permite una selección de parámetros de suavizado basada en los datos, mejorando el rendimiento del riesgo en comparación con reglas fijas.
5. Generalidad: El método funciona tanto cuando el número de fuentes es mayor que el número de variables ($n > p$) como cuando es menor ($n < p$).

4. Resultados

Experimentos de Simulación

Se comparó el método propuesto (ULS - Unbiased Linear Shrinkage y LLS - Local Linear Shrinkage) con métodos existentes:

• UTMOST: Asume estructura de rango reducido/esparcidad.

• ISA (Iterated stable autoencoder): Asume estructura de bajo rango.

• MASH (Multivariate Adaptive Shrinkage): Enfoque Bayes Empírico con mezclas.

• Rendimiento:

  • En escenarios donde la verdad subyacente era de bajo rango, el método ISA funcionó bien, pero el método propuesto (ULS) tuvo un rendimiento muy cercano.
  • En escenarios densos, dispersos aproximados, o de mezcla (Horseshoe/Mixture), el método propuesto superó consistentemente a todos los demás, especialmente en términos de Error Cuadrático Medio (MSE).
  • El método MASH mostró alta variabilidad y problemas de convergencia en escenarios de alta dimensión o con priors de cola pesada (Horseshoe).
  • El estimador local (LLS) mostró un rendimiento superior o comparable a los diseñados específicamente para estructuras complejas, sin necesidad de asumir dicha estructura explícitamente.

Aplicación a Datos Reales (GTEx)

• Datos: Se aplicó al proyecto GTEx (Genotype-Tissue Expression) para predecir la expresión génica en 32 tejidos utilizando 3,201 genes y sus SNPs cis asociados.
• Validación: Se utilizó validación cruzada de 10 pliegues.
• Resultados:
  • El método propuesto (LS) logró el menor Error Cuadrático Medio de Predicción (PMSE) promedio, superando ligeramente a OLS, UTMOST y MASH.
  • En genes con pocos SNPs (hasta 5 o 10), el método propuesto y UTMOST fueron los mejores.
  • Interpretabilidad: Los mapas de calor de los efectos estimados mostraron que el método propuesto permite variabilidad en los efectos entre tejidos (capturando diferencias biológicas reales), mientras que MASH tendía a fusionar efectos a valores idénticos y UTMOST a contraerlos excesivamente a cero.

5. Significancia y Conclusión

Este trabajo ofrece una solución robusta y teóricamente fundamentada para la integración de datos en regresión multi-respuesta, un problema común en bioinformática y otras disciplinas.

• Impacto Científico: Proporciona una alternativa superior a los métodos de "rango reducido" o "sparse" cuando la estructura de los datos es desconocida o compleja. La capacidad de manejar datos de alta dimensión ($p > n$) y la eficiencia computacional lo hacen ideal para estudios genómicos modernos.
• Innovación Metodológica: La conexión explícita entre la estimación de coeficientes de regresión y la estimación óptima de matrices de covarianza bajo pérdida de ahorro relativo es una contribución teórica significativa.
• Aplicabilidad: La aplicación exitosa a datos de GTEx demuestra su utilidad práctica para mejorar la predicción de expresión génica y la identificación de marcadores genéticos, facilitando la comprensión de la arquitectura genética de enfermedades complejas.

En resumen, el método propuesto logra un equilibrio óptimo entre flexibilidad, precisión estadística y eficiencia computacional, superando a los enfoques actuales en una amplia gama de escenarios de simulación y en datos biológicos reales.