Ze Stanem i pakietów frontend rstanarmczy brmsmogę łatwo analizować dane droga Bayesa jak ja zanim z mieszanych modeli takich jak lme. Chociaż na biurku mam większość książek i artykułów Kruschke-Gelman-Wagenmakers itp., Nie mówią mi one, jak podsumować wyniki dla medycznej publiczności, rozdartej między Skylla gniewu Bayesiana a Charybdą recenzentów medycznych ( „chcemy znaczeń, a nie rzeczy rozproszonych”).
Przykład: Częstotliwość żołądka (1 / min) mierzy się w trzech grupach; referencje stanowią zdrowe kontrole. Dla każdego uczestnika jest kilka pomiarów, więc często korzystałem z następującego modelu mieszanego lme:
summary(lme(freq_min~ group, random = ~1|study_id, data = mo))Lekko edytowane wyniki:
Fixed effects: freq_min ~ group
Value Std.Error DF t-value p-value
(Intercept) 2.712 0.0804 70 33.7 0.0000
groupno_symptoms 0.353 0.1180 27 3.0 0.0058
groupwith_symptoms 0.195 0.1174 27 1.7 0.1086
Dla uproszczenia użyję 2 * standardowego błędu jako 95% CI.
W kontekście częstokrzyskim podsumowałbym to następująco:
- W grupie kontrolnej szacowana częstotliwość wyniosła 2,7 / min (może dodać tutaj CI, ale czasami tego unikam z powodu zamieszania wywołanego przez absolut i różnicę CI).
- W grupie no_symptoms częstotliwość była wyższa o 0,4 / min, CI (0,11 do 0,59) / min, p = 0,006 niż kontrola.
- W grupie with_symptoms częstotliwość była wyższa o 0,2 / min, CI (-0,04 do 0,4) / min, p = 0,11 niż kontrola.
Chodzi o maksymalną dopuszczalną złożoność publikacji medycznej, recenzent prawdopodobnie poprosi mnie o dodanie „nieistotnego” w drugim przypadku.
Oto to samo z stan_lmerdomyślnymi priorytetami.
freq_stan = stan_lmer(freq_min~ group + (1|study_id), data = mo)
contrast lower_CredI frequency upper_CredI
(Intercept) 2.58322 2.714 2.846
groupno_symptoms 0.15579 0.346 0.535
groupwith_symptoms -0.00382 0.188 0.384
gdzie CredI to 90% wiarygodnych przedziałów (patrz winieta rstanarm, dlaczego 90% jest używane jako domyślne).
Pytania:
- Jak przetłumaczyć powyższe podsumowanie na świat bayesowski?
- W jakim stopniu wymagana jest wcześniejsza dyskusja? Jestem całkiem pewien, że artykuł powróci ze zwykłym „subiektywnym założeniem”, kiedy wspomnę o aurorze; lub przynajmniej „bez dyskusji technicznej, proszę”. Ale wszystkie władze bayesowskie żądają, aby interpretacja była ważna tylko w kontekście a priori.
- Jak mogę dostarczyć surogat „znaczenia” w sformułowaniu, nie zdradzając pojęć bayesowskich? Coś w rodzaju „wiarygodnie inny” (uuuh ...) lub prawie wiarygodnie inny (buoha ..., brzmi jak „na granicy znaczenia).
Jonah Gabry i Ben Goodrich (2016). rstanarm: Bayesian Applied Regression Modeling via Stan. Wersja pakietu R 2.9.0-3. https://CRAN.R-project.org/package=rstanarm
Stan Development Team (2015). Stan: Biblioteka C ++ dla prawdopodobieństwa i próbkowania, wersja 2.8.0. URL http://mc-stan.org/ .
Paul-Christian Buerkner (2016). brms: modele regresji bayesowskiej przy użyciu Stana. Wersja pakietu R 0.8.0. https://CRAN.R-project.org/package=brms
Pinheiro J, Bates D, DebRoy S, Sarkar D and R Core Team (2016). nlme: Liniowe i nieliniowe modele efektów mieszanych . Wersja pakietu R 3.1-124, http://CRAN.R-project.org/package=nlme>.
źródło
mean(as.matrix(freq_stan)[,"groupwith_symptoms"] < 0).group_nosymptomsa następnie powiedzieć, że prawdopodobieństwo, że będzie ujemne, jest1 / draws. Ale w przypadku przechwytywania łańcuch nigdy nie będzie wędrować na terytorium negatywne dla tych danych, więc myślę, że można powiedzieć, że prawdopodobieństwo jest mniejsze niż1 / draws.Odpowiedzi:
Szybkie przemyślenia:
1) Kluczową kwestią jest pytanie, na które pytanie próbujesz odpowiedzieć swoim odbiorcom, ponieważ to określa, jakie informacje chcesz uzyskać z analizy statystycznej. W tym przypadku wydaje mi się, że chcesz oszacować wielkość różnic między grupami (a może wielkość stosunków grup, jeśli jest to miara bardziej znana odbiorcom). Wielkość różnic nie wynika bezpośrednio z analiz przedstawionych w pytaniu. Ale łatwo jest uzyskać to, czego chcesz od analizy bayesowskiej: chcesz tylnego rozkładu różnic (lub stosunków). Następnie na podstawie tylnego rozkładu różnic (lub stosunków) można wykonać bezpośrednie stwierdzenie prawdopodobieństwa, takie jak to:
„Najbardziej 95% wiarygodne różnice mieszczą się między [dolnym limitem 95% HDI] a [wysokim limitem 95% HDI]” (tutaj używam przedziału 95% największej gęstości [HDI] jako przedziału wiarygodnego, a ponieważ są one definicja najwyższych wartości parametrów gęstości są one określane jako „najbardziej wiarygodne”)
Publiczność w czasopiśmie medycznym intuicyjnie i poprawnie zrozumiałaby to stwierdzenie, ponieważ to, co zwykle uważa, to znaczenie częstego przedziału ufności (nawet jeśli nie jest to znaczenie częstego przedziału ufności).
Jak uzyskać różnice (lub stosunki) od Stana lub JAGS? Tylko przez przetwarzanie końcowe skompletowanego łańcucha MCMC. Na każdym etapie łańcucha obliczyć odpowiednie różnice (lub stosunki), a następnie zbadać tylny rozkład różnic (lub stosunki). Przykłady podano w DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ dla MCMC ogólnie na ryc. 7.9 (s. 177), dla JAGS na rycinie 8.6 (s. 211) i dla Stan w sekcji 16.3 (s. . 468) itd.!
2) Jeśli tradycja zmusza cię do złożenia oświadczenia o odrzuceniu różnicy zero, masz dwie opcje bayesowskie.
2A) Jedną z opcji jest złożenie oświadczeń prawdopodobieństwa dotyczących przedziałów bliskich zeru i ich relacji do HDI. W tym celu ustawiłeś region praktycznej równoważności (LINA) wokół zera, który jest jedynie progiem decyzyjnym odpowiednim dla twojej zastosowanej domeny - jak duża różnica jest banalnie mała? Wyznaczanie takich granic jest rutynowo wykonywane na przykład w klinicznych testach nie-gorszości. Jeśli masz miarę „rozmiaru efektu” na swoim polu, mogą istnieć konwencje dla „małego” rozmiaru efektu, a limity LINY mogą być, powiedzmy, połową małego efektu. Następnie możesz wykonać bezpośrednie stwierdzenia prawdopodobieństwa, takie jak te:
„Tylko 1,2% tylnego rozkładu różnic jest praktycznie równoważne zeru”
i
„Wszystkie 95% najbardziej wiarygodnych różnic nie są praktycznie równe zeru (tzn. 95% HDI i LINIA nie pokrywają się) i dlatego odrzucamy zero”. (zauważ różnicę między stwierdzeniem prawdopodobieństwa z rozkładu a posteriori a późniejszą decyzją opartą na tym stwierdzeniu)
Możesz także zaakceptować różnicę zerową, dla celów praktycznych, jeśli wszystkie 95% najbardziej wiarygodnych wartości są praktycznie równoważne zeru.
2B) Drugą opcją bayesowską jest testowanie Bayesowskiej hipotezy zerowej. (Zauważ, że powyższa metoda nie byłanazywane „testowaniem hipotez”!) Testowanie Bayesowskiej zerowej hipotezy porównuje model Bayesa z wcześniejszym rozkładem, który zakłada, że różnica może wynosić tylko zero, z alternatywnym rozkładem wcześniejszym, który zakłada, że różnica może być pewnym rozproszonym zakresem możliwości. Wynik takiego porównania modelu (zwykle) zależy bardzo silnie od konkretnego wyboru alternatywnego rozkładu, dlatego należy dokonać starannego uzasadnienia przy wyborze alternatywy wcześniej. Najlepiej jest używać co najmniej średnio poinformowanych priorów zarówno dla wartości zerowej, jak i alternatywnej, aby porównanie modelu było naprawdę znaczące. Zauważ, że porównanie modelu dostarcza innych informacji niż oszacowanie różnic między grupami, ponieważ porównanie modelu dotyczy innego pytania. Zatem nawet przy porównaniu modeli
Mogą istnieć sposoby wykonania testu Bayesa na zerową hipotezę z wyjścia Stan / JAGS / MCMC, ale nie wiem w tym przypadku. Na przykład można wypróbować przybliżenie Savage'a-Dickeya do współczynnika Bayesa, ale polegałoby to na znajomości wcześniejszej gęstości różnic, co wymagałoby pewnej analizy matematycznej lub dodatkowego przybliżenia MCMC w stosunku do wcześniejszego.
Dwie metody decydowania o wartościach zerowych omówiono w rozdz. 12 DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ . Ale tak naprawdę nie chcę, aby ta dyskusja była śledzona przez debatę na temat „właściwego” sposobu oceny wartości zerowych; są po prostu różne i dostarczają różnych informacji. Głównym punktem mojej odpowiedzi jest punkt 1 powyżej: Spójrz na tylny rozkład różnic między grupami.
źródło
Zgodnie z etykietą SO powinno to być napisane jako komentarz do @Johna K. Kruschke, ale dłuższe komentarze są trudne do sformułowania. Przepraszam.
lower_CredIaupper_CredIw oryginalnym poście zostały obliczone, jak wspomniałeś, z pełnych łańcuchów MCMC i są tylko nieznacznie sformatowane, aby lepiej porównać zlmewynikami. Chociaż wolisz HDI, są to proste kwantyle; z symetrycznym tylnym w tym przykładzie nie robi to dużej różnicy.Widziałem zastosowania w komisjach etycznych, w których obliczono moc statystyczną bez podania założenia o wielkości efektu. Nawet w przypadku, gdy nie ma możliwości zdefiniowania „klinicznie istotnego efektu”, trudno jest wyjaśnić tę koncepcję naukowcom medycznym. Jest to nieco łatwiejsze w przypadku prób nie-niższości, ale nie są one tak często przedmiotem badań.
Jestem więc całkiem pewien, że wprowadzenie LIN nie będzie możliwe do zaakceptowania - inne założenia, ludzie nie mogą mieć na myśli więcej niż jednej liczby. Czynniki Bayesa mogą działać, ponieważ istnieje tylko jedna liczba do zabrania do domu, podobnie jak wartości p.
Dziwi mnie, że ani @John K. Kruschke, ani @Ben Goodrich z zespołu Stan nie wspominają o priors; większość artykułów na ten temat wymaga szczegółowych dyskusji na temat wcześniejszej wrażliwości podczas prezentacji wyników.
Byłoby miło, gdyby w następnym wydaniu książki - miejmy nadzieję, że ze Stanem - można dodać pola „Jak opublikować to (w artykule niestatystycznym) zawierającym 100 słów” dla wybranych przykładów. Gdybym wziął twój rozdział 23.1 słowem, typowy artykuł z badań medycznych miałby długość 100 stron i cyfr ...
źródło