Correction des valeurs de p pour plusieurs tests où les tests sont corrélés (génétique)

24

J'ai des valeurs de p provenant de nombreux tests et j'aimerais savoir s'il y a réellement quelque chose d'important après correction pour plusieurs tests. La complication: mes tests ne sont pas indépendants. La méthode à laquelle je pense (une variante de la méthode des produits de Fisher, Zaykin et al., Genet Epidemiol , 2002) a besoin de la corrélation entre les valeurs de p.

Afin d'estimer cette corrélation, je pense actuellement à des cas d'amorçage, à l'exécution des analyses et à la corrélation des vecteurs résultants des valeurs de p. Quelqu'un a-t-il une meilleure idée? Ou même une meilleure idée de mon problème d'origine (correction pour plusieurs tests dans des tests corrélés)?

Contexte: Je régresse logistiquement, que mes sujets souffrent ou non d'une maladie particulière sur l'interaction entre leur génotype (AA, Aa ou aa) et une covariable. Cependant, le génotype est en fait beaucoup (30-250) de polymorphismes mononucléotidiques (SNP), qui ne sont certainement pas indépendants mais dans le déséquilibre de liaison.

S. Kolassa - Réintégrer Monica
la source

Réponses:

29

C'est en fait un sujet brûlant dans les études d'analyse à l'échelle du génome (GWAS)! Je ne suis pas sûr que la méthode à laquelle vous pensez soit la plus appropriée dans ce contexte. La mise en commun des valeurs de p a été décrite par certains auteurs, mais dans un contexte différent (études de réplication ou méta-analyse, voir par exemple (1) pour une revue récente). La combinaison de valeurs p SNP par la méthode de Fisher est généralement utilisée lorsque l'on veut dériver une valeur p unique pour un gène donné; cela permet de travailler au niveau du gène et de réduire la dimensionnalité des tests ultérieurs, mais comme vous l'avez dit, la non-indépendance entre les marqueurs (résultant de la colocation spatiale ou du déséquilibre de la liaison, LD) introduit un biais. Des alternatives plus puissantes reposent sur des procédures de rééchantillonnage,

