3 BMM code structure

Adding a new model is straightforward using the use_model_template() function, which will be described in the next section. You do not have to edit any of the files below, but it will be helpful to understand the structure of the package.

3.1 The main workhorse - `bmm()`

The main function for fitting models is bmm(). This function is the main entry point for users to fit models. It is set-up to be independent of the specific models that are implemented in the package.

R/bmm.R

bmm <- function(formula, data, model,
                prior = NULL,
                sort_data = getOption("bmm.sort_data", "check"),
                silent = getOption("bmm.silent", 1),
                backend = getOption("brms.backend", NULL),
                file = NULL, file_compress = TRUE,
                file_refit = getOption("bmm.file_refit", FALSE), ...) {
  deprecated_args(...)
  dots <- list(...)
  local_brms_threads(dots)

  # check if the model has been previously fit and return it if requested
  x <- try_read_bmmfit(file, file_refit)
  if (!is.null(x)) {
    return(x)
  }

  # set temporary global options and return modified arguments for brms
  configure_opts <- nlist(
    sort_data, silent, backend,
    parallel = dots$parallel,
    cores = dots$cores
  )
  opts <- configure_options(configure_opts)
  dots$parallel <- NULL

  # check model, formula and data, and transform data if necessary
  user_formula <- formula
  model <- check_model(model, data, formula)
  data <- check_data(model, data, formula)
  formula <- check_formula(model, data, formula)

  # generate the model specification to pass to brms later
  config_args <- configure_model(model, data, formula)

  # configure the default prior and combine with user-specified prior
  prior <- configure_prior(model, data, config_args$formula, prior)

  # configure initial values if necessary
  config_args$init <- create_initfun(model, data, config_args$formula)

  # estimate the model
  fit_args <- combine_args(nlist(config_args, opts, dots, prior))
  fit <- brms::do_call(brms::brm, fit_args)

  # model post-processing
  fit <- postprocess_brm(
    model, fit,
    fit_args = fit_args,
    user_formula = user_formula,
    configure_opts = configure_opts
  )

  # save the fitted model object if file argument provided and return the object
  try_save_bmmfit(fit, file, compress = file_compress)
}

It calls several subroutines, implemented as generic S3 methods, to:

configure_options() - configure local options for fitting, such as parallel sampling
check_model() - check if the model exists and is valid
check_data() - check whether the data contains all necessary information and transform it
check_formula() - check if the formula is specified correctly
configure_model(model, data, formula) - configure the model specification (formula, family, Stan code) to pass to brms
configure_prior() - set the default priors for the model and combine them with the user prior
create_initfun() - create an initial values function for MCMC sampling
postprocess_brm() - post-process the fitted model (e.g., adjust link functions, add metadata)

In addition, it also tests if the specified bmmodel has already been estimated and saved to a file. This is done via the try_read_bmmfit and try_save_bmmfit functions.

The first line worth a second look is local_brms_threads(dots). threads is a brms::brm() argument and reaches it through ... like any other, so bmm() does not need to know about it — except that a model may need to generate different Stan code when threading is on. local_brms_threads() (R/utils.R:128) normalizes whatever the user passed (a bare number becomes brms::threading(n)) and sets the brms.threads option for the duration of the call, so that a configure_model() method further down the pipeline can read it. Only sdm does this today; see Within-chain threading.

3.2 The `bmmformula` class

Before understanding models, it’s important to understand how users specify what to estimate. The bmmformula (alias: bmf()) is a list of formulas, one per model parameter. Unlike brmsformula, there is no response variable on the left-hand side — the response is specified when constructing the model object.

# Each formula predicts one model parameter
f <- bmmformula(
  kappa ~ 0 + set_size + (1 | id),
  c ~ 1 + (1 | id),
  a ~ 1
)

# Equivalent shorthand
f <- bmf(kappa ~ 0 + set_size + (1 | id), c ~ 1 + (1 | id), a ~ 1)

Parameters can be fixed to a constant using = instead of ~:

f <- bmf(kappa ~ condition, a = 0.5)  # 'a' is fixed to 0.5, not estimated

Parameters not mentioned in the formula are automatically given an intercept-only formula (param ~ 1) and a message is printed.

The bmmformula is converted to a brmsformula during model configuration via the bmf2bf() S3 method. This conversion adds the response variable as the first formula line and translates each parameter formula into either brms::lf() (linear) or brms::nlf() (non-linear) components.

