The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.

Working with Mixtures

Mixture components describe different sets of active accumulators. Two choices determine how a mixture is used:

library(AccumulatR)

A common architecture

Both sections use the same fast-versus-slow architecture. Response R1 comes from a target route and response R2 comes from a competitor. The components differ only in whether the target route is fast or slow on a given trial.

Fixed mixture weights with observed components

In this first example, the mixture weights are fixed and the component label is observed on each trial.

model_fixed <- race_spec() |>
  add_accumulator("target_fast", "lognormal") |>
  add_accumulator("target_slow", "lognormal") |>
  add_accumulator("competitor", "lognormal") |>
  add_pool("TARGET", c("target_fast", "target_slow")) |>
  add_outcome("R1", "TARGET") |>
  add_outcome("R2", "competitor") |>
  add_component("fast", members = c("target_fast", "competitor")) |>
  add_component("slow", members = c("target_slow", "competitor")) |>
  set_parameters(separate = list(m = TRUE, s = TRUE)) |>
  set_mixture(mode = "fixed", weights = c(fast = 0.4, slow = 0.6)) |>
  finalize_model()

true_params_fixed <- c(
  target_fast.m = log(0.25),
  target_fast.s = 0.15,
  target_slow.m = log(0.45),
  target_slow.s = 0.20,
  competitor.m = log(0.35),
  competitor.s = 0.18
)

We now create a trial table that says which component each trial belongs to.

set.seed(123456)

n_trials_fixed <- 1200
trial_component <- sample(
  c("fast", "slow"),
  n_trials_fixed,
  replace = TRUE,
  prob = c(0.4, 0.6)
)
trial_conditions <- data.frame(
  trials = rep(seq_len(n_trials_fixed), each = 3L),
  racer = rep(c("target_fast", "target_slow", "competitor"), n_trials_fixed),
  component = rep(trial_component, each = 3L),
  stringsAsFactors = FALSE
)

params_df_fixed <- build_param_matrix(
  model_fixed,
  true_params_fixed,
  n_trials = n_trials_fixed
)

sim_fixed <- simulate(
  model_fixed,
  params_df_fixed,
  trial_df = trial_conditions
)

data_fixed <- sim_fixed[, c("trials", "R", "rt", "component")]

table(data_fixed$component, data_fixed$R)
##       
##         R1  R2
##   fast 452  30
##   slow 131 587

Fixed mixtures include the selected component label in simulation output by default. For real data, retain a component column when trial type is observed. Omit it to evaluate the data as a latent mixture.

prepared_fixed <- prepare_data(model_fixed, data_fixed)
ctx_fixed <- make_context(model_fixed)

params_df_fixed_true <- build_param_matrix(
  model_fixed,
  true_params_fixed,
  n_trials = n_trials_fixed
)

ll_fixed <- as.numeric(log_likelihood(ctx_fixed, prepared_fixed, params_df_fixed_true))
ll_fixed
## [1] 1554.949

Sampled mixture weights with latent components

Now we use the same architecture as a latent mixture. The component is not observed per trial. Instead, the model samples a hidden component for each trial, with p.fast controlling the proportion of fast trials.

model_sampled <- race_spec() |>
  add_accumulator("target_fast", "lognormal") |>
  add_accumulator("target_slow", "lognormal") |>
  add_accumulator("competitor", "lognormal") |>
  add_pool("TARGET", c("target_fast", "target_slow")) |>
  add_outcome("R1", "TARGET") |>
  add_outcome("R2", "competitor") |>
  add_component("fast", members = c("target_fast", "competitor")) |>
  add_component("slow", members = c("target_slow", "competitor")) |>
  set_parameters(separate = list(m = TRUE, s = TRUE)) |>
  set_mixture(mode = "sample", reference = "slow") |>
  finalize_model()

true_params_sampled <- c(
  target_fast.m = log(0.25),
  target_fast.s = 0.15,
  target_slow.m = log(0.45),
  target_slow.s = 0.20,
  competitor.m = log(0.35),
  competitor.s = 0.18,
  p.fast = 0.35
)

By default, sampled mixtures do not expose the component label in the output, because that label is latent.

set.seed(123456)

n_trials_sampled <- 1500
params_df_sampled <- build_param_matrix(
  model_sampled,
  true_params_sampled,
  n_trials = n_trials_sampled
)

sim_sampled <- simulate(model_sampled, params_df_sampled)
names(sim_sampled)
## [1] "trials" "R"      "rt"

If you want to inspect the hidden generating component, you can ask for it explicitly.

sim_sampled_with_component <- simulate(
  model_sampled,
  params_df_sampled,
  keep_component = TRUE
)

table(sim_sampled_with_component$component)
## 
## fast slow 
##  523  977

For likelihood evaluation, we keep the sampled-mixture data in the same form as real latent-mixture data: only the observed response R and response time rt. Because component is omitted here, the likelihood marginalizes over the possible components. If you do supply a component column, that component is treated as observed on those trials.

data_sampled <- sim_sampled[, c("trials", "R", "rt")]

prepared_sampled <- prepare_data(model_sampled, data_sampled)
ctx_sampled <- make_context(model_sampled)

params_df_sampled_true <- build_param_matrix(
  model_sampled,
  true_params_sampled,
  n_trials = n_trials_sampled
)

ll_sampled_true <- as.numeric(log_likelihood(ctx_sampled, prepared_sampled, params_df_sampled_true))
ll_sampled_true
## [1] 1177.892

Estimate a sampled-mixture weight

Latent weights can also be estimated. To keep the example simple, we estimate only the latent mixture weight p.fast, while holding the accumulator timing parameters fixed.

neg_loglik <- function(theta) {
  est <- true_params_sampled
  est["p.fast"] <- plogis(theta[["logit_p_fast"]])
  params_df <- build_param_matrix(
    model_sampled,
    est,
    n_trials = n_trials_sampled
  )
  ll <- log_likelihood(ctx_sampled, prepared_sampled, params_df)
  -as.numeric(ll)
}

start <- c(logit_p_fast = qlogis(0.60))

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

fit_params <- c(p.fast = plogis(fit$par[["logit_p_fast"]]))
target <- c(p.fast = true_params_sampled[["p.fast"]])

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

Use mode = "fixed" when mixture weights are known and mode = "sample" when mixture weights should be estimated. Independently of that choice, include a component column only when the component is observed on a trial. Otherwise omit it, or set it to NA, so the likelihood marginalizes over components.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.