Pearson VS Desviación Residuales en regresión logística

Question

Sé que los residuos de Pearson estandarizados se obtienen de forma probabilística tradicional:

$$ r_i = \frac{y_i-\hat{\pi}_i}{\sqrt{\hat{\pi}_i(1-\hat{\pi}_i)}}$$

y los residuos de desviación se obtienen de una forma más estadística (la contribución de cada punto a la probabilidad):

$$ d_i = s_i \sqrt{-2[y_i \log \hat{\pi_i} + (1 - y_i)\log(1-\hat{\pi}_i)]} $$

donde $s_i$ = 1 si $y_i$ = 1 y $s_i$ = -1 si $y_i$ = 0.

¿Puedes explicarme, intuitivamente, cómo interpretar la fórmula de los residuos de desviación?

Además, si quiero elegir una, ¿cuál es la más adecuada y por qué?

BTW, algunas referencias afirman que derivamos los residuos de desviación basándonos en el término

$$-\frac{1}{2}{r_i}^2$$

donde $r_i$ se menciona más arriba.

Answer 1

La regresión logística trata de maximizar la función logarítmica de verosimilitud

$LL = \sum^k \ln(P_i) + \sum^r \ln(1-P_i)$

donde $P_i$ es la probabilidad prevista de que el caso i sea $\hat Y=1$ ; $k$ es el número de casos observados como $Y=1$ y $r$ es el número de casos (el resto) observados como $Y=0$ .

Esa expresión es igual a

$LL = ({\sum^k d_i^2} + {\sum^r d_i^2})/-2$

porque un caso desviación residual se define como:

$d_i = \begin{cases} \sqrt{-2\ln(P_i)} &\text{if } Y_i=1\\ -\sqrt{-2\ln(1-P_i)} &\text{if } Y_i=0\\ \end{cases}$

Así, la regresión logística binaria busca directamente minimizar la suma de los residuos de desviación al cuadrado. Son los residuos de desviación los que están implícitos en el algoritmo ML de la regresión.

El estadístico Chi-sq del ajuste del modelo es $2(LL_\text{full model} - LL_\text{reduced model})$ donde el modelo completo contiene predictores y el modelo reducido no.

Answer 2

En respuesta a esta pregunta He añadido un código R para mostrar cómo aplicar manualmente la fórmula para el cálculo de los residuos de desviación

El modelo del código es un modelo logit en el que

$$p_i := Pr(Y_i = 1) = \frac{\exp(b_0 + b_1x_i)}{1+\exp(b_0 + b_1x_i)}.$$

Defino $v_i := b_0 + b_1x_i$ de forma que el modelo pueda escribirse como

$$p_i := Pr(Y_i = 1) = \frac{\exp(v_i)}{1+\exp(v_i)}.$$

Estimando el modelo obtengo estimaciones $\hat b_0$ y $\hat b_1$ . A partir de estas estimaciones, los valores latentes previstos

$$\hat v_i := \hat b_0 + \hat b_1 x_i,$$

y, a continuación, se calculan las probabilidades previstas

$$\hat p_i =\frac{\exp(\hat v_i )}{1+\exp(\hat v_i )}.$$

A partir de estas probabilidades previstas, se aplica la fórmula de los residuos de desviación en la etapa de codificación

sign(y-pred_p) * ifelse(y==1,sqrt(-2*log(pred_p)),sqrt(-2*log(1-pred_p)))

que no es más que una aplicación de la fórmula

$d_i = \begin{cases} \sqrt{-2\ln(\hat p_i)} &\text{if } Y_i=1\\ -\sqrt{-2\ln(1-\hat p_i)} &\text{if } Y_i=0\\ \end{cases}$

# Simulate some data
N <- 1000
b0 <- 0.5
b1 <- 1
x <- rnorm(N)
v <- b0 + b1*x
p <- exp(v)/(1+exp(v))
y <- as.numeric(runif(N)<p)

# Estimate model
model <- glm(y~x,family=binomial)
summary_model <- summary(model)
summary_dev_res <- summary_model$deviance.resid
# This is the output you get:
quantile(summary_dev_res)

# Calculate manually deviance residuals
# First calculate predicted v's
pred_v <- coef(model)[1] + coef(model)[2]*x
# The calculate predicted probabilities
pred_p <- exp(pred_v)/(1+exp(pred_v))
# Apply formula for deviance residuals
dev_res <- sign(y-pred_p) * ifelse(y==1,sqrt(-2*log(pred_p)),sqrt(-2*log(1-pred_p)))
# Check that it is the same as deviance residuals returned from summary
plot(summary_dev_res,dev_res)
points(seq(-3,3,length.out=100),seq(-3,3,length.out=100),type="l",col="red",lwd=2)
# all points should be on the red line 

# Also compare the quantiles ... 
quantile(summary_dev_res)
quantile(dev_res)