R Academy · Lekcja

Posterior predictive checks

Weryfikuj dopasowanie modelu, porównując rozkłady danych symulowanych i obserwowanych.

Lekcja 4 z 413 kroki

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 model

ppc_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 misfit

Wię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ów
  • ppc_scatter_avg(y, yrep) — wykres rozrzutu danych zaobserwowanych względem średniej yrep
  • ppc_intervals(y, yrep) — przedziały niepewności wokół każdej obserwacji
  • ppc_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] matrix

Walidacja 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 reliable

Najlepsze 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-check

Szybkie 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()
Bezpłatny start

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

  1. Wprowadzenie do myślenia bayesowskiego
  2. Pisanie modeli Stan w R
  3. Próbkowanie MCMC i diagnostyka
  4. Posterior predictive checks
← Powrót do R Academy