1
$\begingroup$

I created some data using the following code:

set.seed(1221)
x <- runif(500)
y <- runif(500, 0, 2)
z <- rep(0, 500)
z[-0.8*x + y - 0.75 > 0] <- 1
plot(x, y, col=as.factor(z))

This produces the following plot

enter image description here

The data is linearly separable. Then, I applied the glm function to create a logistic regression model.

df <- data.frame(class = z, x = x, y = y)
model <- glm(z ~ x + y, family = binomial, data = df)

This produces the following output:

summary(model)
Call:
glm(formula = z ~ x + y, family = binomial, data = df)

Deviance Residuals: 
       Min          1Q      Median          3Q         Max  
-8.127e-04  -2.000e-08  -2.000e-08   2.000e-08   7.699e-04  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)
(Intercept)    -1062      52666   -0.02    0.984
x              -1163      57197   -0.02    0.984
y               1433      70408    0.02    0.984

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 6.8274e+02  on 499  degrees of freedom
Residual deviance: 1.3345e-06  on 497  degrees of freedom
AIC: 6

Number of Fisher Scoring iterations: 25

The result surprised me, first because the parameter estimates are huge, and second because I was expecting such estimates to be close to the original decision boundary function, i.e. -0.8x + y - 0.75 = 0.

I then used the glmnet package to see if I could solve this issue. This package creates a penalised logistic regression model in order to deal with the large values in the parameter estimates. The code I used is the following:

library(glmnet)
cvfit <- cv.glmnet(as.matrix(df[,-1]), as.factor(df$class), 
                   family =   "binomial", 
                   type.measure = "class")
plot(cvfit)

enter image description here

And the coefficients for the optimal penalty strength are:

coef(cvfit, s = "lambda.min")
3 x 1 sparse Matrix of class "dgCMatrix"
                    1
(Intercept) -84.01446
x           -91.40983
y           113.18736

Such coefficients are smaller than the ones obtained with the glm function. Still they are not the same as the decision boundary function.

Does anybody know why this is happening? Any help is greatly appreciated.

$\endgroup$
5
  • 3
    $\begingroup$ Review posts bearing the tag Hauck-Donner effect. $\endgroup$
    – Sycorax
    Dec 15, 2015 at 14:37
  • 2
    $\begingroup$ Also note that the coefficients are log odds ratios - not of a decision boundary function. See here. (In fact the decision boundary function derived from your logistic regression is pretty close to what it should be: I make it $0.81x + y - 0.74 = 0$.) $\endgroup$ Dec 15, 2015 at 14:46
  • $\begingroup$ Thanks @user777 and Scortchi, following your advises I found an interesting link that clarifies my question. $\endgroup$ Dec 15, 2015 at 15:58
  • $\begingroup$ Yes, I just noticed that @Scortchi I completely forgot the fact that the coefficients are log odds ratios. I was expecting to get the boundary function right away. Sorry, my mistake... $\endgroup$ Dec 15, 2015 at 16:01
  • 2
    $\begingroup$ @jroberayalas: Don't be sorry. It's perhaps interesting to note here that even though the maximum-likelihood estimates are running away from each other, they're still constrained to give a sensible answer to your question. The Wald standard errors are the only thing that's unequivocally wrong. $\endgroup$ Dec 15, 2015 at 17:01

0

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.