Ma principale préoccupation avec le bootstraping (avec remplacement) serait que vous introduisiez une forme artificielle de parenté, ou en d'autres termes que vous créiez des jumeaux virtuels, altérant ainsi l'équilibre Hardy-Weinberg (mais aussi la fréquence minimale des allèles et le taux d'appel). Ce ne serait pas le cas avec une approche par permutation où vous permutez des étiquettes individuelles et conservez les données de génotypage telles quelles. Habituellement, le logiciel Plink peut vous donner des valeurs p brutes et permutées, bien qu'il utilise (par défaut) une stratégie de test adaptative avec une fenêtre coulissante qui permet d'arrêter l'exécution de toutes les permutations (disons 1000 par SNP) s'il apparaît que le SNP sous la considération n'est pas "intéressante"; il a également une option pour calculer maxT, voir l' aide en ligne .

Mais étant donné le faible nombre de SNP que vous envisagez, je suggérerais de s'appuyer sur des tests basés sur FDR ou maxT tels qu'implémentés dans le package multtest R (voir mt.maxT), mais le guide définitif des stratégies de rééchantillonnage pour l'application génomique est Procédures de test multiples avec applications à Génomique , de Dudoit et van der Laan (Springer, 2008). Voir aussi le livre d'Andrea Foulkes sur la génétique avec R , qui est examiné dans le JSS. Elle a un excellent matériel sur plusieurs procédures de test.

Notes complémentaires

De nombreux auteurs ont souligné le fait que de simples méthodes de correction de tests multiples telles que Bonferroni ou Sidak sont trop strictes pour ajuster les résultats pour les SNP individuels. De plus, aucune de ces méthodes ne prend en compte la corrélation qui existe entre les SNP en raison de LD qui marque la variation génétique entre les régions géniques. D'autres alternatives ont été proposées, comme un dérivé de la méthode de Holm pour la comparaison multiple (3), le modèle de Markov caché (4), le FDR conditionnel ou positif (5) ou un dérivé de celui-ci (6), pour n'en nommer que quelques-uns. Les soi-disant statistiques d'écart ou fenêtre coulissante se sont avérées efficaces dans certains cas, mais vous trouverez une bonne revue dans (7) et (8).

J'ai également entendu parler de méthodes qui utilisent efficacement la structure des haplotypes ou LD, par exemple (9), mais je ne les ai jamais utilisées. Ils semblent cependant plus liés à l'estimation de la corrélation entre les marqueurs, et non à la valeur p comme vous le vouliez. Mais en fait, vous pourriez mieux penser en termes de structure de dépendance entre les statistiques de test successives, qu'entre les valeurs p corrélées.

Les références

  1. Cantor, RM, Lange, K et Sinsheimer, JS. Hiérarchisation des résultats du GWAS: examen des méthodes statistiques et recommandations pour leur application . Suis J Hum Genet. 2010 86 (1): 6–22.
  2. Corley, RP, Zeiger, JS, Crowley, T et al. Association de gènes candidats à la pharmacodépendance antisociale chez les adolescents . Dépendance aux drogues et à l'alcool 2008 96: 90–98.
  3. Dalmasso, C, Génin, E et Trégouet DA. Une procédure pondérée Holm tenant compte des fréquences alléliques dans les études d'association génomique . Génétique 2008 180 (1): 697–702.
  4. Wei, Z, Sun, W, Wang, K et Hakonarson, H. Tests multiples dans les études d'association à l'échelle du génome via des modèles de Markov cachés . Bioinformatics 2009 25 (21): 2802-2808.
  5. Broberg, P. Un examen comparatif des estimations de la proportion de gènes inchangés et du taux de fausses découvertes . BMC Bioinformatics 2005 6: 199.
  6. Besoin, AC, Ge, D, Weale, ME, et a. Une enquête à l'échelle du génome des SNP et CNV dans la schizophrénie . PLoS Genet. 2009 5 (2): e1000373.
  7. Han, B, Kang, HM et Eskin, E. Correction de plusieurs tests rapides et précis et estimation de la puissance pour des millions de marqueurs corrélés . PLoS Genetics 2009
  8. Liang, Y et Kelemen, A. Progrès statistiques et défis pour l'analyse des données corrélées de haute dimension snp dans l'étude génomique des maladies complexes . Enquêtes statistiques 2008 2: 43–60. - le meilleur avis récent de tous les temps
  9. Nyholt, DR. Une correction simple pour les tests multiples des polymorphismes mononucléotidiques dans le déséquilibre de liaison les uns avec les autres . Suis J Hum Genet. 2004 74 (4): 765–769.
  10. Nicodemus, KK, Liu, W, Chase, GA, Tsai, YY et Fallin, MD. Comparaison de l'erreur de type I pour plusieurs corrections de test dans les grandes études de polymorphisme mononucléotidique utilisant les principaux composants par rapport aux algorithmes de blocage des haplotypes . BMC Genetics 2005; 6 (Suppl 1): S78.
  11. Peng, Q, Zhao, J et Xue, F. Tests d'intervalle de confiance bootstrap basés sur PCA pour l'association gène-maladie impliquant plusieurs SNP . BMC Genetics 2010, 11: 6
  12. Li, M, Romero, R, Fu, WJ et Cui, Y (2010). Cartographie des interactions haplotype-haplotype avec LASSO adaptatif . BMC Genetics 2010, 11:79 - bien qu'il ne soit pas directement lié à la question, il couvre l'analyse basée sur l'haplotype / l'effet épistatique
chl
la source
1
Wow, merci pour tous ces ennuis! Je comprends vos scrupules à propos du bootstrap et je suis presque convaincu. Je pense que ma principale complication est la covariable numérique que j'ai, qui sera certainement nécessaire (soit en elle-même, soit en interaction avec le génotype), et qui semble exclure mt.maxT et plink, bien que je doive peut-être examiner à nouveau plink. Mais je vais certainement fouiller dans les références que vous avez fournies!
S.Kolassa - Reinstate Monica
Vous pouvez toujours travailler avec les résidus de votre GLM pour contourner vos covariables, bien que vous ayez perdu du Df qui sera difficile à prendre en compte ou à réintroduire par la suite (par exemple pour calculer la valeur de p).
chl
Hm, résidus de ma régression logistique? Serait-ce légitime?
S.Kolassa - Reinstate Monica
Oui pourquoi pas? Il n'est pas rare de supprimer la variance prise en compte par d'autres covariables puis de passer à l'analyse de 2e niveau avec vos données résiduelles. C'est souvent plus rapide (par exemple, le plink est assez lent avec des covariables catégorielles, alors qu'il est correct avec des variables continues; snpMatrixou il glm()fonctionne tout simplement bien mieux sur ce point mais vous ne pouvez pas intégrer beaucoup de SNP dans glm()...); le problème est que l'obtention de la valeur p corrigée à la fin de votre 2ème analyse est assez délicate (car il faut tenir compte des paramètres déjà estimés).
chl
Pour une illustration de la façon dont les gens travaillent avec les résidus, voir par exemple p. 466 de Heck et al. L'étude de 17 gènes candidats pour les traits de personnalité confirme les effets du gène HTR2A sur la recherche de nouveautés. Gènes, cerveau et comportement (2009) vol. 8 (4) pp. 464-72
chl
2

L'utilisation d'une méthode comme bonferroni est très bien, le problème est que si vous avez de nombreux tests, vous ne trouverez probablement pas beaucoup de "découvertes".

Vous pouvez aller avec l'approche FDR pour les tests dépendants (voir ici pour plus de détails ) le problème est que je ne sais pas si vous pouvez dire à l'avance si vos corrélations sont toutes positives.

Dans R, vous pouvez faire un FDR simple avec p.adjust. Pour des choses plus complexes, je jetterais un coup d'œil à multcomp , mais je ne l'ai pas parcouru pour voir des solutions en cas de dépendances.

Bonne chance.

Tal Galili
la source
1
Salut Tal, merci! Bonferroni ne me semble pas approprié - si l'un de mes SNP est causal et que d'autres y sont associés, il devrait y avoir un signal, et Bonferroni m'a toujours semblé trop conservateur (je préfère généralement la correction pas à pas de Holm). Le FDR auquel vous vous connectez et p.adjust ne prennent pas en compte les preuves combinées (et le FDR m'oblige toujours à comprendre la corrélation de mes tests, la question d'origine). multcomp peut aider, bien qu'à première vue, il semble qu'il traite plus de plusieurs tests au sein d'un même modèle, alors que j'ai plusieurs modèles. Je vais creuser plus profondément ...
S. Kolassa - Reinstate Monica
Bonjour Stephan. Je comprends, désolé de ne pas avoir aidé davantage. Bonne chance! Tal
Tal Galili
Bonjour Stephan, je pense toujours que vous pouvez toujours utiliser la méthode = BY (pour la procédure Benjamini Hochberg Yekuteli) dans p.adjust in R, comme l'a souligné Tal. Certainement, l'utilisation de Bonferroni peut être conservatrice.
suncoolsu
suncoolsu, je pense que cette méthode ne fonctionne que lorsque la corrélation est positive (et non négative) entre les variables. À votre santé.
Tal Galili
2

Je pense que des modèles normaux multivariés sont utilisés pour modéliser les valeurs de p corrélées et pour obtenir le bon type de corrections de tests multiples. Correction rapide et précise de plusieurs tests et estimation de la puissance pour des millions de marqueurs corrélés. PLoS Genet 2009 en parle et donne également d'autres références. Cela ressemble à ce dont vous parliez, mais je pense qu'en plus d'obtenir une correction globale de la valeur de p plus précise, la connaissance de la structure LD devrait également être utilisée pour supprimer les faux positifs provenant de marqueurs corrélés avec des marqueurs causaux.

bande passante élevée
la source
2

Je cherche une solution de travail pour exactement le même problème. Le meilleur que j'ai trouvé est le Bootstrap Null Unrestricted présenté par Foulkes Andrea dans son livre Applied Statistical Genetics with R (2009) . Contrairement à tous les autres articles et livres, il considère spécifiquement les régressions. Outre d'autres méthodes, il conseille le Null Unrestricted Bootstrap, qui convient lorsque l'on ne peut pas facilement calculer des résidus (comme dans mon cas, où je modélise de nombreuses régressions indépendantes (essentiellement des corrélations simples), chacune avec la même variable de réponse et un snip différent). J'ai trouvé que cette méthode s'appelait également la méthode 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
+    }
+ }

Une fois que nous avons toute la TestStatBootmatrice (dans les lignes, nous avons des réplications bootstrap, et dans les colonnes, nous avons bootstrapT^ statistiques) que nous trouvons pour lesquels Tcrit. on observe exactement α=0,05 pour cent de plus significatif T^ statistiques (plus importantes signifient qu'avec une valeur absolue plus grande que Tcrit.).

Nous rapportons je-th composant du modèle significatif, si son Tje^>Tcrit.

La dernière étape peut être accomplie avec ce code

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
la source