8.1 El buen rey Markov y su reino insular
El buen rey Markov y su reino insular proporcionan una metáfora intuitiva para MCMC. El rey visita islas con probabilidades proporcionales a su población, creando una cadena de Markov que, a largo plazo, refleja la distribución de interés.
8.1.1 El algoritmo de Metropolis
El rey Markov gobierna una cadena de 10 islas. Las poblaciones son diferentes en cada isla. El rey quiere pasar tiempo en cada isla proporcional a su población, pero no conoce la población total — solo conoce la población relativa de las islas vecinas.
El algoritmo de Metropolis resuelve esto:
- Proponer: Proponer aleatoriamente moverse a una isla vecina.
- Aceptar/Rechazar: Si la isla propuesta tiene mayor población, aceptar. Si es menor, aceptar con probabilidad igual al cociente de poblaciones.
- Repetir: Tras muchos pasos, la fracción de tiempo pasado en cada isla converge a las proporciones de población.
# Good King Markov: Metropolis algorithm on 10 islands
set.seed(42)
# Island populations (unnormalized)
populations <- c(5, 12, 3, 8, 20, 15, 7, 10, 4, 6)
n_islands <- length(populations)
# Metropolis algorithm
n_steps <- 10000
current <- 1 # start at island 1
visits <- numeric(n_steps)
for (t in 1:n_steps) {
# Propose a neighboring island
proposal <- current + sample(c(-1, 1), 1)
# Wrap around
if (proposal < 1) proposal <- n_islands
if (proposal > n_islands) proposal <- 1
# Accept or reject
accept_prob <- min(1, populations[proposal] / populations[current])
if (runif(1) < accept_prob) {
current <- proposal
}
visits[t] <- current
}
# Compare visited frequencies to true proportions
visited_freq <- table(visits) / n_steps
true_prop <- populations / sum(populations)
cat("Island Visited True\n")
for (i in 1:n_islands) {
cat(sprintf(" %2d %.3f %.3f\n", i, visited_freq[i], true_prop[i]))
}
# Plot
barplot(rbind(visited_freq, true_prop), beside = TRUE,
col = c("steelblue", "firebrick"),
main = "King Markov: visited frequency vs true proportion",
xlab = "Island", ylab = "Proportion")
legend("topright", legend = c("Visited", "True"), fill = c("steelblue", "firebrick"))8.1.2 Por qué funciona
El algoritmo de Metropolis funciona porque satisface el balance detallado:
donde es la distribución objetivo y es la probabilidad de transición. Esto asegura que la cadena sea reversible y converja a la distribución objetivo.
# Verify detailed balance for the Metropolis algorithm
# For two neighboring islands i and j:
# pi(i) * T(i->j) = pi(j) * T(j->i)
populations <- c(5, 12, 3, 8, 20, 15, 7, 10, 4, 6)
# Check islands 4 and 5
i <- 4; j <- 5
pi_i <- populations[i]
pi_j <- populations[j]
# T(i->j) = min(1, pi_j/pi_i) * 0.5 (0.5 = prob of proposing j from i)
# T(j->i) = min(1, pi_i/pi_j) * 0.5
T*ij <- min(1, pi_j / pi_i) * 0.5
T*ji <- min(1, pi_i / pi_j) * 0.5
cat("Detailed balance check (islands 4 and 5):\n")
cat(" pi(4) _ T(4->5) =", pi_i _ T*ij, "\n")
cat(" pi(5) * T(5->4) =", pi*j * T*ji, "\n")
cat(" Equal?", abs(pi_i * T*ij - pi_j * T_ji) < 1e-10, "\n")El rey no necesita conocer la población total. Lo notable del algoritmo de Metropolis es que solo necesita el cociente de poblaciones, no los valores absolutos. Por esto MCMC funciona para inferencia bayesiana: solo necesitamos el cociente de probabilidades a posteriori (verosimilitud por distribución a priori), no la constante de normalización. El algoritmo del rey es la base de todos los métodos MCMC, desde el simple Metropolis hasta el sofisticado Hamiltonian Monte Carlo.
Resumen de la Sección 8.1
- El algoritmo de Metropolis muestrea de una distribución proponiendo movimientos y aceptando/rechazando basándose en el cociente de probabilidades.
- Satisface el balance detallado, asegurando convergencia a la distribución objetivo.
- Solo necesita cocientes de probabilidades, no la constante de normalización.
- Esta es la base de todos los métodos MCMC.