Como os bayesianos verificam seus métodos usando os métodos de simulação de Monte Carlo?

Formação : Tenho doutorado em psicologia social, onde estatística e matemática teóricas mal foram abordadas em meus cursos quantitativos. Durante o curso de graduação e pós-graduação, fui ensinado (como muitos de vocês também nas ciências sociais, provavelmente) através da estrutura freqüentista "clássica". Agora, eu também adoro R e usando métodos de simulação para verificar que os métodos de trabalho faz maneiramais sentido para mim do que provas matemáticas (novamente: experiência em uma ciência social quantitativa, não estatística teórica). Métodos freqüentistas e métodos de simulação juntos fazem muito sentido para mim. Como os freqüentadores veem a probabilidade como probabilidades de longo prazo (por exemplo, se eu fizer isso um número arbitrariamente grande de vezes, e isso acontece 50% das vezes, existe uma probabilidade de 50%). Podemos simular esse longo prazo com os métodos Monte Carlo!

Complicações : Desde a graduação, eu tenho muita consciência dos métodos bayesianos, e sempre houve pessoas na minha vida me chamando para o lado bayesiano, dizendo que os resultados eram mais fáceis de interpretar, que obtemos probabilidade de uma hipótese em vez dos dados dada uma hipótese, etc. Eu realmente participei disso e fiz uma aula bayesiana, li alguns livros e jornais bayesianos e agora estou bastante familiarizada com Stan e seus pacotes R associados.

Digite Mayo : Depois de pensar "Bayesiano é provavelmente o caminho do futuro" por um tempo, li a Inferência Estatística de Deborah Mayo como Teste Severo . Ela diz que não escolhe um lado no começo do livro, mas escolhe: é freqüentadora e grande parte do livro defende metodologias freqüentistas. Eu não quero necessariamente entrar em uma discussão sobre se achamos ou não válidos o modo como ela vê as evidências, mas isso me fez pensar: Bayes é realmente tudo o que é anunciado? Quero dizer, a multidão de Bayes está tão fraturada que nem sequer sei a maneira "certa" de analisar dados em uma estrutura bayesiana frequentemente. Normalmente, eu usaria apenasrstanarme estimativas pontuais atuais e intervalos confiáveis ​​... que frequentemente se alinham estreitamente com estimativas freqüentes e intervalos de confiança. Posso fazer comparações de modelos, mas sempre tenho medo de descrever os fatores de Bayes como comparações de probabilidade posteriores, etc.

Mais reflexão : O que eu continuei pensando no livro de Mayo foi: existe uma maneira de usar computadores para garantir que nossos métodos freqüentadores funcionem, porque probabilidade é o que vemos a longo prazo e podemos simular isso. Os bayesianos nem conseguem concordar com o que realmente é a probabilidade, ao que parece, dependendo da escola bayesiana (padrão, subjetiva, etc.). O que me leva à minha pergunta:

Pergunta : Como os bayesianos verificam se seus métodos definem a incerteza corretamente (isto é, calculam intervalos credíveis válidos e distribuições posteriores) usando métodos de simulação de Monte Carlo, se a probabilidade não é definida como taxas a longo prazo?

Exemplo : eu crio um gerador de dados. Isso simulará apenas a partir de uma distribuição de Bernoulli com uma probabilidade de 0,5:

set.seed(1839)
p <- .50
n <- 100
gen_dat <- function(n, p) {
  rbinom(n, 1, p)
}

Agora, digamos que eu queira ter certeza de que os intervalos de confiança em uma regressão logística são realmente válidos. Posso simular uma regressão várias vezes e garantir que o valor real da população caia no intervalo de confiança de 95% em 95% das vezes. É um modelo somente de interceptação, então só quero ter certeza de que estou calculando pcorretamente:

set.seed(1839)
iter <- 10000
results <- sapply(seq_len(iter), function(zzz) {
  mod <- glm(gen_dat(n, p) ~ 1, binomial)
  conf <- suppressMessages(confint(mod))
  log(p / (1 - p)) < max(conf) & log(p / (1 - p)) > min(conf)
})
mean(results)

