La siguiente explicación es no se limita a la regresión logística pero se aplica igualmente en la regresión lineal normal y en otros MLG. Por lo general, R
excluye un nivel del categórico y los coeficientes denotan la diferencia de cada clase con esta clase de referencia (o a veces llamada clase de referencia) (esto se llama codificación ficticia o contrastes de tratamiento en R
, ver aquí para un excelente resumen de las diferentes opciones de contraste). Para ver los contrastes actuales en R
, tipo options("contrasts")
. Normalmente, R
ordena los niveles de la variable categórica alfabéticamente y toma el primero como clase de referencia. Esto no siempre es óptimo y puede cambiarse escribiendo (aquí, estableceríamos la clase de referencia como "c" en la nueva variable) new.variable <- relevel(old.variable, ref="c")
. Para cada coeficiente de cada nivel de la variable categórica, un Prueba de Wald se realiza para probar si la diferencia por pares entre el coeficiente de la clase de referencia y la otra clase es diferente de cero o no. Esto es lo que el $z$ y $p$ -valores en la tabla de regresión son. Si sólo una clase categórica es significativa, esto hace no implican que toda la variable carece de sentido y debe eliminarse del modelo. Puede comprobar el efecto global de la variable realizando un prueba de razón de verosimilitud : ajustan dos modelos, uno con y otro sin la variable y el tipo anova(model1, model2, test="LRT")
en R
(véase el ejemplo siguiente). He aquí un ejemplo:
mydata <- read.csv("https://stats.idre.ucla.edu/stat/data/binary.csv")
mydata$rank <- factor(mydata$rank)
my.mod <- glm(admit ~ gre + gpa + rank, data = mydata, family = "binomial")
summary(my.mod)
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.989979 1.139951 -3.500 0.000465 ***
gre 0.002264 0.001094 2.070 0.038465 *
gpa 0.804038 0.331819 2.423 0.015388 *
rank2 -0.675443 0.316490 -2.134 0.032829 *
rank3 -1.340204 0.345306 -3.881 0.000104 ***
rank4 -1.551464 0.417832 -3.713 0.000205 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
El nivel rank1
se ha omitido y cada coeficiente de rank
denota la diferencia entre el coeficiente de rank1
y el correspondiente rank
nivel. Así, la diferencia entre el coeficiente de rank1
y rank2
sería $-0.675$ . El coeficiente de rank1
es simplemente la intercepción. Así que el verdadero coeficiente de rank2
sería $-3.99 - 0.675 = -4.67$ . Las pruebas de Wald prueban aquí si la diferencia entre el coeficiente de la clase de referencia (aquí rank1
) y los niveles correspondientes difieren de cero. En este caso, tenemos pruebas de que los coeficientes de todas las clases difieren del coeficiente de rank1
. También se puede ajustar el modelo sin un intercepto añadiendo - 1
a la fórmula del modelo para ver directamente todos los coeficientes:
my.mod2 <- glm(admit ~ gre + gpa + rank - 1, data = mydata, family = "binomial")
summary(my.mod2) # no intercept model
Coefficients:
Estimate Std. Error z value Pr(>|z|)
gre 0.002264 0.001094 2.070 0.038465 *
gpa 0.804038 0.331819 2.423 0.015388 *
rank1 -3.989979 1.139951 -3.500 0.000465 ***
rank2 -4.665422 1.109370 -4.205 2.61e-05 ***
rank3 -5.330183 1.149538 -4.637 3.54e-06 ***
rank4 -5.541443 1.138072 -4.869 1.12e-06 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Obsérvese que el intercepto ha desaparecido y que el coeficiente de rank1
es exactamente el intercepto del primer modelo. En este caso, la prueba de Wald no comprueba la diferencia entre pares de coeficientes, sino la hipótesis de que cada coeficiente individual es cero. De nuevo, tenemos pruebas de que cada coeficiente de rank
difiere de cero. Por último, para comprobar si toda la variable rank
mejora el ajuste del modelo, ajustamos un modelo con ( my.mod1
) y una sin la variable rank
( my.mod2
) y realizar una prueba de razón de verosimilitud. Esto prueba la hipótesis de que todos los coeficientes de rank
son cero:
my.mod1 <- glm(admit ~ gre + gpa + rank, data = mydata, family = "binomial") # with rank
my.mod2 <- glm(admit ~ gre + gpa, data = mydata, family = "binomial") # without rank
anova(my.mod1, my.mod2, test="LRT")
Analysis of Deviance Table
Model 1: admit ~ gre + gpa + rank
Model 2: admit ~ gre + gpa
Resid. Df Resid. Dev Df Deviance Pr(>Chi)
1 394 458.52
2 397 480.34 -3 -21.826 7.088e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
La prueba de la razón de verosimilitud es altamente significativa y concluiríamos que la variable rank
debe permanecer en el modelo.
Este puesto también es muy interesante.