---
title: "Multiple Outcomes"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Multiple Outcomes}
  %\VignetteEncoding{UTF-8}
  %\VignetteEngine{knitr::rmarkdown}
editor_options: 
  markdown: 
    wrap: 72
---

This vignette shows how to work with data in which more than one ordered
response is observed on each trial. Here each trial records the first
response (`R`, `rt`) and the second response (`R2`, `rt2`), and so on
for any number of responses (`R3`, `rt3` etc.).

```{r}
library(AccumulatR)
```

## Define the model

We use a simple two-accumulator race with direct
responses `A` and `B`. Setting `n_outcomes = 2` tells the model to
retain the first and second finishing responses.

```{r}
model <- race_spec(n_outcomes = 2L) |>
  add_accumulator("A", "lognormal") |>
  add_accumulator("B", "lognormal") |>
  add_outcome("A", "A") |>
  add_outcome("B", "B") |>
  set_parameters(separate = list(m = TRUE, s = TRUE)) |>
  finalize_model()

true_params <- c(
  A.m = log(0.30),
  A.s = 0.18,
  B.m = log(0.38),
  B.s = 0.22
)
```

## Simulate data

We generate data and keep both ordered responses for
each trial.

```{r}
set.seed(123456)

n_trials <- 500
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,
  R2 = factor(sim$R2),
  rt2 = sim$rt2,
  stringsAsFactors = FALSE
)

head(data_df)
```

## Estimate parameters with `optim()`

We estimate `A.m`, `A.s`, `B.m`,
and `B.s`. The lognormal standard deviations are optimized on the log scale and
transformed back inside the function.

```{r eval=FALSE}
prepared <- prepare_data(model, data_df)
ctx <- make_context(model)
neg_loglik <- function(theta) {
  est <- true_params
  est[c("A.m", "A.s", "B.m", "B.s")] <- theta[c("A.m", "A.s", "B.m", "B.s")]
  est[c("A.s", "B.s")] <- exp(est[c("A.s", "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.m = log(0.24),
  A.s = log(0.12),
  B.m = log(0.48),
  B.s = log(0.12)
)

fit <- optim(start, neg_loglik, method = "Nelder-Mead")

fit_params <- fit$par
fit_params[c("A.s", "B.s")] <- exp(fit_params[c("A.s", "B.s")])
target <- true_params[c("A.m", "A.s", "B.m", "B.s")]

data.frame(
  true = target,
  recovered = fit_params,
  miss = abs(target - fit_params)
)
```

Ranked observations require direct accumulator or pool outcomes. Chained
onsets are supported, but logical outcome expressions, guessing, remapping,
censoring, and truncation are not. Observed response times must increase
strictly and response labels must be distinct. Later unobserved ranks use
`NA` in both columns; the first response must be observed.