3.3 Models

All models in the package are defined as S3 classes and follow a strict template. This allows us to implement general methods for handling model fitting, data checking, and post-processing. Each model has an internal constructor (.model_<name>()) that defines the model and its parameters, and a user-facing alias. To keep the specification easy to read and edit, the parameters, links, fixed parameters, priors, and initialization ranges live in a separate .<name>_defaults list — or, for models with several versions, in a .<name>_version_table (see Version tables) — and the constructor spells out each field of the model object inline, reading from that list. For a complete example model file and an explanation, see Section 4. The general model template looks like this:

# the model's parameter specification (one such list per version in a version table)
.my_new_model_defaults <- list(
  parameters = list(
    par1 = "Parameter 1 = description of parameter 1",
    par2 = "Parameter 2 = description of parameter 2"
  ),
  links = list(par1 = "identity", par2 = "log"),
  fixed_parameters = list(mu = 0),
  priors = list(
    par1 = list(main = "normal(0, 1)", effects = "normal(0, 0.5)"),
    par2 = list(main = "normal(0, 0.5)", effects = "normal(0, 0.5)")
  ),
  init_ranges = list(par1 = c(-1, 1), par2 = c(0.5, 1.5))
)

.model_my_new_model <- function(resp_var1 = NULL, required_arg1 = NULL,
                                required_arg2 = NULL, links = NULL,
                                call = NULL, ...) {
  out <- structure(
    list(
      resp_vars = nlist(resp_var1),
      other_vars = nlist(required_arg1, required_arg2),
      domain = "",
      task = "",
      name = "",
      citation = "",
      version = "NA",
      requirements = "",
      parameters = .my_new_model_defaults[["parameters"]],
      links = .my_new_model_defaults[["links"]],
      fixed_parameters = .my_new_model_defaults[["fixed_parameters"]],
      default_priors = .my_new_model_defaults[["priors"]],
      init_ranges = .my_new_model_defaults[["init_ranges"]]
    ),
    class = c("bmmodel", "my_new_model"),
    call = call
  )
  out$links[names(links)] <- links
  out
}

The .my_new_model_defaults list holds the parameter specification (one such list per version in a version table). Its fields are:

