## -----------------------------------------------------------------------------
library(AccumulatR)

## -----------------------------------------------------------------------------
model <- race_spec() |>
  add_accumulator("A1", "lognormal") |>
  add_accumulator("A2", "lognormal") |>
  add_accumulator("A3", "lognormal") |>
  add_accumulator("B", "lognormal") |>
  add_pool("A_pool", c("A1", "A2", "A3"), k = 2L) |>
  add_outcome("A", "A_pool") |>
  add_outcome("B", "B") |>
  set_parameters(separate = list(m = TRUE, s = TRUE)) |>
  finalize_model()

true_params <- c(
  A1.m = log(0.28),
  A1.s = 0.16,
  A2.m = log(0.28),
  A2.s = 0.16,
  A3.m = log(0.28),
  A3.s = 0.16,
  B.m = log(0.28),
  B.s = 0.18
)

## -----------------------------------------------------------------------------
set.seed(123456)

n_trials <- 1500
params_df <- build_param_matrix(model, true_params, n_trials = n_trials)

sim <- simulate(model, params_df)

data_df <- data.frame(
  trials = sim$trials,
  R = factor(sim$R),
  rt = sim$rt,
  stringsAsFactors = FALSE
)

table(data_df$R)

## ----eval=FALSE---------------------------------------------------------------
# prepared <- prepare_data(model, data_df)
# ctx <- make_context(model)
# neg_loglik <- function(theta) {
#   est <- true_params
#   est[c("A1.m", "A2.m", "A3.m")] <- theta[["A_pool.m"]]
#   est[c("A1.s", "A2.s", "A3.s")] <- exp(theta[["log_A_pool.s"]])
#   est["B.m"] <- theta[["B.m"]]
#   est["B.s"] <- exp(theta[["log_B.s"]])
#   params_df <- build_param_matrix(
#     model,
#     est,
#     n_trials = n_trials
#   )
#   ll <- log_likelihood(ctx, prepared, params_df)
#   -as.numeric(ll)
# }
# 
# start <- c(
#   A_pool.m = log(0.24),
#   log_A_pool.s = log(0.12),
#   B.m = log(0.34),
#   log_B.s = log(0.12)
# )
# 
# fit <- optim(
#   start,
#   neg_loglik,
#   method = "Nelder-Mead",
#   control = list(maxit = 4000, reltol = 1e-9)
# )
# 
# fit_params <- c(
#   A_pool.m = fit$par[["A_pool.m"]],
#   A_pool.s = exp(fit$par[["log_A_pool.s"]]),
#   B.m = fit$par[["B.m"]],
#   B.s = exp(fit$par[["log_B.s"]])
# )
# target <- c(
#   A_pool.m = true_params[["A1.m"]],
#   A_pool.s = true_params[["A1.s"]],
#   B.m = true_params[["B.m"]],
#   B.s = true_params[["B.s"]]
# )
# 
# data.frame(true = target, recovered = fit_params, miss = abs(target - fit_params))