Извлечение вероятности и SE из логистической регрессии - PullRequest
0 голосов
/ 01 июля 2018

У меня есть набор данных для выбора задачи (1 или 0) с переменной x. Для использования mtcars в качестве примера

#binomial_smooth() from https://ggplot2.tidyverse.org/reference/geom_smooth.html
binomial_smooth <- function(...) {
  geom_smooth(method = "glm", method.args = list(family = "binomial"), ...)
}

#plot
ggplot(data = mtcars, aes(x = disp, y = am)) + geom_point(alpha = 0.5) + binomial_smooth()

#create a model
model <- glm(am ~ disp, family = "binomial", data = mtcars)

enter image description here

где am будет выбор субъекта и переменная x. Я хотел бы определить значение x +/- SE (которое я представляю, каково то, что binomial_smooth строит, хотя я могу ошибаться) для двоичной переменной = 0.5.

Используя mtcars, я бы хотел выяснить, для каких дисплеев +/- SE am = 0,5. Огляделся и просто запутался, поэтому любая помощь будет высоко ценится!

лучший

1 Ответ

0 голосов
/ 05 июля 2018

Хорошо, я понял это после того, как проследил за кроличьей ношей от Романа Луштрика (ура!).

Использует пакет MASS и функцию, используемую для расчета LD50. Также можно вручную выбрать значение p для поиска.

library(ggplot2)
library(MASS)
#binomial_smooth() from https://ggplot2.tidyverse.org/reference/geom_smooth.html
binomial_smooth <- function(...) {
  geom_smooth(method = "glm", method.args = list(family = "binomial"), ...)
}

#create a model
model <- glm(am ~ disp, family = "binomial", data = mtcars)

#get the 'LD50'- the point at which the binomial regression crosses 50%
LD50 <- dose.p(model, p = 0.5)

#print the details
print(LD50)

#replot the figure with the LD50 vlines
ggplot(data = mtcars, aes(x = disp, y = am)) + 
  geom_point(alpha = 0.5) + 
  binomial_smooth() +
  geom_vline(xintercept = LD50[[1]])

enter image description here

...