 # Bernoulli BRMS

Hello everyone,
(sorry in advance for my extensive message… i’m a beginner in stan and BRMS)
Here’s my problem.
I have a dependent variable, denoted REP, binary (1 or 0). I wish to predict the probability of answering 1 with a non-linear model such as, `REP ~ inv_logit ((X/(1 + k*D)) - X2)`. with X and X2 two numerical variables , D a numerical variable and k : the free parameter to estimate.

I would like to estimate the parameter “k” according to the different experimental conditions I have realized.
Independent variable A: “DOMAIN” with two modalities (1: LOSS, 2: GAIN). Each subject is exposed to these 2 modalities of this independent variable.
Independent variable B: “STATUS” with two modalities (1: Control, 2: Threat). A participant is either “control” or “threat” but never both.

my data frame :

I would like to estimate the parameter “k” for each condition :

1. k1_Control_GAIN,
2. k2_Control_LOSS,
3. k3_Threat_GAIN,
4. k4_Threat_LOSS.

My hypothesis:
k1_Control_GAIN < k3_Threat_GAIN
k2_Control_LOSS < k4_Threat_LOSS

I have several ways to code this I think, the way I tried is this: (I’ve simplified the model here to make it more understandable)

``````main_model = bf(
rep ~ inv_logit(
((X / (1 + exp(k1 * dummy_ctrl * dummy_gain + k2 * dummy_gain * dummy_threat + k3 * dummy_loss * dummy_ctrl + k4 * dummy_loss * dummy_threat)*D ) - X2))
,
k1 ~ ( 1| ID ),
k2 ~ ( 1| ID ) ,
k3 ~ ( 1| ID ),
k4 ~ ( 1| ID ),
nl = TRUE)
``````

with `dummy_ctrl = 1` when the subject is in the “control” group and 0 otherwise
with `dummy_threat = 1` when the subject is in the “threat” group and 0 otherwise
with `dummy_GAIN= 1` when the subject is exposed to a trial of the “GAIN” condition and 0 otherwise
with dummy_LOSS= 1 when the subject is exposed to a trial of the “LOSS” condition and 0 otherwise

``````main_fit = brm(main_model,
data = DF, family = bernoulli(link='identity'),
prior = c(prior(normal(-5,3), nlpar = 'k1'),
prior(normal(-5,3), nlpar = 'k2'),
prior(normal(-5,3), nlpar = 'k3'),
prior(normal(-5,3), nlpar = 'k4')),
cores = parallel::detectCores(),
seed = "12345",
inits = "random",
chains = 4,
iter = 4000,
warmup = 2000,
``````

I don’t think my way is right.
I think it’s better to do it like this?

``````main_model = bf(
rep ~ inv_logit(  ( (X / (1 + exp(k)*D ) - X2)) ,

k ~ DOMAIN*STATUS + (1 + DOMAIN | ID),

nl = TRUE)
``````