Corrección de los valores de p para múltiples pruebas donde las pruebas están correlacionadas (genética)

24

Tengo valores de p de muchas pruebas y me gustaría saber si realmente hay algo significativo después de corregir las pruebas múltiples. La complicación: mis pruebas no son independientes. El método en el que estoy pensando (una variante del Método del producto de Fisher, Zaykin et al., Genet Epidemiol , 2002) necesita la correlación entre los valores de p.

Para estimar esta correlación, actualmente estoy pensando en casos de bootstrapping, ejecutando los análisis y correlacionando los vectores resultantes de los valores de p. ¿Alguien tiene una idea mejor? ¿O incluso una mejor idea para mi problema original (corregir las pruebas múltiples en pruebas correlacionadas)?

Antecedentes: estoy retrocediendo logísticamente si mis sujetos sufren o no una enfermedad particular en la interacción entre su genotipo (AA, Aa o aa) y una covariable. Sin embargo, el genotipo es en realidad una gran cantidad (30-250) de polimorfismos de nucleótido único (SNP), que ciertamente no son independientes, sino que están en desequilibrio de ligamiento.

S. Kolassa - Restablece a Monica
fuente

Respuestas:

29

¡Este es realmente un tema candente en los estudios de análisis de Genomewide (GWAS)! No estoy seguro de que el método en el que esté pensando sea el más apropiado en este contexto. Algunos autores describieron la agrupación de los valores de p, pero en un contexto diferente (estudios de replicación o metanálisis, ver, por ejemplo, (1) para una revisión reciente). La combinación de los valores p de SNP por el método de Fisher se usa generalmente cuando se quiere derivar un valor p único para un gen dado; esto permite trabajar a nivel genético y reducir la cantidad de dimensionalidad de las pruebas posteriores, pero como dijiste, la no independencia entre los marcadores (derivados de la colocación espacial o el desequilibrio de enlace, LD) introduce un sesgo. Las alternativas más potentes dependen de los procedimientos de remuestreo,

Mi principal preocupación con bootstraping (con reemplazo) sería que esté introduciendo una forma artificial de relación, o en otras palabras, cree gemelos virtuales, alterando así el equilibrio de Hardy-Weinberg (pero también la frecuencia mínima de alelos y la tasa de llamadas). Este no sería el caso con un enfoque de permutación donde permuta etiquetas individuales y mantiene los datos de genotipado tal como están. Por lo general, el software plink puede proporcionarle valores p sin procesar y permutados, aunque utiliza (por defecto) una estrategia de prueba adaptativa con una ventana deslizante que permite detener la ejecución de todas las permutaciones (digamos 1000 por SNP) si parece que el SNP está bajo la consideración no es "interesante"; También tiene la opción de calcular maxT, consulte la ayuda en línea .

Pero dado el bajo número de SNP que está considerando, sugeriría confiar en las pruebas basadas en FDR o maxT implementadas en el paquete R de prueba múltiple (ver mt.maxT), pero la guía definitiva para las estrategias de remuestreo para la aplicación genómica es Procedimientos de prueba múltiple con aplicaciones para Genomics , de Dudoit y van der Laan (Springer, 2008). Ver también el libro de Andrea Foulkes sobre genética con R , que se revisa en el JSS. Ella tiene un gran material sobre múltiples procedimientos de prueba.

Notas adicionales

Muchos autores han señalado el hecho de que los métodos simples de corrección de pruebas múltiples, como Bonferroni o Sidak, son demasiado estrictos para ajustar los resultados para los SNP individuales. Además, ninguno de estos métodos tiene en cuenta la correlación que existe entre los SNP debido a LD que marca la variación genética en las regiones genéticas. Se han propuesto otras alternativas, como una derivada del método de Holm para comparación múltiple (3), Modelo de Markov Oculto (4), FDR condicional o positivo (5) o derivada del mismo (6), por nombrar algunos. Las llamadas estadísticas de huecos o ventana deslizante han tenido éxito en algunos casos, pero encontrará una buena revisión en (7) y (8).

También he oído hablar de métodos que hacen un uso efectivo de la estructura de haplotipo o LD, por ejemplo (9), pero nunca los usé. Parecen, sin embargo, más relacionados con la estimación de la correlación entre los marcadores, no con el valor p como quisiste decir. Pero, de hecho, es mejor pensar en términos de la estructura de dependencia entre estadísticas de pruebas sucesivas, que entre valores p correlacionados.