parameters - Named list of model parameters with descriptions
links - Named list of link functions for each parameter. The available options are "identity", "log", "softplus", "log1p", "logm1", "inverse", "sqrt", "logit", "probit", "tan_half" and "cloglog" (R/helpers-parameters.R:76). "softplus", added in 1.3.2 (#363), is an opt-in alternative to "log" for positively-bounded parameters: it keeps the parameter positive but grows linearly rather than exponentially, so predictor effects are additive on the natural scale and large values do not blow up. Users can override any of these per parameter through the model alias’s links argument
fixed_parameters - Named list of parameters fixed to constant values (e.g., list(mu = 0))
priors - Named list with main (intercept) and effects (slopes) priors per parameter. The constructor stores this under the model object’s default_priors field
init_ranges - Named list of c(lower, upper) vectors for MCMC initialization (see Section 3.9)

The constructor copies these into the model object (renaming priors to default_priors) alongside the response and auxiliary variables (resp_vars, other_vars), the model metadata (domain, task, name, citation, requirements), and the version. Keeping the fields inline means a new field on the model object is a one-line addition you can see and copy directly; the trade-off is that older models need updating when a field is added (a periodic, mechanical refactor).

Each model is accompanied by a user-facing alias, the documentation of which is generated automatically based on the info list in the model definition.

# user facing alias
# information in the title and details sections will be filled in
# automatically based on the information in the .model_modelname()$info
#' @title `r .model_my_new_model()$name`
#' @name my_new_model
#' @details `r model_info(.model_my_new_model())`
#' @param resp_var1 A description of the response variable
#' @param required_arg1 A description of the required argument
#' @param required_arg2 A description of the required argument
#' @param ... used internally for testing, ignore it
#' @return An object of class `bmmodel`
#' @export
#' @examples
#' \dontrun{
#' # put a full example here (see 'R/model_ddm.R' for an example)
#' }
my_new_model <- function(resp_var1, required_arg1, required_arg2,
                         links = NULL, ...) {
  call <- match.call()
  stop_missing_args()
  .model_my_new_model(resp_var1 = resp_var1, required_arg1 = required_arg1,
                      required_arg2 = required_arg2, links = links,
                      call = call, ...)
}

For a versioned model, the alias instead carries a version = c("v1", "v2") argument validated with version <- match.arg(version), and passes the chosen version on to the constructor (see cswald for a worked example). use_model_template() generates whichever variant you ask for via its versions argument.

Then users can fit the model using the bmm() function, and the model will be automatically recognized and handled by the package:

fit <- bmm(formula = my_bmmformula,
           data = my_data,
           model = my_new_model(resp_var1, required_arg1, required_arg2))

3.4 Model domains and class chains

Models use hierarchical S3 class chains that encode the model taxonomy. The class chain determines which S3 methods are dispatched. The general pattern is:

class = c("bmmodel", <domain/intermediate>, <model>, <version>)

The current model families and their class chains are:

Circular / Visual working memory models:

mixture2p: c("bmmodel", "circular", "mixture2p")
mixture3p: c("bmmodel", "circular", "non_targets", "mixture3p")
imm (full): c("bmmodel", "circular", "non_targets", "imm", "imm_full")
sdm: c("bmmodel", "circular", "sdm", "sdm_simple")

Response time / Decision models:

ddm: c("bmmodel", "ddm")
ezdm (3par): c("bmmodel", "ezdm", "ezdm_3par")
cswald (simple): c("bmmodel", "cswald", "cswald_simple")

Categorical response models:

m3 (ss): c("bmmodel", "m3", "m3_ss")

Intermediate classes like "circular" and "non_targets" allow shared check_data methods (e.g., check_data.circular validates that responses are in radians, check_data.non_targets validates non-target feature columns). S3 methods chain via NextMethod() from general to specific: dispatch walks the class vector left to right, so check_data.bmmodel runs first and the model-specific method runs last. See the worked dispatch chain in Chapter 4.

3.5 S3 methods

The package uses S3 methods to handle different models. This means that the same function can behave differently depending on the class of the object it is called with. For example, the configure_model(model, data, formula) function called by bmm(), is generally defined as:

R/helpers-model.R

configure_model <- function(model, data, formula) {
   UseMethod('configure_model')
}

and it will call a function configure_model.modelname() that is specified for each model. The same is true for other functions, such as check_data(), check_formula(), configure_prior(), create_initfun(), and postprocess_brm(). This allows us to add new models without having to edit the main fitting function, bmm().

The key S3 generics in the pipeline are:

Generic	Signature	Defined in
`check_model`	`(model, data = NULL, formula = NULL)`	`helpers-model.R`
`check_data`	`(model, data, formula)`	`helpers-data.R`
`check_formula`	`(model, data, formula)`	`bmmformula.R`
`configure_model`	`(model, data, formula)`	`helpers-model.R`
`configure_prior`	`(model, data, formula, user_prior, ...)`	`helpers-prior.R`
`create_initfun`	`(model, data, formula)`	`helpers-inits.R`
`postprocess_brm`	`(model, fit, ...)`	`helpers-postprocess.R`
`bmf2bf`	`(model, formula)`	`bmmformula.R`

Match these signatures exactly in your methods, including the ... where the generic has one and the = NULL defaults on check_model. R CMD check compares every method against its generic, and a method that drops an argument the generic declares is reported under “checking S3 generic/method consistency”. Every method in the package follows this: all four configure_prior.* methods carry ..., and both check_model.* methods carry the defaults.

The bmf2bf generic converts a bmmformula to a brmsformula. Most models use the default bmf2bf.bmmodel method, but some models require custom methods to handle response variable wrapping (e.g., bmf2bf.ddm wraps the response as rt | dec(response) ~ 1) or multiple response variables (e.g., bmf2bf.ezdm_3par packs variables via vreal(), vint(), trials()). Currently, model-specific bmf2bf methods exist for: ddm, cswald, ezdm_3par, ezdm_4par, and m3.

The configure_model method must return a named list via nlist(). The required and optional fields are:

formula (required) — the brmsformula object, with formula$family set to the appropriate family
data (required) — the data frame (possibly with additional columns added during configuration)
stanvars (optional) — custom Stan code, needed only for models with brms::custom_family(). Models using built-in brms families (e.g., brms::mixture() for the IMM) do not need stanvars.

Initial values are not returned by configure_model — they are handled by create_initfun() in a separate pipeline step.

3.6 Data flow: the attribute bridge pattern

A key architectural pattern is how information flows between check_data() and configure_model(). Since these are separate S3 methods, they cannot share local variables. Instead, check_data stores computed values as attributes on the data frame, and configure_model retrieves them.

# In check_data.non_targets (R/helpers-data.R):
attr(data, "max_set_size") <- max_set_size
attr(data, "lure_idx_vars") <- lure_idx_vars

# Later, in configure_model.imm_full (R/model_imm.R):
max_set_size <- attr(data, "max_set_size")
lure_idx <- attr(data, "lure_idx_vars")

This is the way to pass model-specific metadata between pipeline steps. If your check_data method computes something that configure_model needs (e.g., the maximum set size, index variables, or transformed columns), store it as a data attribute.

Two qualifications. First, the infrastructure uses the same channel for its own bookkeeping: check_data.bmmodel sets data_name and checked on the data frame (R/helpers-data.R:51-52), so do not assume the attributes you find there are all yours. Second, an attribute is not a guarantee. A configure_model method can be called directly, bypassing check_data entirely, and then the attribute is absent. configure_model.sdm handles this by recomputing rather than failing:

run_metadata <- attr(data, "sdm_run_metadata")
if (is.null(run_metadata)) {
  # Guard direct configure_model.sdm() calls that bypass check_data.sdm().
  run_metadata <- sdm_run_metadata(data, formula, model)
}

Whether to guard this way is a judgement call. It costs a little duplicated computation and it means the expensive part has to live in a function both methods can call, which is good structure anyway. If you do not guard, make the failure loud rather than silent.

3.7 Default priors structure

The default_priors field in model definitions has a specific structure that the configure_prior() method interprets:

default_priors = list(
  kappa = list(
    main = "student_t(5, 1.75, 0.75)",  # prior for Intercept
    effects = "normal(0, 1)",            # prior for regression coefficients
    sd = "exponential(1)"               # prior for random effect SDs (optional)
  ),
  c = list(
    main = "student_t(5, 2, 0.75)",
    effects = "normal(0, 1)"
  )
)

main: Applied to the Intercept of the parameter. If the user suppresses the intercept (e.g., kappa ~ 0 + condition), this prior is applied to each level of the factor instead.
effects: Applied to regression coefficients (slopes). If omitted, brms default priors are used for effects.
sd: Applied to all random effect standard deviations for this parameter. This is optional — if omitted, brms default priors apply.

The priors use the same string format as brms::set_prior(). These defaults are automatically applied by bmm() unless the user explicitly overrides them via the prior argument.

For parameters listed in fixed_parameters that also appear in parameters (e.g., mu1 in the IMM, which is fixed by default but can be estimated), the fixed value takes precedence unless the user provides a formula for that parameter.

3.8 Suppressing `mu` for models without a location parameter

The brms package requires every custom family to have a mu parameter (the location parameter). However, some models – particularly response time and decision models like DDM, EZDM, and cswald – do not have a natural location/mu parameter. These models suppress mu rather than arbitrarily promoting one of their parameters to it:

A dummy mu parameter is included in the family’s dpars with an identity link
The mu intercept is fixed to 0 by adding mu = 0 to fixed_parameters
Users do not specify a formula for mu in their bmmformula
The model’s bmf2bf method constructs the first formula line differently (e.g., rt | dec(response) ~ 1 for DDM, rather than using mu as a response)

This is purely a technical workaround to satisfy brms requirements; the fixed mu = 0 ensures it has no effect on estimation.

Note

Earlier versions tracked this with a void_mu flag in the model definition. The flag was never read by the infrastructure – suppression is handled entirely by fixed_parameters and the family’s dpars – and it has been removed. If you see void_mu in older code or notes, it is dead and can be deleted.

3.9 The `init_ranges` field and `create_initfun()`

Models with custom families often require carefully chosen initial values for stable MCMC sampling. When a model’s likelihood has hard boundaries or regions of zero density, the default brms initialization – which draws from a wide symmetric range – will often place chains in impossible territory, causing immediate divergences or chain failures.

This is particularly relevant for response time models, where parameters have hard constraints relative to the observed data. For example, the non-decision time parameter ndt must be positive and smaller than the fastest observed RT – otherwise the implied decision time would be negative and the likelihood undefined. Similarly, boundary separation (bound) must be strictly positive and within a plausible range, and extreme drift rate values can cause numerical overflow in the Wiener likelihood. Circular models (mixture2p, mixture3p) typically do not have these issues and work fine with default initialization.

The init_ranges field in the model definition addresses this by specifying narrow, safe windows on the native (pre-link-transform) scale:

init_ranges = list(
  mu    = c(-0.1, 0.1),
  drift = c(-1, 1),
  bound = c(1, 2),
  ndt   = c(0.01, 0.05)
)

Each entry maps a parameter name to a c(lower, upper) range. The create_initfun() S3 generic (defined in R/helpers-inits.R) uses these ranges to generate an initialization function that:

Parses the Stan code to determine parameter dimensions and types
Generates initial values for intercepts within the specified ranges (applying link transformations)
Generates small random values for regression coefficients, random effect SDs, and correlation matrices

If init_ranges is NULL (the default), create_initfun returns 1, which tells brms to use its default initialization. Most circular models (mixture2p, mixture3p) do not need custom initialization, but RT models (DDM, EZDM, cswald) and SDM do.

The conversion between native and sampling scales is handled by link_transform() (defined in R/helpers-parameters.R). It supports all link functions used by bmm models (identity, log, softplus, log1p, logm1, inverse, sqrt, logit, probit, tan_half, cloglog) and can apply transformations in either direction – from the scale on which parameters are interpretable (native) to the scale on which Stan samples (sampling), and back. Model developers rarely need to call it directly; it is used automatically by create_initfun() when generating initial values, and by the post-processing infrastructure when converting estimates back to the native scale.

3.10 Within-chain threading

Passing threads to bmm() asks brms to split each chain’s likelihood evaluation across cores. For most models nothing is needed: brms slices the data and calls the family’s likelihood inside partial_log_lik, and a custom family written the usual way comes along for the ride. A model that injects Stan code through stanvar() is the exception, because that code has to be written differently depending on whether it ends up inside partial_log_lik or in the model block.

SDM is the one model in 1.3.2 that does this; the NEWS.md entry reports it as reducing fitting time by up to ~45% in benchmarks (#374). It is worth reading if your model has a likelihood with a normalizing constant computed over groups of rows, because the same three pieces will apply.

1. Two versions of the chunk. inst/stan_chunks/ holds both sdm_simple_likelihood.stan and sdm_simple_likelihood_threaded.stan. The threaded one calls a _slice variant of the same function and passes start and end, the row bounds brms makes available inside partial_log_lik:

// serial
target += sdm_simple_run_ldenom(c, kappa, COSN, G_sdm_runs,
                                sdm_run_start, sdm_run_count);
// threaded
target += sdm_simple_run_ldenom_slice(c, kappa, COSN, start, end,
                                      G_sdm_runs, sdm_run_start,
                                      sdm_run_count);

2. A predicate that decides which one to read. configure_model.sdm selects the file at configuration time:

likelihood_file <- if (sdm_use_threaded_likelihood()) {
  "sdm_simple_likelihood_threaded.stan"
} else {
  "sdm_simple_likelihood.stan"
}

The predicate reads the brms.threads option that local_brms_threads() set at the top of bmm() (R/model_sdm.R:234-246). It has to normalize the option the same way brms does, because brms accepts a bare number where brms::threading() is expected; if it did not, brm() would thread while bmm emitted the serial chunk.

3. threading(force = TRUE) is a separate case. With force = TRUE, brms compiles with threading support but emits unsliced code — there is no partial_log_lik, so start and end do not exist and the sliced chunk fails to compile. The predicate therefore ends with && !isTRUE(threads$force), and SDM emits its serial likelihood in that case. This was a bug fix in 1.3.2; if you write a threaded chunk, test threading(2, force = TRUE) explicitly, because it is not the case you will think to try.

Data that the threaded function needs must also be declared for the partial-log-likelihood signature, via the pll_args argument to stanvar(). SDM passes "data matrix COSN" for its precomputed cosine matrix and, through the helper sdm_stanvar_int_array(), "data array[] int <name>" for the run-index arrays. A variable that is not declared in pll_args is simply not in scope inside partial_log_lik, and the model will not compile.

3.11 Version tables

Some models support multiple versions with different parameter sets. Two patterns are used:

If/else pattern (older, used by IMM): Version-specific parameters are selected via conditional logic within the model constructor.

Version table pattern (newer, used by EZDM, cswald and M3): Version-specific parameters, links, priors, and init_ranges are stored in a lookup table:

R/model_ezdm.R

.ezdm_version_table <- list(
  "3par" = list(
    parameters = list(...),
    links = list(...),
    fixed_parameters = list(...),
    priors = list(...),
    init_ranges = list(...)
  ),
  "4par" = list(
    parameters = list(...),
    links = list(...),
    ...
  )
)

The model constructor then indexes into this table for each field, selecting the entry for the requested version:

parameters = .ezdm_version_table[[version]][["parameters"]],
links = .ezdm_version_table[[version]][["links"]],
...

The table is the single source of truth for the version-specific specification. This pattern is cleaner and easier to extend than if/else logic, and new models with multiple versions should prefer it — use_model_template("mymodel", versions = c("v1", "v2")) generates it for you.

M3 nests the table one level deeper

The shape above is not universal, and M3 is the exception you will hit if you build a model whose parameterization depends on something other than the version. In .m3_version_table (R/model_m3.R:6-49), links and priors are indexed by the choice rule before the parameter:

links = list(
  simple  = list(c = "log", a = "log"),
  softmax = list(c = "identity", a = "identity")
),
priors = list(
  simple  = list(a = list(...), c = list(...)),
  softmax = list(a = list(...), c = list(...))
)

so the constructor indexes twice: .m3_version_table[[version]][["links"]][[choice_rule]] (R/model_m3.R:81-82). Three further departures in the same constructor are worth knowing before you copy it as a template:

fixed_parameters is not in the table at all. It is computed inline from the choice rule: list(b = if (choice_rule == "softmax") 0 else 0.1) (:76-78).
parameters is the table entry with a background-activation parameter b prepended in the constructor, because every M3 version has it (:73-75).
The table holds no init_ranges, and it has entries for ss and cs only. M3’s default version is "custom", which is not in the table — the indexing returns NULL and the user supplies the specification. That is deliberate, not an omission.

If your model needs a second axis of variation, this is what it looks like. If it does not, keep the single-level shape the template generates.

Not every model uses either pattern

Four of the eight released models — sdm, mixture2p, mixture3p and imm — carry no .{model}_defaults object and no version table. They write parameters, links, fixed_parameters, default_priors and init_ranges as literal lists inside the constructor (see R/model_sdm.R:5-44). The separated specification is the pattern for new models and the one use_model_template() generates; the inline models are simply older. Do not conclude the book is wrong when you open R/model_sdm.R and find no defaults list.

3.12 File organization

The bmm package is organized into several files. The main files are:

`R/bmm.R`

It contains the main function for fitting models, bmm(). This function is the main entry point for users to fit models. It is set-up to be independent of the specific models that are implemented in the package.

To add new models, you do not have to edit this file. The functions above are generic S3 methods, and they will automatically recognize new models if you add appropriate methods for them (see section Adding new models).

`R/helpers-*.R`

R/helpers-data.R, R/helpers-inits.R, R/helpers-model.R, R/helpers-parameters.R, R/helpers-postprocess.R, and R/helpers-prior.R

These files define the main generic S3 methods for checking data, creating initial value functions, configuring the model, extracting parameters, post-processing the fitted model, and combining priors. They contain the default methods for these functions, which are called by bmm() if no specific method is defined for a model. If you want to add a new model, you will need to add specific methods for these functions for your model. You do not need to edit these files to add a new model.

R/helpers-data.R also contains several user-facing utility functions for response time data preparation. These follow a deliberate design philosophy: each function performs one well-defined task without combining too many features. They are composable – users can pipe them together into a workflow but also inspect intermediate results at each step:

flag_contaminant_rts() – identifies likely contaminant trials by fitting a mixture of the RT distribution and a uniform contaminant process, returning per-trial posterior probabilities of contamination
validate_fast_guesses() – takes flagged contaminants and tests whether fast responses show chance-level accuracy consistent with guessing, using a Bayesian test
ezdm_summary_stats() – computes the summary statistics (mean RT, variance, accuracy) that the EZDM model requires as input, with options for robust estimation and contaminant handling

`R/bmmformula.R`

This file contains the definition of the bmmformula class, which is used to represent the formula for the model. It contains the bmmformula() function and its alias bmf(), which is used to create a new formula object.

In addition, it contains the definition of the bmf2bf S3 method that is used to convert a bmmformula object into a brmsformula object. This is necessary, as brmsformula objects are required to include the response variable in the first formula line. In contrast bmmformula objects only contain formulas predicting the parameters of a bmmodel. The bmf2bf S3 method is used to perform this conversion and add the first formula line including the response variable in the brmsformula created during model configuration. Some models override this with custom methods (e.g., bmf2bf.ddm, bmf2bf.ezdm_3par).

`R/model_*.R`

Each model and its methods are defined in a separate file. The current model files are:

model_mixture2p.R - Two-parameter mixture model
model_mixture3p.R - Three-parameter mixture model
model_imm.R - Interference Measurement Model (versions: full, bsc, abc)
model_sdm.R - Signal Discrimination Model
model_ddm.R - Diffusion Decision Model
model_ezdm.R - EZ-Diffusion Model (versions: 3par, 4par)
model_cswald.R - Censored-Shifted Wald Model (versions: simple, crisk)
model_m3.R - Multinomial Memory Measurement Model (versions: ss, cs, custom)

Each file contains the internal model constructor, the user-facing alias, and model-specific S3 methods. Your new model will exist in a file like this. The name of the file should be model_name_of_your_model.R. You don’t have to add this file manually - see section Adding new models.

`R/distributions.R`

This file contains the definition of the custom distributions used in the package. It specifies the density (d*), random number generation (r*), cumulative probability (p*), and quantile (q*) functions for each model’s distribution. Currently implemented distributions include: sdm, mixture2p, mixture3p, imm, m3, ddm, ezdm, and cswald. If your model requires a custom distribution, you will need to add it to this file. These are not used during model fitting, but are used to generate data from the model and for posterior predictive checks.

`R/conditional_effects.R`, `R/emmeans.R`, `R/pp_check.R`

Post-processing infrastructure for fitted bmmfit objects:

conditional_effects.R - conditional_effects() method for visualizing marginal effects
emmeans.R - Integration with the emmeans package for estimated marginal means
pp_check.R - Posterior predictive checking, including custom support for multinomial models

These work automatically for most models and do not need to be edited when adding a new model.

`R/utils.R`, `R/brms-misc.R`, `R/restructure.R`, `R/summary.R`, `R/update.R`

Various utility functions. Note that R/brms-misc.R is a read-only vendored copy from brms and should not be edited.

`inst/stan_chunks/`

This directory contains the Stan code chunks that are passed to the brms::stanvar() function. These define the custom likelihoods and helper functions for models that require custom brms families. Current Stan files:

sdm_simple_funs.stan, sdm_simple_tdata.stan, sdm_simple_likelihood.stan, sdm_simple_likelihood_threaded.stan - SDM model
ddm_functions.stan - DDM model (Wiener process likelihood)
ezdm_3par_functions.stan, ezdm_4par_functions.stan - EZDM model
cswald_simple_functions.stan, cswald_crisk_functions.stan, cswald_helper_functions.stan - cswald model

When writing Stan code for bmm, keep in mind that Stan’s HMC sampler requires continuous gradients to function correctly. Non-differentiable operations – such as fabs() or fmax() – create discontinuities in the log-density surface that can cause divergent transitions or break sampling entirely. Where such operations are mathematically necessary, bmm uses smooth approximations instead. For example, sqrt(x^2 + epsilon) serves as a differentiable stand-in for |x|. This is a standard Stan modelling technique and should be applied whenever your model’s likelihood involves non-smooth functions.

If you add a new model with a custom family, you will need to add Stan chunks to this directory. The naming convention is {model_name}_{block}.stan (e.g., my_model_functions.stan).

3.1 The main workhorse - bmm()

3.2 The bmmformula class

3.3 Models

3.4 Model domains and class chains

3.5 S3 methods

3.6 Data flow: the attribute bridge pattern

3.7 Default priors structure

3.8 Suppressing mu for models without a location parameter

3.9 The init_ranges field and create_initfun()

3.10 Within-chain threading

3.11 Version tables

M3 nests the table one level deeper

Not every model uses either pattern

3.12 File organization

R/bmm.R

R/helpers-*.R

R/bmmformula.R

R/model_*.R

R/distributions.R

R/conditional_effects.R, R/emmeans.R, R/pp_check.R

R/utils.R, R/brms-misc.R, R/restructure.R, R/summary.R, R/update.R

inst/stan_chunks/