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.).
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.
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
)We generate data and keep both ordered responses for each trial.
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)## trials R rt R2 rt2
## 1 1 A 0.3394130 B 0.3514512
## 2 2 A 0.2373401 B 0.4565748
## 3 3 A 0.3709420 B 0.4913625
## 4 4 B 0.2974072 A 0.3441526
## 5 5 A 0.3297834 B 0.4790856
## 6 6 A 0.3671582 B 0.4663272optim()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.
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.