Com os pacotes Stan e frontend, rstanarm
ou brms
posso analisar facilmente os dados da maneira bayesiana, como fiz antes com modelos mistos como lme
. Embora eu tenha a maioria dos livros e artigos de Kruschke-Gelman-Wagenmakers-etc em minha mesa, eles não me dizem como resumir os resultados para uma audiência médica, divididos entre a ira de Skylla da Bayesiana e os Charybdis de revisores médicos ( "queremos significados, não coisas difusas").
Um exemplo: a frequência gástrica (1 / min) é medida em três grupos; controles saudáveis são a referência. Existem várias medidas para cada participante, então, como freqüentador, usei o seguinte modelo misto lme
:
summary(lme(freq_min~ group, random = ~1|study_id, data = mo))
Resultados ligeiramente editados:
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
Para simplificar, usarei o erro 2 * std como IC de 95%.
No contexto frequentista, eu teria resumido isso como:
- No grupo de controle, a frequência estimada foi de 2,7 / min (talvez adicione IC aqui, mas às vezes evito isso devido à confusão criada pelo IC absoluto e pela diferença).
- No grupo no_symptoms, a frequência foi maior em 0,4 / min, IC (0,11 a 0,59) / min, p = 0,006 que o controle.
- No grupo com sintomas, a frequência foi maior em 0,2 / min, IC (-0,04 a 0,4) / min, p = 0,11 que o controle.
Trata-se da complexidade máxima aceitável para uma publicação médica, o revisor provavelmente me pedirá para adicionar "não significativo" no segundo caso.
Aqui está o mesmo com stan_lmer
e anteriores padrão.
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
onde CredI são intervalos de 90% credíveis (veja a vinheta do rstanarm por que 90% é usado como padrão).
Questões:
- Como traduzir o resumo acima para o mundo bayesiano?
- Até que ponto é necessária uma discussão prévia? Tenho certeza de que o artigo voltará com a "suposição subjetiva" usual quando menciono os anteriores; ou pelo menos com "nenhuma discussão técnica, por favor". Mas todas as autoridades bayesianas solicitam que a interpretação seja válida apenas no contexto de anteriores.
- Como posso fornecer algum substituto de "significado" na formulação, sem trair os conceitos bayesianos? Algo como "credivelmente diferente" (uuuh ...) ou quase credivelmente diferente (buoha ..., soa como "à beira do significado).
Jonah Gabry e Ben Goodrich (2016). rstanarm: Modelagem de regressão aplicada bayesiana via Stan. Pacote R versão 2.9.0-3. https://CRAN.R-project.org/package=rstanarm
Equipe de Desenvolvimento Stan (2015). Stan: Uma biblioteca C ++ para probabilidade e amostragem, versão 2.8.0. URL http://mc-stan.org/ .
Paul-Christian Buerkner (2016). brms: Modelos de regressão bayesiana usando Stan. Pacote R versão 0.8.0. https://CRAN.R-project.org/package=brms
Pinheiro J, Bates D, DebRoy S, Sarkar D e R Core Team (2016). nlme: Modelos de efeitos mistos lineares e não lineares . Pacote R versão 3.1-124, http://CRAN.R-project.org/package=nlme>.
fonte
mean(as.matrix(freq_stan)[,"groupwith_symptoms"] < 0)
.group_nosymptoms
e, em seguida, dizer que a probabilidade de ser negativo é1 / draws
. Mas para a interceptação, a cadeia nunca vai entrar em território negativo para esses dados, então acho que você poderia dizer que a probabilidade é menor que1 / draws
.Respostas:
Pensamentos rápidos:
1) A questão principal é qual pergunta aplicada você está tentando responder para o seu público, porque isso determina quais informações você deseja da sua análise estatística. Nesse caso, parece-me que você deseja estimar a magnitude das diferenças entre os grupos (ou talvez a magnitude das proporções dos grupos, se essa for a medida mais familiar para o seu público). A magnitude das diferenças não é fornecida diretamente pelas análises que você apresentou na pergunta. Mas é fácil obter o que você quer da análise bayesiana: você quer a distribuição posterior das diferenças (ou proporções). Então, a partir da distribuição posterior das diferenças (ou proporções), você pode fazer uma declaração de probabilidade direta como esta:
"As diferenças de 95% com maior credibilidade caem entre [limite baixo de 95% do IDH] e [limite alto de 95% do IDH]" (aqui estou usando o intervalo de densidade de 95% mais alto [IDH] como o intervalo de credibilidade, e porque são definição dos valores dos parâmetros de densidade mais alta, eles são considerados "mais credíveis")
Um público de revistas médicas entenderia intuitivamente e corretamente essa afirmação, porque é o que a audiência normalmente pensa que é o significado de um intervalo de confiança freqüentista (mesmo que isso não signifique um intervalo de confiança freqüentista).
Como você obtém as diferenças (ou proporções) de Stan ou JAGS? Apenas pelo pós-processamento da cadeia MCMC concluída. Em cada etapa da cadeia, calcule as diferenças (ou proporções) relevantes e examine a distribuição posterior das diferenças (ou proporções). Exemplos são dados no DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ para o MCMC geralmente na Figura 7.9 (p. 177), para o JAGS na Figura 8.6 (p. 211) e para Stan na Seção 16.3 (p 468), etc.!
2) Se você é compelido pela tradição a fazer uma declaração sobre se uma diferença de zero é ou não rejeitada, você tem duas opções bayesianas.
2A) Uma opção é fazer declarações de probabilidade em relação a intervalos próximos de zero e sua relação com o IDH. Para isso, você configura uma região de equivalência prática (ROPE) em torno de zero, que é apenas um limite de decisão apropriado para seu domínio aplicado - qual o tamanho da diferença é trivialmente pequeno? O estabelecimento desses limites é rotineiramente realizado em testes clínicos de não inferioridade, por exemplo. Se você tiver uma medida de 'tamanho do efeito' em seu campo, pode haver convenções para o tamanho do efeito 'pequeno', e os limites do ROPE podem ser, digamos, metade do efeito pequeno. Em seguida, você pode fazer declarações de probabilidade direta como estas:
"Apenas 1,2% da distribuição posterior das diferenças é praticamente equivalente a zero"
e
"As 95% de diferenças mais credíveis nem todas são praticamente equivalentes a zero (ou seja, o IDH e o ROPE de 95% não se sobrepõem) e, portanto, rejeitamos zero". (observe a distinção entre a declaração de probabilidade da distribuição posterior e a decisão subsequente com base nessa declaração)
Você também pode aceitar uma diferença de zero, para fins práticos, se os valores com mais de 95% de credibilidade forem praticamente equivalentes a zero.
2B) Uma segunda opção bayesiana é o teste de hipótese nula bayesiana. (Observe que o método acima não foichamado "teste de hipótese"!) O teste de hipótese nula bayesiana faz uma comparação do modelo bayesiano de uma distribuição anterior que assume a diferença só pode ser zero em relação a uma distribuição anterior alternativa que assume que a diferença pode ser uma gama difusa de possibilidades. O resultado dessa comparação de modelo (geralmente) depende muito fortemente da escolha particular da distribuição alternativa, e, portanto, deve ser feita uma justificativa cuidadosa para a escolha da alternativa anterior. É melhor usar priors pelo menos levemente informados, tanto para o nulo quanto para a alternativa, para que a comparação do modelo seja genuinamente significativa. Observe que a comparação do modelo fornece informações diferentes da estimativa das diferenças entre os grupos, porque a comparação do modelo está abordando uma questão diferente. Assim, mesmo com uma comparação de modelos,
Pode haver maneiras de fazer um teste de hipótese nula bayesiana a partir da saída Stan / JAGS / MCMC, mas eu não sei neste caso. Por exemplo, pode-se tentar uma aproximação de Savage-Dickey a um fator de Bayes, mas isso dependeria do conhecimento da densidade anterior das diferenças, o que exigiria alguma análise matemática ou alguma aproximação adicional do MCMC do anterior.
Os dois métodos para decidir sobre valores nulos são discutidos no cap. 12 do DBDA2E https://sites.google.com/site/doingbayesiandataanalysis/ . Mas eu realmente não quero que essa discussão seja desviada por um debate sobre a maneira "adequada" de avaliar valores nulos; eles são diferentes e fornecem informações diferentes. O ponto principal da minha resposta é o ponto 1 acima: Observe a distribuição posterior das diferenças entre os grupos.
fonte
Seguindo a etiqueta do SO, isso deveria ter sido escrito como um comentário para @ John K. Kruschke, mas comentários mais longos são difíceis de estruturar. Desculpe.
lower_CredI
eupper_CredI
no post original foram computados como você mencionou nas cadeias completas do MCMC e são reformatados apenas ligeiramente para melhor comparação com alme
saída. Enquanto você é a favor do IDH, esses são quantis simples; com o posterior simétrico neste exemplo, não faz grande diferença.Eu já vi aplicações em comitês de ética em que o poder estatístico foi calculado sem afirmar a suposição sobre o tamanho do efeito. Mesmo no caso em que não há como definir um "efeito clinicamente relevante", é difícil explicar o conceito aos pesquisadores médicos. É um pouco mais fácil para estudos de não inferioridade, mas estes não são frequentemente objeto de um estudo.
Portanto, tenho certeza de que a introdução do ROPES não será aceitável - outras suposições, as pessoas não podem ter mais de um número em mente. Os fatores de Bayes podem funcionar, porque há apenas um número para levar para casa, como os valores de p antes.
Estou surpreso que nem John K. Kruschke nem Ben Goodrich, da equipe de Stan, mencionem anteriores; a maioria dos trabalhos sobre o assunto solicita uma discussão detalhada da sensibilidade anterior ao apresentar resultados.
Seria bom que, na próxima edição do seu livro - espero que com Stan - você pudesse adicionar as caixas "Como publicar isso (em um artigo não estatístico) com 100 palavras" para exemplos selecionados. Quando eu pegava seu capítulo 23.1 por palavra, um trabalho típico de pesquisa médica teria 100 páginas e figuras ...
fonte