I am getting different results (close but not exact the same) from R GLM and manual solving logistic regression optimization. Could anyone tell me where is the problem?
BFGS does not converge?
Numerical problem with finite precision?
Thanks
# logistic regression without intercept
fit=glm(factor(vs) ~ hp+wt-1, mtcars, family=binomial())
# manually write logistic loss and use BFGS to solve
x=as.matrix(mtcars[,c(4,6)])
y=ifelse(mtcars$vs==1,1,-1)
lossLogistic <- function(w){
L=log(1+exp(-y*(x %*% w)))
return(sum(L))
}
opt=optim(c(1,1),lossLogistic, method="BFGS")
Best Answer
Short answer: Optimise harder.
Your loss function is fine, no numeric issues there. For instance you can easily check that:
What happens is that you assume that
optim
's BFGS implementation can get as good as an routine that use actual gradient information - remember Fisher scoring is essentially a Newton-Raphson routine. BFGS converged (opt$convergence
equals0
) but the best of BFGS was not the best you could get because as you did not provide a gradient function the routine had to numerically approximate the gradient. If you used a better optimisation procedure that could use more gradient-like information you would get the same results. Here, because the log-likelihood is actually a very well behaved function I can even use a quadratic approximation procedure to "fake" gradient information.It works.