Monte Carlo via Cadeias de Markov

Hoje é muito difícil usar Estatística Bayesiana sem usar métodos de Monte Carlo via Cadeias de Markov (Markov Chain Monte Carlo, MCMC). Isto ocorre devido a complexidade dos modelos que os dados reais exigem. Quando começamos a estudar Estatística Bayesiana, os parâmetros de interesse são costumeiramente, a média e/ou variância de uma população normalmente distribuída, o parâmetro de proporção de uma distribuição Binomial, ou até mesmo a taxa lambda de uma distribuição Poisson. Nestas situações, caso as distribuições a priori sejam escolhidas adequadamente, temos distribuições a posteriori conhecidas e as inferências são realizadas facilmente. Para estes modelos citados anteriormente, tem-se as chamadas famílias conjugadas de distribuições. Por exemplo, se escolhemos uma distribuição a priori gama para a taxa lambda de uma distribuição Poisson, a distribuição a posteriori também terá distribuição gama. Isto também ocorre com a distribuição Normal quando escolhemos uma distribuição Normal para média, o que retorna uma distribuição posteriori a Normal para a média.

Toda distribuição pertencente a família exponencial tem uma família conjugada de distribuições a priori. Entretanto, caso tivéssemos indícios de que o modelo mais adequado aos nossos dados é o modelo t-Student, distribuição não pertencente à família exponencial e que não conheço um família conjugada para a média desta distribuição. Uma alternativa é o uso de métodos de MCMC.

Monte Carlo via Cadeias de Markov são métodos de simulação baseados em Cadeias de Markov ergódicas onde a distribuição estacionária do processo estocástico é a distribuição a posteriori de interesse. Existem diversos tipos de métodos MCMC, dentre estes, o Algoritmo de Metropolis-Hastings e Algoritmo (Amostrador) de Gibbs são os mais utilizados. O Algoritmo de Metropolis-Hastings é muito utilizado por ser pouco restritivo com relação a distribuição a posteriori. Pois para o uso deste algoritmo, é suficiente e necessário ter apenas a distribuição a posteriori a menos de uma constante de proporcionalidade e escolher uma distribuição proposta adequada, o que nos trará uma taxa de rejeição, pois alguns valores das distribuição proposta serão rejeitados. O algoritmo de Gibbs é mais restritivo, pois para seu uso, é necessário conhecer as distribuições condicionais completas. Entretanto, não temos que escolher uma distribuição proposta e com isso não existe taxa de rejeição. Em cada um dos casos torna-se imprescindível verificar se a distribuição a posteriori é realmente uma distribuição de probabilidade (distribuição própria). É sugerido o uso de distribuições a priori próprias nestes caso.

Dado que temos um modelo estatístico, por exemplo um modelo linear generalizado (MLG) onde nossos parâmetros de interesse são os coeficientes regressores. Infelizmente neste caso não temos uma família conjugada e uma alternativa é utilizar MCMC. Uma pergunta comum é: Gibbs ou Metropolis-Hastings? Costuma-se tentar usar Gibbs na primeira tentativa, por ter uma implementação mais fácil. Entretanto, é difícil encontrar as distribuições condicionais completas em MLG. Então uma opção é utilizar Metropolis-Hastings. Outra opção, caso você goste de usar o algoritmo de Gibbs é tentar aumentar os dados, e com isso obter as distribuições condicionais completas.

Assumindo que o método MCMC foi escolhido e que já temos uma cadeia, o próximo passo é encontrar o ponto cujo a cadeia atingiu a distribuição estacionária. Encontrado este ponto, todos os valores obtidos a partir deste ponto serão assumidos uma amostra da distribuição a posteriori. Entretanto, neste momento entra um estudo muito importante que é a análise de convergência da Cadeia. Muitas vezes não é fácil encontrar este ponto. Existem vários meios de testar se a cadeia atingiu um comportamento estacionário; métodos baseados em representações gráficas, monitoração das médias ergódicas, método de Gelman e Rubin, método dos quantis de Raftey e Lewis, método de Geweke, além do uso de cadeias paralelas.

Assumindo agora que a cadeia atingiu o comportamento estacionário após um período de aquecimento, devemos notar que a Cadeia de Markov nos retorna valores correlacionados da distribuição a posteriori, principalmente quando utilizamos o Amostrador de Gibbs. Neste caso, devemos saltar observações para tentar reduzir a correlação da cadeia, caso contrário a variância da distribuição a posteriori será superestimada. O tamanho do salto é escolhido de acordo com correlação existente na cadeia. Uma forma de ver esta correlação é visualizar os gráficos de autocorrelações seriais das cadeiais, tanto os individuais como os cruzados, caso forem gerados cadeias simultâneas para várias quantidades de interesse. Após estes passos, temos uma amostra não correlacionada da distribuições a posteriori dos nossos parâmetros de interesse e as inferências podem ser realizadas.

Um breve resumo: Estamos interessados em uma certa quantidade (parâmetro) de uma população de interesse. Assume-se um modelo probabilístico adequado para os dados e deseja-se fazer inferência Bayesiana.

Primeira pergunta: Posso e é interessante utilizar uma família conjugada?

SIM -> Ótimo, só olhar num livro quais os parâmetros da distribuição a posteriori. Os quais dependem dos parâmetros da priori (hiperparâmetros) e dos dados.
NÃO -> Tentarei usar um método baseado em MCMC.

Segunda pergunta: Qual método MCMC devo utilizar?

Esta é uma pergunta difícil de ser respondida, cada um tem seu gosto por um método. Para o amostrador de Gibbs deve-se ter as distribuições condicionais completas, o que torna o algoritmo bastante fácil de se implementar. Entretanto, no algoritmo de Metropolis-Hastings estas distribuições não são necessárias, mas tem-se que escolher uma distribuição proposta adequada.

Terceira Pergunta: Qual deve ser meu período de aquecimento?

O interessante é utilizar vários métodos, desde gráficos à testes mais sofisticados.

Quarta Pergunta: De quanto em quanto devo saltar para assumir que minha amostra da distribuição a posteriori é não correlacionada?

Para isso é necessário ver os correlogramas, a partir deles, pode-se ter um idéia do salto.

Após todos estes passos, o que nos resta é escolher as melhores medidas para descrever nosso conhecimento sobre o parâmetro de interesse. Alguns casos, onde queremos uma estimativa pontual, pode-se usar a média a posteriori, outro casos, quando a distribuição é bastante assimétrica, a mediana pode ser a medida mais interessante.

Neste texto quis mostrar superficialmente o uso dos métodos de Monte Carlo via Cadeias de Markov em Estatística Bayesiana. Tentar mostrar de forma fácil a metodologia. Em um novo post irei escrever um pouco mais sobre análise de convergência e comparação de algoritmos, pois são assuntos que tenho interesse. Espero melhorar este texto com o tempo, e caso você tenha alguma sugestão ou crítica, será muito bem vinda.

Em breve colocarei hiperlinks nos assuntos.