Isso leva alguns minutos para ser executado, mas terminamos com a mean(results)ligação 0.9416. Isso é cerca de 95%, e estou confiante em dizer que o glmcomando está descrevendo incertezas de uma maneira válida. Tenho certeza de que teria chegado mais perto do nariz em 95% se eu aumentasse itere quisesse esperar aqui no meu laptop por mais tempo.

Por outro lado, vamos ajustar um modelo bayesiano para a mesma coisa:

library(rstanarm)
set.seed(1839)
dat <- data.frame(y = gen_dat(n, p))
stan_mod <- stan_glm(y ~ 1, binomial, dat)
summary(stan_mod)

Em parte, isso me dá:

Estimates:
                mean   sd    2.5%   25%   50%   75%   97.5%
(Intercept)    -0.1    0.2  -0.5   -0.2  -0.1   0.0   0.3  
mean_PPD        0.5    0.1   0.3    0.4   0.5   0.5   0.6  
log-posterior -73.0    0.7 -75.1  -73.1 -72.7 -72.5 -72.5  

Como os bayesianos não definem probabilidade como o que vemos a longo prazo, como posso usar métodos de simulação para verificar se stan_glmestá capturando com precisão a incerteza? Ou seja, como eu poderia confiar que esses intervalos credíveis são válidos, usando métodos de simulação? E, no momento, nem estou definindo um prior - como a inclusão de priors entra em jogo aqui, pois isso afetará nossas medidas de incerteza?

Quando eu estava tentando escrever uma regressão beta com um componente de modelo de obstáculo no Stan uma vez, pedi a alguém que me recomendasse: "Simule dados. Faça isso várias vezes, e as estimativas verdadeiras devem estar no intervalo credível de 95 % do tempo." Mas para mim, isso vai contra o que os bayesianos acreditam! Isso depende de entendimentos freqüentes de probabilidade! Então, como um bayesiano me convenceria de que o intervalo credível que estou recebendo da summary()chamada para o meu modelo está descrevendo com precisão a incerteza, usando métodos de simulação?

Objetivo da pergunta : Este é um exemplo trivial, mas muitas vezes os clientes me fornecem problemas difíceis. E eu tento coisas com as quais não estou familiarizado, então, frequentemente, faço um estudo de simulação para garantir que o que estou fazendo seja válido. Se eu escrevesse um modelo personalizado em Stan, como saberia que o que estou fazendo é legítimo? Como eu poderia usar métodos de simulação para verificar se o que estou fazendo em Stan realmente me diz o que quero saber?

Eu acho que você vê o problema lógico na sua pergunta. No paradigma freqüentista, é bom presumir uma verdade da população, gerar dados e ver se as estimativas têm uma boa cobertura, porque é isso que eles devem fazer. No paradigma bayesiano, no entanto, não existe uma verdade básica para gerar dados! Os bayesianos perguntam a probabilidade de tais verdades fornecerem dados; portanto, na simulação, precisamos de verdades diferentes que gerem dados e depois condicionem os dados. Na prática, acaba-se simulando a lei da probabilidade condicional, que, felizmente, sempre é válida por definição. Abordo esse problema exatamente em Rouder, 2014, Psychonomic Bulletin and Review. https://dx.doi.org/10.3758/s13423-014-0595-4

Como os bayesianos verificam se seus métodos definem a incerteza adequadamente (isto é, calculam intervalos credíveis válidos e distribuições posteriores) usando métodos de simulação de Monte Carlo, se a probabilidade não é definida como taxas a longo prazo?

Acredito que a confusão aqui seja sobre o propósito dos métodos de simulação nas estatísticas bayesianas. O único objetivo dos métodos Monte Carlo da Cadeia de Markov, como Gibbs Sampling ou Hamiltonian Monte Carlo, é calcular o denominador da regra de Bayes.

Obviamente, existem outros métodos disponíveis que tornariam desnecessário o MCMC. Alguns modelos podem ser expressos usando conjugação, outros através da aplicação de uma grade fina sobre o espaço do parâmetro; outros ainda podem ser resolvidos com testes de aceitação-rejeição. O MCMC é útil quando a integral é mal comportada.

