Wskaźnik prawdopodobieństwa (inaczej dewiacja) Statystyka i test braku dopasowania (lub dobroci dopasowania) jest dość prosty do uzyskania dla modelu regresji logistycznej (dopasowanie przy użyciu funkcji) w R. Jednak może być łatwe jest, aby niektóre liczby komórek były wystarczająco niskie, aby test był niewiarygodny. Jednym ze sposobów weryfikacji wiarygodności testu współczynnika wiarygodności dla braku dopasowania jest porównanie jego statystyki testowej i wartości P z wynikami testu chi-Pearsona (lub ) braku dopasowania.glm(..., family = binomial)
Ani glm
obiekt, ani jego summary()
metoda nie podają statystyki testu testu chi kwadrat Pearsona z powodu braku dopasowania. W moich poszukiwaniach jedyną rzeczą, jaką wymyśliłem, jest chisq.test()
funkcja (w stats
pakiecie): w jej dokumentacji jest napisane: „ chisq.test
wykonuje testy tabeli kontyngencji chi-kwadrat i testy dopasowania”. Jednak dokumentacja jest ograniczona, jak przeprowadzić takie testy:
Jeśli
x
jest to macierz z jednym rzędem lub kolumną, lub jeślix
jest wektorem iy
nie jest podana, wówczas przeprowadzany jest test zgodności pod względem dopasowania (x
jest traktowany jako jednowymiarowa tabela kontyngencji). Wpisyx
muszą być liczbami całkowitymi nieujemnymi. W tym przypadku testowana hipoteza sprawdza, czy prawdopodobieństwo populacji jest równe tymp
, czy są równe, jeślip
nie zostanie podane.
Wyobrażam sobie, że możesz użyć y
komponentu glm
obiektu jako x
argumentu chisq.test
. Nie można jednak użyć fitted.values
komponentu glm
obiektu jako p
argumentu chisq.test
, ponieważ pojawi się błąd: „ probabilities must sum to 1.
”
Jak mogę (w R) przynajmniej obliczyć statystykę testową Pearson kątem braku dopasowania bez konieczności ręcznego wykonywania kroków?
źródło
Statystyka Pearsona ma rozkład zdegenerowany, więc nie jest ogólnie zalecana dla dopasowania modelu logistycznego. Wolę testy strukturalne (liniowość, addytywność). Jeśli chcesz testu zbiorczego, zobacz pojedynczy stopień swobody le Cessie - van Houwelingen - Copas - Hosmer test nieważonej sumy kwadratów zaimplementowany w funkcji
rms
pakietu R.residuals.lrm
źródło
ResourceSelection
pakiecie, a jego wynik różni się od tego, co otrzymujęresid(lrm_object, 'gof')
po uruchomieniu mojego modelu regresji logistycznej jakolrm_object <- lrm(...)
. Jeśli rzeczywiście są różne, czy możesz skomentować, w jaki sposób test HL ma się do tego, o którym tu wspominasz? Dziękuję Ci!Dzięki, nie zdawałem sobie sprawy, że to było tak proste jak: suma (reszty (f1, type = "pearson") ^ 2) Jednak zauważ, że reszta Pearsona różni się w zależności od tego, czy jest obliczana według grupy współzmiennej, czy indywidualnie. Prosty przykład:
m1 jest matrycą (ta jest głową większej matrycy):
Gdzie x1-3 są predyktorami, obs jest nie. obserwacje w każdej grupie, pi jest prawdopodobieństwem przynależności do grupy (przewidywanym na podstawie równania regresji), lew jest dźwignią, przekątną macierzy kapelusza, czyli przewidywaną nie. (y = 1) w grupie iy faktyczna liczba nie.
To da ci Pearsona według grup. Zauważ, że jest inaczej, jeśli y == 0: ' 'fun1 <- function(j){
if (m1[j,"y"] ==0){ # y=0 for this covariate pattern
Pr1 <- sqrt( m1[i,"pi"] / (1-m1[i,"pi"]))
Pr2 <- -sqrt (m1[i,"obs"])
res <- round( Pr1 * Pr2, 3)
return(res)
} else {
Pr1 <- m1[j,"y"] - m1[j,"yhat"]
Pr2 <- sqrt( m1[j,"yhat"] * ( 1-(m1[j,"pi"]) ) )
res <- round( Pr1/Pr2, 3)
return(res)
}
}
A zatem
Jeśli istnieje duża liczba osobników o wzorcach współzmiennych y = 0, wówczas resztkowa wartość Pearons będzie znacznie większa, jeśli zostanie obliczona przy użyciu metody „według grupy”, a nie „według osoby”.
Patrz np. Hosmer & Lemeshow „Applied Logistic Regression”, Wiley, 200.
źródło
Możesz także użyć
c_hat(mod)
tego, co da ten sam wynik, cosum(residuals(mod, type = "pearson")^2)
.źródło
c_hat
się?