Избегание цикла for в R в попытке оценить процент истинных положительных / отрицательных значений при использовании логистической регрессии

Question

Что я получил: матрица, в которой я получил прогнозируемую вероятность результата (из модели логистической регрессии) и известный результат. Для этих любопытных я фактически получил две регрессионные модели и независимый набор тестовых данных, в котором я хотел бы сравнить эти две модели, выполнив это.

> head(matrixComb)
      probComb outComb
[1,] 0.9999902       1
[2,] 0.9921736       0
[3,] 0.9901175       1
[4,] 0.9815581       0
[5,] 0.7692992       0
[6,] 0.7369990       0

Что я хочу: график, на котором я могу построить график того, как часто моя модель прогнозирования дает правильные результаты (одна строка для положительных значений и одна строка для отрицательных значений) в зависимости от предельного значения вероятности. Моя проблема в том, что я не могу понять, как это сделать, не переключаясь на Perl и не используя цикл For для перебора матрицы.

В Perl я бы просто начал с вероятности 0,1, а при прогоне цикла for увеличивал значение на 0,1. На первой итерации я буду считать все вероятности <0,1 и исход = 0 как истинные отрицания, вероятность <0,1 и результат 1 как вероятность ложных отрицаний> 0,1, а результат = 0 - как ложные срабатывания, а вероятность> 0,1 и результат = 1 - как истинные положительные результаты

Затем процесс будет повторяться, и результаты каждой итерации будут печататься как [вероятность, истинные положительные результаты / общее количество положительных результатов, истинные отрицательные значения / общее количество отрицательных значений]. Таким образом, мне будет легко распечатать его в open office calc.

Причина, по которой я спрашиваю это, состоит в том, что операция слишком сложна для меня, чтобы найти подобный случай здесь в stackoverflow или в учебном пособии. Но мне бы очень хотелось узнать, как сделать это эффективно в среде R.

Answer 1

Вот способ сделать это вручную:

#Create some sample data
dat <- data.frame(x=runif(100),y=sample(0:1,100,replace=TRUE))

#Function to compute tp and tn
myFun <- function(x){
    tbl <- table(dat$x > x,dat$y)
    marg <- margin.table(tbl,2)
    tn <- tbl[1,1]/marg[1]
    tp <- tbl[2,2]/marg[2]
    rs <- c(tp,tn)
    names(rs) <- c('truePos','trueNeg')
    return(rs)
}


#Decision thresholds
thresh <- seq(0.1,0.9, by = 0.1)
#Loop using lapply
temp <- as.data.frame(do.call(rbind,lapply(thresh,myFun)))
temp$thresh <- thresh

#Melt and plot using ggplot
tempMelt <- melt(temp,id.vars="thresh")
ggplot(tempMelt,aes(x=thresh,y=value)) + 
    geom_line(aes(group=variable,colour=variable))

В качестве альтернативы, как упомянуто выше в комментариях, существует множество функций ROC или ROC, которые могутнайдено с помощью ??ROC.Например, используя roc из пакета caret:

temp <- as.data.frame(roc(dat$x,factor(dat$y)))
tempMelt <- melt(temp,id.vars="cutoff")
ggplot(tempMelt,aes(x=cutoff,y=value)) + 
    geom_line(aes(group=variable,colour=variable))

Answer 2

Вы можете получить R, чтобы нарисовать кривые, основанные на анализе ROC. Это грубая версия, использующая пакет ROCR, и ее можно легко сделать красивее

ss <- 1000   # sample size
mydf <- data.frame(probComb = runif(ss)) # predictions illustration
mydf$outComb <- 0 + (runif(ss) < mydf$probComb) # actuals illustration

library(ROCR)
pred <- prediction(mydf$probComb, mydf$outComb)
perfp <- performance(pred, "tpr")
perfn <- performance(pred, "tnr")
plot(perfp, col="green", ylab="True positive (green) and true negative (red) rates")
plot(perfn, col="red", ylab="True negative rate", add=TRUE)

для производства

Если нужно, вы можете найти данные в perfp и perfn.

Answer 3

Может быть, что-то вроде этого:

# A function for counting outcomes for a certain probability
f <- function(d, p) {
  lp <- d$prob < p
  c(TNeg=sum(lp & d$out==0), TPos=sum(!lp & d$out==1))
}

# Make it accept a vector of probabilities
vf <- Vectorize(f, 'p')

# Sample data
n <- 100
d <- data.frame(prob=runif(n), out=round(runif(n)))
# Probabilities to plot
p <- seq(0,1, len=20)

res <- vf(d, p)
colnames(res) <- paste('p(', p, ')', sep='')
matplot(p, t(res), type='l', xlab='prob', ylab='count')