С помощью пакетов Stan и frontend, rstanarm
или brms
я могу легко анализировать данные байесовским способом, как я делал раньше со смешанными моделями, такими как lme
. Хотя у меня на столе лежит большая часть книг и статей Крушке-Гельмана-Вагенмейкера и т. Д., Они не говорят мне, как суммировать результаты для медицинской аудитории, разрываясь между Скиллой гнева Байесиана и Харибдой медицинских рецензентов ( «мы хотим значимости, а не этого рассеянного материала»).
Пример: желудочная частота (1 / мин) измеряется в трех группах; здоровые контроли являются эталоном. Есть несколько измерений для каждого участника, поэтому я часто использовал следующую смешанную модель lme
:
summary(lme(freq_min~ group, random = ~1|study_id, data = mo))
Немного отредактированные результаты:
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
Для простоты я буду использовать 2 * стандартную ошибку как 95% CI.
В частом контексте я бы суммировал это как:
- В контрольной группе предполагаемая частота составляла 2,7 / мин (возможно, добавьте здесь CI, но иногда я избегаю этого из-за путаницы, создаваемой абсолютной и разностной CI).
- В группе без симптомов частота была выше на 0,4 / мин, ДИ (от 0,11 до 0,59) / мин, р = 0,006 по сравнению с контролем.
- В группе с симптомами частота была выше на 0,2 / мин, ДИ (от -0,04 до 0,4) / мин, р = 0,11, чем в контрольной группе.
Речь идет о максимально приемлемой сложности для медицинской публикации, рецензент, вероятно, попросит меня добавить «незначительный» во втором случае.
Здесь то же самое с stan_lmer
и по умолчанию приоры.
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
где CredI - 90% вероятные интервалы (см. виньетку rstanarm, почему 90% используется по умолчанию.)
Вопросов:
- Как перевести вышеприведенное резюме в байесовский мир?
- В какой степени требуется предварительное обсуждение? Я совершенно уверен, что статья вернется с обычным «субъективным допущением», когда я упомяну приоры; или, по крайней мере, с "без технической дискуссии, пожалуйста". Но все Байесовские власти требуют, чтобы толкование действовало только в контексте приоры.
- Как я могу дать некоторую «значимость» суррогата в формулировке, не изменяя байесовским концепциям? Что-то вроде "достоверно другой" (уууу ...) или почти достоверно другой (буоха ..., звучит как "на краю значимости").
Джона Габри и Бен Гудрич (2016). rstanarm: Байесовское моделирование прикладной регрессии через Стэн. Версия пакета R 2.9.0-3. https://CRAN.R-project.org/package=rstanarm
Команда разработчиков Stan (2015). Stan: библиотека C ++ для вероятности и выборки, версия 2.8.0. URL http://mc-stan.org/ .
Пауль-Кристиан Бюркнер (2016). brms: модели байесовской регрессии с использованием Stan. Версия пакета R 0.8.0. https://CRAN.R-project.org/package=brms
Pinheiro J, Bates D, DebRoy S, Sarkar D и R Core Team (2016). nlme: линейные и нелинейные модели смешанных эффектов . Версия пакета R 3.1-124, http://CRAN.R-project.org/package=nlme>.
mean(as.matrix(freq_stan)[,"groupwith_symptoms"] < 0)
.group_nosymptoms
а затем сказать, что вероятность того, что он будет отрицательным, равна1 / draws
. Но для перехвата цепочка никогда не собирается уходить в отрицательную территорию для этих данных, поэтому я думаю, вы могли бы сказать, что вероятность меньше, чем1 / draws
.Ответы:
Быстрые мысли:
1) Ключевой вопрос заключается в том, какой прикладной вопрос вы пытаетесь ответить для своей аудитории, потому что это определяет, какую информацию вы хотите получить из своего статистического анализа. В этом случае мне кажется, что вы хотите оценить величину различий между группами (или, возможно, величину соотношений групп, если эта мера более знакома вашей аудитории). Величина различий напрямую не указывается в анализе, который вы представили в вопросе. Но получить байесовский анализ просто так: вы хотите последующее распределение различий (или отношений). Затем из апостериорного распределения различий (или отношений) вы можете сделать прямое утверждение вероятности, такое как это:
«Наиболее вероятные различия в 95% находятся между [низким пределом ИЧР в 95%] и [высоким пределом ИЧР в 95%]» (здесь я использую интервал с высокой плотностью 95% [ИЧР] в качестве вероятного интервала, и потому что определение значений параметров наивысшей плотности, которые они обозначают как «наиболее достоверные»)
Аудитория медицинского журнала интуитивно и правильно поняла бы это утверждение, потому что это то, о чем обычно думает аудитория, это значение доверительного интервала (даже если это не означает доверительный интервал).
Как вы получаете различия (или отношения) от Стэна или JAGS? Только путем постобработки завершенной цепочки MCMC. На каждом этапе цепочки вычислите соответствующие различия (или соотношения), затем изучите последующее распределение различий (или соотношений). Примеры приведены в DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ для MCMC, как правило, на рисунке 7.9 (с. 177), для JAGS на рисунке 8.6 (с. 211) и для Stan в разделе 16.3 (с. 468) и т. Д.!
2) Если вы по традиции вынуждены делать заявления о том, отклоняется или нет нулевая разница, у вас есть два байесовских варианта.
2А) Один из вариантов - сделать вероятностные заявления относительно интервалов, близких к нулю, и их связи с ИЧР. Для этого вы устанавливаете область практической эквивалентности (ROPE) около нуля, которая является просто порогом принятия решения, подходящим для вашей прикладной области - насколько большая разница тривиально мала? Установление таких границ обычно выполняется, например, в клинических тестах на неполноценность. Если у вас есть мера «размера эффекта» в вашей области, могут быть соглашения для «маленького» размера эффекта, и пределы ROPE могут быть, скажем, половиной небольшого эффекта. Затем вы можете сделать прямые вероятностные заявления, такие как эти:
«Только 1,2% заднего распределения различий практически эквивалентно нулю»
и
«95% наиболее вероятных различий практически не эквивалентны нулю (т. Е. 95% ИЧР и ВЕРЕВКА не перекрываются), и поэтому мы отвергаем ноль». (обратите внимание на различие между утверждением вероятности и последующим распределением и последующим решением, основанным на этом утверждении)
Вы также можете принять нулевую разницу для практических целей, если все 95% наиболее вероятных значений практически эквивалентны нулю.
2B) Второй байесовский вариант - это проверка байесовской нулевой гипотезы. (Обратите внимание, что метод выше не былназывается «тестирование гипотез»!) Байесовское тестирование нулевых гипотез делает сравнение байесовской модели предыдущего распределения, которое предполагает, что разница может быть только нулевой, по сравнению с альтернативным предварительным распределением, которое предполагает, что разница может быть некоторой разбросанной областью возможностей. Результат такого сравнения моделей (обычно) очень сильно зависит от конкретного выбора альтернативного распределения, и поэтому необходимо сделать тщательное обоснование для выбора альтернативы. Лучше использовать как минимум умеренно информированные априоры как для нулевого, так и для альтернативного, чтобы сравнение моделей было действительно значимым. Обратите внимание, что сравнение моделей дает иную информацию, нежели оценка различий между группами, потому что сравнение моделей касается другого вопроса. Таким образом, даже при сравнении моделей,
Могут быть способы сделать тест байесовской нулевой гипотезы из вывода Stan / JAGS / MCMC, но я не знаю в этом случае. Например, можно попробовать приближение Сэвиджа-Дики к коэффициенту Байеса, но это будет зависеть от знания предыдущей плотности различий, что потребует некоторого математического анализа или некоторого дополнительного приближения MCMC из предыдущего.
Два метода принятия решения о нулевых значениях обсуждаются в гл. 12 из DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ . Но я действительно не хочу, чтобы это обсуждение было отвлечено дискуссией о «правильном» способе оценки нулевых значений; они просто разные и предоставляют разную информацию. Основным пунктом моего ответа является пункт 1 выше: посмотрите на последующее распределение различий между группами.
источник
Следуя SO этикету, это должно было быть написано как комментарий к @John K. Kruschke, но более длинные комментарии трудно структурировать. Сожалею.
lower_CredI
иupper_CredI
в исходном посте были вычислены, как вы упомянули из полных цепочек MCMC, и только слегка переформатированы для лучшего сравнения сlme
выходными данными. Хотя вы предпочитаете ИЧР, это простые квантили; с симметричным апостериором в этом примере это не имеет большого значения.Я видел заявления в комитеты по этике, в которых статистическая мощность была рассчитана без указания предположения о величине эффекта. Даже для случая, когда нет никакого способа определить «клинически значимый эффект», трудно объяснить эту концепцию медицинским исследователям. Для испытаний, не связанных с неполноценностью, это немного проще, но они не так часто становятся предметом исследования.
Поэтому я совершенно уверен, что введение ROPES не будет приемлемым - другие предположения, люди не могут иметь в виду более одного числа. Байесовские факторы могут сработать, потому что есть только одно число, которое можно вернуть домой, как p-значения ранее
Я удивлен, что ни @Джон К. Крушке, ни @Бен Гудрич из команды Стэна не упоминают приоры; большинство работ по этому вопросу требуют подробного обсуждения предварительной чувствительности при представлении результатов.
Было бы неплохо, если бы в следующем издании вашей книги - надеюсь, со Стэном - вы могли бы добавить блоки «Как опубликовать это (в нестатистическом документе) со 100 словами» для выбранных примеров. Когда я взял бы вашу главу 23.1 словом, типичный медицинский исследовательский документ был бы длиной в 100 страниц и цифрами ...
источник