Regresión logística: prueba de chi-cuadrado de anova versus significancia de coeficientes (anova () versus resumen () en R)

35

Tengo un modelo logístico GLM con 8 variables. Realicé una prueba de chi-cuadrado en R anova(glm.model,test='Chisq')y 2 de las variables resultan ser predictivas cuando se ordenan en la parte superior de la prueba y no tanto cuando se ordenan en la parte inferior. El summary(glm.model)sugiere que sus coeficientes son insignificantes (alto valor de p). En este caso parece que las variables no son significativas.

Quería preguntar cuál es una mejor prueba de significación de las variables: la significación del coeficiente en el resumen del modelo o la prueba de ji al cuadrado anova(). Además, ¿cuándo es uno mejor que el otro?

Supongo que es una pregunta amplia, pero cualquier sugerencia sobre qué considerar será apreciada.

StreetHawk
fuente
44
Esto es análogo a la distinción entre las sumas de cuadrados de tipo I y tipo III para probar coeficientes en modelos lineales. Podría ayudarlo a leer mi respuesta aquí: cómo interpretar ANOVA secuencial tipo I y MANOVA .
gung - Restablece a Monica

Respuestas:

61

Además de la respuesta de @ gung, intentaré proporcionar un ejemplo de lo anovaque realmente prueba la función. Espero que esto le permita decidir qué pruebas son apropiadas para las hipótesis que le interesan probar.

Supongamos que tiene un resultado y 3 variables predictoras: x 1 , x 2 y x 3 . Ahora, si tu modelo de regresión logística fuera . Cuando ejecuta , la función compara los siguientes modelos en orden secuencial:yX1X2X3my.mod <- glm(y~x1+x2+x3, family="binomial")anova(my.mod, test="Chisq")

  1. glm(y~1, family="binomial") vs. glm(y~x1, family="binomial")
  2. glm(y~x1, family="binomial") vs. glm(y~x1+x2, family="binomial")
  3. glm(y~x1+x2, family="binomial") vs. glm(y~x1+x2+x3, family="binomial")

Por lo tanto, compara secuencialmente el modelo más pequeño con el siguiente modelo más complejo al agregar una variable en cada paso. Cada una de esas comparaciones se realiza mediante una prueba de razón de verosimilitud (prueba LR; ver el ejemplo a continuación). Que yo sepa, estas hipótesis rara vez son de interés, pero esto debe decidirlo usted.

Aquí hay un ejemplo en R:

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 

# The sequential analysis
anova(my.mod, test="Chisq")

Terms added sequentially (first to last)    

     Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
NULL                   399     499.98              
gre   1  13.9204       398     486.06 0.0001907 ***
gpa   1   5.7122       397     480.34 0.0168478 *  
rank  3  21.8265       394     458.52 7.088e-05 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

# We can make the comparisons by hand (adding a variable in each step)

  # model only the intercept
mod1 <- glm(admit ~ 1,                data = mydata, family = "binomial") 
  # model with intercept + gre
mod2 <- glm(admit ~ gre,              data = mydata, family = "binomial") 
  # model with intercept + gre + gpa
mod3 <- glm(admit ~ gre + gpa,        data = mydata, family = "binomial") 
  # model containing all variables (full model)
mod4 <- glm(admit ~ gre + gpa + rank, data = mydata, family = "binomial") 

anova(mod1, mod2, test="LRT")

Model 1: admit ~ 1
Model 2: admit ~ gre
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
1       399     499.98                          
2       398     486.06  1    13.92 0.0001907 ***

anova(mod2, mod3, test="LRT")

Model 1: admit ~ gre
Model 2: admit ~ gre + gpa
  Resid. Df Resid. Dev Df Deviance Pr(>Chi)  
1       398     486.06                       
2       397     480.34  1   5.7122  0.01685 *

anova(mod3, mod4, test="LRT")

Model 1: admit ~ gre + gpa
Model 2: admit ~ gre + gpa + rank
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
1       397     480.34                          
2       394     458.52  3   21.826 7.088e-05 ***

pagssummary(my.mod)

  • Para coeficiente de x1: glm(y~x2+x3, family="binomial")vs. glm(y~x1+x2+x3, family="binomial")
  • Para coeficiente de x2: glm(y~x1+x3, family="binomial")vs.glm(y~x1+x2+x3, family="binomial")
  • Para coeficiente de x3: glm(y~x1+x2, family="binomial")vs.glm(y~x1+x2+x3, family="binomial")

Entonces cada coeficiente contra el modelo completo que contiene todos los coeficientes. Las pruebas de Wald son una aproximación de la prueba de razón de probabilidad. También podríamos hacer las pruebas de razón de probabilidad (prueba LR). Aquí es cómo:

mod1.2 <- glm(admit ~ gre + gpa,  data = mydata, family = "binomial")
mod2.2 <- glm(admit ~ gre + rank, data = mydata, family = "binomial")
mod3.2 <- glm(admit ~ gpa + rank, data = mydata, family = "binomial")

anova(mod1.2, my.mod, test="LRT") # joint LR test for rank

Model 1: admit ~ gre + gpa
Model 2: admit ~ gre + gpa + rank
  Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
1       397     480.34                          
2       394     458.52  3   21.826 7.088e-05 ***

anova(mod2.2, my.mod, test="LRT") # LR test for gpa

Model 1: admit ~ gre + rank
Model 2: admit ~ gre + gpa + rank
  Resid. Df Resid. Dev Df Deviance Pr(>Chi)  
1       395     464.53                       
2       394     458.52  1   6.0143  0.01419 *

anova(mod3.2, my.mod, test="LRT") # LR test for gre

Model 1: admit ~ gpa + rank
Model 2: admit ~ gre + gpa + rank
  Resid. Df Resid. Dev Df Deviance Pr(>Chi)  
1       395     462.88                       
2       394     458.52  1   4.3578  0.03684 *

pagssummary(my.mod)

rankanova(my.mod, test="Chisq")rankanova(mod1.2, my.mod, test="Chisq")pags7.08810-5 5rank

COOLSerdash
fuente
1
+1, esta es una explicación buena y completa. 1 pequeño punto: creo que cuando test="Chisq"no estás ejecutando una prueba de razón de probabilidad, debes configurar test="LRT"eso, ¿ ves ? Anova.glm .
gung - Restablece a Monica
66
@gung Gracias por el cumplido. test="LRT"y test="Chisq"son sinónimos (lo dice en la página que ha vinculado).
COOLSerdash
2
No hay problema, pero creo que en realidad es un buen punto. test="LRT"es mejor ya que está claro de inmediato que es una prueba de razón de probabilidad. Lo cambie. Gracias.
COOLSerdash
44
+1 Estoy impresionado con su rápido progreso aquí en solo un mes y su capacidad para proporcionar una explicación clara y bien trabajada. ¡Gracias por tus esfuerzos!
whuber
1
Gran respuesta. ¿Puedo preguntar cómo 7.088e-05, 0.01419, 00.03684deben interpretarse los valores p ( )?
TheSimpliFire