Extracting Model Coefficients in brms - Output Unclear

Interfaces brms
1
  • Operating System: Windows 10
  • brms Version: 2.8.0

I am currently fitting a generalized logistic function model.
The code of the fitted model is given below:

fit<- brm(
bf(LDL ~ 15 + (alpha-15) / (1 + exp ((gamma-Time) / delta)),
alpha ~ 1 + (1|Cat/Sample),
gamma ~ 1 + (1|Cat/Sample),
delta ~ 1 + (1|Cat),
nl = TRUE),
prior = c(
prior(normal(3800, 50), nlpar = “alpha”),
prior(normal(3,0.4),lb=0,ub=5, nlpar = “gamma”),
prior(normal(0.4,0.04),lb=0,ub=1, nlpar = “delta”)),
data = ldlData, family = gaussian(),
chains = 1,
iter=6000,
warmup = 1000,
cores = getOption(“mc.cores”,1L),
control = list(adapt_delta = 0.999)
)

Since I am fitting a multilevel model, I am curious in the parameter values calculated (fixef + ranef) for the model corresponding to each individual. I am using

coef(fit)

to obtain the individuals’ parameter values, however, I have some doubts that the output generated is not doing exactly what I want. Here are extracts of the output:

coef(fit)
$Cat
, , alpha_Intercept

Estimate Est.Error Q2.5 Q97.5
0 3394.274 651.5656 2060.987 4388.694
1 3407.948 662.1964 2081.623 4491.990

, , gamma_Intercept

Estimate Est.Error Q2.5 Q97.5
0 -367.9757 385.5372 -1305.909 213.6622
1 -368.5811 388.0063 -1260.914 205.1106

, , delta_Intercept

Estimate Est.Error Q2.5 Q97.5
0 -781.3041 706.3485 -2450.643 -157.0636
1 -782.9636 819.0308 -2443.098 -143.9825

$Cat:Sample
, , alpha_Intercept

 Estimate Est.Error     Q2.5    Q97.5

0_1 3735.186 320.0096 2945.746 4327.521
0_11 3675.675 345.7927 2774.851 4238.050
1_12 3766.073 312.7897 3068.951 4415.879
1_14 3796.811 321.7548 3134.178 4507.588

, , gamma_Intercept

    Estimate Est.Error      Q2.5    Q97.5

0_1 -27.7285129 153.9734 -396.9019 236.0779
0_11 -64.7150180 165.5760 -506.8677 165.7535
0_13 -38.1695172 151.3383 -390.7412 216.7660
0_15 15.8747282 149.8950 -265.1621 337.6533
0_17 -7.0927562 138.9808 -319.4078 275.6221
0_47 -11.6233737 144.0304 -328.6369 276.9613
1_60 20.8967926 145.8143 -252.7480 348.3330
1_8 -40.6679938 155.7832 -410.7497 222.7164

, , delta_Intercept

  Estimate  Est.Error      Q2.5     Q97.5

0_15 0.3997916 0.03944224 0.3239464 0.4761564
1_60 0.3997916 0.03944224 0.3239464 0.4761564
1_8 0.3997916 0.03944224 0.3239464 0.4761564

According to this, the output of coef() is supposed to be “model coefficients, which are the sum of population-level effects and corresponding group-level effects”.
However, I am not sure if this is the case - the output provided seems to be the group effects wrt Cat and Cat:Sample (where an intercept for delta is given even for Cat:Sample, which is odd, as by definition of delta it is not affected by this group effect).
I also read the Value section of the link provided, however I find it very confusing.

So, to sum up, how do I extract the values of the model coefficients for the individual models (fixef+ranef, not just ran)?

2

In your case, I would recommend extracting the fixed effects via fixef(..., summary = FALSE) and the random effects via ranef(..., summary = FALSE) and then add together whatever you want to add. Afterwards, you can summarize the posterior via posterior_summary.

3

I have a very similiar problem. I have a term that is similar alpha ~ 0 + (1|Cat/Sample). So it has two levels on the random effects (Cat and Sample) and non on the fixed effect size. With ranef() and coef() I get the same values once for the levels of cat and then for the levels of cat:samples.

Is there a build in function that gives me a summary of the combined effect? So a function that adds the right cat effect to the cat:sample effect? If not, how do I do it most efficient tidy style?

4

Did you find a solution to this?
Many thanks

5

I found some a workarounds. First a approach in base R (pritty cumbersome):

Since my experience with R grows I was able to use the dataframes where the relation of samples and categories are specified (e.g. by using unique() on the both columns of the longformat dataset or making a wide representation of it) and then usind this relation to reference the correct categorie ranefs to add it to the sample ranefs.

A more convinient way I found later was the package tidybayes where I spread the draws of sample and cat and then juste add the two columns.

1 Like