5 Adding a new model

If you have read Section 3 and Section 4, you should have a pretty good idea of how the bmm package functions. Now it’s time to add your new model.

You don’t have to create the files by hand. use_model_template() generates the model file with stubs for every method a model needs, so you can focus on filling them in. The chapter ends with a checklist of everything a new model has to touch outside its own file, because the test suite fails on some of those and review catches the rest.

What `use_model_template()` does

Run it from the root of the package. It writes R/model_<model_name>.R containing the parameter specification, the model constructor, the user-facing alias, and stub methods for check_data, bmf2bf, configure_model and postprocess_brm. With custom_family = TRUE it also creates one empty inst/stan_chunks/<model_name>_<block>.stan file per entry in stanvar_blocks, wires the code that reads them into configure_model, and adds empty log_lik, posterior_predict and posterior_epred functions that the family already points to.

The full argument list is on the reference page for use_model_template(). Three things about the call are worth knowing:

versions is the second positional argument. use_model_template("gcm", TRUE) means versions = TRUE, not custom_family = TRUE. Name the arguments.
The stan files are created before the R file is written. If the call fails partway through, empty .stan files are already on disk and you have to remove them by hand before you can call the function again.
An unversioned model gets version = "NA" — the literal string, not NA.

What the template leaves for you

The generated file follows the package’s conventions. Two things in it are placeholders that check() will not flag:

The prediction functions are empty. With custom_family = TRUE, the family passes log_lik_<model>(), posterior_predict_<model>() and posterior_epred_<model>() to brms, but their bodies are empty, so loo(), pp_check() and posterior_epred() cannot work on your model until you fill them (step 8 below).
The settable_links stub is commented out. Uncomment it, so that its tag line reads #' @exportS3Method, if configure_model() builds links into the family or the formulas; otherwise delete it (step 2 below).

There are no stubs for check_formula(), configure_prior(), create_initfun(), revert_check_data() or the optional post-fit methods; add those only if your model needs them.

5.1 Example

Let’s add a new model called gcm. Assume you have tested the model in Stan and have the Stan code ready. We want a custom family with two Stan blocks, likelihood and functions (see ?brms::stanvar for the blocks).

First set up your system and git environment as described in Section 1. Then load the development version with devtools::load_all() and run, from the root directory of the bmm package:

use_model_template("gcm", custom_family = TRUE, stanvar_blocks = c("likelihood", "functions"))

This creates R/model_gcm.R and the empty files inst/stan_chunks/gcm_likelihood.stan and inst/stan_chunks/gcm_functions.stan, and opens them in your editor. The top of the R file holds the specification and the constructor:

.gcm_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
  ),
  # pick the sd rate by meaning: 1 for sensitivity, strength, concentration
  # (kappa), mixing weights and identity-linked drift; 2 for criteria,
  # thresholds, boundary, ndt, start point, log-linked drift, diffusion
  # constant, log ratios, log SDs and correlations; 4 for circular bias
  priors = list(
    par1 = list(main = "normal(0, 1)", effects = "normal(0, 0.5)", sd = "exponential(1)"),
    par2 = list(main = "normal(0, 0.5)", effects = "normal(0, 0.5)", sd = "exponential(2)")
  ),
  init_ranges = list(
    par1 = c(-1, 1),
    par2 = c(0.5, 1.5)
  )
)


.model_gcm <- 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 = character(),
      version = "NA",
      requirements = "",
      parameters = .gcm_defaults[["parameters"]],
      links = .gcm_defaults[["links"]],
      fixed_parameters = .gcm_defaults[["fixed_parameters"]],
      default_priors = .gcm_defaults[["priors"]],
      init_ranges = .gcm_defaults[["init_ranges"]]
    ),
    class = c("bmmodel", "gcm"),
    call = call
  )
  out <- set_links(out, links)
  out
}

# uncomment if configure_model() builds the links into the family or into
# the non-linear formulas rather than reading them from the list above,
# so that a link set by the user is refused instead of silently ignored:
# #' @exportS3Method
# settable_links.gcm <- function(model) character(0)

