# Sharing parameters across distributions in a mixture model

Hi,
I was wondering if I could share parameters across the distributions of a mixture in brms. Imagine a situation where a school grade depends on a predictor, but a school teacher might have favorite students who they don’t want to fail and so another distribution (of favorites) is truncated at 59.5:

``````library(extraDistr)
library(brms)
set.seed(123)
# number of students
N <- 400
alpha <- 55
beta <- 2
sigma <- 8
predictor <- rnorm(N)
scores_1 <- rtnorm(N, alpha + beta * predictor,sigma, a =0,b=100)
# scores for grades that are always a pass
scores_2 <- rtnorm(N, alpha + beta * predictor,sigma, a =59.5,b=100)
# prob that the student is not a favorite:
p <- .7
# which students are a favorite
z <- rbern(N, p)
# distribution of scores:
scores <- ifelse(z, scores_1, scores_2)
hist(scores)
``````

Created on 2023-05-03 with reprex v2.0.2

Now how do I code this in brms? I tried the following:

``````df <- data.frame(scores, predictor)

brm(bf(scores ~ predictor ,
mu1 | rep_trunc(lb=0,ub=100)~ 0,
mu2 | rep_trunc(lb=59.5,ub=100)~ 0
)  , data= df,
family = mixture(gaussian, gaussian, order = 'none'),
prior = ...
)
``````

But I can see with `make_stancode` that the model is ignoring `predictor`. But if I move it to the lines of mu1 and mu2:

``````mu1 | rep_trunc(lb=0,ub=100)~ predictor,
mu2 | rep_trunc(lb=59.5,ub=100)~ predictor
``````

I end up with two betas. (Ideally I would like to share sigma across models, but it’s less critical).

(I know how to do it in Stan, and I know how to edit the brms model, but I wonder if there’s a brms solution that I can use for teaching. And by the way this is a very realistic situation, most datasets with grades have a bump around 60).

• Operating System: Ubuntu 20.04.5 LTS
• brms Version: brms_2.19.0

I’m not totally sure I understand the model you’re trying to fit, or what sort of priors might be reasonable here, but I think you can do this in brms by using a non-linear model. For example:

``````mix <- mixture(gaussian, nmix=2, order = "none")
formula <- bf(scores | resp_trunc(lb=0, ub=100) ~ alpha + beta*predictor,
nlf(mu2 | resp_trunc(lb=59.5, ub=100) ~ alpha + beta*predictor),
nlf(sigma1 ~ x),
nlf(sigma2 ~ x),
alpha ~ 1,
beta ~ 1,
x ~ 1,
nl = TRUE,
family=mix)
get_prior(formula, data=df)
priors = c(set_prior("normal(50, 20)", nlpar="alpha", class="b"),
set_prior("normal(0, 5)", nlpar="beta", class="b"),
set_prior("normal(0, 10)", nlpar="x", class="b"))
fit <- brm(formula, data=df, prior = priors)
``````

This seems to recover the parameter values okay, though not the mixing proportion.

1 Like