# Bmm: R package for easy and flexible Bayesian Measurement Modeling (v1.0.0)

**URL:** <https://discourse.mc-stan.org/t/bmm-r-package-for-easy-and-flexible-bayesian-measurement-modeling-v1-0-0/35483>\
**Category:** Publicity\
**Tags:** hierarchical-model, r, cognitive-science, brms\
**Created:** [June 13, 2024, 8:46pm UTC](https://discourse.mc-stan.org/t/bmm-r-package-for-easy-and-flexible-bayesian-measurement-modeling-v1-0-0/35483 "2024-06-13T20:46:18Z")\
**Posts on this page:** 1\
**Showing post:** 3

<div class="post-metadata">

**Author:** ![mitzimorris](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/mitzimorris/32/29_2.png) [@mitzimorris](https://discourse.mc-stan.org/u/mitzimorris)\
**Post date:** [June 19, 2024, 10:44am UTC](https://discourse.mc-stan.org/t/bmm-r-package-for-easy-and-flexible-bayesian-measurement-modeling-v1-0-0/35483/3 "2024-06-19T10:44:13Z")

</div>

> [@Bob\_Carpenter](#):
>
> For example, @mitzimorris just extended brms to handle clinical measurement errors for diagnostic tests with differing sensitivity and specificity. This requires a custom link function and couldn’t be done directly in brms. And that’s a really really simple thing to do in Stan code.

here’s the brms code that adds a custom link function to a binomial model:

```no-highlight
library(brms)

# alldata is a dataframe of per-category binomial outcomes: number of tests, number of positive tests, categorical covariates sex, eth, age, time

alldata$spec = 0.95
alldata$sens = 0.7

hyper_data <- list(
  intercept_prior_mean = 0,
  intercept_prior_scale = 5,
  coef_prior_scale = 2.5
)

# define a *vectorized* custom family (no loop over observations)
binomial_sens_spec_vec <- custom_family(
"binomial_sens_spec", dpars = c("mu"),
links = c("logit"), lb = c(NA),
type="int", vars= c("trials", "vreal1", "vreal2"), loop = FALSE
)

# define the corresponding Stan density function
stan_density_vec <- "
real binomial_sens_spec_lpmf(int[] y, vector mu, int[] N, real[] vreal1, real[] vreal2) {
  return binomial_lpmf(y | N, mu * vreal1[1] + (1 - mu) * (1 - vreal2[1]));
}
"

stanvars_vec <- stanvar(scode = stan_density_vec, block = "functions")

priors <- c(
  set_prior(paste("normal(", hyper_data$intercept_prior_mean, ",", hyper_data$intercept_prior_scale, ")"), class = "Intercept"),
  set_prior(paste("normal(0,", hyper_data$coef_prior_scale, ")"), class = "b", coef = "sex"),
  set_prior(paste("normal(0,", hyper_data$coef_prior_scale, ")"), class = "sd", group = "eth"),
  set_prior(paste("normal(0,", hyper_data$coef_prior_scale, ")"), class = "sd", group = "age"),
  set_prior(paste("normal(0,", hyper_data$coef_prior_scale, ")"), class = "sd", group = "time")
)

fit = brm(cases | trials(tests) + vreal(sens, spec) ~ sex + (1 | eth) + (1 | age) + (1 | time),
          data = alldata,
          family = binomial_sens_spec_vec,
          prior = priors,
          stanvars = stanvars_vec)

```

More code was required to do `proj_pred` checks. The project team included @jonah, and it took a lot of digging on both our parts to get things running.

---

_[View the full topic](https://discourse.mc-stan.org/t/bmm-r-package-for-easy-and-flexible-bayesian-measurement-modeling-v1-0-0/35483)._