Referencias

  1. Cantor, RM, Lange, K y Sinsheimer, JS. Priorizar los resultados de GWAS: una revisión de métodos estadísticos y recomendaciones para su aplicación . Soy J Hum Genet. 2010 86 (1): 6–22.
  2. Corley, RP, Zeiger, JS, Crowley, T et al. Asociación de genes candidatos con drogodependencia antisocial en adolescentes . Dependencia de drogas y alcohol 2008 96: 90–98.
  3. Dalmasso, C, Génin, E y Trégouet DA. Un procedimiento de Holm ponderado que explica las frecuencias alélicas en los estudios de asociación de Genomewide . Genética 2008 180 (1): 697–702.
  4. Wei, Z, Sun, W, Wang, K y Hakonarson, H. Pruebas múltiples en estudios de asociación de genoma completo a través de modelos ocultos de Markov . Bioinformática 2009 25 (21): 2802-2808.
  5. Broberg, P. Una revisión comparativa de las estimaciones de la proporción de genes sin cambios y la tasa de descubrimiento falso . BMC Bioinformatics 2005 6: 199.
  6. Need, AC, Ge, D, Weale, ME, y otros. Una investigación de todo el genoma de SNP y CNV en la esquizofrenia . PLoS Genet. 2009 5 (2): e1000373.
  7. Han, B, Kang, HM y Eskin, E. Corrección de pruebas múltiples rápida y precisa y estimación de potencia para millones de marcadores correlacionados . PLoS Genetics 2009
  8. Liang, Y y Kelemen, A. Avances y desafíos estadísticos para analizar datos de snp correlacionados de alta dimensión en el estudio genómico de enfermedades complejas . Encuestas estadísticas 2008 2: 43–60. - la mejor reseña reciente
  9. Nyholt, DR. Una corrección simple para múltiples pruebas de polimorfismos de un solo nucleótido en el desequilibrio de ligamiento entre sí . Soy J Hum Genet. 2004 74 (4): 765–769.
  10. Nicodemus, KK, Liu, W, Chase, GA, Tsai, YY y Fallin, MD. Comparación del error tipo I para múltiples correcciones de prueba en grandes estudios de polimorfismo de un solo nucleótido utilizando componentes principales versus algoritmos de bloqueo de haplotipos . BMC Genetics 2005; 6 (Supl. 1): S78.
  11. Peng, Q, Zhao, J y Xue, F. Pruebas de intervalo de confianza de arranque basadas en PCA para la asociación de enfermedades genéticas que involucran múltiples SNP . BMC Genetics 2010, 11: 6
  12. Li, M, Romero, R, Fu, WJ y Cui, Y (2010). Mapeo de interacciones haplotipo-haplotipo con LASSO adaptativo . BMC Genetics 2010, 11:79 - aunque no está directamente relacionado con la pregunta, cubre el análisis basado en haplotipos / efecto epistático
chl
fuente
1
¡Guau, gracias por pasar por todos estos problemas! Entiendo tus reparos sobre bootstrapping, y estoy casi convencido. Creo que mi principal complicación es la covariable numérica que tengo, que sin duda será necesaria (ya sea por sí misma o en interacción con el genotipo), y eso parece descartar mt.maxT y plink, aunque es posible que deba buscar en plink nuevamente. ¡Pero sin duda profundizaré en las referencias que proporcionó!
S. Kolassa - Restablece a Monica el
Siempre puede trabajar con los residuos de su GLM para obtener el valor de sus covariables, aunque perdió algo de Df que será difícil de contabilizar o reintroducir después (por ejemplo, para calcular el valor p).
chl
Hm, ¿residuos de mi regresión logística? ¿Sería eso legítimo?
S. Kolassa - Restablece a Monica el
¿Si por qué no? No es raro eliminar la varianza explicada por otras covariables y luego pasar al análisis de segundo nivel con sus datos residuales. A menudo es más rápido (por ejemplo, el plink es bastante lento con las covariables categóricas, mientras que está bien con las continuas; snpMatrixo simplemente glm()funciona bastante mejor en este punto, pero no puede incorporar muchos SNP dentro de glm()...); el problema es que obtener el valor p corregido al final de su segundo análisis es bastante complicado (porque debe tener en cuenta los parámetros ya estimados).
chl
Para ver una ilustración de cómo las personas trabajan con los residuos, ver por ejemplo p. 466 de Heck et al. La investigación de 17 genes candidatos para los rasgos de personalidad confirma los efectos del gen HTR2A en la búsqueda de novedades. Genes, cerebro y comportamiento (2009) vol. 8 (4) págs. 464-72
chl
2

