[R-sig-ME] Is there a way of getting “marginal effects” from a `glmer` object
Henrik Singmann
henrik.singmann at psychologie.uni-freiburg.de
Mon Jun 16 22:14:43 CEST 2014
An alternatvie to the effects package could be either the lsmeans package or even lsmeans in combination with multcomp (see the great multcomp book for more examples).
Note however that this approach to my knowledge completely ignores the random effects. You most likely will need to use precit.merMod and marginalize on your own given the random effects structure you desire if you want to incorporate the random effects.
Hope that helps,
Henrik
require(lsmeans)
(l1 <- lsmeans(cfelr, ~ rank | ran, at = list(ran = quantile(mydata$ran, c(0.25, 0.5, 0.75)))))
## ran = 0.3027645:
## rank lsmean SE df asymp.LCL asymp.UCL
## 1 -0.04255557 0.2966504 NA -0.6240500 0.53893886
## 2 -0.72272008 0.2089092 NA -1.1322242 -0.31321594
## 3 -1.46011460 0.2828895 NA -2.0146349 -0.90559430
## 4 -1.62713014 0.3568840 NA -2.3266946 -0.92756565
##
## ran = 0.5592498:
## rank lsmean SE df asymp.LCL asymp.UCL
## 1 0.06329396 0.2809776 NA -0.4874786 0.61406651
## 2 -0.61687056 0.1872073 NA -0.9838346 -0.24990654
## 3 -1.35426508 0.2619736 NA -1.8677861 -0.84074405
## 4 -1.52128062 0.3490954 NA -2.2055778 -0.83698340
##
## ran = 0.7915128:
## rank lsmean SE df asymp.LCL asymp.UCL
## 1 0.15914709 0.3000271 NA -0.4289664 0.74726057
## 2 -0.52101743 0.2157033 NA -0.9438392 -0.09819563
## 3 -1.25841195 0.2785249 NA -1.8043768 -0.71244708
## 4 -1.42542749 0.3689908 NA -2.1487238 -0.70213122
##
## Confidence level used: 0.95
# on probability scale
binomial()$linkinv(summary(l1)$lsmean)
## [1] 0.4893627 0.3267943 0.1884498 0.1642239 0.5158182 0.3504935
## 0.2051740 0.1792730 0.5397030 0.3726144 0.2212474 0.1938121
require(multcomp)
lapply(as.glht(l1, df = 0), summary)
## $`ran = 0.302764488209505`
##
## Simultaneous Tests for General Linear Hypotheses
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 1 == 0 -0.04256 0.29665 -0.143 0.99982
## 2 == 0 -0.72272 0.20891 -3.459 0.00217 **
## 3 == 0 -1.46011 0.28289 -5.161 < 1e-04 ***
## 4 == 0 -1.62713 0.35688 -4.559 < 1e-04 ***
## ---
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
## (Adjusted p values reported -- single-step method)
##
##
## $`ran = 0.559249840560369`
##
## Simultaneous Tests for General Linear Hypotheses
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 1 == 0 0.06329 0.28098 0.225 0.99898
## 2 == 0 -0.61687 0.18721 -3.295 0.00393 **
## 3 == 0 -1.35427 0.26197 -5.169 < 1e-05 ***
## 4 == 0 -1.52128 0.34910 -4.358 5.27e-05 ***
## ---
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
## (Adjusted p values reported -- single-step method)
##
##
## $`ran = 0.791512808995321`
##
## Simultaneous Tests for General Linear Hypotheses
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 1 == 0 0.1591 0.3000 0.530 0.970969
## 2 == 0 -0.5210 0.2157 -2.415 0.060143 .
## 3 == 0 -1.2584 0.2785 -4.518 < 1e-04 ***
## 4 == 0 -1.4254 0.3690 -3.863 0.000449 ***
## ---
## Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
## (Adjusted p values reported -- single-step method)
Am 13.06.2014 10:54, schrieb Prof. Karthik D.:
> Hi
> I have asked this question on StackOverflow and was asked by Ben Bolker to
> post it here.
>
> thanks for your help.
>
> Karthik
>
>
>
> *reproducing the SO
> question http://stackoverflow.com/questions/24177197/is-there-a-way-of-getting-marginal-effects-from-a-glmer-object
> <http://stackoverflow.com/questions/24177197/is-there-a-way-of-getting-marginal-effects-from-a-glmer-object>*
>
> I am estimating random effects logit model using glmer and I would like to
> report Marginal Effects for the independent variables. For glm models,
> package mfx helps compute marginal effects. Is there any package or
> function for glmer objects?
>
> Thanks for your help.
>
> A reproducible example is given below
>
> mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv")
> mydata$rank <- factor(mydata$rank) #creating ranks
> id <- rep(1:ceiling(nrow(mydata)/2), times=c(2)) #creating ID variable
> mydata <- cbind(mydata,data.frame(id,stringsAsFactors=FALSE))
> set.seed(12345)
> mydata$ran <- runif(nrow(mydata),0,1) #creating a random variable
>
> library(lme4)
> cfelr <- glmer(admit ~ (1 | id) + rank + gpa + ran + gre, data=mydata
> ,family = binomial)
> summary(cfelr)
>
> [[alternative HTML version deleted]]
>
--
Dr. Henrik Singmann
Albert-Ludwigs-Universität Freiburg, Germany
http://www.psychologie.uni-freiburg.de/Members/singmann
More information about the R-sig-mixed-models
mailing list