Above the specification, a header comment lists the registries outside the file that a new model must update and the format citation must follow. One line of that header is out of date in bmm 1.4.1: the stored-frame test case no longer goes into stored_frame_cases() in test-update.R but into your model’s test-case file (see Testing). The other entries, data_column_roles, revert_check_data(), response_annotations() and the citation format, are correct. Below it come the user-facing alias gcm() with its roxygen block, and the stubs check_data.gcm(), bmf2bf.gcm(), configure_model.gcm() (with a custom_family() skeleton and the lines that read your two Stan files), log_lik_gcm(), posterior_predict_gcm(), posterior_epred_gcm() and postprocess_brm.gcm(). The check_data, bmf2bf and postprocess_brm stubs carry a header comment saying when you may delete them.

For a complete model written to the current conventions, read R/model_sdt_yn.R: a custom family on aggregated counts, with its own log_lik, posterior_predict and posterior_epred functions and a settable_links method. Most of its length is documentation.

5.2 Filling in the template

The specification and the constructor. Most of your work is in .gcm_defaults: the parameters (with descriptions), their links, which parameters are fixed and to what value, the default priors, and the initialization ranges. Set default priors (main, effects and sd) for every parameter, informed by knowledge in the field; .ddm_defaults in R/model_ddm.R is a worked example. Give init_ranges only to parameters whose likelihood is undefined in part of their range, such as a non-decision time that must stay below the fastest RT (see Section A.3). Prefer parameter names of more than one letter that do not match a likely data column (see Failure modes to expect).

In .model_gcm(), rename the response and required arguments and fill in the metadata (domain, task, name, citation, requirements). citation holds one reference per element, each on a single line, without a - bullet, and ending in . or a https://doi.org/ URL; model_citation() hands them to users for their reference lists. Two more things about the constructor are checked by tests: it must run with no arguments (the tests call every .model_*() bare, which is why every argument defaults to NULL), and its fields must stay in the order the template emits them. For domain, reuse the string of the models your model sits with, so that bmm_models() and the website list it under the right task group: "Visual working memory", "Working Memory (categorical), Categorical Decision Making", "Categorical decision making, memory, and reasoning" (MPT), "Perception & Recognition Memory" (or "Recognition Memory", which lands in the same group) or "Decision Making / Response times". A new string becomes a group of its own.

For a model with several versions, generate the scaffold with use_model_template("gcm", versions = c("v1", "v2")). You get a .gcm_version_table with one entry per version, a constructor that reads each field from it and appends paste0("gcm_", version) to the class, and an alias with version = c("v1", "v2") validated by match.arg(). Keep the default as a literal c(...): model_versions() reads it to build the per-version documentation tables.
Links. The constructor passes the user’s links through set_links(), which refuses names that are not parameters, repairs one-letter typos with a warning, and warns when a link widens the parameter’s range (say identity instead of log). By default every link in your links list can be changed by the user. If configure_model() writes a link into the family, the Stan code or the non-linear formulas instead of reading it from the list, a user’s choice would be silently ignored, so declare which links can be set:
```
#' @exportS3Method
settable_links.gcm <- function(model) character(0)   # or c("par1")
```
Name the method for your model class, not for bmmodel, which comes first in every class vector and would shadow it. A model that substitutes the inverse link into its formulas itself also overrides settable_link_functions(); the multinomial SDT models are the examples.
The user-facing alias. Rename the required arguments and fill in @examples with a full example. Keep the generated @keywords bmmodel: the pkgdown test needs it and is skipped inside check() (see Chapter 6). Title and details are generated from the constructor.
check_data.gcm() checks the data and returns it, ending with NextMethod("check_data"). It may compute transformations or store values as data attributes for configure_model (see Section 3.6). If it adds, renames or packs columns, update() needs to undo that before it can re-run it on a fitted model’s data; see Section 3.7.
bmf2bf.gcm(), if needed. You need a custom method if your model has several response variables (DDM: rt | dec(response) ~ 1), uses aggregate data through vreal(), vint() or trials() (EZDM, sdt_yn), or needs other response wrapping. brms reads this formula as the linear model for the family’s mu. If your model has no natural mu, suppress it: put mu = 0 in fixed_parameters, add a dummy mu (identity link) to the family’s dpars, and write the first formula line yourself (see Section A.2). Leave the dummy mu out of parameters, so that prior_info() and parameter_info() do not report it.
configure_model.gcm() builds the formula, the family and the stanvars and returns nlist(formula, data, stanvars). It does not return initial values. brms connects the family to your Stan code by name: a custom family "gcm" with dpars = c("mu", "par1", "par2") calls a Stan function gcm_lpdf(y | mu, par1, par2) (gcm_lpmf for an integer response) that you define in gcm_functions.stan. The DDM method in Chapter 4 shows the whole pattern, and the brms custom families vignette explains it. Add further Stan files to inst/stan_chunks/ and read them here. If your Stan code has to differ when the user asks for within-chain threading, see Section A.5.
postprocess_brm.gcm() changes the fitted object (renaming parameters and the like). Delete the section if you need none. If you keep it, also write revert_postprocess_brm.gcm(), which undoes the change: update() calls it to get the fit back into the shape brms built before refitting. R/model_sdm.R has both.
log_lik, posterior_predict and posterior_epred, for a custom family. These are plain R functions passed to brms::custom_family(); the template generates them empty and already passes them to the family, and their contracts are in Chapter 6. They power LOO and WAIC, posterior predictive checks, and posterior_epred()/fitted(). Build them on d*() and r*() functions for your model in R/distributions.R, which you write first (see Distribution functions). If the mean of your response is not a useful quantity, pass posterior_epred_undefined("gcm"), as SDM does, so the user gets an error that points to native_parameters(). Models built on a native brms family (brms::mixture(), brms::multinomial()) get all three from brms.

