model-diagnostics

# Model Diagnostics

## Key Convergence Metrics

| Metric | Good Value | Concern |
|--------|------------|---------|
| Rhat | < 1.01 | > 1.1 indicates non-convergence |
| ESS bulk | > 400 | < 100 unreliable estimates |
| ESS tail | > 400 | < 100 unreliable intervals |
| Divergences | 0 | Any indicates geometry issues |
| Max treedepth | 0 hits | Hitting limit = slow exploration |

## Stan Diagnostics (cmdstanr)

```r
library(cmdstanr)

fit <- mod$sample(data = stan_data, ...)

# Quick check
fit$cmdstan_diagnose()

# Summary with diagnostics
fit$summary()

# Detailed diagnostics
fit$diagnostic_summary()

# Extract specific metrics
draws <- fit$draws()
rhat <- posterior::rhat(draws)
ess_bulk <- posterior::ess_bulk(draws)
ess_tail <- posterior::ess_tail(draws)

# Divergences
np <- fit$sampler_diagnostics()
sum(np[,,"divergent__"])

# Treedepth
sum(np[,,"treedepth__"] == 10)  # Default max
```

## JAGS Diagnostics (R2jags)

```r
library(R2jags)
library(coda)

fit <- jags(...)

# Summary (includes Rhat, n.eff)
print(fit)
fit$BUGSoutput$summary

# Rhat
max(fit$BUGSoutput$summary[,"Rhat"])

# Effective sample size
min(fit$BUGSoutput$summary[,"n.eff"])

# Convert to coda
mcmc_obj <- as.mcmc(fit)

# Gelman-Rubin
gelman.diag(mcmc_obj)

# Autocorrelation
autocorr.diag(mcmc_obj)
autocorr.plot(mcmc_obj)

# Geweke diagnostic
geweke.diag(mcmc_obj)
```

## Visual Diagnostics

### Trace Plots
```r
# Stan (bayesplot)
library(bayesplot)
mcmc_trace(fit$draws(), pars = c("mu", "sigma"))

# JAGS
traceplot(fit)
```

### Rank Histograms
```r
# Should be uniform if chains mixed well
mcmc_rank_hist(fit$draws(), pars = "mu")
```

### Pairs Plot (Detect Correlations)
```r
mcmc_pairs(fit$draws(), pars = c("mu", "sigma", "tau"))
```

## Divergence Diagnosis (Stan)

```r
# Identify divergent transitions
np <- nuts_params(fit)
divergent <- np[np$Parameter == "divergent__" & np$Value == 1, ]

# Pairs plot highlighting divergences
mcmc_pairs(fit$draws(), np = np,
           pars = c("mu", "tau"),