Manual calculation of logistic regression residuals

Question

I believe I understand how logistic regression works, and what it approximates. However, I do not understand how residuals of a logistic regression are calculated, especially when there are no predictors. Consider this simple example using R syntax:

> residuals(glm(cbind(c(1, 2, 4), c(0, 0, 2)) ~ 0, family=binomial))
        1         2         3 
1.1774100 1.6651092 0.8243762 

The first two values make sense: of course the evidence in the 2:0 case is stronger than in the 1:0 case, which should lead to more positive log odds estimation. But how do I get to those particular values? Why don't we get infinity in both 0 cases if we're calculating log(1/0)? The residuals for the 4:2 case are not log(2), they're log(2.28), as I would wrongly predict on using the same logic.

Code showing what's going on would be greatly appreciated.

Answer 1

The particular type of residual you're asking for here are deviance residuals, which are the signed square root of the deviance contribution from the observation.

The deviance is the -2log-likelihood at the fitted value (modelling $y\sim Binom(n,\hat\mu))$ minus the -2log-likelihood for a saturated model (modelling $y\sim Binom(n, y/n)$)

So:

dev.r<-function(y,n,muhat){
    diff <- dbinom(y,n,muhat,log=TRUE)-dbinom(y,n, y/n,log=TRUE)
    sign(-2*diff)*sqrt(-2*diff)
}

And $\hat\mu=0.5$ for a logistic regression with no parameters, so

> dev.r(y=c(1,2,4), n=c(1,2,6), muhat=0.5)
[1] 1.1774100 1.6651092 0.8243762

Nothing bad happens with the division by zero, because there really isn't one. The worst you get is $y\log y$ with $y=0$ or the equivalent, and that's not a problem.