R Academy · Lección

Muestreo MCMC y diagnósticos

Ejecute el muestreo, inspeccione las cadenas e interprete los diagnósticos Rhat y ESS.

Muestreo MCMC y diagnósticos es una lección gratuita de R Academy en CoddyKit. Esta es la lección 3 de 4. Puedes leer la lección completa abajo gratuitamente — luego la practicas en el navegador con un editor de código integrado y un tutor de IA 24/7. Forma parte de la ruta de aprendizaje de R Academy, y tu progreso se sincroniza en la web y la app de CoddyKit. El curso de R Academy incluye 4 lecciones en total.

¿Qué es MCMC?

Monte Carlo mediante cadenas de Markov (MCMC) es una familia de algoritmos para extraer muestras de una distribución de probabilidad cuando el muestreo directo es imposible. En estadística bayesiana, MCMC extrae muestras de la distribución posterior P(parameters | data).

RStan implementa el muestreador No-U-Turn (NUTS), un algoritmo MCMC de vanguardia.

Un modelo sencillo de Stan

Un modelo de Stan es un bloque de texto que define tipos de datos, parámetros y la distribución posterior logarítmica. El modelo más sencillo estima la media de una distribución normal con varianza conocida.

# 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
# }
# '

Llamar a stan() para muestrear

stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000) compila el modelo (una vez) y después extrae muestras. Con iter=2000 y warmup=1000, cada cadena produce 1000 muestras posteriores al calentamiento; en total, son 4000 muestras entre las 4 cadenas.

# 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): la tabla de resumen

print(fit) muestra una tabla de resumen para cada parámetro con la media posterior, la desviación estándar, los cuantiles, Rhat y n_eff. Estos dos diagnósticos son los primeros que debe comprobar.

# 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    1

Criterio de convergencia Rhat

Rhat (factor de reducción de la escala potencial) compara la varianza dentro de las cadenas con la varianza entre cadenas. Los valores cercanos a 1.0 indican que todas las cadenas han convergido a la misma distribución.

Rhat < 1.01: convergencia (estándar actual)
Rhat > 1.01: las cadenas no se han mezclado; ejecute más iteraciones
Rhat > 1.1: problema grave de convergencia

# 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: tamaño de muestra efectivo

n_eff (tamaño de muestra efectivo) tiene en cuenta la autocorrelación entre muestras MCMC sucesivas. Las muestras correlacionadas contienen menos información que las independientes.

n_eff cercano al total de iteraciones: muestras casi independientes, excelente
n_eff / total_samples > 0.1: generalmente aceptable
n_eff muy bajo: autocorrelación alta; considere reparametrizar el 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 la convergencia visual

traceplot(fit, pars = 'mu') representa los valores muestreados de mu a lo largo de las iteraciones para cada cadena. Las cadenas convergentes tienen el aspecto de una oruga difusa: todas se superponen, sin tendencias ni desplazamientos. Las cadenas divergentes deambulan o permanecen 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 las correlaciones posteriores

pairs(fit, pars = c('mu', 'sigma')) muestra una matriz de diagramas de dispersión de las muestras posteriores. Revela las correlaciones entre parámetros y destaca las transiciones divergentes (representadas en rojo), que indican regiones con las que el muestreador tiene dificultades.

# 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!
')

Extraer muestras posteriores

extract(fit, pars = 'mu')$mu devuelve un vector numérico con todas las muestras posteriores al calentamiento de mu. Use estas muestras para calcular cualquier resumen posterior: media, intervalos de credibilidad o probabilidad de una condición.

# 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')

Iniciar ShinyStan para obtener diagnósticos interactivos

shinystan::launch_shinystan(fit) abre una aplicación interactiva de Shiny con gráficos de trazas, distribuciones posteriores, gráficos de pares y diagnósticos de NUTS en un mismo lugar. Es la herramienta más completa para explorar un ajuste de 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 chain

Soluciones habituales para la convergencia

Cuando Rhat > 1.01 o se observan divergencias:

Aumente iter y warmup
Aumente adapt_delta hacia 1.0 (por ejemplo, a 0.95) en control
Reparametrice: use una parametrización no centrada para modelos jerárquicos
Establezca distribuciones previas más restrictivas si son demasiado difusas
Compruebe los datos: los valores atípicos o las diferencias de escala causan problemas al muestreador

# 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
# )

Comprobación rápida: umbral de Rhat

¿Cuál es el umbral estándar actual de Rhat que indica que un modelo de Stan ha convergido?

Repaso del muestreo y los diagnósticos MCMC

Flujo de trabajo MCMC principal de RStan:

stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000) ajusta el modelo
print(fit) muestra Rhat y n_eff, los principales diagnósticos de convergencia
Rhat < 1.01 y n_eff / total > 0.1 indican una muestra adecuada
traceplot(): comprobación visual de la mezcla; pairs(): revela problemas en la geometría posterior
extract(fit, pars='mu')$mu: permite acceder a las muestras posteriores sin procesar
shinystan::launch_shinystan(fit): diagnósticos interactivos completos

Preguntas frecuentes

¿La lección «Muestreo MCMC y diagnósticos» es gratis?

Sí — el texto completo de «Muestreo MCMC y diagnósticos» es gratis para leer aquí en la web. Para practicarla de forma interactiva (editor de código integrado y tutor de IA 24/7) y desbloquear el resto del curso de R Academy, actualiza a CoddyKit PRO. El curso de R Academy incluye 4 lecciones en total.

¿Qué aprenderé en «Muestreo MCMC y diagnósticos»?

Ejecute el muestreo, inspeccione las cadenas e interprete los diagnósticos Rhat y ESS. Practicas R Academy con código real que ejecutas directamente en el navegador, y un tutor de IA 24/7 responde tus preguntas mientras trabajas en la lección.

¿Necesito experiencia previa para empezar R Academy?

No se requiere experiencia previa. R Academy en CoddyKit está estructurado para principiantes hasta estudiantes avanzados, así que puedes empezar aquí o desde el inicio y avanzar a tu ritmo. Esta es la lección 3 de 4.

¿Cuánto tiempo toma la lección «Muestreo MCMC y diagnósticos»?

La mayoría de las lecciones de CoddyKit toman alrededor de 5–10 minutos. Cada una es compacta e interactiva, así que avanzas constantemente y retomas exactamente por donde dejaste en la web y la app.

¿Puedo escribir y ejecutar código en esta lección de R Academy?

Sí. Cada lección de R Academy incluye un editor de código integrado, así que escribes y ejecutas código real directamente en tu navegador y obtienes retroalimentación instantánea de IA — sin configuración local necesaria.

Todas las lecciones de este curso

← Volver a R Academy