Échantillonnage MCMC et diagnostics
Lancez l’échantillonnage, examinez les chaînes et interprétez les diagnostics Rhat et ESS.
Échantillonnage MCMC et diagnostics est une leçon R Academy gratuite sur CoddyKit. Ceci est la leçon 3 sur 4. Tu peux lire la leçon complète ci-dessous gratuitement — puis la pratiquer en direct dans le navigateur avec un éditeur de code intégré et un tuteur IA 24/7. Elle fait partie du parcours d'apprentissage R Academy, et ta progression se synchronise sur le web et l'application CoddyKit. Le cours R Academy comprend 4 leçons au total.
Qu’est-ce que le MCMC ?
La méthode de Monte-Carlo par chaînes de Markov (MCMC) regroupe des algorithmes permettant de tirer des échantillons d’une distribution de probabilités lorsque l’échantillonnage direct est impossible. En statistique bayésienne, le MCMC échantillonne à partir de la distribution a posteriori P(paramètres | données).
RStan implémente l’échantillonneur No-U-Turn (NUTS), un algorithme MCMC de pointe.
Un modèle Stan simple
Un modèle Stan est un bloc de texte qui définit les types de données, les paramètres et la log-probabilité a posteriori. Le modèle le plus simple estime la moyenne d’une distribution normale dont la variance est connue.
# 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
# }
# 'Appeler stan() pour échantillonner
stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000) compile le modèle (une seule fois), puis tire des échantillons. Avec iter=2000 et warmup=1000, chaque chaîne produit 1000 échantillons après la phase d’échauffement, soit 4000 au total pour les 4 chaînes.
# 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) — Le tableau récapitulatif
print(fit) affiche un tableau récapitulatif pour chaque paramètre, avec la moyenne a posteriori, l’écart type, les quantiles, Rhat et n_eff. Ces deux diagnostics sont les premiers éléments à vérifier.
# 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 1Le critère de convergence Rhat
Rhat (facteur potentiel de réduction d’échelle) compare la variance au sein des chaînes à la variance entre les chaînes. Des valeurs proches de 1,0 indiquent que toutes les chaînes ont convergé vers la même distribution.
- Rhat < 1.01 — convergence atteinte (norme actuelle)
- Rhat > 1.01 — les chaînes ne se sont pas mélangées ; exécutez davantage d’itérations
- Rhat > 1.1 — problème de convergence sérieux
# 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 — Taille d’échantillon effective
n_eff (taille d’échantillon effective) tient compte de l’autocorrélation entre les échantillons MCMC successifs. Les échantillons corrélés contiennent moins d’informations que les échantillons indépendants.
- n_eff proche du nombre total d’itérations — échantillons presque indépendants, excellent résultat
- n_eff / total_samples > 0.1 — généralement acceptable
- n_eff très faible — forte autocorrélation ; envisagez de reparamétrer le modèle
# 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() pour visualiser la convergence
traceplot(fit, pars = 'mu') représente les valeurs échantillonnées de mu au fil des itérations pour chaque chaîne. Des chaînes ayant convergé ressemblent à une chenille floue : elles se recouvrent toutes, sans tendance ni dérive. Des chaînes divergentes errent ou restent séparées.
# 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() pour les corrélations a posteriori
pairs(fit, pars = c('mu', 'sigma')) affiche une matrice de nuages de points des échantillons a posteriori. Elle révèle les corrélations entre les paramètres et met en évidence les transitions divergentes (représentées en rouge), qui indiquent les régions posant problème à l’échantillonneur.
# 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!
')Extraire les échantillons a posteriori
extract(fit, pars = 'mu')$mu renvoie un vecteur numérique contenant tous les échantillons de mu obtenus après la phase d’échauffement. Utilisez ces échantillons pour calculer n’importe quel résumé a posteriori : moyenne, intervalles de crédibilité ou probabilité d’une condition.
# 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')Lancer ShinyStan pour des diagnostics interactifs
shinystan::launch_shinystan(fit) ouvre une application Shiny interactive contenant les graphiques de traces, les distributions a posteriori, les graphiques pairs et les diagnostics NUTS au même endroit. C’est l’outil le plus complet pour explorer un ajustement 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 chainSolutions courantes aux problèmes de convergence
Lorsque Rhat > 1.01 ou que vous observez des divergences :
- Augmentez
iteretwarmup - Augmentez
adapt_deltavers 1,0 (par exemple, 0,95) danscontrol - Reparamétrez le modèle — utilisez une paramétrisation non centrée pour les modèles hiérarchiques
- Resserrez les lois a priori si elles sont trop diffuses
- Vérifiez les données — les valeurs aberrantes ou les différences d’échelle peuvent provoquer des problèmes pour l’échantillonneur
# 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
# )Vérification rapide : seuil de Rhat
Quel est le seuil actuellement admis pour Rhat indiquant qu’un modèle Stan a convergé ?
Récapitulatif de l’échantillonnage et des diagnostics MCMC
Flux de travail MCMC avec RStan :
stan(model_code=..., data=..., chains=4, iter=2000, warmup=1000)ajuste le modèleprint(fit)affiche Rhat et n_eff — les principaux diagnostics de convergence- Rhat < 1.01 et n_eff / total > 0.1 indiquent un échantillon de bonne qualité
traceplot()— vérification visuelle du mélange ;pairs()— révèle les problèmes de géométrie a posterioriextract(fit, pars='mu')$mu— accès aux échantillons a posteriori brutsshinystan::launch_shinystan(fit)— diagnostics interactifs complets
Questions Fréquemment Posées
La leçon « Échantillonnage MCMC et diagnostics » est-elle gratuite ?
Oui — le texte complet de « Échantillonnage MCMC et diagnostics » est gratuit à lire ici sur le web. Pour la pratiquer de manière interactive (un éditeur de code intégré et un tuteur IA 24/7) et déverrouiller le reste du cours R Academy, passe à CoddyKit PRO. Le cours R Academy comprend 4 leçons au total.
Qu'est-ce que j'apprendrai dans « Échantillonnage MCMC et diagnostics » ?
Lancez l’échantillonnage, examinez les chaînes et interprétez les diagnostics Rhat et ESS. Tu pratiques R Academy avec du code pratique que tu exécutes directement dans le navigateur, et un tuteur IA 24/7 répond à tes questions au fur et à mesure que tu avances dans la leçon.
Dois-je avoir de l'expérience pour commencer R Academy ?
Aucune expérience préalable n'est requise. R Academy sur CoddyKit est structuré pour les débutants jusqu'aux apprenants avancés, donc tu peux commencer ici ou depuis le début et avancer à ton rythme. Ceci est la leçon 3 sur 4.
Combien de temps prend la leçon « Échantillonnage MCMC et diagnostics » ?
La plupart des leçons CoddyKit prennent environ 5–10 minutes. Chacune est courte et interactive, tu progresses régulièrement et tu repiques exactement où tu t'es arrêté sur le web et l'app.
Peux-tu écrire et exécuter du code dans cette leçon R Academy ?
Oui. Chaque leçon R Academy inclut un éditeur de code intégré, tu écris et exécutes du vrai code directement dans ton navigateur et tu reçois des retours IA instantanés — aucune configuration locale requise.
Toutes les leçons de ce cours
- Introduction au raisonnement bayésien
- Écrire des modèles Stan dans R
- Échantillonnage MCMC et diagnostics
- Vérifications prédictives a posteriori