Embora eu adorasse evitar a matemática, isso realmente não pode ser evitado. Ao examinar a regra de Bayes o numerador é composto de e . é uma probabilidade e não uma probabilidade; portanto, não soma a um, exceto por acaso. O denominador garante que igual a um. O objetivo do MCMC é determinar o número inferior. Observe que o número inferior é uma constante. É a probabilidade esperada.

A precisão desse número determina algumas estimativas de parâmetros, mas não todas. Se você estava usando o estimador a posteriori máximo, o MCMC é uma etapa desnecessária. Você deve criar um algoritmo de escalada. Por outro lado, é necessário determinar a média posterior ou um intervalo. Isso ocorre porque o intervalo de 95% deve ser 95% de alguma coisa e o denominador determina qual é a escala dessa coisa.

O objetivo do MCMC nas metodologias bayesianas é fazer com que as cadeias de Markov convergam para a densidade posterior. É isso. Não testa a validade de nada. É apenas uma tentativa de determinar um valor de ponto fixo. É uma forma de integração numérica. Como não há como saber, sem deixar que o algoritmo corra até o infinito, se todas as regiões densas foram cobertas, existe algum julgamento humano. O algoritmo terá um ponto de corte quando acreditar que está pronto, mas isso não significa que esteja realmente pronto.

Nas metodologias freqüentistas, o MCMC é frequentemente usado para testar a razoabilidade de um modelo ou aproximar numericamente uma solução quando não há uma analítica disponível. Não serve a nenhum propósito semelhante aqui.

Se eu escrevesse um modelo personalizado em Stan, como saberia que o que estou fazendo é legítimo? Como eu poderia usar métodos de simulação para verificar se o que estou fazendo em Stan realmente me diz o que quero saber?

Esta questão é muito mais difícil. Stan é um algoritmo rápido, o que significa que negocia velocidade para aumentar o risco de imprecisão. Stan, por construção, estará mais frequentemente correto do que incorreto. Existem outros algoritmos projetados para pesquisar amplamente no espaço de parâmetros os máximos locais que podem ser mais precisos, mas que serão muito lentos.

O que você deve fazer, antes de usar um algoritmo específico, é ler a literatura sobre esse algoritmo e examinar suas limitações funcionais. Infelizmente, esse geralmente é um trabalho matemático, pois o único objetivo real de qualquer método não conjugado é estimar

A segunda coisa que você pode fazer é validá-lo com um algoritmo alternativo. Os números nunca coincidirão, mas se você os considerar próximos o suficiente, estará bem.

Terceiro, a maioria dos pacotes pré-criados fornece avisos de que algo pode estar errado. Se um aviso aparecer, use outra coisa depois de investigar a origem do problema, para não recriá-lo em outro algoritmo.

Quarto, observe sua densidade anterior. Imagine que você tivesse uma densidade anterior de com conhecida apenas para simplificá-lo e uma probabilidade de . No mínimo, você deveria estar uau, ou eu estava errado, a amostra estava ruim ou há algo mais que eu deveria investigar.$\Pr(\mu)=\mathcal{N}(7,2^2)$$\sigma^2$$\mathcal{N}(25,.1^2)$

Quinto, e você deve fazer isso antes de iniciar Stan, em primeiro lugar, faça um gráfico das suas probabilidades marginais em uma ou duas dimensões. Há surpresas em algum lugar que possam interferir no algoritmo?

Como os bayesianos não definem probabilidade como o que vemos a longo prazo, como posso usar métodos de simulação para verificar se stan_glm está capturando com precisão a incerteza? Ou seja, como eu poderia confiar que esses intervalos credíveis são válidos, usando métodos de simulação? E, no momento, nem estou definindo um prior - como a inclusão de priors entra em jogo aqui, pois isso afetará nossas medidas de incerteza?

Se você não definir um prior, seu modelo não será válido. Se você não está definindo uma densidade anterior razoável, por que você usaria um modelo bayesiano? Os modelos freqüentistas minimizam o risco de perda máxima que poderia acontecer ao coletar uma amostra ruim. Eles são muito pessimistas e geralmente são necessárias mais informações para produzir o mesmo resultado que um método bayesiano produziria.

