Korekta wartości p dla wielu testów, w których testy są skorelowane (genetyka)

Mam wartości p z wielu testów i chciałbym wiedzieć, czy rzeczywiście jest coś znaczącego po poprawieniu wielu testów. Komplikacja: moje testy nie są niezależne. Metoda, o której myślę (wariant metody produktu Fishera, Zaykin i wsp., Genet Epidemiol , 2002) wymaga korelacji między wartościami p.

Aby oszacować tę korelację, obecnie myślę o przypadkach ładowania początkowego, przeprowadzaniu analiz i korelowaniu powstałych wektorów wartości p. Czy ktoś ma lepszy pomysł? A może nawet lepszy pomysł na mój pierwotny problem (poprawianie wielokrotnych testów w testach skorelowanych)?

Tło: loguję się logicznie, czy moi pacjenci cierpią na określoną chorobę w związku z interakcją między ich genotypem (AA, Aa lub aa) a współzmienną. Jednak genotyp jest w rzeczywistości dużą ilością (30–250) polimorfizmów pojedynczego nukleotydu (SNP), które z pewnością nie są niezależne, ale znajdują się w zaburzeniu równowagi sprzężenia.

To właściwie gorący temat w badaniach analizy genomewide (GWAS)! Nie jestem pewien, czy metoda, o której myślisz, jest najbardziej odpowiednia w tym kontekście. Łączenie wartości p zostało opisane przez niektórych autorów, ale w innym kontekście (badania replikacji lub metaanaliza, patrz np. (1) w ostatnim przeglądzie). Łączenie wartości p SNP metodą Fishera jest ogólnie stosowane, gdy chce się uzyskać unikalną wartość p dla danego genu; pozwala to na pracę na poziomie genów i zmniejsza wymiar wymiarowości kolejnych testów, ale jak już powiedziałeś, brak niezależności między markerami (wynikającymi z kolokacji przestrzennej lub nierównowagi połączeń, LD) wprowadza błąd systematyczny. Bardziej zaawansowane alternatywy polegają na procedurach ponownego próbkowania,

Moje główne obawy związane z bootstrapowaniem (z wymianą) polegałyby na tym, że wprowadzasz sztuczną formę powiązania, lub innymi słowy tworzysz wirtualne bliźnięta, zmieniając w ten sposób równowagę Hardy'ego-Weinberga (ale także minimalną częstotliwość alleli i szybkość połączeń). Nie byłoby tak w przypadku podejścia permutacyjnego, w którym permutujesz poszczególne etykiety i zachowujesz dane genotypowania w obecnej postaci. Zazwyczaj oprogramowanie Plink może dać surowe i permutowane wartości p, chociaż wykorzystuje (domyślnie) strategię testowania adaptacyjnego z przesuwanym oknem, które pozwala zatrzymać wszystkie permutacje (powiedzmy 1000 na SNP), jeśli wydaje się, że SNP pod rozważanie nie jest „interesujące”; ma również opcję obliczania maxT, patrz pomoc online .

Biorąc jednak pod uwagę małą liczbę SNP, które rozważasz, sugerowałbym poleganie na testach opartych na FDR lub maxT, jak zaimplementowano w pakiecie wielokrotnego testu R (patrz mt.maxT), ale ostatecznym przewodnikiem po strategiach ponownego próbkowania dla aplikacji genomowej jest wiele procedur testowych z aplikacjami do Genomics , od Dudoit & van der Laan (Springer, 2008). Zobacz także książkę Andrei Foulkes na temat genetyki z R , która jest recenzowana w JSS. Ma świetny materiał na temat wielu procedur testowych.

Dalsze uwagi

Wielu autorów wskazało na fakt, że proste metody wielokrotnego korygowania testów, takie jak Bonferroni lub Sidak, są zbyt rygorystyczne, aby dostosować wyniki dla poszczególnych SNP. Ponadto żadna z tych metod nie bierze pod uwagę korelacji między SNP z powodu LD, który oznacza zmienność genetyczną między regionami genowymi. Inne alternatywy zostały zaproponowane, na przykład pochodna metody Holma do wielokrotnego porównania (3), ukryty model Markowa (4), warunkowy lub dodatni FDR (5) lub jego pochodna (6). Tak zwane statystyki luk lub przesuwane okna okazały się w niektórych przypadkach skuteczne, ale dobrą recenzję znajdziesz w (7) i (8).

Słyszałem również o metodach, które skutecznie wykorzystują strukturę haplotypów lub LD, np. (9), ale nigdy ich nie używałem. Wydaje się jednak, że są one bardziej związane z oszacowaniem korelacji między markerami, a nie wartością p, jak zamierzałeś. Ale w rzeczywistości lepiej jest pomyśleć o strukturze zależności między kolejnymi statystykami testów, niż między skorelowanymi wartościami p.

Referencje

