0

I have created a logistic regression model and a corresponding ROC using pROC. I got the threshold for the "best" value that maximizes sensitivity and specificity. The predictor is a score that goes from 4 to 13 points. The predicted variable is survival. I need to know which value in my score (from 4 to 13) is represented by the threshold value (e.g. 0.043). Will appreaciate your help.

Code looks like this>

#MULTIVARIATE ANALYSIS

summary(glm((VIVO_AL_ALTA==0)~ CALCULADORA_CALL_SCORE + 
SOPORTE_VENT_AL_INGRESO + ETE_DURANTE_HOSP + 
SOBREINFECC_BACT_DURANTE_HOSP + COP_DURANTE_HOSP + TOCILIZUMAB + 
CORTICOIDES_HOSP, family = binomial, data = work_data))

exp(coef(glm((VIVO_AL_ALTA==0)~ CALCULADORA_CALL_SCORE + 
SOPORTE_VENT_AL_INGRESO + ETE_DURANTE_HOSP + 
SOBREINFECC_BACT_DURANTE_HOSP + COP_DURANTE_HOSP + TOCILIZUMAB + 
CORTICOIDES_HOSP , family = binomial, data = work_data)))

exp(confint.default(glm((VIVO_AL_ALTA==0) ~ 
CALCULADORA_CALL_SCORE + SOPORTE_VENT_AL_INGRESO + 
ETE_DURANTE_HOSP + SOBREINFECC_BACT_DURANTE_HOSP + 
COP_DURANTE_HOSP + TOCILIZUMAB + CORTICOIDES_HOSP , family= 
binomial, data= work_data), level = .95))

mod_vivo_alta_multi<-glm((VIVO_AL_ALTA==0)~ 
CALCULADORA_CALL_SCORE + SOPORTE_VENT_AL_INGRESO + 
ETE_DURANTE_HOSP + SOBREINFECC_BACT_DURANTE_HOSP + 
COP_DURANTE_HOSP + TOCILIZUMAB + CORTICOIDES_HOSP, family = 
binomial, data = work_data, na.action = "na.exclude")

#ROC Curves
library(pROC)

#ROC VIVO ALTA MULTI
work_data$pred_vivo_alta_multi<-predict(mod_vivo_alta_multi, type 
= "response", na.action = "na.omit")
pROC_obj_pred_vivo_alta_multi <- 
roc((work_data$VIVO_AL_ALTA==0),work_data$pred_vivo_alta_multi,
                               smoothed = TRUE, direction="<",
                               # arguments for ci
                               ci=TRUE, ci.alpha=0.95, 
                               stratified=FALSE,
                               # arguments for plot
                               plot=TRUE, auc.polygon=F, 
                               max.auc.polygon=TRUE, grid=TRUE, 
                               print.thres=T,
                               print.auc=TRUE, show.thres=TRUE)

coords(pROC_obj_pred_vivo_alta_multi,x= "best", 
input="threshold", ret=c("threshold", "specificity", 
"sensitivity", "npv", "ppv","youden"), as.list=FALSE, drop=TRUE, 
best.method=c("youden"), best.weights=c(1, 0.5), transpose = 
FALSE, as.matrix=FALSE)

Sorry for some of the variables are in Spanish. Basically, these are my variables for the model: CALCULADORA_CALL_SCORE + SOPORTE_VENT_AL_INGRESO + ETE_DURANTE_HOSP + SOBREINFECC_BACT_DURANTE_HOSP + COP_DURANTE_HOSP + TOCILIZUMAB + CORTICOIDES_HOSP

And my predicted variable is VIVO_AL_ALTA

  • It's easier to help you if you provide a [reproducible example](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) with sample input and desired output that can be used to test and verify possible solutions. – MrFlick May 20 '22 at 18:42

1 Answers1

0

You can use coords(my_roc, "best").

Here's a demonstration using a built in example from pROC

library(pROC)

ROC <- roc(aSAH$outcome, aSAH$s100b, levels=c("Good", "Poor"))
#> Setting direction: controls < cases

coords(ROC, "best")
#>   threshold specificity sensitivity
#> 1     0.205   0.8055556   0.6341463

Created on 2022-05-20 by the reprex package (v2.0.1)

Allan Cameron
  • 147,086
  • 7
  • 49
  • 87
  • Thanks for your response. It does work for a univariate logistic regression. However, for the next analysis I created a multivariable logistic regression with 7 predictor variables and the predicted variable is still survival. I used this model with the function 'predict' to get the predicted probability for each of my scores. Then using 'pROC', I created a 'ROC object'. Using 'coords' I got the "best" probability (e.g. 0.067). How can I transform this probability to my score that goes from 4 to 13? Thanks again for your help. – Juan Cristóbal Pedemonte May 21 '22 at 11:24
  • Why not use the score directly for the roc instead of the output of your logistic regression probabilities? That would be the usual way to do it – Allan Cameron May 21 '22 at 11:32
  • How can I include the score plus my 6 covariables (SCORE + X1 + X2 + X3 ... X6) without including the predicted probabilities? I don't know how to do this. Thanks again. – Juan Cristóbal Pedemonte May 21 '22 at 20:10
  • Oh I thought your score was _made_ from your co-variables. Perhaps you could share your data / model – Allan Cameron May 21 '22 at 20:18
  • Don't attempt in the comments - edit your question to include please – Allan Cameron May 21 '22 at 20:44
  • I edited my question with the corresponding code. I hope it is easier to follow know. Thanks again! – Juan Cristóbal Pedemonte May 21 '22 at 21:12