## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")

## ----setup--------------------------------------------------------------------
library(intraclass)

## ----ci-bootstrap, eval = requireNamespace("glmmTMB", quietly = TRUE)---------
mc <- tidy(icc(ratings, score, subject, rater, seed = 1))
bs <- tidy(icc(ratings, score, subject, rater,
  ci_method = "bootstrap", boot_samples = 999, seed = 1
))
data.frame(
  term = mc$term,
  estimate = round(mc$estimate, 3),
  mc = sprintf("[%.2f, %.2f]", mc$conf.low, mc$conf.high),
  bootstrap = sprintf("[%.2f, %.2f]", bs$conf.low, bs$conf.high)
)

## ----ci-oneway-optin, eval = requireNamespace("glmmTMB", quietly = TRUE)------
mc <- tidy(icc(ratings, score, subject, rater, model = "oneway", seed = 1))
se <- tidy(icc(ratings, score, subject, rater,
  model = "oneway", ci_method = "searle"
))
bu <- tidy(icc(ratings, score, subject, rater,
  model = "oneway", ci_method = "burch"
))
np <- tidy(icc(ratings, score, subject, rater,
  model = "oneway", ci_method = "npbootstrap", boot_samples = 199, seed = 1
))
data.frame(
  term = mc$term,
  estimate = round(mc$estimate, 3),
  montecarlo = sprintf("[%.2f, %.2f]", mc$conf.low, mc$conf.high),
  searle = sprintf("[%.2f, %.2f]", se$conf.low, se$conf.high),
  burch = sprintf("[%.2f, %.2f]", bu$conf.low, bu$conf.high),
  npbootstrap = sprintf("[%.2f, %.2f]", np$conf.low, np$conf.high)
)

## ----ci-mpl, eval = requireNamespace("glmmTMB", quietly = TRUE)---------------
set.seed(88)
n_s <- 20
n_r <- 4
subj_eff <- rnorm(n_s, sd = sqrt(0.6))
rater_eff <- rnorm(n_r, sd = sqrt(0.1))
noise <- matrix(rnorm(n_s * n_r, sd = sqrt(0.2)), n_s, n_r)
sim <- data.frame(
  subject = factor(rep(seq_len(n_s), times = n_r)),
  rater = factor(rep(seq_len(n_r), each = n_s)),
  score = as.numeric(outer(subj_eff, rep(1, n_r)) +
    outer(rep(1, n_s), rater_eff) + noise)
)
mc2 <- tidy(icc(sim, score, subject, rater, type = "agreement", seed = 1))
ml <- tidy(icc(sim, score, subject, rater, type = "agreement", ci_method = "mpl"))
data.frame(
  term = mc2$term,
  estimate = round(mc2$estimate, 3),
  montecarlo = sprintf("[%.2f, %.2f]", mc2$conf.low, mc2$conf.high),
  mpl = sprintf("[%.2f, %.2f]", ml$conf.low, ml$conf.high)
)

## ----posterior, eval = FALSE--------------------------------------------------
# icc(ratings, score, subject, rater, engine = "brms", type = "agreement", seed = 1)

## ----posterior-hpdi, eval = FALSE---------------------------------------------
# icc(ratings, score, subject, rater, engine = "brms",
#   type = "agreement", posterior_summary = "hpdi", seed = 1)