1. Cantor, RM, Lange, K i Sinsheimer, JS. Priorytetyzacja wyników GWAS: przegląd metod statystycznych i zaleceń dotyczących ich stosowania . Am J Hum Genet. 2010 86 (1): 6–22.
2. Corley, RP, Zeiger, JS, Crowley, T i in. Związek genów kandydujących z uzależnieniem od narkotyków aspołecznych u młodzieży . Uzależnienie od narkotyków i alkoholu 2008 96: 90–98.
3. Dalmasso, C, Génin, E i Trégouet DA. Procedura ważonego holma uwzględniająca częstotliwości alleli w badaniach stowarzyszenia genomewide . Genetics 2008 180 (1): 697–702.
4. Wei, Z, Sun, W, Wang, K i Hakonarson, H. Wielokrotne testy w badaniach asocjacyjnych całego genomu za pomocą ukrytych modeli Markowa . Bioinformatics 2009 25 (21): 2802–2808.
5. Broberg, P. Porównawczy przegląd szacunków odsetka niezmienionych genów i odsetka fałszywych odkryć . BMC Bioinformatics 2005 6: 199.
6. Need, AC, Ge, D, Weale, ME, i in. Badanie całego genomu SNP i CNV w schizofrenii . PLoS Genet. 2009 5 (2): e1000373.
7. Han, B, Kang, HM i Eskin, E. Szybka i dokładna wielokrotna korekcja testowa i oszacowanie mocy dla milionów skorelowanych markerów . PLoS Genetics 2009
8. Liang, Y i Kelemen, A. Postępy i wyzwania statystyczne w analizie skorelowanych danych snp w badaniach genomowych złożonych chorób . Ankiety statystyczne 2008 2: 43–60. - najnowsza najnowsza recenzja
9. Nyholt, DR. Prosta poprawka do wielokrotnego testowania polimorfizmów pojedynczego nukleotydu w nierównowagach sprzężonych ze sobą . Am J Hum Genet. 2004 74 (4): 765–769.
10. Nikodem, KK, Liu, W, Chase, GA, Tsai, YY i Fallin, MD. Porównanie błędu typu I dla wielu poprawek testowych w dużych badaniach polimorfizmu pojedynczego nukleotydu przy użyciu głównych składników w porównaniu z algorytmami blokującymi haplotyp . BMC Genetics 2005; 6 (suplement 1): S78.
11. Peng, Q, Zhao, J i Xue, F. Testy przedziału ufności oparte na PCA bootstrap dla powiązania choroby genowej z udziałem wielu SNP . BMC Genetics 2010, 11: 6
12. Li, M, Romero, R, Fu, WJ i Cui, Y (2010). Mapowanie Haplotyp-haplotyp Interakcje z adaptacyjnym LASSO . BMC Genetics 2010, 11:79 - chociaż nie jest bezpośrednio związany z pytaniem, obejmuje analizę opartą na haplotypie / efekt epistatyczny

Używanie metody takiej jak bonferroni jest w porządku, problem polega na tym, że jeśli masz wiele testów, prawdopodobnie nie znajdziesz wielu „odkryć”.

Możesz zastosować podejście FDR do testów zależnych (szczegóły tutaj ), problem polega na tym, że nie jestem pewien, czy możesz powiedzieć z góry, czy wszystkie korelacje są pozytywne.

W R możesz zrobić prosty FDR za pomocą p.adjust. W przypadku bardziej skomplikowanych rzeczy rzuciłbym okiem na multcomp , ale nie przejrzałem go, aby znaleźć rozwiązania w przypadkach zależności.

Powodzenia.

Myślę, że wielowymiarowe modele normalne są używane do modelowania skorelowanych wartości p i uzyskania odpowiedniego typu wielu poprawek testowych. Szybka i dokładna wielokrotna korekcja testowa i oszacowanie mocy dla milionów skorelowanych markerów. PLoS Genet 2009 mówi o nich, a także podaje inne referencje. Brzmi podobnie do tego, o czym mówiłeś, ale myślę, że oprócz uzyskania dokładniejszej globalnej korekty wartości p, znajomość struktury LD powinna również być użyta do usuwania fałszywych wyników dodatnich wynikających ze znaczników skorelowanych ze znacznikami przyczynowymi.

Szukam rozwiązania dla dokładnie tego samego problemu. Najlepsze, co znalazłem, to Null Unrestricted Bootstrap wprowadzony przez Foulkesa Andreę w jego książce Applied Statistics Genetics with R (2009) . W przeciwieństwie do wszystkich innych artykułów i książek, rozważa w szczególności regresje. Oprócz innych metod radzi zerowy nieograniczony bootstrap, który jest odpowiedni tam, gdzie nie można łatwo obliczyć reszt (jak w moim przypadku, gdzie modeluję wiele niezależnych regresji (w zasadzie prostych korelacji), każda z tą samą zmienną odpowiedzi i innym wycinaniem). Odkryłem, że ta metoda nazywa się również metodą 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
+    }
+ }

TestStatBoot$\hat{\vec{T}^*}$$T_{\text{crit.}}$$\alpha=0.05$$\hat{\vec{T}^*}$$T_{\text{crit.}}$

$i$$\hat{\vec{T}_i} > T_{\text{crit.}}$

Ostatni krok można wykonać za pomocą tego kodu

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.