Amostragem MCMC e diagnósticos
Execute a amostragem, inspecione cadeias e interprete os diagnósticos Rhat e ESS.
Amostragem MCMC e diagnósticos é uma aula grátis de R Academy no CoddyKit. Esta é a aula 3 de 4. Você pode ler a aula completa abaixo gratuitamente — depois pratica ao vivo no navegador com um editor de código integrado e um tutor de IA 24/7. Faz parte do caminho de aprendizado de R Academy, e seu progresso é sincronizado entre a web e o app CoddyKit. O curso de R Academy inclui 4 aulas no total.
O que é MCMC?
Markov Chain Monte Carlo (MCMC) é uma família de algoritmos para extrair amostras de uma distribuição de probabilidade quando a amostragem direta é impossível. Na estatística bayesiana, MCMC extrai amostras da distribuição a posteriori P(parameters | data).
O RStan implementa o No-U-Turn Sampler (NUTS), um algoritmo MCMC de última geração.
Um modelo Stan simples
Um modelo Stan é um bloco de texto que define tipos de dados, parâmetros e o logaritmo da distribuição a posteriori. O modelo mais simples estima a média de uma distribuição normal com variância conhecida.
# library(rstan)
#
# stan_code <- '
# data {
# int<lower=0> N;
# vector[N] y;
# }
# parameters {
# real mu;
# real<lower=0> sigma;
# }
# model {
# mu ~ normal(0, 10); // prior
# sigma ~ exponential(1); // prior
# y ~ normal(mu, sigma); // likelihood
# }
# 'Chamando stan() para amostrar
stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000) compila o modelo (uma vez) e, em seguida, extrai amostras. Com iter=2000 e warmup=1000, cada cadeia produz 1000 amostras após o aquecimento — 4000 no total em 4 cadeias.
# library(rstan)
# options(mc.cores = parallel::detectCores())
#
# y <- c(2.1, 1.8, 2.4, 1.9, 2.3, 2.0, 1.7, 2.2)
# stan_data <- list(N = length(y), y = y)
#
# fit <- stan(
# model_code = stan_code,
# data = stan_data,
# chains = 4,
# iter = 2000,
# warmup = 1000,
# seed = 42
# )print(fit) — A tabela de resumo
print(fit) mostra uma tabela de resumo para cada parâmetro, com média a posteriori, desvio-padrão, quantis, Rhat e n_eff. Esses dois diagnósticos são os primeiros itens a verificar.
# print(fit)
#
# Example output:
# mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff Rhat
# mu 2.05 0.00 0.15 1.76 1.95 2.05 2.15 2.34 3842 1
# sigma 0.22 0.00 0.06 0.13 0.18 0.21 0.25 0.37 3521 1
# lp__ 4.38 0.02 1.01 1.60 3.91 4.71 5.18 5.50 2148 1O critério de convergência Rhat
Rhat (fator de redução potencial da escala) compara a variância dentro das cadeias com a variância entre as cadeias. Valores próximos de 1,0 indicam que todas as cadeias convergiram para a mesma distribuição.
- Rhat < 1.01 — convergência alcançada (padrão atual)
- Rhat > 1.01 — as cadeias não se misturaram; execute mais iterações
- Rhat > 1.1 — problema sério de convergência
# Check Rhat for all parameters:
# s <- summary(fit)$summary
# rhat_vals <- s[, 'Rhat']
# cat('Max Rhat:', max(rhat_vals, na.rm = TRUE), '
')
# if (any(rhat_vals > 1.01, na.rm = TRUE)) {
# warning('Convergence issue detected!')
# } else {
# cat('All Rhat < 1.01 — chains converged
')
# }n_eff — Tamanho efetivo da amostra
n_eff (tamanho efetivo da amostra) leva em conta a autocorrelação entre amostras MCMC sucessivas. Amostras correlacionadas carregam menos informação que amostras independentes.
- n_eff próximo do total de iterações — amostras quase independentes, excelente
- n_eff / total_samples > 0.1 — geralmente aceitável
- n_eff muito baixo — autocorrelação alta; considere reparametrizar o modelo
# s <- summary(fit)$summary
# n_eff_vals <- s[, 'n_eff']
# total_samples <- 4 * 1000 # chains * post-warmup iter
# ratio <- n_eff_vals / total_samples
# cat('n_eff ratio (mu) :', round(ratio['mu'], 2), '
')
# cat('n_eff ratio (sigma):', round(ratio['sigma'], 2), '
')traceplot() para convergência visual
traceplot(fit, pars = 'mu') plota os valores amostrados de mu ao longo das iterações para cada cadeia. Cadeias convergentes se parecem com uma lagarta difusa — todas sobrepostas, sem tendências ou desvios. Cadeias divergentes vagueiam ou permanecem separadas.
# library(rstan)
#
# traceplot(fit, pars = c('mu', 'sigma'), inc_warmup = FALSE)
#
# Good traceplot characteristics:
# - All 4 chains overlapping completely (same range)
# - No visible drift or trend
# - Rapid mixing (values jump around quickly)
# - No flat regions (stuck sampler)
cat('A healthy traceplot looks like a fuzzy caterpillar
')pairs() para correlações a posteriori
pairs(fit, pars = c('mu', 'sigma')) mostra uma matriz de gráficos de dispersão das amostras a posteriori. Ela revela correlações entre parâmetros e destaca transições divergentes (plotadas em vermelho), que indicam regiões com as quais o amostrador tem dificuldade.
# pairs(fit, pars = c('mu', 'sigma'))
#
# What to look for:
# - Elliptical clouds: mild correlation (OK)
# - Banana / funnel shapes: reparameterization needed
# - Red dots (divergences): geometry problem in posterior
# => increase adapt_delta: stan(..., control=list(adapt_delta=0.95))
cat('Red dots in pairs() indicate divergent transitions — investigate!
')Extraindo amostras a posteriori
extract(fit, pars = 'mu')$mu retorna um vetor numérico com todas as amostras de mu obtidas após o aquecimento. Use essas amostras para calcular qualquer resumo da distribuição a posteriori: média, intervalos de credibilidade ou probabilidade de uma condição.
# mu_samples <- extract(fit, pars = 'mu')$mu
# cat('Posterior mean :', mean(mu_samples), '
')
# cat('95% CI:', quantile(mu_samples, c(0.025, 0.975)), '
')
# cat('P(mu > 2):', mean(mu_samples > 2), '
')
# hist(mu_samples, main = 'Posterior of mu', xlab = 'mu', col = 'steelblue')Iniciando o ShinyStan para diagnósticos interativos
shinystan::launch_shinystan(fit) abre um aplicativo Shiny interativo com gráficos de traços, distribuições a posteriori, gráficos pairs e diagnósticos NUTS, tudo em um só lugar. É a ferramenta mais abrangente para explorar um ajuste do RStan.
# install.packages('shinystan')
# library(shinystan)
#
# shinystan::launch_shinystan(fit)
#
# ShinyStan tabs:
# - Diagnose: Rhat, n_eff, divergences, energy
# - Explore: marginal posteriors, scatter plots
# - Model: Stan code, data
# - NUTS: step size, tree depth per chainCorreções comuns de convergência
Quando Rhat > 1.01 ou quando houver divergências:
- Aumente
iterewarmup - Aumente
adapt_deltaem direção a 1,0 (por exemplo, 0,95) emcontrol - Reparametrize — use uma parametrização não centralizada para modelos hierárquicos
- Restrinja as distribuições a priori se forem muito difusas
- Verifique os dados — valores discrepantes ou diferenças de escala causam problemas para o amostrador
# Re-run with higher adapt_delta to reduce divergences:
# fit2 <- stan(
# model_code = stan_code,
# data = stan_data,
# chains = 4,
# iter = 4000,
# warmup = 2000,
# control = list(adapt_delta = 0.95, max_treedepth = 12),
# seed = 42
# )Verificação rápida: limite de Rhat
Qual é o limite padrão atual de Rhat que indica que um modelo Stan convergiu?
Recapitulação da amostragem e dos diagnósticos MCMC
Fluxo de trabalho MCMC principal do RStan:
stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000)ajusta o modeloprint(fit)mostra Rhat e n_eff — os principais diagnósticos de convergência- Rhat < 1.01 e n_eff / total > 0.1 indicam uma amostra com bom comportamento
traceplot()— verificação visual da mistura;pairs()— revela problemas na geometria da distribuição a posterioriextract(fit, pars='mu')$mu— acessa as amostras brutas da distribuição a posteriorishinystan::launch_shinystan(fit)— diagnósticos interativos abrangentes
Perguntas Frequentes
A aula “Amostragem MCMC e diagnósticos” é grátis?
Sim — o texto completo de “Amostragem MCMC e diagnósticos” é grátis para ler aqui na web. Para praticá-la interativamente (um editor de código integrado e um tutor de IA 24/7) e desbloquear o restante do curso de R Academy, atualize para CoddyKit PRO. O curso de R Academy inclui 4 aulas no total.
O que vou aprender em “Amostragem MCMC e diagnósticos”?
Execute a amostragem, inspecione cadeias e interprete os diagnósticos Rhat e ESS. Você pratica R Academy com código prático que executa diretamente no navegador, e um tutor de IA 24/7 responde suas dúvidas enquanto trabalha na aula.
Preciso ter experiência prévia para começar R Academy?
Nenhuma experiência prévia é necessária. R Academy no CoddyKit é estruturado para alunos iniciantes até avançados, então você pode começar aqui ou desde o início e aprender no seu ritmo. Esta é a aula 3 de 4.
Quanto tempo leva a aula “Amostragem MCMC e diagnósticos”?
A maioria das aulas CoddyKit leva cerca de 5–10 minutos. Cada uma é compacta e interativa, então você faz progresso constante e retoma exatamente de onde parou entre web e app.
Posso escrever e executar código nesta aula de R Academy?
Sim. Cada aula de R Academy inclui um editor de código integrado, então você escreve e executa código real direto no navegador e recebe feedback de IA instantaneamente — nenhuma configuração local necessária.
Todas as aulas deste curso
- Introdução ao pensamento bayesiano
- Escrevendo modelos Stan em R
- Amostragem MCMC e diagnósticos
- Verificações preditivas a posteriori