[R-sig-ME] [Lme4-authors] Error in lme4 package

[David Federer]

Hi Ben, Steve,

See attached the dataset I performed the analysis on.

Employee_Id: Employee Number
Dpmt_Code: Department Code
HadIncident: Flags whether the employee has had an incident
Yrs_Service: The number of service years with the company

First I ran a regular logit, see below, which works out fine. The estimate of the intercept equals -0.9301 which equals the natural logarithm of the odd ratio:
p = 430/1520 = 0.2829; ln(0.2829/(1-0.2829)) = -0.9301

setwd("C:\\Users\\david.federer\\Documents")

mydata <- read.csv("SafetyAnalytics.csv", header=TRUE)

fixedmodel <- glm(HadIncident ~ 1, data = mydata, family = binomial("logit"))

summary(fixedmodel)



Call:

glm(formula = HadIncident ~ 1, family = binomial("logit"), data = mydata)



Deviance Residuals:

    Min       1Q   Median       3Q      Max

-0.8155  -0.8155  -0.8155   1.5891   1.5891



Coefficients:

            Estimate Std. Error z value Pr(>|z|)

(Intercept) -0.93015    0.05695  -16.33   <2e-16 ***

---

Signif. codes:

0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1



(Dispersion parameter for binomial family taken to be 1)



    Null deviance: 1810.8  on 1519  degrees of freedom

Residual deviance: 1810.8  on 1519  degrees of freedom

AIC: 1812.8



Number of Fisher Scoring iterations: 4

The I went ahead and ran a multilevel logit (intercept-only) which gave an estimated intercept of logistic of -1.2589 (see below). This figure roughly equates to the natural logarithm of the proportion 0.2829 (-1.2627) which is remarkable. If I would take -1.2627 to be the ln of the odds ratio the proportion would work out to be 0.2205 which is a 22% decrease relative to the fixed effect model.
It might be that ‘Dpmt_Code’ adds a lot of variance. I just wanted to ascertain that that is the case and glmer comes back with the ln of the odds ratio, rather than the proportion. Thanks guys for making time to look into this. David


library(lme4)

fit <- glmer(HadIncident ~ (1|Dpmt_Code), family=binomial("logit"), data=mydata)

summary(fit)

Generalized linear mixed model fit by maximum likelihood

  (Laplace Approximation) [glmerMod]

 Family: binomial ( logit )

Formula: HadIncident ~ (1 | Dpmt_Code)

   Data: mydata



     AIC      BIC   logLik deviance df.resid

  1728.6   1739.2   -862.3   1724.6     1518



Scaled residuals:

    Min      1Q  Median      3Q     Max

-1.2337 -0.6240 -0.4391  0.9691  2.8970



Random effects:

 Groups    Name        Variance Std.Dev.

 Dpmt_Code (Intercept) 0.7986   0.8936

Number of obs: 1520, groups: Dpmt_Code, 145



Fixed effects:

            Estimate Std. Error z value Pr(>|z|)

(Intercept)  -1.2589     0.1164  -10.81   <2e-16 ***

---

Signif. codes:

0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


From: Ben Bolker [mailto:bbolker at gmail.com]
Sent: Friday, 4 July 2014 7:46 AM
To: Steve Walker
Cc: David Federer; r-sig-mixed-models at r-project.org
Subject: Re: [Lme4-authors] Error in lme4 package

Agree with what Steve said.  It's quite possible that accounting for the random effects is giving you results that are correct but initially surprising. See for example this thread: http://article.gmane.org/gmane.comp.lang.r.lme4.devel/11951
  [taking the liberty of cc'ing to r-sig-mixed-models]

On Thu, Jul 3, 2014 at 11:14 AM, Steve Walker <steve.walker at utoronto.ca<mailto:steve.walker at utoronto.ca>> wrote:
Hi David,

I explored this issue with the build-in cbpp example in lme4 (attached) but didn't find any problems.  We would need a reproducible example to pursue this any farther.

Cheers,
Steve


On 2014-07-03, 2:44 AM, David Federer wrote:
Hi Ben,

I need to run a multilevel logistic regression on a dataset that embeds a nested structure.  I used the 'glmer' function to run a multilevel logit, setting the 'family' parameter to 'binomial("logit")'. I expected the function to regress the natural logarithm of the odds ratio ( ln(p/1-p) )against the independent variable that I had specified. Instead, it returned the natural logarithm of the proportion ( ln(p) ). Here is a numerical example (for the intercept only model):

In my dataset, p equals 0.2829. ln(p/p-1) should lie in the vicinity of -0.9301. Instead glmer returned -1.2589 which approximates ln(0.2829)

Is this an error in the 'glmer' function?


Kind Regards,



David

This message is for the designated recipient only and may contain privileged, proprietary, or otherwise private information. If you have received it in error, please notify the sender immediately and delete the original. Any other use of the email by you is prohibited.


_______________________________________________
Lme4-authors mailing list
Lme4-authors at lists.r-forge.r-project.org<mailto:Lme4-authors at lists.r-forge.r-project.org>
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/lme4-authors



This message is for the designated recipient only and may contain privileged, proprietary, or otherwise private information. If you have received it in error, please notify the sender immediately and delete the original. Any other use of the email by you is prohibited.