Comment appliquer la correction de tests multiples pour le chevauchement de la liste de gènes à l'aide de R

J'ai 2 études examinant la réponse du patient au même médicament. L'étude 1 a trouvé 10 000 gènes exprimés au-dessus du fond et 500 d'entre eux sont différentiellement exprimés et appelés signature de réponse au médicament. L'étude 2 a trouvé 1 000 gènes représentant la signature de la réponse au médicament. Le chevauchement entre les deux signatures est de 100 gènes.

Je veux calculer la signification statistique du chevauchement entre les signatures. Si je comprends bien, une façon de le faire (sur la base des articles ici: Calculer la probabilité de chevauchement de la liste de gènes entre un seq d'ARN et un ensemble de données de puce ChLP et ici: Utiliser le phyper de R pour obtenir la probabilité de chevauchement de liste ) est via phyper():

> overlap  <- 100
> list1    <- 500
> totalPop <- 10000
> list2    <- 1000
> 
> 1-phyper(overlap-1, list1, totalPop-list1, list2)
[1] 4.103051e-12

1. Cela vous semble-t-il raisonnable?

2. Si je voulais appliquer la correction de Bonferroni, je devrais multiplier cette valeur de p par le nombre de comparaisons. À quoi correspondrait le nombre de comparaisons dans ce cas? List2? Sinon, quel serait un moyen rapide de faire une correction moins conservatrice (par exemple, Benjamini-Hochberg)?

Je ne connais rien aux études sur l'expression des gènes, mais j'ai un certain intérêt pour les inférences multiples, je vais donc risquer une réponse sur cette partie de la question de toute façon.

Personnellement, je n'aborderais pas le problème de cette façon. Je voudrais ajuster le niveau d'erreur dans les études originales, calculer le nouveau chevauchement et laisser le test à la fin seul. Si le nombre de gènes exprimés différentiellement (et tout autre résultat que vous utilisez) est déjà basé sur des tests ajustés, je dirais que vous n'avez rien à faire.

Si vous ne pouvez pas revenir aux données d'origine et que vous voulez vraiment ajuster la valeur de p , vous pouvez en effet la multiplier par le nombre de tests mais je ne vois pas pourquoi cela devrait avoir quelque chose à voir avec la taille de list2. Il serait plus logique d'ajuster le nombre total de tests effectués dans les deux études (c'est-à-dire deux fois la population). Cela va être brutal, cependant.

Pour ajuster les valeurs de p dans R, vous pouvez utiliser p.adjust(p), où pest un vecteur de valeurs de p .

p.adjust(p, method="bonferroni") # Bonferroni method, simple multiplication
p.adjust(p, method="holm") # Holm-Bonferroni method, more powerful than Bonferroni
p.adjust(p, method="BH") # Benjamini-Hochberg

Comme indiqué dans le fichier d'aide, il n'y a aucune raison de ne pas utiliser Holm-Bonferroni sur Bonferroni car il fournit également un contrôle fort du taux d'erreur par famille dans tous les cas mais est plus puissant. Benjamini-Hochberg contrôle le taux de fausses découvertes, qui est un critère moins strict.

Modifié après le commentaire ci-dessous:

Plus je pense au problème, plus je pense qu'une correction pour les comparaisons multiples est inutile et inappropriée dans cette situation. C'est là que la notion de «famille» d'hypothèses entre en jeu. Votre dernier test n'est pas tout à fait comparable à tous les tests précédents, il n'y a aucun risque de «capitaliser sur le hasard» ou de sélectionner des résultats significatifs, il n'y en a qu'un test d'intérêt et il est légitime d'utiliser le niveau d'erreur ordinaire pour celui-ci.

Même si vous corrigez agressivement les nombreux tests effectués auparavant, vous ne répondriez pas directement à la principale préoccupation, à savoir le fait que certains des gènes des deux listes auraient pu être détectés à tort comme exprimés différentiellement. Les résultats des tests antérieurs restent «valables» et si vous souhaitez interpréter ces résultats tout en contrôlant le taux d'erreur familial, vous devez également tous les corriger.

Mais si l'hypothèse nulle est vraiment vraie pour tous les gènes, tout résultat significatif serait un faux positif et vous ne vous attendriez pas à ce que le même gène soit à nouveau signalé dans l'échantillon suivant. Le chevauchement entre les deux listes ne se produirait donc que par hasard et c'est exactement ce que le test basé sur la distribution hypergéométrique teste. Donc, même si les listes de gènes sont inutiles, le résultat de ce dernier test est sûr. Intuitivement, il semble que tout ce qui se trouve entre les deux (un mélange d'hypothèses vraies et fausses) devrait également convenir.

Peut-être que quelqu'un ayant plus d'expérience dans ce domaine pourrait peser, mais je pense qu'un ajustement ne deviendrait nécessaire que si vous voulez comparer le nombre total de gènes détectés ou savoir lesquels sont exprimés de manière différentielle, c'est-à-dire si vous voulez interpréter les milliers d'individus tests effectués dans chaque étude.

Vous n'avez pas besoin de corriger la valeur p pour votre seul test de chevauchement. Cependant, supposons que vous vouliez déterminer si le médicament affecte les gènes de la même voie. Comment détermineriez-vous quelle voie avait le plus de chevauchement? Disons que vous avez 500 ensembles de gènes de voie. Vous exécutez le test de chevauchement des ensembles hypergéométriques 500 fois et les avez classés par valeur p. Puisque vous avez effectué ce test 500 fois (ou même plus selon la quantité de données dont vous disposez), il est possible que vous obteniez un bon score juste par hasard (faux positif). Alors, vous devez corriger cela et effectuer un ajustement de la valeur ... soit Bonferroni (le plus conservateur) ou Benjamini Hochberg.