No entanto, isso é inútil sem o uso de uma boa densidade anterior. A densidade anterior permite que o método bayesiano minimize a perda média ao escolher uma amostra ruim. As informações anteriores agem como um esquema de ponderação, de modo que, se alguma amostra extrema for escolhida por infeliz acaso, a prévia enfraquece o papel que os dados desempenham.

EDIT Eu percebi que não forneci uma resposta específica. Foi à pergunta

Como eu poderia usar métodos de simulação para verificar se o que estou fazendo em Stan realmente me diz o que quero saber?

O que torna essa questão desafiadora é que, no paradigma bayesiano, os pontos fixos são, , a amostra. Nos métodos freqüentistas, os parâmetros são fixos e milhares de amostras invisíveis são criadas. No lado bayesiano da moeda, é a amostra que é fixada. Você precisa simular milhares de universos paralelos.$X$

Para ver o que pode ser como, imagine todas as funções de densidade de um sorteio com uma probabilidade desconhecida de serem cabeças e de serem caudas. Você observa seis cabeças e duas caudas. Imagine um pequeno espaço de parâmetro em que . Sua simulação consideraria todos os casos em que seis cabeças poderiam ser obtidas nas três distribuições binomiais objetivas. O posterior seria a média ponderada de cada parâmetro, sendo o valor verdadeiro. Sua distribuição preditiva seria a soma das distribuições binomiais ponderadas.$p$$1-p$$p\in\{1/3,1/2,2/3\}$

De importância para você, é impossível que a previsão bayesiana seja a verdadeira distribuição. Uma das três distribuições é a verdadeira distribuição. Os métodos bayesianos ponderam sua probabilidade com base no valor observado e no anterior. O posterior nunca pode ser a verdadeira distribuição, nem a densidade preditiva.

Ele está perguntando "qual é a probabilidade de ver seis caras e duas caudas no conjunto de todas as explicações possíveis (parâmetros, modelos, etc.)".

O freqüentista afirmava que uma das três opções era o valor verdadeiro, tornando-o nulo. Seis caras e duas caudas falsificariam mas não as outras. Se, por acaso, você escolher a correta das três distribuições, estará perfeitamente correta. Caso contrário, você estará errado.$H_0:p=1/3,$

Se você usasse simulações para manter uma amostra fixa, descobriria que Stan teria um desempenho admirável, pois o teorema de Bayes é um teorema matemático. É ex-post ideal. Tudo o que você descobriria é que o algoritmo implementou corretamente o teorema de Bayes até o nível de erro natural na estimativa do denominador.

Há três coisas que você pode fazer. Primeiro, você pode usar os métodos de pontuação do modelo para dados fora da amostra. Segundo, você pode usar uma seleção de modelo bayesiano ou processo de média de modelo. Terceiro, você pode tratá-lo como um problema freqüentista e construir a distribuição amostral dos estimadores.

Para o primeiro, os métodos de pontuação são uma literatura inteira em si. Você deve pesquisá-los. A seleção do modelo bayesiano e a média do modelo tratam os modelos como parâmetros. Para a seleção de modelos, é calculada a probabilidade dos modelos serem verdadeiros. Para a média do modelo, calcula-se a probabilidade de cada modelo ser verdadeiro e isso serve como ponderação sobre o espaço do modelo. Finalmente, você pode tratá-lo como um modelo freqüentista.

O último será um problema em muitos casos padrão por causa do anterior. Para modelos com três ou mais dimensões e uma distribuição normal, a densidade posterior não será integrada à unidade se a densidade anterior não for uma densidade adequada. Em outras palavras, você precisa morder o marcador e escolher um prior para qualquer modelo com qualquer complexidade real.

A presença de um anterior adequado corretamente centralizado força o caso em que o método bayesiano será superior ao método frequentista correspondente devido à melhoria da informação. O método bayesiano vencerá sob qualquer padrão razoável. Isso não se deve a uma falha no método Frequentist, mas o método Bayesiano assume informações externas. O método Frequentist, considerando apenas as informações da amostra, terá menos informações se você tiver um real anterior.

Novamente, se você não tem um real anterior, então por que você está usando um método bayesiano?