[R] Different goodness of fit tests leads to contradictory conclusions

[Kiyoshi Sasaki]

I am trying to test goodness of fit for my legalistic regression using several 
options as shown below.  Hosmer-Lemeshow test (whose function I borrowed from 
a previous post), Hosmer–le Cessie omnibus lack of fit test (also borrowed 
from a previous post), Pearson chi-square test, and deviance test.  All the 
tests, except the deviance tests, produced p-values well above 0.05.  Would 
anyone please teach me why deviance test produce p-values such a small value 
(0.001687886) that suggest lack of fit, while other tests suggest good fit? Did 
I do something wrong?
 
Thank you for your time and help!
 
Kiyoshi
 
 
mod.fit <- glm(formula = no.NA$repcnd ~  no.NA$svl, family = binomial(link = 
logit), data =  no.NA, na.action = na.exclude, control = list(epsilon = 
0.0001, maxit = 50, trace = F))
 
> # Option 1: Hosmer-Lemeshow test
> mod.fit <- glm(formula = no.NA$repcnd ~  no.NA$svl, family = binomial(link = 
> logit), data =  no.NA, na.action = na.exclude, control = list(epsilon = 
> 0.0001, maxit = 50, trace = F))
             
 >  hosmerlem <- function (y, yhat, g = 10) 
{
cutyhat <- cut(yhat, breaks = quantile(yhat, probs = seq(0, 1, 1/g)), 
include.lowest = T)
obs <- xtabs(cbind(1 - y, y) ~ cutyhat)
expect <- xtabs(cbind(1 - yhat, yhat) ~ cutyhat)
chisq <- sum((obs - expect)^2/expect)
P <- 1 - pchisq(chisq, g - 2)
c("X^2" = chisq, Df = g - 2, "P(>Chi)" = P)
}
 
> hosmerlem(no.NA$repcnd, fitted(mod.fit))
 X^2                       Df   
                        P(>Chi) 
7.8320107            8.0000000            0.4500497
 
 
> # Option 2 - Hosmer–le Cessie omnibus lack of fit test:
> library(Design)
> lrm.GOF <- lrm(formula = no.NA$repcnd ~  no.NA$svl, data =  no.NA, y = T, x 
> = T)
> resid(lrm.GOF,type = "gof")
Sum of squared errors     Expected value|H0        
SD                     
Z                     P 
                 48.3487115           
               48.3017399            
              0.1060826     0.4427829     0.6579228
 
> # Option 3 - Pearson chi-square p-value:
> pp <- sum(resid(lrm.GOF,typ="pearson")^2)
> 1-pchisq(pp, mod.fit$df.resid)
[1] 0.4623282
 
 
> # Option 4 - Deviance (G^2) test:
> 1-pchisq(mod.fit$deviance, mod.fit$df.resid)
[1] 0.001687886


      
        [[alternative HTML version deleted]]


______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.