Is unit deviance (statistics) equivalent to the loss function (machine learning)

Question

In this page from scikit learn, about GLM, the notion of unit deviance is introduced as loss function (from the machine learning perspective).

I want to know if there is equivalence between these two notions: unit devance vs. loss function.

For linear regression, there is equivalence, since the loss function is the squared error and it is equivalent to maximizing the Gaussian distribution in MLE (Maximum Likelihood Estimation).

For logistic regression, I struggled to understand, since the deviance is:

$$2({y}\log \frac{y}{\hat{y}}+({1}-{y})\log \frac{{1}-{y}}{{1}-\hat{y}})$$

Whereas the classic log loss is:

$$y \log(p) + (1-y) \log(1 - p)$$

And this post seems to confirm that the deviance is the same as log loss.

As for Gamma regression, I found this post discussing the loss function for gamma in XGBoost and the result is different from the one seen in scikit-learn.

Answer 1

The total deviance is twice the difference in log-likelihood between a "saturated model" and the fitted model. Using notation from Dobson and Barnett's "An Introduction to Generalized Linear Models", the total deviance is equal to

$$ D = 2 \left[ l(\mathbf{b}_{\max}; \mathbf{y}) - l(\mathbf{b}; \mathbf{y}) \right] \>.$$

Here, $\mathbf{b}_{\max}$ is a maximum likelihood estimator for the saturated model with true parameter $\boldsymbol{\beta}_{\max}$, and $\mathbf{b}$ is a maximum likelihood estimator for the fitted model. The easiest way to get $\mathbf{b}_{\max}$ is to have 1 parameter per observation, which then just makes each prediction equal to the observation.

You can verify all this in R. The glm function will compute the deviance for you (listed under residual deviance)

library(tidyverse)

set.seed(0)
N <- 1000
x <- rnorm(N)
p <- plogis(-2 + 0.8*x)
y <- rbinom(N, 1, p)
fit <- glm(y~x,family = binomial())
p_est <- predict(fit, type = 'response')

# Deviance
D = sum(2*(dbinom(y, 1, y, log = T) - dbinom(y, 1, p_est, log = T)))
print(D)
#> [1] 739.2842

s <- summary(fit)
print(s$deviance)
#> [1] 739.2842

^{Created on 2023-03-24 with reprex v2.0.2}

Note that for a given $\mathbf{y}$, $l(\mathbf{b}_{\max}; \mathbf{y})$ is constant, so optimizing $D$ is the same as optimizing the log likelihood. In this sense, the deviance is the same as the loss. Again, easy to verify in R

loss<-function(b, X, y){
  p <- plogis(X %*% b)
  ll <- 2*sum(dbinom(y, 1, y, log = T) - dbinom(y, 1, p, log = T))
  
  return(ll)
}


b <- c(0,0)

r <- optim(b, fn=loss, X=model.matrix(~x), y=y)
round(r$par, 4)
#> [1] -2.0182  0.7623

round(coef(fit), 4)
#> (Intercept)           x 
#>     -2.0182      0.7622

# Coefficients from optimizing the log-likelihood are the same
# As coefficients from optimizing the deviance.

^{Created on 2023-03-24 with reprex v2.0.2}

A little algebra goes a long way. First, let's note that $\hat{y}$ is the estimated risk from the logistic regression model (it is not a 0 or 1). Let's re-write half the binomial deviance using log rules as

$$ y\log(y) - (1-y)\log(1-y) - y \log(\hat{y}) + (1-y)\log(1-\hat{y}) $$

Now using what we've discussed

$$ \underbrace{y\log(y) - (1-y)\log(1-y)}_{l(\mathbf{b}_{\max}; \mathbf{y})} - \overbrace{y \log(\hat{y}) + (1-y)\log(1-\hat{y})}^{l(\mathbf{b}; \mathbf{y})} $$

OK so, what have we know now:

• The log-likelihood is the loss for generalized linear models
• The deviance is proportional to the log likelihood minus a constant
• We showed this with the binomial deviance (the poisson is perhaps the next easiest, I suggest you do that as practice).
• This means optimizing the deviance is equivalent to maximizing the optimizing log likelihood.
• Hence, deviance and loss should be equivalent in so far as optimizing both leads to the same model. This might be why XGboost has a different loss function for gamma regression (it uses the loglikelihood as opposed to the deviance).

Personally, I find it a bit strange to optimize the deviance as opposed to the log likelihood directly, but sklearn has made it very clear they are NOT a statistics package so I'm ultimately not surprised.