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.).
##
## Attaching package: 'AccumulatR'
## The following object is masked from 'package:stats':
##
## simulate
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.
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.
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.4663272
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.
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.
