Ho to make the Receiver operating characteristic (ROC) and identifying the Youden index

user_9323156 0
Viewed 60

I would like to make a ROC curve and identify the Youden-Index. I have a subsample as below. Where I need to find cut point for the column "val" based on the reference column "ref". How can I make the ROC plot, identifying the cutoff with specificities and sensitivities as well as Youden Index, in Stata or R?

cut_= structure(list(val = c("2.5", "3.5", "1.5", "1.1", "2.4", "1.6", 
    "1.9", "2.7", "1.2", "1.5", "2.1", "1.4", "1.8", "3.5", "2.5", 
    "2.4"), ref = c(1, 0, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 0, 
    1)), class = c("tbl_df", "tbl", "data.frame"), row.names = c(NA, 
    -16L))
roc stata r
2 Answers

As mentioned in the comments, it is pretty easy to obtain in R with pROC. You need to convert val to a numeric vector first, then you can create the ROC curve, let's call it cut_roc:

cut_$val <- as.numeric(cut_$val)
library(pROC)
cut_roc <- roc(cut_$ref, cut_$val)

Then it's as simple as a call to coords with x="best" to get the best threshold(s) (youden is the default so the last argument is optional):

coords(cut_roc, x="best", best.method="youden")
#   threshold specificity sensitivity
# 1      1.15       0.125       1.000
# 2      1.45       0.250       0.875
# 3      2.25       0.625       0.500
# 4      2.45       0.750       0.375
# 5      2.60       0.875       0.250

Note that in this specific ROC curve, multiple points of the curve maximize the Youden index.

user_3335990

Consider this a comment, not an answer; see the answer by @Calimo instead.

I thought it might be interesting to plot the ROC curve with the "best" thresholds. Unfortunately, as the figure shows, the binary classifier behind this ROC curve doesn't have much discriminative ability at any threshold.

Here is some R code. Run it after @Calimo's code which shows how to create the cut_roc object.

par(pty = "s") # Make this ROC plot square
plot(cut_roc, print.thres = "best", print.thres.best.method = "youden")


par(pty = "m") # Reset the `pty` parameter

Created on 2022-07-28 by the reprex package (v2.0.1)

user_95664840
Related