Recherche d’indices d’exposition pouvant expliquer l’apparition de nouveau diabète post-transplantation
##Chargement des packages nécessaires
load("NODAT.RData")
library(rmarkdown)
library(knitr)
library(kableExtra)
library(dplyr)
library(plyr)
library(purrr)
library(readr)
library(tidyr)
library(ggplot2)
library(pastecs)
library(AUC)
library(stringr)
library(reshape)
library(survival)
library(ggfortify)
library(survminer)
library(broom)
library(GGally, quietly = TRUE)
library(pROC)
library(rstatix)
library(tidymodels)
library(ranger)
library(randomForest)
library(magick)
library(funModeling)Création d’une fonction permettant de calculer la durée au dessus d’un vecteur de seuil
La fonction s’appelle “Tps_sup”, elle permet de calculer le temps passé supérieur à une valeur seuil que nous avons déterminé. Grâce à la fonction map_df du package purrr, nous pouvons appliquer un vecteur à cette fonction pour déterminer plusieurs temps cumulés à différents seuils.
# Seuil cible
Threshold<-c(seq(10,25,2.5))
## fonction
Tps_sup<-function(Threshold){
IntervalRatio<-ifelse(Threshold==0,1,0.5)
nom<-unique(NODATvubft$id)
rest<-data.frame()
for (i in 1:length(nom)){
idpat<-nom[i]
patient<-subset(NODATvubft, id == idpat)
previousVisit<-NULL
cumulatedTime<-0
# ItÈration sur les visites
for (j in 1:nrow(patient)) {
visite<-patient[j,]
delai<-visite$delai_C0_Y_cc
delta<-ifelse(is.null(previousVisit), 0, delai-previousVisit$delai_C0_Y_cc)
intervalOfInterest<-delta*IntervalRatio
cumulatedTime<- cumulatedTime + ifelse(visite$CO_cc >= Threshold, intervalOfInterest, 0) # valeur cible ‡ dÈfinir essayer sinon ==
previousVisit<-visite
}
rest<-data.frame(rbind(rest, data.frame(id=idpat, CumulatedTime=cumulatedTime, Threshold=Threshold)))
print(rest)
}
return(rest)}
reg<-map_df(Threshold,Tps_sup)## préparation données pour courbe de survie et modèle de COX
reg1<-pivot_wider(reg, names_from = Threshold, values_from = CumulatedTime)
survivalnodat15d<-merge(survivalnodat15c,reg1, by="id")
survivalnodat15d<-survivalnodat15d%>% select(-"15")
write_csv(survivalnodat15d, path="survivalnodat15d.csv")
survivalnodat15d$`10_semaine`<-survivalnodat15d$`10`*52
survivalnodat15d$`12.5_semaine`<-survivalnodat15d$`12.5`*52
survivalnodat15d$`17.5_semaine`<-survivalnodat15d$`17.5`*52
survivalnodat15d$`20_semaine`<-survivalnodat15d$`20`*52
survivalnodat15d$`22.5_semaine`<-survivalnodat15d$`22.5`*52
survivalnodat15d$`25_semaine`<-survivalnodat15d$`25`*52
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_10<-survivalnodat15d$`10_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_10<-survivalnodat15d$`10_semaine`/survivalnodat15d$delai_C0_Y_cc
}
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_12.5<-survivalnodat15d$`12.5_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_12.5<-survivalnodat15d$`12.5_semaine`/survivalnodat15d$delai_C0_Y_cc
}
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_17.5<-survivalnodat15d$`17.5_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_17.5<-survivalnodat15d$`17.5_semaine`/survivalnodat15d$delai_C0_Y_cc
}
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_20<-survivalnodat15d$`20_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_20<-survivalnodat15d$`20_semaine`/survivalnodat15d$delai_C0_Y_cc
}
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_22.5<-survivalnodat15d$`22.5_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_22.5<-survivalnodat15d$`22.5_semaine`/survivalnodat15d$delai_C0_Y_cc
}
if(survivalnodat15d$nodat=1){
survivalnodat15d$ratio_semaine_25<-survivalnodat15d$`25_semaine`/survivalnodat15d$delai.y
} else (survivalnodat15d$nodat=0){
survivalnodat15d$ratio_semaine_25<-survivalnodat15d$`25_semaine`/survivalnodat15d$delai_C0_Y_cc
}Comparaison des différents seuils sur différent datent de censures (Création de modèle de Cox)
Test de Mann-Whitney entre les ratios obtenus dans les 2 populations (nodat et non nodat) aux différents thresholds
w10<-survivalnodat15d %>% wilcox_test(ratio_semaine_10~nodat)
w12.5<-survivalnodat15d %>% wilcox_test(ratio_semaine_12.5~nodat)
w15<-survivalnodat15d %>% wilcox_test(ratio_semaine_15~nodat)
w17.5<-survivalnodat15d %>% wilcox_test(ratio_semaine_17.5~nodat)
w20<-survivalnodat15d %>% wilcox_test(ratio_semaine_20~nodat)
w22.5<-survivalnodat15d %>% wilcox_test(ratio_semaine_22.5~nodat)
w25<-survivalnodat15d %>% wilcox_test(ratio_semaine_25~nodat)
W_all<-bind_rows(w10,w12.5,w15,w17.5,w20, w22.5, w25)
names(W_all)[1] <- "Ratio"Résultats des modèles de COX aux seuils significativement différents entre les 2 populations
Dans un premier temps, nous testons le paramètres BMIr (BMI >25).
test20<-survfit(Surv(delai,nodat)~BMIr,data=survivalnodat15d)
ggsurvplot(test20, surv.scale="percent", pval=T, fun="event")
Nous pouvons voir que ce critère permet de bien différencier les 2 populations.
- Modèle de COX pour un seuil à 15 microg/L.
coxfit20<-coxph(Surv(delai,nodat)~ratio_semaine_15+BMIr, data=survivalnodat15d,method=c("breslow"))
summary(coxfit20)Call:
coxph(formula = Surv(delai, nodat) ~ ratio_semaine_15 + BMIr,
data = survivalnodat15d, method = c("breslow"))
n= 404, number of events= 30
(92 observations deleted due to missingness)
coef exp(coef) se(coef) z Pr(>|z|)
ratio_semaine_15 0.14303 1.15376 0.02977 4.805 1.55e-06 ***
BMIr> 25 0.94064 2.56161 0.39070 2.408 0.0161 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
ratio_semaine_15 1.154 0.8667 1.088 1.223
BMIr> 25 2.562 0.3904 1.191 5.509
Concordance= 0.691 (se = 0.055 )
Likelihood ratio test= 19.13 on 2 df, p=7e-05
Wald test = 26.33 on 2 df, p=2e-06
Score (logrank) test = 33.21 on 2 df, p=6e-08
Dans ce modèle, on peut observer un Hazard Ratio de 1,2 concernant le paramètre “ratio semaine”. Soit pour une semaine/an de suivi au dessus de ce seuil, 1,2 fois plus de chance de développer un NODAT. On peut également observer que la proportionnalité des risques est respectée.
- Modèle de COX pour un seuil à 10 microg/L.
coxfit21<-coxph(Surv(delai,nodat)~ratio_semaine_10+BMIr, data=survivalnodat15d,method=c("breslow"))
summary(coxfit21)Call:
coxph(formula = Surv(delai, nodat) ~ ratio_semaine_10 + BMIr,
data = survivalnodat15d, method = c("breslow"))
n= 404, number of events= 30
(92 observations deleted due to missingness)
coef exp(coef) se(coef) z Pr(>|z|)
ratio_semaine_10 0.10411 1.10973 0.02606 3.995 6.47e-05 ***
BMIr> 25 0.83214 2.29823 0.38740 2.148 0.0317 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
ratio_semaine_10 1.110 0.9011 1.054 1.168
BMIr> 25 2.298 0.4351 1.076 4.911
Concordance= 0.7 (se = 0.049 )
Likelihood ratio test= 19.37 on 2 df, p=6e-05
Wald test = 20.5 on 2 df, p=4e-05
Score (logrank) test = 21.96 on 2 df, p=2e-05
Dans ce modèle, on peut observer un Hazard Ratio de 0,90. Soit pour une semaine/an de suivi au dessus de ce seuil, 0,91 fois plus de chance de développer un NODAT. Mais la proportionnalité des risques n’est pas respectée donc on ne peut pas utiliser ce modèle.
Résultats des modèles de COX aux seuils significativement différents entre les 2 populations
Dans un premier temps, nous testons le paramètres BMIr (BMI >25). 
Nous pouvons voir que ce critère permet de bien différencier les 2 populations.
- Modèle de COX pour un seuil à 10 microg/L.
Call:
coxph(formula = Surv(delai, nodat) ~ ratio_semaine_10 + BMIr,
data = survivalnodat15d6m, method = c("breslow"))
n= 393, number of events= 30
(89 observations deleted due to missingness)
coef exp(coef) se(coef) z Pr(>|z|)
ratio_semaine_10 0.07822 1.08137 0.02893 2.704 0.00686 **
BMIr> 25 0.82333 2.27808 0.38747 2.125 0.03359 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
ratio_semaine_10 1.081 0.9248 1.022 1.144
BMIr> 25 2.278 0.4390 1.066 4.868
Concordance= 0.655 (se = 0.05 )
Likelihood ratio test= 12.4 on 2 df, p=0.002
Wald test = 12.15 on 2 df, p=0.002
Score (logrank) test = 12.79 on 2 df, p=0.002Warning in .get_data(model, data = data): The `data` argument is not provided.
Data will be extracted from model fit.
Dans ce modèle, on peut observer un Hazard Ratio de 1,10. Soit pour une semaine/an de suivi au dessus de ce seuil, 1,10 fois plus de chance de développer un NODAT. Mais la proportionnalité des risques n’est pas respectée donc on ne peut pas utiliser ce modèle.
Calcul du temps passé au dessus de la valeur supérieur de la fenêtre thérapeutique au cours du temps (>15 ug/L avant 3 mois et >10 ug/L après 3 mois)
IntervalRatio<-ifelse(Threshold==0,1,0.5)
NODATvubft03<- NODATvubft %>% filter(delai_C0_Y_cc<0.25)
NODATvubft312<- NODATvubft %>% filter(delai_C0_Y_cc>=0.25)
nom<-unique(NODATvubft03$id)
rest_15<-data.frame()
nom1<-unique(NODATvubft312$id)
rest_10<-data.frame()
for (i in 1:length(nom)){
idpat<-nom[i]
patient<-subset(NODATvubft03, id == idpat)
previousVisit<-NULL
cumulatedTime<-0
# ItÈration sur les visites
for (j in 1:nrow(patient)) {
visite<-patient[j,]
delai<-visite$delai_C0_Y_cc
delta<-ifelse(is.null(previousVisit), 0, delai-previousVisit$delai_C0_Y_cc)
intervalOfInterest<-delta*IntervalRatio
cumulatedTime<- cumulatedTime + ifelse(visite$CO_cc >= threshold, intervalOfInterest, 0) # valeur cible ‡ dÈfinir essayer sinon ==
previousVisit<-visite
}
rest_15<-data.frame(rbind(rest_15, data.frame(id=idpat, CumulatedTime=cumulatedTime)))
print(rest_15)
}
##
for (i in 1:length(nom1)){
idpat<-nom1[i]
patient<-subset(NODATvubft312, id == idpat)
previousVisit<-NULL
cumulatedTime<-0
# ItÈration sur les visites
for (j in 1:nrow(patient)) {
visite<-patient[j,]
delai<-visite$delai_C0_Y_cc
delta<-ifelse(is.null(previousVisit), 0, delai-previousVisit$delai_C0_Y_cc)
intervalOfInterest<-delta*IntervalRatio
cumulatedTime<- cumulatedTime + ifelse(visite$CO_cc >= Threshold_v, intervalOfInterest, 0) # valeur cible ‡ dÈfinir essayer sinon ==
previousVisit<-visite
}
rest_10<-data.frame(rbind(rest_10, data.frame(id=idpat, CumulatedTime=cumulatedTime)))
print(rest_15)
}
###
rest_10_15<-merge(rest_10, rest_15, by="id", all = T)
rest_10_15$CumulatedTime.y<- rest_10_15$CumulatedTime.y %>% replace_na(0)
rest_10_15$CumulatedTime.x<- rest_10_15$CumulatedTime.x %>% replace_na(0)
rest_10_15$CumulatedTime<-rest_10_15$CumulatedTime.x+rest_10_15$CumulatedTime.y
NODATvubftid<-survivalnodat15d %>% select(id, nodat, delai.y, delai_C0_Y_cc) %>% group_by(id) %>%slice(1)
survivalnodat15d3m<-merge(rest_10_15, NODATvubftid, by="id")
survivalnodat15d3m$CumulatedTime<-ifelse(survivalnodat15d3m$id=="AMAF087", 1/365,survivalnodat15d3m$CumulatedTime)
survivalnodat15d3m$CumulatedTime_semaine<-survivalnodat15d3m$CumulatedTime*52
if (survivalnodat15d3m$nodat==1){
survivalnodat15d3m$ratio<-survivalnodat15d3m$CumulatedTime_semaine/survivalnodat15d3m$delai.y
} else (survivalnodat15d3m$nodat==0){
survivalnodat15d3m$ratio<-survivalnodat15d3m$CumulatedTime_semaine/survivalnodat15d3m$delai_C0_Y_cc
} Comparaison des ratios des patients non diabétiques entre 0-6 mois & 6- 12 mois
ratio_semaine.x= 0-6 mois
ratio_semaine.y= 6-12 mois
ratio_semaine = 0-12 mois
ggplot(test_comp3) +
aes(x = "", y = value, fill = name) +
geom_boxplot() +
scale_fill_hue() +
labs(x = " ") +
theme_minimal() +
facet_wrap(vars(name))
A priori on peut en conclure à une non-proportionalité des ratios au cours du temps. Il serait donc pertinent de censurer à une date avant 1 an.
Description apparition NODAT
On peut observer que 95% des diabètes sont apparues à 0.28 année soit 3 mois et 11 jours environ.
Censure à 6 mois (soit 0.5 année)
## temps passé 6mois
Tps_sup6m<-function(Threshold){
IntervalRatio<-ifelse(Threshold==0,1,0.5)
nom<-unique(NODATvubft6m$id)
resta<-data.frame()
for (i in 1:length(nom)){
idpat<-nom[i]
patient<-subset(NODATvubft6m, id == idpat)
previousVisit<-NULL
cumulatedTime<-0
# ItÈration sur les visites
for (j in 1:nrow(patient)) {
visite<-patient[j,]
delai<-visite$delai_C0_Y_cc
delta<-ifelse(is.null(previousVisit), 0, delai-previousVisit$delai_C0_Y_cc)
intervalOfInterest<-delta*IntervalRatio
cumulatedTime<- cumulatedTime + ifelse(visite$CO_cc >= Threshold, intervalOfInterest, 0) # valeur cible ‡ dÈfinir essayer sinon ==
previousVisit<-visite
}
resta<-data.frame(rbind(resta, data.frame(id=idpat, CumulatedTime=cumulatedTime, Threshold=Threshold)))
print(resta)
}
return(resta)}
reg6m<-map_df(Threshold,Tps_sup6m)survivalnodat15d6m$nodat<-as.factor(survivalnodat15d6m$nodat)
w106m<-survivalnodat15d6m %>%wilcox_test(ratio_semaine_10~nodat)
w12.56m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_12.5~nodat)
w156m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_15~nodat)
w17.56m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_17.5~nodat)
w206m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_20~nodat)
w22.56m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_22.5~nodat)
w256m<-survivalnodat15d6m %>% wilcox_test(ratio_semaine_25~nodat)
W_all_6m<-bind_rows(w106m,w12.56m,w156m,w17.56m,w206m, w22.56m, w256m)
names(W_all_6m)[1] <- "Ratio"Censure à 3 mois (soit 0.25 année)
## CRéation de la base à 3 mois
NODAT1vubft3m<-NODAT1vubf %>%
filter(delai.y<0.25)
NODAT0vubft3m<-NODAT0vubf %>%
filter(delai_C0_Y_cc<0.25)
NODATvubft3m<-bind_rows(NODAT0vubft3m,NODAT1vubft3m)
NODATvubft3m<-NODATvubft3m%>% group_by(id)%>% arrange(id, delai_C0_Y_cc)
NODATvubft3m$nodat<-as.factor(NODATvubft3m$nodat)
nom2<-unique(NODATvubft3m$id)
rest3m<-data.frame()
for (i in 1:length(nom2)){
idpat<-nom2[i]
patient<-subset(NODATvubft3m, id == idpat)
previousVisit<-NULL
cumulatedTime<-0
# ItÈration sur les visites
for (j in 1:nrow(patient)) {
visite<-patient[j,]
delai<-visite$delai_C0_Y_cc
delta<-ifelse(is.null(previousVisit), 0, delai-previousVisit$delai_C0_Y_cc)
intervalOfInterest<-delta*IntervalRatio
cumulatedTime<- cumulatedTime + ifelse(visite$CO_cc >= threshold, intervalOfInterest, 0) # valeur cible ‡ dÈfinir essayer sinon ==
previousVisit<-visite
}
rest3m<-data.frame(rbind(rest3m, data.frame(id=idpat, CumulatedTime=cumulatedTime, Threshold=threshold)))
print(rest3m)
}
NODATvubft3m<-NODATvubft3m%>% group_by(id)%>% arrange(id, desc(delai_C0_Y_cc))
reg13m<-pivot_wider(rest3m, names_from = Threshold, values_from = CumulatedTime)
reg13m$`15`<-ifelse (reg13m$id=="AMAF087", 1/365,reg13m$`15`)
survivalnodat3m<-merge(NODATvubft3m,reg13m, by="id")
survivalnodat15d3m<-survivalnodat3m%>%group_by(id)%>%slice(1)
survivalnodat15d3m$`15_semaine`<-survivalnodat15d3m$`15`*52
if(survivalnodat15d3m$nodat=1){
survivalnodat15d3m$ratio_semaine_15<-survivalnodat15d3m$`15_semaine`/survivalnodat15d3m$delai.y
} else (survivalnodat15d3m$nodat=0){
survivalnodat15d3m$ratio_semaine_15<-survivalnodat15d3m$`15_semaine`/survivalnodat15d3m$delai_C0_Y_cc
}
survivalnodat15d3m$IMC=(survivalnodat15d3m$poids)/((survivalnodat15d3m$taille/100)*(survivalnodat15d3m$taille/100))
survivalnodat15d3m$BMIr<-as.factor(ifelse(survivalnodat15d3m$IMC>25, " > 25", " < 25" ))
survivalnodat15d3m$nodat<-as.factor(survivalnodat15d3m$nodat)
write_csv(survivalnodat15d3m, path = "survivalnodat15d3m.csv")
w153m<-survivalnodat15d3m %>% wilcox_test(ratio_semaine_15~nodat)
names(w153m)[1] <- "Ratio"Essai au Machine Learning
Utilisation de Tidymodels
Essai du Machine Learning pour estimer les patients ayant un NODAT.
Partage de la population totale en 0.6 (Training)/0.4 (Testing) soit 298 patients pour le training et 198 patients pour le testing.
Package Ranger ntrees=100, mode=“classification”
survivalnodat15dLM<-survivalnodat15d%>% select(nodat, ratio_semaine_15,ratio_semaine_10, ratio_semaine_12.5, ratio_semaine_17.5 )
survivalnodat15dLM$nodat<-as.factor(survivalnodat15dLM$nodat)
set.seed(100)
survivalnodat15d_split <- initial_split(survivalnodat15dLM, prop = 0.6)
survivalnodat15d_split #<Training/Validation/Total> <298/198/496>
survivalnodat15d_split %>%
training() %>%
glimpse()
##Recipe
survivalnodat15d_recipe <- training(survivalnodat15d_split) %>%
recipe(nodat ~.) %>%
step_corr(all_predictors()) %>%
step_center(all_predictors(), -all_outcomes()) %>%
step_scale(all_predictors(), -all_outcomes()) %>%
prep()
survivalnodat15d_recipe
##Bake
survivalnodat15d_testing <- survivalnodat15d_recipe %>%
bake(testing(survivalnodat15d_split))
glimpse(survivalnodat15d_testing)
##Juice
survivalnodat15d_training <- juice(survivalnodat15d_recipe)
glimpse(survivalnodat15d_training)
##Rand Forest
survivalnodat15d_ranger <- rand_forest(trees = 100, mode = "classification") %>%
set_engine("ranger") %>%
fit(nodat ~ ., data = survivalnodat15d_testing)
survivalnodat15d_rf <- rand_forest(trees = 100, mode = "classification") %>%
set_engine("randomForest") %>%
fit(nodat ~ ., data = survivalnodat15d_training)
##Predictions
predict(survivalnodat15d_ranger, survivalnodat15d_testing)
survivalnodat15d_ranger %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
glimpse()
## Model Validation
survivalnodat15d_ranger %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
metrics(truth = nodat, estimate = .pred_class)
survivalnodat15d_rf %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
metrics(truth = nodat, estimate = .pred_class)
## Per lassifier metrics
survivalnodat15d_ranger %>%
predict(survivalnodat15d_testing, type = "prob") %>%
glimpse()
survivalnodat15d_probs <- survivalnodat15d_ranger %>%
predict(survivalnodat15d_testing, type = "prob") %>%
bind_cols(survivalnodat15d_testing)
glimpse(survivalnodat15d_probs)
survivalnodat15d_probs%>%
glimpse()
testpred<-predict(survivalnodat15d_ranger, survivalnodat15d_testing, type = "prob") %>%
bind_cols(predict(survivalnodat15d_ranger, survivalnodat15d_testing)) %>%
bind_cols(select(survivalnodat15d_testing, nodat)) %>%
glimpse()
testpred1<-predict(survivalnodat15d_ranger, survivalnodat15d_testing, type = "prob") %>%
bind_cols(predict(survivalnodat15d_ranger, survivalnodat15d_testing)) %>%
bind_cols(select(survivalnodat15d_testing, nodat)) %>%
metrics(nodat, estimate = .pred_class)
mat <- as.table.confusionMatrix(data=testpred$.pred_class,reference=testpred$nodat, positive="1")Confusion Matrix and Statistics
Reference
Prediction 0 1
0 183 5
1 0 10
Accuracy : 0.9747
95% CI : (0.9421, 0.9918)
No Information Rate : 0.9242
P-Value [Acc > NIR] : 0.002113
Kappa : 0.7871
Mcnemar's Test P-Value : 0.073638
Sensitivity : 0.66667
Specificity : 1.00000
Pos Pred Value : 1.00000
Neg Pred Value : 0.97340
Prevalence : 0.07576
Detection Rate : 0.05051
Detection Prevalence : 0.05051
Balanced Accuracy : 0.83333
'Positive' Class : 1
Package Ranger ntrees=200, mode=“classification”
##Rand Forest
survivalnodat15d_ranger2 <- rand_forest(trees = 200, mode = "classification") %>%
set_engine("ranger") %>%
fit(nodat ~ ., data = survivalnodat15d_testing)
##Predictions
predict(survivalnodat15d_ranger2, survivalnodat15d_testing)
survivalnodat15d_ranger2 %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
glimpse()
## Model Validation
survivalnodat15d_ranger2 %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
metrics(truth = nodat, estimate = .pred_class)
survivalnodat15d_rf %>%
predict(survivalnodat15d_testing) %>%
bind_cols(survivalnodat15d_testing) %>%
metrics(truth = nodat, estimate = .pred_class)
## Per lassifier metrics
survivalnodat15d_ranger2 %>%
predict(survivalnodat15d_testing, type = "prob") %>%
glimpse()
survivalnodat15d_probs2 <- survivalnodat15d_ranger2 %>%
predict(survivalnodat15d_testing, type = "prob") %>%
bind_cols(survivalnodat15d_testing)
glimpse(survivalnodat15d_probs2)
survivalnodat15d_probs2%>%
glimpse()
testpred2<-predict(survivalnodat15d_ranger2, survivalnodat15d_testing, type = "prob") %>%
bind_cols(predict(survivalnodat15d_ranger2, survivalnodat15d_testing)) %>%
bind_cols(select(survivalnodat15d_testing, nodat)) %>%
glimpse()
testpred12<-predict(survivalnodat15d_ranger2, survivalnodat15d_testing, type = "prob") %>%
bind_cols(predict(survivalnodat15d_ranger2, survivalnodat15d_testing)) %>%
bind_cols(select(survivalnodat15d_testing, nodat)) %>%
metrics(nodat, estimate = .pred_class)
mat2 <- confusionMatrix(data=testpred2$.pred_class,reference=testpred2$nodat, positive="1")Confusion Matrix and Statistics
Reference
Prediction 0 1
0 183 7
1 0 8
Accuracy : 0.9646
95% CI : (0.9285, 0.9857)
No Information Rate : 0.9242
P-Value [Acc > NIR] : 0.01495
Kappa : 0.6787
Mcnemar's Test P-Value : 0.02334
Sensitivity : 0.53333
Specificity : 1.00000
Pos Pred Value : 1.00000
Neg Pred Value : 0.96316
Prevalence : 0.07576
Detection Rate : 0.04040
Detection Prevalence : 0.04040
Balanced Accuracy : 0.76667
'Positive' Class : 1