Zbudowałem regresję logistyczną, w której zmienna wynikowa jest leczona po otrzymaniu leczenia ( Cure
vs. No Cure
). Wszyscy pacjenci w tym badaniu zostali poddani leczeniu. Interesuje mnie, czy cukrzyca jest związana z tym wynikiem.
W R mój wynik regresji logistycznej wygląda następująco:
Call:
glm(formula = Cure ~ Diabetes, family = binomial(link = "logit"), data = All_patients)
...
Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.2735 0.1306 9.749 <2e-16 ***
Diabetes -0.5597 0.2813 -1.990 0.0466 *
...
Null deviance: 456.55 on 415 degrees of freedom
Residual deviance: 452.75 on 414 degrees of freedom
(2 observations deleted due to missingness)
AIC: 456.75
Jednak przedział ufności dla ilorazu szans obejmuje 1 :
OR 2.5 % 97.5 %
(Intercept) 3.5733333 2.7822031 4.646366
Diabetes 0.5713619 0.3316513 1.003167
Kiedy wykonuję test chi-kwadrat na tych danych, otrzymuję:
data: check
X-squared = 3.4397, df = 1, p-value = 0.06365
Jeśli chcesz to obliczyć samodzielnie, rozkład cukrzycy w wyleczonych i nieutwardzonych grupach przedstawia się następująco:
Diabetic cure rate: 49 / 73 (67%)
Non-diabetic cure rate: 268 / 343 (78%)
Moje pytanie brzmi: dlaczego wartości p i przedział ufności, w tym 1, nie są zgodne?
r
hypothesis-testing
logistic
generalized-linear-model
odds-ratio
SniperBro2000
źródło
źródło
confint()
. Tj. Prawdopodobieństwo zostało sprofilowane. W ten sposób otrzymujesz CI, które są analogiczne do LRT. Twoje obliczenia są prawidłowe, ale zamiast tego stanowią Wald CI. Więcej informacji znajduje się w mojej odpowiedzi poniżej.Odpowiedzi:
Dzięki uogólnionym modelom liniowym można uruchomić trzy różne typy testów statystycznych. Są to: testy Walda, testy współczynnika prawdopodobieństwa i testy punktowe. Doskonała strona pomocy na temat statystyk UCLA omawia je tutaj . Poniższy rysunek (skopiowany z ich witryny) pomaga je zilustrować:
Testy, które są dołączone,1,96 χ2)
summary.glm()
to testy Walda. Nie mówisz, jak uzyskałeś przedziały ufności, ale zakładam, że wykorzystałeśconfint()
, co z kolei wymagaprofile()
. Mówiąc dokładniej, te przedziały ufności są obliczane przez profilowanie prawdopodobieństwa (co jest lepszym podejściem niż pomnożenie SE przez ). Oznacza to, że są one analogiczne do testu współczynnika wiarygodności, a nie testu Walda. Z kolei test jest testem punktacji.Gdy twoje staje się nieskończenie duże, trzy różne powinny zbiegać się na tej samej wartości, ale mogą się nieznacznie różnić, gdy nie masz nieskończonych danych. Warto zauważyć, że wartość (Wald) w początkowym wyniku jest ledwie znacząca i istnieje niewielka prawdziwa różnica między nieco powyżej i poniżej ( cytat ). Ta linia nie jest „magią”. Biorąc pod uwagę, że dwa bardziej wiarygodne testy są nieco ponad , powiedziałbym, że twoje dane nie są „znaczące” według konwencjonalnych kryteriów.N. p p α = 0,05 .05.05
Poniżej profiluję współczynniki na skali predyktora liniowego i przeprowadzam test prawdopodobieństwa jawnie (via
anova.glm()
). Otrzymuję takie same wyniki jak Ty:Jak zauważył @JWilliman w komentarzu (teraz usunięty), w
R
można również uzyskać wartość p opartą na wynikach, używającanova.glm(model, test="Rao")
. W poniższym przykładzie należy pamiętać, że wartość p nie jest całkiem takie same jak w teście chi-kwadrat powyżej, ponieważ domyślnieR
jestchisq.test()
zastosowanie korekty ciągłości. Jeśli zmienimy to ustawienie, wartości p będą zgodne:źródło