5.3 Testing

A new model adds two files to tests/testthat/:

test-model_gcm.R tests your model’s own contract: constructor validation, check_data(), configure_model(), the content of the generated Stan code, and your d*/r* functions. tests/testthat/test-model_ddm.R and test-model_sdt_yn.R are examples.

helper-model-gcm.R registers your model with the tests that run over every model (default priors, initial values, update(), links, posterior predictive checks). It defines one function, model_cases_gcm(), which returns small cases those tests use. A stored_frame, an sd_default, an init and a links case are required; check_data, observables and fixture cases apply only to some models. Copy tests/testthat/helper-model-ddm.R and adapt it. For a model without random-effect quirks or extra observables, it looks like this:

model_cases_gcm <- function() {
  model <- gcm("resp", "exemplar", "category")
  data <- toy_gcm_data()  # write a small dataset next to those in helper-model-data.R
  re_formula <- bmf(par1 ~ 1 + (1 | id), par2 ~ 1 + (1 | id))
  list(
    stored_frame = list(gcm = list(model = model, formula = bmf(par1 ~ 1, par2 ~ 1), data = data)),
    sd_default = list(gcm = list(
      model = model, formula = re_formula, data = data,
      sd = c(par1 = "exponential(1)", par2 = "exponential(2)")
    )),
    init = list(gcm = list(model = model, formula = re_formula, data = data)),
    links = list(gcm = list(model = model, settable = c("par1", "par2")))
  )
}

sd repeats your default sd priors and settable what settable_links() returns. The header of helper-model-cases.R describes the optional cases.

Tests must not sample, so pipeline tests use brms’ mock backend. backend = "mock" is a brms::brm() argument, and whatever you pass as mock_fit is returned in place of the fitted Stan object. Everything before the fit (data checks, formula, family, Stan code, priors, initial values) runs for real, which is what you want to test:

test_that("gcm model runs without errors", {
  withr::local_options("bmm.silent" = 2)
  skip_on_cran()

  dat <- data.frame(
    resp = rnorm(10),
    exemplar = rep(c("A", "B"), 5),
    category = rep(c(1, 2), each = 5)
  )

  f <- bmmformula(par1 ~ 1, par2 ~ 1)
  mock_fit <- bmm(
    f, dat,
    gcm(resp_var1 = "resp", required_arg1 = "exemplar", required_arg2 = "category"),
    backend = "mock", mock_fit = 1, rename = FALSE
  )
  expect_equal(mock_fit$fit, 1)
})

Every argument of the alias without a default has to be supplied, because the alias calls stop_missing_args(). Older test files write mock = 1, which works only because R partial-matches it to mock_fit.

The rest of the rules are numbered R1 to R10 in bmm’s AGENTS.md, with the reasoning in ADR 0004. The ones a new model meets:

Seed with withr::local_seed(), never set.seed().
A test that compiles Stan code opens with skip_unless_stan_tests(). A test that takes more than 0.4 s with NOT_CRAN unset on your machine opens with skip_on_cran(), so that the package’s check stays within CRAN’s time limit.
The suite has a time budget. Check your tests with devtools::test(reporter = "slow").
tests/testthat/test-conventions.R checks several of these rules and fails with the rule’s number.

