Логистическая регрессия: критерий хи-квадрат anova против значимости коэффициентов (anova () против суммарного () в R)

35

У меня есть логистическая модель GLM с 8 переменными. Я anova(glm.model,test='Chisq')выполнил тест хи-квадрат в R, и 2 переменные оказываются прогнозирующими, если их упорядочивать в верхней части теста, и не так сильно, когда их упорядочивают в нижней части. Предполагается, summary(glm.model)что их коэффициенты незначительны (высокое значение p). В этом случае кажется, что переменные не являются значимыми.

Я хотел спросить , что является лучшим тестом переменных значимости - коэффициент значимости в кратком изложении модели или хи-квадрат теста с anova(). Кроме того - когда один из них лучше другого?

Я предполагаю, что это широкий вопрос, но любые советы о том, что следует учитывать, будут признательны.

StreetHawk
источник
4
Это аналогично различию между суммами квадратов типа I и типа III для проверки коэффициентов в линейных моделях. Это может помочь вам прочитать мой ответ здесь: как интерпретировать тип I Последовательная ANOVA и MANOVA .
gung - Восстановить Монику

Ответы:

61

В дополнение к ответу @ gung я постараюсь привести пример того, что anova самом деле тестирует функция. Я надеюсь, что это позволит вам решить, какие тесты подходят для гипотез, которые вы интересуетесь тестированием.

yx1x2x3my.mod <- glm(y~x1+x2+x3, family="binomial")anova(my.mod, test="Chisq")

  1. glm(y~1, family="binomial") против glm(y~x1, family="binomial")
  2. glm(y~x1, family="binomial") против glm(y~x1+x2, family="binomial")
  3. glm(y~x1+x2, family="binomial") против glm(y~x1+x2+x3, family="binomial")

Таким образом, он последовательно сравнивает меньшую модель со следующей более сложной моделью, добавляя одну переменную на каждом шаге. Каждое из этих сравнений делается через теста отношения правдоподобия (тест LR; см. Пример ниже). Насколько мне известно, эти гипотезы редко представляют интерес, но вы должны решить это.

Вот пример в 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 ***

п-значениями в выходных данных summary(my.mod)являются тесты Вальда, которые проверяют следующие гипотезы (обратите внимание, что они взаимозаменяемы и порядок тестов не имеет значения ):

  • Для коэффициента x1:glm(y~x2+x3, family="binomial") против glm(y~x1+x2+x3, family="binomial")
  • Для коэффициента x2:glm(y~x1+x3, family="binomial") противglm(y~x1+x2+x3, family="binomial")
  • Для коэффициента x3:glm(y~x1+x2, family="binomial") противglm(y~x1+x2+x3, family="binomial")

Таким образом, каждый коэффициент против полной модели, содержащей все коэффициенты. Тесты Вальда являются приближением теста отношения правдоподобия. Мы также могли бы сделать тесты отношения правдоподобия (тест LR). Вот как:

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 *

п-значения из тестов отношения правдоподобия очень похожи на те, которые получены тестами Вальда summary(my.mod)выше.

Примечание. Сравнение третьей модели для rankof anova(my.mod, test="Chisq")- то же самое, что и сравнение для rankпримера ниже ( anova(mod1.2, my.mod, test="Chisq")). Каждый разп-значение одинаковое, 7,08810-5, Каждый раз проводится сравнение между моделью без модели rankи моделью, содержащей ее.

COOLSerdash
источник
1
+1, это хорошее, исчерпывающее объяснение. Небольшое замечание: я считаю, что когда test="Chisq"вы не проводите тест отношения правдоподобия, вам нужно установить test="LRT"это, см. ? Anova.glm .
gung - Восстановить Монику
6
@ Gung Спасибо за комплимент. test="LRT"и test="Chisq"являются синонимами (это говорится на странице, на которую вы ссылаетесь).
COOLSerdash
2
Нет проблем, но я думаю, что это хороший момент. test="LRT"лучше, поскольку сразу становится ясно, что это критерий отношения правдоподобия. Я изменил это. Спасибо.
COOLSerdash
4
+1 Я впечатлен вашим быстрым прогрессом здесь всего за один месяц и вашей способностью дать хорошо проработанное, четкое объяснение. Спасибо за ваши старания!
whuber
1
Отличный ответ. Могу ли я спросить, как 7.088e-05, 0.01419, 00.03684следует интерпретировать p-значения ( )?
TheSimpliFire