# Unexpected behaviour of logit regression with glm in R

I recently was puzzled by the behaviour of R's glm when trying to compute a logistic regression

data <- data.frame(
response = rbinom(600, 1, prob=rep(c(1,0.5,0),each=200)),
predictor = rep(c("A","B","C"), each = 200)
)

llm1 <- glm(response ~ predictor, data = data, family=binomial(link='logit'))
summary(llm1)


The output tells me that there is apparently no difference in my predictors

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept)    20.57    1253.73   0.016    0.987
predictorB    -20.81    1253.73  -0.017    0.987
predictorC    -41.13    1773.04  -0.023    0.981



When I compute glm with repsonse values where not all values are ones or zeros in predictors "A" and "C" respectively, I receive results as I had expected them.

data2 <- data.frame(
response = rbinom(600, 1, prob=rep(c(0.95,0.5,0.05),each=200)),
predictor = rep(c("A","B","C"), each = 200)
)

llm2 <- glm(response ~ predictor, data = data2, family=binomial(link='logit'))
summary(llm2)

# output:
...
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept)   2.9444     0.3244   9.075  < 2e-16 ***
predictorB   -2.8243     0.3540  -7.978 1.49e-15 ***
predictorC   -6.2612     0.5033 -12.441  < 2e-16 ***
...


I suspect it has something to do with the link function, but I can't really pinpoint it. Has anyone experienced this and can offer an explanation for the observed behaviour?

This is due to perfect separation of the response variable, where the covariates perfectly predict the outcome.

This is the case in your example when predictor = "A"(always predict 1), and since A goes into the intercept you get an intercept of 20, pushing the probability toward 1. It's also the case when prediction = "C" but now the opposite way, since you now get $$P(y = 1) = \exp(20 - 40) / (1 + \exp(20 - 40))$$ giving you a prediction of 0.

The maximum likelihood estimate is not valid when you have perfect separation and the estimates should, in theory, be $$\infty$$, since the bigger the coefficient, the closer the probability is to 1 (or 0) and the better your prediction. Why is the estimate not infinity? Probably due to some stoping rule inside glm. The likelihood becomes very flat (see figure 2 here) and there is not a big enough change in the likelihood when doing an additional iteration in the optimization, even though you are not at a global maximum.

The flatness of the likelihood is also the reason behind the big standard errors, since these are based on the inverse Hessian (matrix of the second derivatives of the likelihood wrt the parameters), so a small curvature (flat likelihood) leads to large elements in the inverse Hessian.

As you can see when you don't have perfect separation, the parameter estimates are more sensible.

EDIT

One way to solve the problem of perfect separation is to add a prior to your parameters/penalize big estimates. This can be done using Ridge or Lasso, where you penalize the likelihood when $$\beta$$ becomes big. For your example:

data <- data.frame(
response = rbinom(600, 1, prob=rep(c(1,0.5,0),each=200)),
predictor = rep(c("A","B","C"), each = 200)
)

mod3 = ridge::logisticRidge(response ~ predictor, data = data)

summary(mod)

Call:
ridge::logisticRidge(formula = response ~ predictor, data = data)

Coefficients:
Estimate Scaled estimate Std. Error (scaled) t value (scaled) Pr(>|t|)
(Intercept)    2.438              NA                  NA               NA       NA
predictorB    -2.471         -28.531               3.349            -8.52   <2e-16 ***
predictorC    -4.446         -51.339               3.835           -13.39   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Ridge paramter: 0.0005737279, chosen automatically, computed using 2 PCs

Degrees of freedom: model 2.97 , variance 2.94


where now the parameter estimates are a lot closer to zero, but you still get problems with the standard errors.

• (+1) Nice answer. Using the curvature of the likelihood surface to relate to the estimation of both the coefficients and standard errors really helped me understand this better.
– Noah
Commented Aug 3, 2020 at 23:55
• Thanks for the great answer. it addressed all my questions and even some more! Commented Aug 4, 2020 at 7:14