Skip to content

lme4 categorical outcome probability values at each time point

2 messages · Igor Yakovenko, Evan Palmer-Young

Original
Igor Yakovenko
# 
Hello Everyone,

 

Does anyone know how to graph the results of a lme4 longitudinal multilevel
model as a bar graph of Time (3 points) on x-axis with probability of
getting the 1 value on the outcome after taking into account the predictors
as the y-axis? Alternatively, if there is a way to derive the mean
probability of the outcome being 1 (vs. 0) during an intervention at time 1,
2 and 3 divided by Condition (two levels), I can graph those values
separately. The model is based on data from a two-group RCT evaluating an
intervention over 3 time points. Here is a sample model with a dichotomous
outcome and fixed effects of time (continuous), condition (2 levels) and
whether an additional workshop was completed (2 levels).

 

workshop<-glmer(CurrentAsst~Time*Condition + Workshop_completion_reordered +
(1|id), data = vse, family = binomial, nAGQ = 25, control =
glmerControl(optimizer = "bobyqa"))

 

Thank you,

Igor Yakovenko
Evan Palmer-Young
# 
Dear Igor,
Am new to the list, so defer by all means to more experienced opinions.
My suggestion would be to use the "lsmeans" package and the inverse logit
transformation
(to get back to the scale of probabilities from the hard-to-interpret
logits)

If you have time as a factor this will be no problem.
So before you run the model, specify
vse$Time<-as.factor(vse$Time)
#then use lsmeans by R Lenth
#see
https://cran.r-project.org/web/packages/lsmeans/vignettes/using-lsmeans.pdf
#example
# lsmeans(oranges.lm1, "day", at = list(price1 = 50, price2 = c(40,60), day
= c("2","3","4")) )
Newdf<- summary(lsmeans(workshop, ~Time*Condition))
#or if Time is numeric:
Newdf<- summary(lsmeans(workshop, ~Time*Condition, at=list(Time=c(1, 2, 3))

#if you need to graph standard errors, create a new column with raw std
errors first,
#then transform the mean +/- SE
Newdf$sehi<-Newdf$lsmean+ Newdf$SE
Newdf$selo<-Newdf$lsmeans - Newdf$SE
#now apply inverse logit
library(boot)
Newdf$meanprop<-inv.logit (Newdf$lsmean)
Newdf$prophi<-inv.logit(Newdf$sehi)
Newdf$proplo<-inv.logit(Newdf$selo)

You should be able to plot the resulting variables using ggplot2 or your
package of choice.
You could also use inv.logit on the 95% confidence limits, though the error
bars tend to get pretty enormous with these reverse transformations and may
not seem to reflect any significant differences in your model.
Hope this helps!
Evan


On Fri, Mar 25, 2016 at 12:20 AM, Igor Yakovenko <iyakoven at ucalgary.ca>
wrote: