Intervalles crédibles a posteriori

Concentrons-nous sur le paramètre de pente $b$, soit le taux de variation du poids en fonction de la taille. La moyenne a posteriori de $b$ traduit la tendance du modèle a posteriori pour la pente. À l'inverse, un intervalle crédible a posteriori fournit un intervalle de valeurs plausibles de la pente a posteriori, ce qui reflète l'incertitude a posteriori autour de $b$. Par exemple, l'intervalle crédible à 95 % pour $b$ s'étend du 2,5e au 97,5e quantile de la postérieure de $b$. Il y a donc 95 % de probabilité (a posteriori) que $b$ se trouve dans cet intervalle.

Vous utiliserez la sortie de simulation de RJAGS pour approcher les intervalles crédibles de $b$. La simulation RJAGS de 100 000 itérations de la postérieure, weight_sim_big, est dans votre espace de travail ainsi qu'un cadre de données des sorties de la chaîne de Markov, weight_chains.

Cette activité fait partie du cours

Modélisation bayésienne avec RJAGS

Voir le cours

Instructions de l’exercice

Obtenez les statistiques de summary() pour les chaînes de weight_sim_big.
Les quantiles a posteriori 2.5% et 97.5% pour $b$ sont présentés au tableau 2 du summary(). Appliquez quantile() aux valeurs brutes de weight_chains pour vérifier ces calculs. Enregistrez le tout dans ci_95 et affichez-le.
De la même façon, utilisez les données de weight_chains pour construire un intervalle crédible à 90 % pour $b$. Enregistrez-le dans ci_90 et affichez-le.
Construisez un graphique de densité des valeurs de la chaîne de Markov de $b$. Superposez des lignes verticales représentant l'intervalle crédible à 90 % pour $bà l'aide degeomvline()avecxintercept = ci90`.

Exercice interactif pratique

Essayez cet exercice en complétant ce code d’exemple.

# Summarize the posterior Markov chains


# Calculate the 95% posterior credible interval for b
ci_95 <- quantile(___, probs = c(___, ___))
ci_95

# Calculate the 90% posterior credible interval for b
ci_90 <- ___
ci_90

# Mark the 90% credible interval 
ggplot(___, aes(x = ___)) + 
    geom_density() + 
    geom_vline(xintercept = ___, color = "red")

Modifier et exécuter le code