Коэффициент отношения правдоподобия (он же отклонение) и критерий несоответствия (или качества соответствия) довольно просто получить для модели логистической регрессии (подгонка с использованием функции) в R. Однако это может быть легко подсчитать количество клеток в конечном итоге достаточно низко, чтобы тест был ненадежным. Один из способов проверить надежность теста отношения правдоподобия при отсутствии соответствия - это сравнить его статистику теста и P- значение с данными теста Пирсона хи-квадрат (или ) на отсутствие соответствия.glm(..., family = binomial)
Ни glm
объект, ни его summary()
метод не сообщают статистику теста для критерия хи-квадрат Пирсона на отсутствие соответствия. В моем поиске единственной вещью, которую я придумал, является chisq.test()
функция (в stats
пакете): ее документация гласит: « chisq.test
выполняет тесты таблицы непредвиденных обстоятельств по критерию хи-квадрат и проверки на соответствие». Тем не менее, документации о том, как выполнять такие тесты, немного:
Если
x
это матрица с одной строкой или столбцом, или еслиx
это вектор иy
он не задан, то выполняется проверка на соответствие критерия соответствия (x
рассматривается как одномерная таблица сопряженности). Записиx
должны быть неотрицательными целыми числами. В этом случае проверяется гипотеза о том, равны ли вероятности совокупности с вероятностямиp
или все они равны, еслиp
не даны.
Я полагаю, что вы можете использовать y
компонент glm
объекта для x
аргумента chisq.test
. Однако вы не можете использовать fitted.values
компонент glm
объекта для p
аргумента chisq.test
, потому что вы получите ошибку: " probabilities must sum to 1.
"
Как я могу (в R) по крайней мере рассчитать статистику теста Пирсона на отсутствие соответствия без необходимости выполнять шаги вручную?
источник
Статистика Пирсона имеет вырожденное распределение, поэтому вообще не рекомендуется для соответствия логистической модели. Я предпочитаю структурированные тесты (линейность, аддитивность). Если вы хотите провести комплексный тест, посмотрите единую степень свободы теста невзвешенной суммы квадратов Ле Сесси - ван Хоуилинген - Копас - Хосмер, реализованную в функции
rms
пакета R.residuals.lrm
источник
ResourceSelection
пакете, и его результат отличается от того, что я получаю при запускеresid(lrm_object, 'gof')
после подгонки моей модели логистической регрессии какlrm_object <- lrm(...)
. Если они действительно разные, можете ли вы прокомментировать, как тест HL соотносится с тем, который вы упомянули здесь? Спасибо!Спасибо, я не осознавал, что это было так просто: сумма (остатки (f1, type = "pearson") ^ 2) Однако, обратите внимание, что остаток Пирсона варьируется в зависимости от того, рассчитывается ли он по ковариантной группе или по отдельным. Простой пример:
m1 - это матрица (это голова более крупной матрицы):
Где x1-3 - предикторы, obs - нет. наблюдения в каждой группе, pi - вероятность членства в группе (предсказано по уравнению регрессии), lev - плечо, диагональ шляпной матрицы, yhat - предсказанное нет. (из у = 1) в группе и у фактического нет.
Это даст вам Пирсона группой. Обратите внимание, как это отличается, если 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)
}
}
таким образом
Если имеется большое количество субъектов с ковариатными шаблонами y = 0, то остаток Пирона будет намного больше при расчете с использованием метода «по группам», а не «по индивидуальным».
См., Например, Hosmer & Lemeshow "Прикладная логистическая регрессия", Wiley, 200.
источник
Вы также можете использовать,
c_hat(mod)
что даст тот же результат, что иsum(residuals(mod, type = "pearson")^2)
.источник
c_hat
?