2016-07-22 7 views
1

Im R Paket {Epi} die ROC() Funktion ein Diagramm aus dem Datensatz aSAH in im {pROC} Paket wie folgt erzeugen:Wie werden die Grenzwerte oder Grenzwerte im Paket {Epi} R ausgewählt?

enter image description here

mit den folgenden Befehlen:

require(Epi) 
require(pROC) 
data(aSAH) 
rock = ROC(form = outcome ~ s100b, data=aSAH, plot = "ROC", MX = T) 

Die Empfindlichkeit und Spezifität wurden für 51 Punkte berechnet, die in dem Objekt nrow(rock$res) enthalten sind. Beachten Sie in diesem Zusammenhang, dass nrow(aSAH) stattdessen 113 ist.

Welche Punkte wurden verwendet, um rock$res zu generieren?

Wenn wir stattdessen die Funktion roc() im Paket {pROC} verwenden, können wir dies über: roc(aSAH$outcome, aSAH$s100b)$threshold. Aber da es sich um verschiedene Pakete handelt, sind sie wahrscheinlich anders.

Antwort

2

Die Antwort ist ... natürlich ... im package documentation:

res Datenrahmen mit Variablen sens, spec, pvp, PVN und den Namen der Test Variable. Letzteres ist der eindeutige Wert des Tests oder linearen Prädiktors aus der logistischen Regression in aufsteigender Reihenfolge mit -Inf vorangestellt. So

was sind die einzigartigen Werte:

points = unique(aSAH$s100b); length(points) [1] 50 plus die Pre-ended -Inf!

Nizza Ahnung davon, aber wir können es beweisen ... Ich denke so:

require(Epi) 
require(pROC) 
data(aSAH) 
rock = ROC(form = outcome ~ s100b, data=aSAH, plot = "ROC", MX = T) 

d = aSAH 
points = sort(unique(d$s100b)) 
length(points) 

## Logistic regression coefficients: 
beta.0 = as.numeric(rock$lr$coefficients[1]) 
beta.1 = as.numeric(rock$lr$coefficients[2]) 

## Sigmoid function: 
sigmoid = 1/(1 + exp(-(beta.0 + beta.1 * points))) 
sigmoid = as.numeric(c("-Inf", sigmoid)) 

lr.eta = rock$res$lr.eta 
length(lr.eta) 
head(lr.eta) 
head(sigmoid) 

> head(lr.eta) 
[1]  -Inf 0.1663429 0.1732556 0.1803934 0.1877585 0.1953526 
> head(sigmoid) 
[1]  -Inf 0.1663429 0.1732556 0.1803934 0.1877585 0.1953526 

## Trying to get the lr.eta number 0.304 on the plot: 
which.max(rowSums(rock$res[, c("sens", "spec")])) 0.30426295405785  18 

## Excluding the "-Inf" introduced by the ROC function (position 17 as opposed to 18): 
max.sens.sp.cut = points[17] 
1/(1 + exp(-(beta.0 + beta.1 * max.sens.sp.cut))) [1] 0.304263 !!! 

Fertig!

+1

@gung Danke für die Bearbeitung! – Toni