Posterior predictive checks
Weryfikuj dopasowanie modelu, porównując rozkłady danych symulowanych i obserwowanych.
Posterior predictive checks to bezpłatna lekcja R Academy na CoddyKit. To lekcja 4 z 4. Możesz przeczytać całą lekcję poniżej za darmo — a potem ćwiczyć ją interaktywnie w przeglądarce z wbudowanym edytorem kodu i tutorem AI dostępnym 24/7. To część ścieżki edukacyjnej R Academy, a Twój postęp synchronizuje się między webem a aplikacją CoddyKit. Kurs R Academy zawiera 4 lekcji w sumie.
Czym są predykcyjne sprawdzenia a posteriori?
Po dopasowaniu modelu bayesowskiego należy zadać pytanie: czy model generuje dane podobne do zaobserwowanych? Predykcyjne sprawdzenia a posteriori (PPC) odpowiadają na to pytanie, symulując powtórzone zbiory danych yrep z rozkładu a posteriori i porównując je graficznie z zaobserwowanymi danymi y.
Wyodrębnianie próbek a posteriori
extract(fit, pars='mu') zwraca nazwaną listę; $mu to wektor liczbowy wszystkich próbek tego parametru po fazie rozgrzewki. Przy 4 łańcuchach i 1000 iteracji po fazie rozgrzewki otrzymuje się 4000 próbek.
# library(rstan)
# mu_samples <- extract(fit, pars = 'mu')$mu
# sigma_samples <- extract(fit, pars = 'sigma')$sigma
#
# cat('Samples drawn:', length(mu_samples), '
')
# cat('Posterior mean of mu:', mean(mu_samples), '
')
# cat('90% CI:', quantile(mu_samples, c(0.05, 0.95)), '
')Generowanie yrep z próbek a posteriori
Dla każdego losowania a posteriori (mu_s, sigma_s) należy zasymulować powtórzony zbiór danych o takim samym rozmiarze jak oryginalne dane. Zapisz te zbiory w macierzy yrep, w której każdy wiersz reprezentuje jeden zasymulowany zbiór danych.
# y <- c(2.1, 1.8, 2.4, 1.9, 2.3, 2.0, 1.7, 2.2)
# n_obs <- length(y)
# S <- length(mu_samples) # 4000 posterior draws
#
# yrep <- matrix(NA, nrow = S, ncol = n_obs)
# for (s in seq_len(S)) {
# yrep[s, ] <- rnorm(n_obs, mean = mu_samples[s], sd = sigma_samples[s])
# }
# dim(yrep) # [4000, 8]ppc_dens_overlay() — porównanie gęstości
bayesplot::ppc_dens_overlay(y, yrep[1:50,]) nakłada gęstość jądrową zaobserwowanych danych (ciemna linia) na gęstości 50 losowo wybranych zasymulowanych zbiorów danych (jasne linie). Dobre dopasowanie modelu oznacza, że ciemna linia znajduje się w obszarze wyznaczonym przez chmurę jasnych linii.
# library(bayesplot)
#
# ppc_dens_overlay(y, yrep[1:50, ])
#
# Interpretation:
# - Dark line (y_obs) surrounded by light lines (yrep): good fit
# - Dark line systematically outside the cloud: model misfit
# - Light lines much wider than dark: overdispersed model
# - Light lines much narrower than dark: underdispersed modelppc_stat() — sprawdzanie statystyki testowej
ppc_stat(y, yrep, stat = 'mean') pokazuje histogram statystyki testowej (np. średniej) obliczonej dla każdego zasymulowanego zbioru danych, z pionową linią oznaczającą zaobserwowaną wartość statystyki. Jeśli zaobserwowana wartość znajduje się w głównej części histogramu, model poprawnie odwzorowuje ten aspekt danych.
# library(bayesplot)
#
# ppc_stat(y, yrep, stat = 'mean') # does model capture the mean?
# ppc_stat(y, yrep, stat = 'sd') # does model capture spread?
# ppc_stat(y, yrep, stat = 'max') # does model capture extremes?
#
# If observed stat is in the tail of the histogram,
# the model fails to reproduce that statistic.Wartość p bayesowska
Wartość p bayesowska (predykcyjna wartość p a posteriori) to odsetek zasymulowanych zbiorów danych, w których statystyka testowa jest bardziej ekstremalna niż wartość zaobserwowana. Wartości bliskie 0,5 wskazują na dobre skalibrowanie, a wartości bliskie 0 lub 1 — na niedopasowanie modelu dla danej statystyki.
# Bayesian p-value for the mean:
# obs_mean <- mean(y)
# rep_means <- apply(yrep, 1, mean)
# pval <- mean(rep_means >= obs_mean)
# cat('Bayesian p-value (mean):', round(pval, 3), '
')
# # 0.5 is perfect; < 0.05 or > 0.95 suggests misfitWięcej funkcji PPC w bayesplot
bayesplot oferuje wiele wizualizacji PPC wykraczających poza nakładanie gęstości:
ppc_hist(y, yrep[1:8,])— siatka histogramówppc_scatter_avg(y, yrep)— wykres rozrzutu danych zaobserwowanych względem średniej yrepppc_intervals(y, yrep)— przedziały niepewności wokół każdej obserwacjippc_rootogram(y, yrep)— dla danych zliczeniowych
# library(bayesplot)
#
# # Grid of 8 simulated histograms vs the observed
# ppc_hist(y, yrep[1:8, ])
#
# # Scatter: y_obs (x) vs mean of yrep (y) — should hug diagonal
# ppc_scatter_avg(y, yrep)
#
# # 50% and 90% posterior predictive intervals around each y_i
# ppc_intervals(y, yrep)Interpretacja wykresów PPC — niedopasowanie modelu
Typowe wzorce niedopasowania i ich przyczyny:
- yrep jest zbyt szeroki — rozkład a priori jest zbyt rozproszony lub model ma nadmierną dyspersję
- yrep jest przesunięty — niewłaściwa rodzina rozkładów wiarygodności (np. rozkład normalny dla danych skośnych)
- yrep nie odwzorowuje wielomodalności — potrzebny jest model mieszaninowy
- yrep nie odwzorowuje wartości ekstremalnych — potrzebny jest rozkład z ciężkimi ogonami
# Example: if data has a long right tail but yrep does not,
# consider switching:
# y ~ normal(mu, sigma) => y ~ student_t(nu, mu, sigma)
#
# Or for count data:
# y ~ poisson(lambda) => y ~ neg_binomial_2(mu, phi) (overdispersion)
cat('PPCs guide model improvement by revealing specific failure modes
')PPC w bloku modelu Stan
Można generować yrep bezpośrednio w Stanie, używając bloku generated quantities. Pozwala to uniknąć ponownego wyodrębniania parametrów w R i jest równoważne obliczeniowo.
# Stan model with generated quantities:
# '
# generated quantities {
# array[N] real y_rep;
# for (n in 1:N) {
# y_rep[n] = normal_rng(mu, sigma);
# }
# }
# '
# Then extract in R:
# yrep <- extract(fit, pars = 'y_rep')$y_rep # [S, N] matrixWalidacja krzyżowa leave-one-out
Oprócz PPC funkcja loo::loo(fit) oblicza walidację krzyżową leave-one-out, umożliwiając porównywanie konkurencyjnych modeli. Preferowany jest model o wyższym ELPD (oczekiwanej logarytmicznej gęstości predykcyjnej). Użyj loo::loo_compare(loo1, loo2), aby uszeregować modele.
# library(loo)
# loo1 <- loo(fit1) # normal model
# loo2 <- loo(fit2) # student-t model
#
# comparison <- loo_compare(loo1, loo2)
# print(comparison)
#
# Model with elpd_diff > 0 is preferred
# se_diff > |elpd_diff| means difference is not reliableNajlepsze praktyki dotyczące PPC
Przestrzegaj poniższych zasad, aby przeprowadzać rygorystyczne predykcyjne sprawdzanie a posteriori:
- Zawsze zaczynaj od
ppc_dens_overlay()jako ogólnej kontroli poprawności - Następnie używaj statystyk specyficznych dla danej dziedziny (
ppc_stat()), odpowiednich do celów analizy - Do kontroli wizualnych używaj co najmniej 50 losowań yrep, a do bayesowskich wartości p — wszystkich 4000
- Nieudane PPC wskazują kierunek ulepszania modelu — są narzędziem diagnostycznym, a nie oznaką porażki
# Workflow:
# 1. Fit model -> extract() -> generate yrep matrix
# 2. ppc_dens_overlay(y, yrep[1:50,]) -- visual global check
# 3. ppc_stat(y, yrep, stat='mean') -- check mean
# 4. ppc_stat(y, yrep, stat='sd') -- check spread
# 5. ppc_stat(y, yrep, stat='max') -- check tails
# 6. If misfit found -> revise model -> refit -> re-checkSzybkie sprawdzenie: interpretacja bayesowskiej wartości p
Co oznacza dla modelu bayesowska wartość p równa 0,03 dla statystyki maksimum?
Podsumowanie: predykcyjne sprawdzenia a posteriori
Przebieg pracy z PPC w RStan i bayesplot:
- Wyodrębnij próbki:
extract(fit, pars='mu')$mu - Wygeneruj yrep: wykonaj pętlę po losowaniach a posteriori, wywołując
rnorm(n, mu_s, sigma_s) - Kontrola ogólna:
ppc_dens_overlay(y, yrep[1:50,]) - Sprawdzanie statystyki:
ppc_stat(y, yrep, stat='mean') - Bayesowska wartość p:
mean(apply(yrep,1,stat) >= stat(y))— wartość bliska 0,5 jest dobra - Do porównywania modeli używaj
loo::loo_compare()
Ucz się R dzięki korepetycjom AI — za darmo
Pisz i uruchamiaj kod w przeglądarce, otrzymuj natychmiastową pomoc od korepetytora AI dostępnego 24/7 i kontynuuj naukę w sieci lub w aplikacji.
- Kursy
- 43
- Lekcje
- 159
Często zadawane pytania
Czy lekcja „Posterior predictive checks” jest bezpłatna?
Tak — pełny tekst „Posterior predictive checks” jest dostępny za darmo tutaj w sieci. Aby ćwiczyć ją interaktywnie (wbudowany edytor kodu i tutor AI dostępny 24/7) i odblokować resztę kursu R Academy, przejdź na CoddyKit PRO. Kurs R Academy zawiera 4 lekcji w sumie.
Co nauczysz się w „Posterior predictive checks”?
Weryfikuj dopasowanie modelu, porównując rozkłady danych symulowanych i obserwowanych. Ćwiczysz R Academy z praktycznym kodem, który uruchamiasz bezpośrednio w przeglądarce, a tutor AI dostępny 24/7 odpowiada na Twoje pytania podczas pracy nad lekcją.
Czy potrzebuję doświadczenia, aby zacząć R Academy?
Nie wymagamy żadnego doświadczenia. R Academy w CoddyKit jest strukturyzowany dla początkujących i zaawansowanych użytkowników, więc możesz zacząć tutaj lub od początku i uczyć się w swoim tempie. To lekcja 4 z 4.
Ile czasu zajmuje lekcja „Posterior predictive checks”?
Większość lekcji CoddyKit trwa około 5–10 minut. Każda lekcja to mały, interaktywny krok, dzięki czemu robisz systematyczne postępy i zawsze wracasz dokładnie do tego samego miejsca — na webie i w aplikacji.
Czy mogę pisać i uruchamiać kod w tej lekcji R Academy?
Tak. Każda lekcja R Academy zawiera wbudowany edytor kodu, więc piszesz i uruchamiasz prawdziwy kod bezpośrednio w przeglądarce i od razu otrzymujesz sprzężenie zwrotne od AI — bez konfiguracji na komputerze.
Wszystkie lekcje w tym kursie
- Wprowadzenie do myślenia bayesowskiego
- Pisanie modeli Stan w R
- Próbkowanie MCMC i diagnostyka
- Posterior predictive checks