Uzyskiwanie przewidywanych wartości (Y = 1 lub 0) z dopasowania modelu regresji logistycznej

Powiedzmy, że mam obiekt klasy glm(odpowiadający modelowi regresji logistycznej) i chciałbym zamienić przewidywane prawdopodobieństwa podane za predict.glmpomocą argumentu type="response"na odpowiedzi binarne, tj. lub . Jaki jest najszybszy i najbardziej kanoniczny sposób to zrobić w R?$Y=1$$Y=0$

Chociaż znowu jestem tego świadomy predict.glm, nie wiem, gdzie dokładnie znajduje się wartość odcięcia - i myślę, że to jest mój główny problem.$P(Y_i=1|\hat X_{i})$

Kiedy masz przewidywane prawdopodobieństwo, od Ciebie zależy, jakiego progu chciałbyś użyć. Możesz wybrać próg, aby zoptymalizować czułość, specyficzność lub jakąkolwiek miarę, która jest najważniejsza w kontekście aplikacji (niektóre dodatkowe informacje byłyby tutaj pomocne w celu uzyskania bardziej szczegółowej odpowiedzi). Możesz przyjrzeć się krzywym ROC i innym miernikom związanym z optymalną klasyfikacją.

Edycja: Aby nieco wyjaśnić tę odpowiedź, podam przykład. Prawdziwa odpowiedź jest taka, że ​​optymalne odcięcie zależy od tego, jakie właściwości klasyfikatora są ważne w kontekście aplikacji. Niech będzie prawdziwą wartością dla obserwacji , a będzie przewidywaną klasą. Niektóre typowe miary wydajności to$Y_{i}$$i$$\hat{Y}_{i}$

(1) Czułość: - proporcja „1” poprawnie zidentyfikowanych jako tak.$P(\hat{Y}_i=1 | Y_i=1)$

(2) Swoistość: - proporcja , które są poprawnie zidentyfikowane jako tak$P(\hat{Y}_i=0 | Y_i=0)$

(3) (Prawidłowo) Współczynnik klasyfikacji: - odsetek poprawnych prognoz.$P(Y_i = \hat{Y}_i)$

(1) jest również nazywany prawdziwą częstością dodatnią, (2) jest również nazywany prawdziwą częstością ujemną.

Na przykład, jeśli twój klasyfikator dążył do oceny testu diagnostycznego dla poważnej choroby, która ma stosunkowo bezpieczne wyleczenie, czułość jest znacznie ważniejsza niż specyficzność. W innym przypadku, jeśli choroba byłaby stosunkowo niewielka, a leczenie ryzykowne, ważniejsza byłaby specyficzność. W przypadku ogólnych problemów z klasyfikacją uważa się za „dobre” wspólne optymalizowanie czułości i specyfikacji - na przykład można użyć klasyfikatora, który minimalizuje ich odległość euklidesową od punktu :$(1,1)$

$\delta$ może być ważona lub modyfikowana w inny sposób, aby odzwierciedlić bardziej rozsądną miarę odległości od w kontekście aplikacji - odległość euklidesowa od (1,1) została tu wybrana arbitralnie w celach ilustracyjnych. W każdym razie wszystkie te cztery środki mogą być najbardziej odpowiednie, w zależności od zastosowania.$(1,1)$

Poniżej znajduje się symulowany przykład wykorzystujący prognozę z modelu regresji logistycznej do klasyfikacji. Wartość graniczna jest zróżnicowana, aby zobaczyć, która wartość graniczna daje „najlepszy” klasyfikator w ramach każdej z tych trzech miar. W tym przykładzie dane pochodzą z modelu regresji logistycznej z trzema predyktorami (patrz wykres R poniżej wykresu). Jak widać z tego przykładu, „optymalne” odcięcie zależy od tego, który z tych środków jest najważniejszy - jest to całkowicie zależne od aplikacji.

Edycja 2: i , dodatnia wartość predykcyjna i ujemna wartość predykcyjna (zauważ, że NIE są one takie same jako czułość i swoistość) mogą być również użytecznymi miernikami wydajności.$P(Y_i = 1 | \hat{Y}_i = 1)$$P(Y_i = 0 | \hat{Y}_i = 0)$

# data y simulated from a logistic regression model 
# with with three predictors, n=10000
x = matrix(rnorm(30000),10000,3)
lp = 0 + x[,1] - 1.42*x[2] + .67*x[,3] + 1.1*x[,1]*x[,2] - 1.5*x[,1]*x[,3] +2.2*x[,2]*x[,3] + x[,1]*x[,2]*x[,3]
p = 1/(1+exp(-lp))
y = runif(10000)<p

# fit a logistic regression model
mod = glm(y~x[,1]*x[,2]*x[,3],family="binomial")

# using a cutoff of cut, calculate sensitivity, specificity, and classification rate
perf = function(cut, mod, y)
{
   yhat = (mod$fit>cut)
   w = which(y==1)
   sensitivity = mean( yhat[w] == 1 ) 
   specificity = mean( yhat[-w] == 0 ) 
   c.rate = mean( y==yhat ) 
   d = cbind(sensitivity,specificity)-c(1,1)
   d = sqrt( d[1]^2 + d[2]^2 ) 
   out = t(as.matrix(c(sensitivity, specificity, c.rate,d)))
   colnames(out) = c("sensitivity", "specificity", "c.rate", "distance")
   return(out)
}

s = seq(.01,.99,length=1000)
OUT = matrix(0,1000,4)
for(i in 1:1000) OUT[i,]=perf(s[i],mod,y)
plot(s,OUT[,1],xlab="Cutoff",ylab="Value",cex.lab=1.5,cex.axis=1.5,ylim=c(0,1),type="l",lwd=2,axes=FALSE,col=2)
axis(1,seq(0,1,length=5),seq(0,1,length=5),cex.lab=1.5)
axis(2,seq(0,1,length=5),seq(0,1,length=5),cex.lab=1.5)
lines(s,OUT[,2],col="darkgreen",lwd=2)
lines(s,OUT[,3],col=4,lwd=2)
lines(s,OUT[,4],col="darkred",lwd=2)
box()
legend(0,.25,col=c(2,"darkgreen",4,"darkred"),lwd=c(2,2,2,2),c("Sensitivity","Specificity","Classification Rate","Distance"))