## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment  = "#>",
  fig.width  = 7,
  fig.height = 4.5,
  dpi = 110
)

## ----message = FALSE, warning = FALSE-----------------------------------------
library(bayesqm)
library(ggplot2)

## ----setup-cmdstan, eval = FALSE----------------------------------------------
# install.packages("cmdstanr", repos = c("https://stan-dev.r-universe.dev", getOption("repos")))
# cmdstanr::check_cmdstan_toolchain(fix = TRUE)
# cmdstanr::install_cmdstan()
# cmdstanr::cmdstan_version()   # should print a version

## ----eval = FALSE-------------------------------------------------------------
# qdata <- read_qsort("mystudy.csv")   # CSV / Excel / HTMLQ / FlashQ
# qdata <- read_qsort("mystudy.DAT")   # PQMethod
# qdata <- read_qsort("mystudy.zip")   # KADE project

## -----------------------------------------------------------------------------
sim   <- generate_data(N = 33, J = 42, K = 3, seed = 1)
qdata <- qsort_data(sim$Y, distribution = sim$distribution)
qdata

## ----eval = FALSE-------------------------------------------------------------
# fit <- fit_bayesian(qdata, K = 3, chains = 4, iter = 2000,
#                     warmup = 1000, seed = 1)

## -----------------------------------------------------------------------------
fit <- demo_fit(N = 33, J = 42, K = 3, seed = 1)
fit

## -----------------------------------------------------------------------------
fit$diagnostics

## ----tucker, fig.cap = "MatchAlign alignment quality by factor."--------------
plot_tucker(fit)

## ----loadings, fig.cap = "Loadings with nested 50% and 95% credible intervals."----
plot_loading_posterior(fit)

## -----------------------------------------------------------------------------
head(compute_loadings(fit$Lambda_draws, prob = 0.95))

## ----zscores, fig.cap = "Posterior z-scores across factors with 50% and 95% credible intervals."----
plot(fit)

## ----zscore-one, fig.cap = "Posterior z-score for a single statement across factors."----
plot_zscore_posterior(fit, statement = 1)

## ----eval = FALSE-------------------------------------------------------------
# run <- run_bayes(qdata, K_max = 4, seed = 1,
#                  chains = 2, iter = 1500, warmup = 800)

## -----------------------------------------------------------------------------
run <- demo_run(K_max = 5, k_peak = 3, k_sivula = 1, case = "gap")
run

## ----elpd, fig.width = 7, fig.height = 4.2, fig.cap = "ΔELPD across K with the Sivula non-promotion band at |ΔELPD| < 4. Triangle marks the Sivula K, square marks the ELPD peak (the adopted K)."----
make_elpd_diff(run)

## -----------------------------------------------------------------------------
critical_delta(fit$Lambda_draws)
head(fit$qdc)

## ----dist-cons, fig.cap = "Posterior viewpoint divergence D_j with 95% credible intervals; statements ordered by P(D_j > delta | Y)."----
plot_dist_cons(fit)

## -----------------------------------------------------------------------------
head(classify_membership(fit$Lambda_draws))

## ----membership, fig.width = 6.5, fig.height = 5.5, fig.cap = "Blue tiles: posterior probability that each factor is dominant for a given participant (values printed in each cell). The right column shows the overall assignment verdict on an orange-red scale."----
make_dominant_panel(fit)

## ----ppc, fig.cap = "PPC on the by-person correlation matrix."----------------
plot_ppc(fit)

## ----hyper, fig.cap = "Posterior densities of nu, sigma, and tau."------------
plot_hyper(fit)

## ----eval = FALSE-------------------------------------------------------------
# save_bayesqm_plot("fig_loadings.pdf", plot_loading_posterior(fit))
# save_bayesqm_plot("fig_elpd.pdf",      plot_elpd(run),
#                   width = 3.5, height = 3)

## -----------------------------------------------------------------------------
caption_bayesqm(fit)

## ----eval = FALSE-------------------------------------------------------------
# bayesqm_set_colors("teal")   # "blue" (default), "red", "purple", "grey"
# plot(fit)
# bayesqm_set_colors("blue")

