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

## ----pywitness-availability, include=FALSE------------------------------------
library(r4lineups)

# The model fitting runs on an isolated pyWitness environment, which is not
# present on every machine and is never installed by package checks. Chunks that
# need the engine are guarded on this flag rather than being switched off
# outright, so that they run wherever the engine is available.
pywitness_ready <- FALSE
try({
  status <- check_pywitness_deps(initialize = TRUE)
  pywitness_ready <- isTRUE(status$installed) &&
                     isTRUE(status$audited_revision_matches)
}, silent = TRUE)

## ----engine-status------------------------------------------------------------
check_pywitness_deps(initialize = FALSE)

## ----not-available, echo=FALSE, results='asis', eval=!pywitness_ready---------
# cat("> **Note.** The pyWitness engine is not available in the environment where",
#     "this vignette was built, so the fitting output below is not shown. The code",
#     "is exactly what would run. Install the engine as above to reproduce it.\n")

## ----data---------------------------------------------------------------------
condition_a <- as.data.frame(simulate_lineup_data(
  n_tp        = 500,
  n_ta        = 500,
  d_prime     = 1.6,
  c_criterion = 0.9,
  lineup_size = 6,
  conf_levels = 9,
  foil_bias   = c(0.4, 0.2, 0.0, -0.4, -1.0),
  seed        = 404
))

# Put confidence on the 0-100 scale used elsewhere in the package
condition_a$confidence <- round(seq(20, 100, length.out = 9))[condition_a$confidence]

str(condition_a[c("target_present", "identification", "confidence")])

## ----response_table-----------------------------------------------------------
knitr::kable(
  as.data.frame.matrix(table(condition_a$target_present,
                             condition_a$identification)),
  caption = "Response counts by lineup type, simulated data (500 target-present and 500 target-absent lineups of size 6). Rows: FALSE = target-absent, TRUE = target-present."
)

## ----fitting, eval=pywitness_ready--------------------------------------------
fits <- fit_lineup_models(
  condition_a,
  models          = c("independent", "ensemble", "integration"),
  lineup_size     = 6,
  confidence_bins = c(0, 60, 80, 100),
  variance        = "equal",
  shared_variance = "estimated",
  control         = list(maxiter = 5000L, method = "Nelder-Mead")
)

knitr::kable(
  as.data.frame(fits$comparison)[, c("model", "converged", "n_parameters",
                                     "pearson_chisq", "df", "p_value")],
  digits = 4, row.names = FALSE,
  caption = "Model fits, simulated data generated by the MAX (Independent Observations) rule."
)

## ----starts, eval=pywitness_ready---------------------------------------------
fit_started <- fit_lineup_models(
  condition_a,
  models = c("independent", "ensemble"),
  confidence_bins = c(0, 60, 80, 100),
  starts = list(
    independent = c(targetMean = 1.5, c1 = 0.8, c2 = 1.4, c3 = 2.0),
    ensemble    = c(targetMean = 1.5, c1 = 0.7, c2 = 1.2, c3 = 1.7)
  )
)

fit_started$comparison$converged

## ----recovery, eval=pywitness_ready-------------------------------------------
simulate_under <- function(rule, seed = 404) {
  d <- as.data.frame(simulate_lineup_data(
    n_tp = 500, n_ta = 500, d_prime = 1.6, c_criterion = 0.9,
    lineup_size = 6, conf_levels = 9, decision_rule = rule,
    foil_bias = c(0.4, 0.2, 0.0, -0.4, -1.0), seed = seed
  ))
  d$confidence <- round(seq(20, 100, length.out = 9))[d$confidence]
  d
}

recovery <- do.call(rbind, lapply(c("max", "ensemble", "integration"), function(rule) {
  f <- fit_lineup_models(simulate_under(rule),
                         models = c("independent", "ensemble", "integration"),
                         lineup_size = 6, confidence_bins = c(0, 60, 80, 100))
  cc <- as.data.frame(f$comparison)
  get <- function(m, col) cc[[col]][cc$model == m]
  data.frame(
    `Generated by`     = rule,
    `Independent chi2` = round(get("independent", "pearson_chisq"), 1),
    `Independent p`    = signif(get("independent", "p_value"), 3),
    `Ensemble chi2`    = round(get("ensemble", "pearson_chisq"), 1),
    `Ensemble p`       = signif(get("ensemble", "p_value"), 3),
    `Integration chi2` = round(get("integration", "pearson_chisq"), 1),
    `Integration p`    = signif(get("integration", "p_value"), 3),
    check.names = FALSE
  )
}))

knitr::kable(
  recovery, row.names = FALSE,
  caption = "Pearson chi-square and its p-value by fitted model, for data generated under each decision rule. Simulated, 500 target-present and 500 target-absent lineups per row."
)

## ----results, eval=pywitness_ready--------------------------------------------
knitr::kable(
  as.data.frame(fits$parameters), digits = 4, row.names = FALSE,
  caption = "Parameter estimates, simulated data."
)

## ----diagnostics, eval=pywitness_ready----------------------------------------
fits$diagnostics
fits$engine

## ----plots, eval=pywitness_ready, fig.width=7, fig.height=5-------------------
plot(fits, type = "comparison")

## ----plot_fit, eval=pywitness_ready, fig.width=7, fig.height=5----------------
plot(fits, type = "fit")

## ----plot_roc, eval=pywitness_ready, fig.width=7, fig.height=5----------------
plot(fits, type = "roc")