El uso de un método como bonferroni está bien, el problema es que si tiene muchas pruebas, es probable que no encuentre muchos "descubrimientos".

Puede utilizar el enfoque FDR para las pruebas dependientes (consulte aquí para más detalles ). El problema es que no estoy seguro si puede decir por adelantado si sus correlaciones son todas positivas.

En R puede hacer FDR simple con p.adjust. Para cosas más complejas, echaría un vistazo a multcomp , pero no lo revisé para buscar soluciones en casos de dependencias.

Buena suerte.

Tal Galili
fuente
1
Hola tal gracias Bonferroni no me parece apropiado: si uno de mis SNP es causal y otros están asociados con él, debería haber una señal, y Bonferroni siempre me ha parecido demasiado conservador (generalmente prefiero la corrección gradual de Holm). El FDR con el que se vincula y no se ajustan a la evidencia combinada (y el FDR aún requiere que comprenda la correlación de mis pruebas, la pregunta original). multcomp puede ayudar, aunque a primera vista parece que trata más con múltiples pruebas dentro de un solo modelo, mientras que tengo varios modelos.
Profundizaré
Hola stephan Entiendo, lo siento por no ayudar más. ¡Buena suerte! Tal
Tal Galili
Hola Stephan, sigo pensando que todavía puedes usar el método = BY (para el procedimiento Benjamini Hochberg Yekuteli) en p.ajustar en R, como lo señala Tal. Definitivamente, usar Bonferroni puede ser conservador.
suncoolsu
suncoolsu, creo que este método solo funciona cuando la correlación es positiva (no negativa) entre las variables. Aclamaciones.
Tal Galili
2

Creo que los modelos normales multivariados se están utilizando para modelar los valores p correlacionados y para obtener el tipo correcto de múltiples correcciones de prueba. Corrección de pruebas múltiples rápida y precisa y estimación de potencia para millones de marcadores correlacionados. PLoS Genet 2009 habla de ellos y también ofrece otras referencias. Suena similar a lo que estaba hablando, pero creo que, además de obtener una corrección global del valor p más precisa, el conocimiento de la estructura LD también debería usarse para eliminar los falsos positivos que surgen de los marcadores correlacionados con los marcadores causales.

highBandWidth
fuente
2

Estoy buscando una solución de trabajo para exactamente el mismo problema. Lo mejor que encontré es el Bootstrap sin restricciones nulo introducido por Foulkes Andrea en su libro Genética estadística aplicada con R (2009) . Al contrario de todos los demás artículos y libros, considera específicamente las regresiones. Además de otros métodos, aconseja el Bootstrap no restringido nulo, que es adecuado cuando no se pueden calcular fácilmente los residuos (como en mi caso, donde modelo muchas regresiones independientes (básicamente correlaciones simples), cada una con la misma variable de respuesta y un recorte diferente). Encontré que este método también se llama método maxT .

> attach(fms)
> Actn3Bin <- > data.frame(actn3_r577x!="TT",actn3_rs540874!="AA",actn3_rs1815739!="TT",actn3_1671064!="GG")
> Mod <- summary(lm(NDRM.CH~.,data=Actn3Bin))
> CoefObs <- as.vector(Mod$coefficients[-1,1]) 
> B <-1000
> TestStatBoot <- matrix(nrow=B,ncol=NSnps)
> for (i in 1:B){
+    SampID <- sample(1:Nobs,size=Nobs, replace=T)
+    Ynew <- NDRM.CH[!MissDat][SampID]
+    Xnew <- Actn3BinC[SampID,]
+    CoefBoot <- summary(lm(Ynew~.,data=Xnew))$coefficients[-1,1]
+    SEBoot <- summary(lm(Ynew~.,data=Xnew))$coefficients[-1,2]
+    if (length(CoefBoot)==length(CoefObs)){
+       TestStatBoot[i,] <- (CoefBoot-CoefObs)/SEBoot
+    }
+ }

TestStatBootT^Tcrit.α=0,05T^Tcrit.

yoTyo^>Tcrit.

El último paso se puede lograr con este código

p.value<-0.05 # The target alpha threshold
digits<-1000000
library(gtools) # for binsearch

pValueFun<-function(cj)
{
   mean(apply(abs(TestStatBoot)>cj/digits,1,sum)>=1,na.rm=T)
}
ans<-binsearch(pValueFun,c(0.5*digits,100*digits),target=p.value)
p.level<-(1-pnorm(q=ans$where[[1]]/digits))*2 #two-sided.
Adam Ryczkowski
fuente