tests/testthat/helper-epred.R has epred_vs_predict() and expect_epred_by_cell() for testing a posterior_epred function against posterior_predict. If your Stan code changes with threading, test the generated code with and without it (Section A.5).

5.3.1 Parameter recovery

A new model needs a parameter recovery study on simulated data. The pull request description shows true against recovered values, as a table or a plot, from a hierarchical design with several subjects whose parameters vary, on data simulated with an independent generator rather than your own r* function. Attach the script as a gist or in a collapsed block. Chapter 6 explains why the generator matters, and CONTRIBUTING.md has the requirements.

5.3.2 Distribution functions

A custom family also needs distribution functions in R/distributions.R: d* and r*, optionally p* and q*. They simulate data, back log_lik, posterior_predict and posterior_epred, and drive the plots. Write them before the prediction functions.

5.4 Before you open the pull request

Everything below lives outside your model file.

The test suite fails without these:

The two test files above, with all four required cases in your model_cases_gcm().
An entry in data_column_roles (R/helpers-model.R) describing in a few words what the column named by each alias argument without a default holds. Use <model>_<version> keys where versions differ, and NA for an argument that names no column. The Get started article builds its model tables from it.
A revert_check_data() method if your check_data() adds, renames or packs columns (Section 3.7). The stored_frame case is what tests it.
@keywords bmmodel on the alias, and a constructor that runs without arguments and keeps the template’s field order.
At least one reference in citation, in the format of step 1.

Review expects these:

A branch in response_annotations() (R/helpers-model.R) giving the coding each response variable needs (“radians in [-pi, pi]”, “seconds”). print() shows it next to each response.
A domain that places the model in the right task group (step 1).
A parameter recovery study, an example dataset, an article, and a NEWS.md bullet.

Optional extension points. The defaults cover most models; write a method only when its condition applies.

Method	Write it when	Example
`pp_observables()`, `pp_simulate()`	the likelihood has more than one observable (RT and choice), so `pp_check()` should check each. Add an `observables` case to your test cases.	`R/model_ddm.R`
`native_transform()`	the map from sampling to native scale is not elementwise `links` plus at most one softmax group	`native_transform.non_targets()`
`data_check_findings()`	users commonly make a data mistake that `check_data()` tolerates; `bmm_data_check()` reports it	`data_check_findings.circular()`
`create_initfun()`, `configure_prior()`	`init_ranges` and the default priors cannot express what the model needs	`create_initfun.m3()`

A model fitted to counts or summary statistics has a few more things to set up (Section A.4).

5.5 Failure modes to expect

configure_prior() errors when fixed_parameters names a parameter the formula never gets. Every fixed parameter has to be wired into the formula by configure_model(). The usual cause is mu = 0 without a dummy mu in the family (Section A.2).
Method signatures must match the generics. create_initfun() receives the prior as its fourth argument, and configure_prior() methods receive data2, knots or drop_unused_levels through ... when the user passes them; a method without ... fails only in that case. R CMD check reports both.
Short parameter names collide with Stan parameter names and data columns. Prefer names of more than one letter (Section A.8).

5.6 Post-processing

Most post-fit tools read the model object and need no code from you. native_parameters(), parameter_info() and prior_info() read the specification; conditional_effects() and emmeans() work through the bmmfit class; posterior_epred() and fitted() call your family’s posterior_epred function; pp_check() checks the response, and further observables if you declared them; bmm_data_check() runs your check_data() without compiling.

5.7 Add an example dataset

Every new model needs an example dataset that users can load and that the examples use: a new dataset added to the package, one already in it, or one from another package loaded with data(). The IMM and the SDT models use datasets shipped with bmm; the mixture2p and mixture3p articles use data('bays2009_full', package = 'mixtur'). R Packages explains how to add one.

5.8 Add an article

Every new model needs an article explaining the model and how to fit it. Create it with usethis::use_article(); the existing ones are in vignettes/articles/, and the published versions are here. Articles load their fits from cached files, which are not committed; CONTRIBUTING.md says how to hand over the fit. The README, the home page and the Get started tables list your model automatically, from its constructor and its data_column_roles entry.

What use_model_template() does