4
$\begingroup$

I used the brms package to carry out a mixed-effects logistic regression analysis with random intercepts. I'm having trouble interpreting the coefficients, in particular the transformation of the coefficients from log odds to probabilities with the plogis() function, which is equal to exp() / 1 + exp().

Here are the results:

 Family: bernoulli 
  Links: mu = logit 
Formula: OUTCOME ~ PREDICTOR + (1 | WORD) 
   Data: lang.data (Number of observations: 584) 
Samples: 6 chains, each with iter = 5000; warmup = 1000; thin = 1;
         total post-warmup samples = 24000

Group-Level Effects: 
~WORD (Number of levels: 196) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     3.57      0.65     2.45     4.99 1.00     7968    13158

Population-Level Effects: 
                 Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept          -5.49      0.79    -7.19    -4.09 1.00    13082    15203
PREDICTOR:LEVEL1    6.45      0.80     5.00     8.14 1.00    18012    16581
PREDICTOR:LEVEL2    0.80      0.74    -0.64     2.27 1.00    26702    20032

Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

The PREDICTOR variable has three levels and the OUTCOME is a dichotomous variable, with the values 0 for absence and 1 for presence.

Here are my questions:

  1. When I transform the estimate of the intercept with plogis(), I get the following:
plogis(-5.49)
[1] 0.004110875

Does this mean that when PREDICTOR is at the reference level, the probability of OUTCOME being 1 is 0.41 percent?

  1. When I add the transformed intercept and the estimate of PREDICTOR:LEVEL1, I get the following result:
plogis(-5.49 + 6.45)
[1] 0.7231218

I understand this to mean that the probability of OUTCOME being 1 is about 72 percent when PREDICTOR has the value of LEVEL1.

When I add the transformed intercept and the estimate of PREDICTOR:LEVEL2, I get the following result:

plogis(-5.49 + 0.80)
[1] 0.009103059

What I don't understand is why the three values above (0.004110875, 0.7231218, 0.009103059) do not sum to 1. If these values exhaust the space of predictor-variable values, the outcome has to occur in one of these three conditions, so why then do they not sum to 1?

  1. How do interpret the value of sd(Intercept)? I know that it refers to the standard deviation around the value of the intercepts, but what do I make of the value 3.57? Is it high? Is it low?
$\endgroup$

1 Answer 1

5
$\begingroup$

Example 1

Think of a simple example first. You have a single school from which you select a random sample of 100 students. For each student, you measure the values of an OUTCOME value and a PREDICTOR value, where:

OUTCOME = whether or not the student passed a particular exam (1 = passed; 0 = did not pass)

PREDICTOR = number of hours the student spent preparing for the exam, categorized into "1 day or less", "2-4 days" and "5-7 days"

When you fit a binary logistic regression model of the form brms(OUTCOME ~ PREDICTOR, family = bernoulli, etc.), what is it that you are ultimately estimating with that model? You are estimating the proportion of students in that school who passed the exam given they spent "1 day or less", "2-4 days" or "5-7 days", respectively, preparing for the exam.

There is no reason to believe that the 3 proportions should add up to 1 - in fact they can add up to either less than 1 or more than 1. For example, the 3 proportions of interest could be 0.72, 0.42 and 0.14, say.

What you can expect to add up to 1 would be:

  1. The proportion of students in the school who passed the exam given they studied for "1 day or less" and the proportion of students who did NOT pass the exam given they studied for "1 day or less";

  2. The proportion of students in the school who passed the exam given they studied for "2 - 4 days" and the proportion of students who did NOT pass the exam given they studied for "2 - 4 days";

  3. The proportion of students in the school who passed the exam given they studied for "5 - 7 days" and the proportion of students who did NOT pass the exam given they studied for "5 - 7 days".

Example 2

Now, assume you expand on the first example and sample 20 schools at random and then within each school you sample 100 students at random. For each student, you measure the OUTCOME and PREDICTOR variables defined above; you then use the data collected on these variables to fit a mixed effects binary logistic regression model of the form brms(OUTCOME ~ PREDICTOR + (1|SCHOOL), family = bernoulli, etc.).

With this more complicated set up, the parameters reported in your brms model summary refer to the typical school; if you plug them into the plogis() function as you did in your post, you'll be able to estimate the proportion of students at the "typical" school who passed the exam given they studied a specified number of days (e.g., "5 - 7 days"). Again, there is no reason to believe that these 3 proportions should add up to 1! They are whatever they are - for example, these proportions could increase in as a function of the number of days spent studying and their total could be less than 1 or greater than 1. But, given a number of days spent studying, the proportion of students who passed the exam and the proportion of students who did NOT pass the exam will add up to 1.

sd(Intercept)

To get a handle on how big/small the SD(Intercep) value (SD for short) is, why not compute plogis(-5.49 - 1.96 * SD) and plogis(-5.49 + 1.96 * SD) and see how variable the resulting probabilities are? If very variable, then SD is "large"; if not very variable, then SD is "small".

$\endgroup$
2
  • 1
    $\begingroup$ Thank you! Could you just explain how calculating plogis(-5.49 - 1.96 * SD) and plogis(-5.49 + 1.96 * SD) provides insight into how big/small the SD is? Where did the number 1.96 come from? $\endgroup$
    – Namenlos
    Dec 23, 2020 at 22:15
  • $\begingroup$ If you fit your mixed effects model with lmer or glmer, you can look into the concept of 95% random effects confidence intervals. (Not sure what the equivalent of this would be in the Bayesian setting for a model fitted with brms.) See Lesa Hoffman’s slide on Quantification of Random Effects Variances from this presentation: lesahoffman.com/CLDP944/CLDP944_Lecture05_Random_Effects.pdf. The slide refers to a random effect in a linear mixed effects model. $\endgroup$ Dec 29, 2020 at 2:01

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.