INTRODUCCIÓN

Un hospital ha contactado para que se haga un estudio sobre la probabilidad de que sus pacientes sufran un infarto.

Para ello, nos ha facilitado un archivo donde recoge a sus pacientes y los enfermedades-características que poseen en relación a los infartos.

Se van a realizar distintos modelos empleando Machine Learning para conocer la raiz de sufrir infartos y poder mostrar un camino más certero para futuras investigaciones.

1.- IMPORTACIÓN Y MUESTREO

Preparación del entorno

Eliminación de la notación cientifica

options(scipen=999)

Se cargan e instalan las librerias

paquetes<- c('data.table',
             'dplyr',
             'tidyr',
             'ggplot2',
             'randomForest',
             'ROCR',
             'purrr',
             'smbinning',
             'rpart',
             'rpart.plot',
             'h2o',
             'faraway')
instalados<-paquetes%in%installed.packages()
if(sum(instalados==FALSE)>0){
  install.packages(paquetes[!instalados])
}
lapply(paquetes,require,character.only=TRUE)

Se cargan los datos en un df

df <- read.csv("healthcare-dataset-stroke-data.csv")

2.- ANÁLISIS EXPLORATORIO

Análisis exploratorio general

Se van a arealizar las siguientes acciones:

  • Conocer el nombre de las variables

  • Tipología de las variables y una muestra de los valores que la conforman

  • Comprobar si la variable objetivo Churn, está balanceada

names(df)
##  [1] "id"                "gender"            "age"              
##  [4] "hypertension"      "heart_disease"     "ever_married"     
##  [7] "work_type"         "Residence_type"    "avg_glucose_level"
## [10] "bmi"               "smoking_status"    "stroke"
glimpse(df)
## Rows: 5,110
## Columns: 12
## $ id                <int> 9046, 51676, 31112, 60182, 1665, 56669, 53882, 10...
## $ gender            <chr> "Male", "Female", "Male", "Female", "Female", "Ma...
## $ age               <dbl> 67, 61, 80, 49, 79, 81, 74, 69, 59, 78, 81, 61, 5...
## $ hypertension      <int> 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0...
## $ heart_disease     <int> 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1, 0, 1...
## $ ever_married      <chr> "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", ...
## $ work_type         <chr> "Private", "Self-employed", "Private", "Private",...
## $ Residence_type    <chr> "Urban", "Rural", "Rural", "Urban", "Rural", "Urb...
## $ avg_glucose_level <dbl> 228.69, 202.21, 105.92, 171.23, 174.12, 186.21, 7...
## $ bmi               <chr> "36.6", "N/A", "32.5", "34.4", "24", "29", "27.4"...
## $ smoking_status    <chr> "formerly smoked", "never smoked", "never smoked"...
## $ stroke            <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1...
table(df$stroke)#Se comprueba si la variable target está bien balanceada
## 
##    0    1 
## 4861  249

Gracias a este análisis inicial, se extraen las siguientes conclusiones:

  • Existen variables que tienen que transformarse de character a factor, hay que dividirlas en niveles

  • La variable bmi (body mass index) presenta valores nulos NA’s

  • La variable stroke (target) está poco balanceada ya que presenta un 4% de 1. Esto puede suponer que los modelos sobreajusten demasiado. Para no falsear los datos, no se van a introducir artificialmente 1 para balancear la target. Así los valores son reales

Calidad de datos

Se realizarán una serie de análisis:

  • Estadísticos básicos

  • Análisis de nulos

  • Análisis de ceros

  • Análisis de atípicos

Estadísticos básicos

lapply(df,summary)
## $id
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##      67   17741   36932   36518   54682   72940 
## 
## $gender
##    Length     Class      Mode 
##      5110 character character 
## 
## $age
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.08   25.00   45.00   43.23   61.00   82.00 
## 
## $hypertension
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## 0.00000 0.00000 0.00000 0.09746 0.00000 1.00000 
## 
## $heart_disease
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## 0.00000 0.00000 0.00000 0.05401 0.00000 1.00000 
## 
## $ever_married
##    Length     Class      Mode 
##      5110 character character 
## 
## $work_type
##    Length     Class      Mode 
##      5110 character character 
## 
## $Residence_type
##    Length     Class      Mode 
##      5110 character character 
## 
## $avg_glucose_level
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   55.12   77.25   91.89  106.15  114.09  271.74 
## 
## $bmi
##    Length     Class      Mode 
##      5110 character character 
## 
## $smoking_status
##    Length     Class      Mode 
##      5110 character character 
## 
## $stroke
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## 0.00000 0.00000 0.00000 0.04873 0.00000 1.00000

Se observa lo siguiente:

  • Existen muchas variables categóricas en formato carácter que conviene transformar para ver su distribución

  • Se observa que la variable age presente valores de 0.8 años como edad mínima. Por eso se va a realizar un visualización de su distribución

Análisis de nulos

data.frame(colSums(is.na(df)))
##                   colSums.is.na.df..
## id                                 0
## gender                             0
## age                                0
## hypertension                       0
## heart_disease                      0
## ever_married                       0
## work_type                          0
## Residence_type                     0
## avg_glucose_level                  0
## bmi                                0
## smoking_status                     0
## stroke                             0

No existen valores nulos en formato numérico, pero como se ha visto anteriormente bmi son números en formato character y posee valores nulos. Por eso, se va a transformar bmi y se observará cuantos son NA´s

df$bmi <- as.numeric(df$bmi)
data.frame(colSums(is.na(df)))
##                   colSums.is.na.df..
## id                                 0
## gender                             0
## age                                0
## hypertension                       0
## heart_disease                      0
## ever_married                       0
## work_type                          0
## Residence_type                     0
## avg_glucose_level                  0
## bmi                              201
## smoking_status                     0
## stroke                             0
#Esxisten 201 Nulos de 5110 personas. Lo que se va hacer es eliminarlos
df<-na.omit(df)

Como son solo 201 pacientes frente a los 5110 totales, se decide eliminarlos del df

Análisis de ceros

#Se va a crear la función para contar los ceros de las variables
contar_ceros <- function(variable){
  temp<-transmute(df,if_else(variable==0,1,0))
  sum(temp)
}
num_ceros<-sapply(df,contar_ceros)#Con sapply se aplica al df completamente
num_ceros<-data.frame(VARIABLE=names(num_ceros),CEROS=as.numeric(num_ceros),stringsAsFactors=F)
num_ceros
##             VARIABLE CEROS
## 1                 id     0
## 2             gender     0
## 3                age     0
## 4       hypertension  4458
## 5      heart_disease  4666
## 6       ever_married     0
## 7          work_type     0
## 8     Residence_type     0
## 9  avg_glucose_level     0
## 10               bmi     0
## 11    smoking_status     0
## 12            stroke  4700

Existen 3 variables con una alta cantidad de ceros:

  • hypertension –> diferencia entre pacientes con hipertensión (1) y los que no la poseen (0). Es normal que haya más personas sin hipertensión que con hipertensión

  • hearth_disease –> personas que presentan enfermedades del corazón (1) y las que no (0). Se ve que hay más sin enfermedades que con enfermedades

  • stroke –> personas que han sufrido un infarto.

Como el estudio es general, es decir, no solo trata a pacientes enfermos, también hay personas sanas o que no padecen enfermedades, no se hará ninguna acción sobre ellos.

Análisis de atípicos

Se estudiarán dos tipos de variables:

  • Tipología double

  • Tipología integer

Variables double

Se creará una función que se encargará de mostrar los 35 valores mayores de aquellas variables doubles y comprobar que el valor decrece gradualmente y no existen grandes saltos

out<-function(variable){
  t(t(head(sort(variable,decreasing=T),35)))
}
lapply(df,function(x){
  if(is.double(x))out(x)
})
## $id
## NULL
## 
## $gender
## NULL
## 
## $age
##       [,1]
##  [1,]   82
##  [2,]   82
##  [3,]   82
##  [4,]   82
##  [5,]   82
##  [6,]   82
##  [7,]   82
##  [8,]   82
##  [9,]   82
## [10,]   82
## [11,]   82
## [12,]   82
## [13,]   82
## [14,]   82
## [15,]   82
## [16,]   82
## [17,]   82
## [18,]   82
## [19,]   82
## [20,]   82
## [21,]   82
## [22,]   82
## [23,]   82
## [24,]   82
## [25,]   82
## [26,]   82
## [27,]   82
## [28,]   82
## [29,]   82
## [30,]   82
## [31,]   82
## [32,]   82
## [33,]   82
## [34,]   82
## [35,]   82
## 
## $hypertension
## NULL
## 
## $heart_disease
## NULL
## 
## $ever_married
## NULL
## 
## $work_type
## NULL
## 
## $Residence_type
## NULL
## 
## $avg_glucose_level
##         [,1]
##  [1,] 271.74
##  [2,] 267.76
##  [3,] 267.61
##  [4,] 267.60
##  [5,] 266.59
##  [6,] 263.56
##  [7,] 263.32
##  [8,] 261.67
##  [9,] 259.63
## [10,] 256.74
## [11,] 255.17
## [12,] 254.95
## [13,] 254.63
## [14,] 254.60
## [15,] 253.86
## [16,] 253.16
## [17,] 252.72
## [18,] 251.99
## [19,] 251.60
## [20,] 251.46
## [21,] 250.89
## [22,] 250.80
## [23,] 250.20
## [24,] 249.31
## [25,] 249.29
## [26,] 248.37
## [27,] 248.24
## [28,] 247.97
## [29,] 247.87
## [30,] 247.69
## [31,] 247.51
## [32,] 247.48
## [33,] 246.53
## [34,] 246.34
## [35,] 244.30
## 
## $bmi
##       [,1]
##  [1,] 97.6
##  [2,] 92.0
##  [3,] 78.0
##  [4,] 71.9
##  [5,] 66.8
##  [6,] 64.8
##  [7,] 64.4
##  [8,] 63.3
##  [9,] 61.6
## [10,] 61.2
## [11,] 60.9
## [12,] 60.9
## [13,] 60.2
## [14,] 59.7
## [15,] 58.1
## [16,] 57.9
## [17,] 57.7
## [18,] 57.5
## [19,] 57.3
## [20,] 57.2
## [21,] 57.2
## [22,] 56.6
## [23,] 56.6
## [24,] 56.1
## [25,] 56.0
## [26,] 55.9
## [27,] 55.9
## [28,] 55.7
## [29,] 55.7
## [30,] 55.7
## [31,] 55.7
## [32,] 55.2
## [33,] 55.1
## [34,] 55.0
## [35,] 55.0
## 
## $smoking_status
## NULL
## 
## $stroke
## NULL

Se va a ver gráficamente como se distribuyen:

Age
hist(df$age,breaks = 9)

Se observa que existen pocas personas con más de 80 años, cabe mencionar que la esperanza de vida ronda los 80 años.

Avg glucose level
hist(df$avg_glucose_level, breaks=10)

Existe una gran concentración de pacientes con una glucosa media de 75, un valor dentro de la media como persona sana.

Avg glucose level
hist(df$bmi)

La distribución del índice corporal es alto en valores normales y bajos en valores muy pequeños o muy altos.

Variables integer

Como solo hay 1 variable en formato integer, que es la identificación del paciente, no es necesario estudiar este tipo de variables

Preselección de variables

Se genera una lista donde se recogen las variables que se van a convertir a factor para continuar con la creación del modelo

a_factor<-(c('gender','ever_married','work_type','Residence_type','smoking_status'))

Análisis longitudinal

longi<-df%>%
  summarise_all(mean)%>%
  t %>%
  as.data.frame()
data.frame(variable=rownames(longi),media=longi$V1)%>%
  arrange(desc(variable))
##             variable          media
## 1          work_type             NA
## 2             stroke     0.04257486
## 3     smoking_status             NA
## 4     Residence_type             NA
## 5                 id 37064.31350581
## 6       hypertension     0.09187207
## 7      heart_disease     0.04950092
## 8             gender             NA
## 9       ever_married             NA
## 10               bmi    28.89323691
## 11 avg_glucose_level   105.30514972
## 12               age    42.86537380

Conclusiones finales:

  • El estudio está enfocado a una población con una edad media de 40 años, con un nivel de glucosa en sangre algo superior a la media y un índice de masa corporal por encima de lo recomendable

  • Se crea una copia del df ahora para emplearlo posteriormente en el Information Value

df1<-df
  • Se pasan las variables preseleccionadas a factor y a parte se pasa Senior Citizen a factor también. Senior Citizen se hace por separado porque su transformación a factor será diferente que las de a_factor de la lista:
df<-df%>%
  mutate_at(a_factor,.funs=factor)
df$hypertension<- as.factor(as.character(df$hypertension))
df$heart_disease<- as.factor(as.character(df$heart_disease))
df$stroke <- as.factor(as.character(df$stroke))
glimpse(df)
## Rows: 4,909
## Columns: 12
## $ id                <int> 9046, 31112, 60182, 1665, 56669, 53882, 10434, 60...
## $ gender            <fct> Male, Male, Female, Female, Male, Male, Female, F...
## $ age               <dbl> 67, 80, 49, 79, 81, 74, 69, 78, 81, 61, 54, 79, 5...
## $ hypertension      <fct> 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0...
## $ heart_disease     <fct> 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0...
## $ ever_married      <fct> Yes, Yes, Yes, Yes, Yes, Yes, No, Yes, Yes, Yes, ...
## $ work_type         <fct> Private, Private, Private, Self-employed, Private...
## $ Residence_type    <fct> Urban, Rural, Urban, Rural, Urban, Rural, Urban, ...
## $ avg_glucose_level <dbl> 228.69, 105.92, 171.23, 174.12, 186.21, 70.09, 94...
## $ bmi               <dbl> 36.6, 32.5, 34.4, 24.0, 29.0, 27.4, 22.8, 24.2, 2...
## $ smoking_status    <fct> formerly smoked, never smoked, smokes, never smok...
## $ stroke            <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1...

3.- TRANSFORMACIÓN DE DATOS

En esta fase se van a realizar tres funciones

Preparación de las variables independientes

Las variables independientes son aquellas distintas a la variable target (stroke) y a la identidad del paciente id. Esta lista se empleará en Random Forest

ind_larga<-names(df) #Lista con los nombres de todas las variables
no_usar<- c('id','stroke') #No se van a usar ni codigo de paciente ni la target
ind_larga<-setdiff(ind_larga,no_usar)
ind_larga
##  [1] "gender"            "age"               "hypertension"     
##  [4] "heart_disease"     "ever_married"      "work_type"        
##  [7] "Residence_type"    "avg_glucose_level" "bmi"              
## [10] "smoking_status"

Se crea una copia del entorno para poder empezar desde este punto

save.image(file='sesion1.RData')

Para poder cargar la copia

load(file='sesion1.RData')

Preselección con RandomForest

El número de variables que cogerá de todas las disponibles será de 1, para poder obtener la mejor correlación entre RandomForest e Information Value

#Es necesario que la variable target esté en formato factor
pre_rf<-randomForest(formula=reformulate(ind_larga,'stroke'),data=df,mtry=1,ntree=400,importance=T)
imp_rf<-importance(pre_rf)[,4]
imp_rf<-data.frame(VARIABLE=names(imp_rf),IMP_RF=imp_rf)
imp_rf<-imp_rf %>% arrange(desc(IMP_RF))%>%mutate(RANKING_RF=1:nrow(imp_rf))
imp_rf
##             VARIABLE    IMP_RF RANKING_RF
## 1                age 19.409105          1
## 2  avg_glucose_level 15.752642          2
## 3                bmi 10.815882          3
## 4       hypertension  4.571845          4
## 5     smoking_status  4.455400          5
## 6      heart_disease  4.002062          6
## 7          work_type  3.578263          7
## 8       ever_married  2.432185          8
## 9     Residence_type  1.419773          9
## 10            gender  1.328978         10

Preselección con Information Value (IV)

Tal y como está el df actualmente, no es posible la realización del estudio sobre él.

Como estos cambios se emplearán exclusivamente para IV, se utilizará la copia del data frame llamada df1, asegurando que el df original seguirá intacto.

Transformación de carácter a número

Como se expuso antes, existen muchas variables categóricas en formato texto. Para poder realizar el estudio de IV, es necesario asignar un valor numérico a aquellos valores que conforman dichas variables y transformar estos valores a factor.

#Transformación de variables character con dos niveles
df1<-df
df1<-df1 %>%
  mutate(gender=ifelse((gender=="Male"),1,0),
         Residence_type=ifelse((Residence_type=="Urban"),1,0),
         ever_married=ifelse((ever_married=="Yes"),1,0))%>%
  mutate(gender=as.factor(as.character(gender)),
         Residence_type=as.factor(as.character(Residence_type)),
         ever_married=as.factor(as.character(ever_married)))

#Transformación de variables character con tres o más niveles
df1<-df1 %>%
  mutate(work_type=ifelse((work_type=="Never_worked"),0,ifelse((work_type=="children"),1,ifelse((work_type=="Private"),2,ifelse((work_type=="Self-employed"),3,4)))),
         smoking_status=ifelse((smoking_status=="never smoked"),0,ifelse((smoking_status=="smokes"),1,ifelse((smoking_status=="formerly smoked"),2,3))))%>%
  mutate(work_type=as.factor(as.character(work_type)),
        smoking_status=as.factor(as.character(smoking_status)))
glimpse(df1)
## Rows: 4,909
## Columns: 12
## $ id                <int> 9046, 31112, 60182, 1665, 56669, 53882, 10434, 60...
## $ gender            <fct> 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0...
## $ age               <dbl> 67, 80, 49, 79, 81, 74, 69, 78, 81, 61, 54, 79, 5...
## $ hypertension      <fct> 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0...
## $ heart_disease     <fct> 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0...
## $ ever_married      <fct> 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1...
## $ work_type         <fct> 2, 2, 2, 3, 2, 2, 2, 2, 2, 4, 2, 2, 3, 2, 2, 2, 4...
## $ Residence_type    <fct> 1, 0, 1, 0, 1, 0, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0...
## $ avg_glucose_level <dbl> 228.69, 105.92, 171.23, 174.12, 186.21, 70.09, 94...
## $ bmi               <dbl> 36.6, 32.5, 34.4, 24.0, 29.0, 27.4, 22.8, 24.2, 2...
## $ smoking_status    <fct> 2, 0, 1, 0, 2, 0, 0, 3, 0, 1, 1, 0, 0, 1, 1, 0, 1...
## $ stroke            <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1...

Una vez transformado df1, se aplica IV

#Para poder aplicar el IV es necesario pasar la target a numeric. Antes de transformarla en numeric, se transformará en character, para que mantenga el formato 0 y 1 y no se pase como 1 y 2. Esto se hace porque para calcular el IV con smbinning.sumiv, es necesario que esté en numérico
temp <- mutate(df1,stroke = as.numeric(as.character(stroke))) %>% as.data.frame()#as.data.frame es por estética
imp_iv<-smbinning.sumiv(temp[c(ind_larga,'stroke')],y="stroke")
imp_iv <- imp_iv %>% mutate(Ranking = 1:nrow(imp_iv), IV = ifelse(is.na(.$IV),0,IV)) %>% select(-Process)
names(imp_iv)<-c('VARIABLE','IMP_IV','RANKING_IV')
imp_iv
##             VARIABLE IMP_IV RANKING_IV
## 1                age 2.3391          1
## 2          work_type 0.5203          2
## 3  avg_glucose_level 0.3808          3
## 4       ever_married 0.3732          4
## 5                bmi 0.3289          5
## 6       hypertension 0.3038          6
## 7      heart_disease 0.2455          7
## 8     smoking_status 0.2002          8
## 9             gender 0.0012          9
## 10    Residence_type 0.0009         10

Ranking Final

Para poder hacer la selección final, se van a enfrentar los resultados de ambos estudios en un ranking final.

imp_final<-inner_join(imp_rf,imp_iv,by='VARIABLE')%>%
  select(VARIABLE,IMP_RF,IMP_IV,RANKING_RF,RANKING_IV)%>%
  mutate(RANKING_TOT=RANKING_RF,RANKING_IV)%>%
  arrange(RANKING_TOT)
imp_final
##             VARIABLE    IMP_RF IMP_IV RANKING_RF RANKING_IV RANKING_TOT
## 1                age 19.409105 2.3391          1          1           1
## 2  avg_glucose_level 15.752642 0.3808          2          3           2
## 3                bmi 10.815882 0.3289          3          5           3
## 4       hypertension  4.571845 0.3038          4          6           4
## 5     smoking_status  4.455400 0.2002          5          8           5
## 6      heart_disease  4.002062 0.2455          6          7           6
## 7          work_type  3.578263 0.5203          7          2           7
## 8       ever_married  2.432185 0.3732          8          4           8
## 9     Residence_type  1.419773 0.0009          9         10           9
## 10            gender  1.328978 0.0012         10          9          10

Como el df solo cuenta con 10 variables, se van a seguir empleando, ya que en los distintos modelos se hará una selección de variables antes de obtener el modelo final.

  • Además, la correlación entre ambos métodos es media: 0.7554114

Limpieza del entorno

Definido el nuevo estado del df, se va a proceder a la limpieza del entorno

ls()
rm(list=setdiff(ls(),'df'))
target<-'stroke'
indep<-setdiff(names(df),c('stroke','id'))

Y se va a guardar una caché, cache1 para poder partir desde este punto.

saveRDS(df,'cache1.rds')

Para poder cargar dicha caché:

df <- readRDS(file = 'cache1.rds')
df1 <- readRDS(file='cache1.rds')

4.- CREACIÓN DE VARIABLES SINTÉTICAS

En este modelo, las variables carecen de histórico, es decir, no se diferencian por meses, solo representan el ahora.

Debido a eso, no se podrán crear ni tenencias, ni tendencias, ni cancelaciones, ni contrataciones ni medias.

La única tipología de variable sintética que se puede aplicar es la discretización.

Discretización

La discretización consiste en dividir cada variable en distintos escalones en función de la penetración de cada escalón en la variable target.

Dicha penetración se debe caracterizar por dos motivos:

  • Que sea similar a la penetración de la variable original

  • Que presente un comportamiento monotónico, es decir, a medida que se avanza en la variable, esta siempre crece o siempre decrece.

Existen dos tipos de discretizaciones:

  • Discretización automática –> age, avg_glucose_level, bmi

  • Discretización manual –> gender, work_type, smoking_status y Residence_type

Discretización automática

Para la discretización automática, se va a crear una función

discretizar <- function (vi,target){ #vi es variuable independiente
  df_temp <- data.frame(vi=vi,target=target)
  df_temp$target <- as.numeric(as.character((df_temp$target)))#smbinning necesita que la target esté en formato numérico
  disc <- smbinning(df_temp,y='target',x='vi')
  return(disc)
}

Se aplica dicha función

#Se discretiza age
disc_temp_age <- discretizar(df$age,df$stroke)
df_temp <- select(df,age,stroke)
df_temp <- smbinning.gen(df_temp,disc_temp_age,chrname='age_disc')
df<- cbind(df,df_temp[3]) %>% select(-age)

#Se discretiza average_glucose_level
disc_temp_avg_glucose_level <- discretizar(df$avg_glucose_level,df$stroke)
df_temp <- select(df,avg_glucose_level,stroke)
df_temp <- smbinning.gen(df_temp,disc_temp_avg_glucose_level,chrname='avg_glucose_level_disc')
df <- cbind(df,df_temp[3]) %>% select(-avg_glucose_level)

#Se discretiza bmi
disc_temp_bmi <- discretizar(df$bmi,df$stroke)
df_temp <- select(df,bmi,stroke)
df_temp <- smbinning.gen(df_temp,disc_temp_bmi,chrname='bmi_disc')
df<- cbind(df,df_temp[3]) %>% select(-bmi)

Discretización manual

Conforme se discretizan las variables, se van borrando las originales

#Gender
df <- df%>%
  mutate(gender_disc=as.factor(case_when(
    gender=='Male'~'01_Male',
    gender=='Female'~'02_Female',
    gender=='Other'~'03_Other',
    TRUE ~ '00_ERROR'
  )))
df$gender<-NULL

#Work_type
df <- df%>%
  mutate(work_type_disc=as.factor(case_when(
    work_type=='Never_worked'~'01_Never_worked',
    work_type=='children'~'02_children',
    work_type=='Govt_job'~'03_Govt_job',
    work_type=='Private'~'04_Private',
    work_type=='Self-employed'~'05_Self_employed',
    TRUE ~ '00_ERROR'
  )))
df$work_type<-NULL

#Smoking_status
df <- df%>%
  mutate(smoking_status_disc=as.factor(case_when(
    smoking_status=='Unknown'~'01_Unknown',
    smoking_status=='never smoked'~'02_never_smoked',
    smoking_status=='smokes'~'03_smokes',
    smoking_status=='formerly smoked'~'04_formerly_smoked',
    TRUE ~ '00_ERROR'
  )))
df$smoking_status<-NULL

#Residence_type
df <- df%>%
  mutate(Residence_type_disc=as.factor(case_when(
    Residence_type=='Urban'~'01_Urban',
    Residence_type=='Rural'~'02_Rural',
    TRUE ~ '00_ERROR'
  )))
df$Residence_type<-NULL

Estudio de las discretizaciones

Se crea una función para la representación de la penetración sobre la variable target

a <- function(var1,var2) {
  df_temp <- data.frame(var1 = df[[var1]],var2 = df[[var2]])
  df_temp %>% 
    group_by(var1) %>% 
    summarise(Conteo = n(), Porc = mean(as.numeric(as.character(var2)))) %>% 
  ggplot(aes(var1,Porc)) + geom_bar(stat='identity') + xlab(var1) + ggtitle("                                                   Penetración")
}
df2_nombres <- df %>% select_if(is.factor) %>% names()
graficos<-lapply(df2_nombres,function(x){a(x,'stroke')})

Se van a visualizar las distintas variables

Hypertension

ggplot(df,aes(hypertension))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[1]]

heart_disease

ggplot(df,aes(heart_disease))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[2]]

ever_married

ggplot(df,aes(ever_married))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[3]]

age_disc

ggplot(df,aes(age_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[5]]

avg_glucose_level_disc

ggplot(df,aes(avg_glucose_level_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[6]]

bmi_disc

ggplot(df,aes(bmi_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[7]]

gender_disc

ggplot(df,aes(gender_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[8]]

work_type_disc

ggplot(df,aes(work_type_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[9]]

smoking_status_disc

ggplot(df,aes(smoking_status_disc))+geom_bar()+theme(axis.text=element_text(size=10))+ggtitle("                                                 Distribución")

graficos[[10]]

Conclusiones

Todas las variables discretizadas se comportan de forma monotónica en la penetración.

Discretización manual 2ª etapa.

Limpieza del entorno y preparación de la caché 2

Se observa el estado del df

glimpse(df)
## Rows: 4,909
## Columns: 12
## $ id                     <int> 9046, 31112, 60182, 1665, 56669, 53882, 1043...
## $ hypertension           <fct> 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1,...
## $ heart_disease          <fct> 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0,...
## $ ever_married           <fct> Yes, Yes, Yes, Yes, Yes, Yes, No, Yes, Yes, ...
## $ stroke                 <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ age_disc               <fct> 03 <= 67, 05 > 76, 03 <= 67, 05 > 76, 05 > 7...
## $ avg_glucose_level_disc <fct> 02 > 170.95, 01 <= 170.95, 02 > 170.95, 02 >...
## $ bmi_disc               <fct> 03 > 23.7, 03 > 23.7, 03 > 23.7, 03 > 23.7, ...
## $ gender_disc            <fct> 01_Male, 01_Male, 02_Female, 02_Female, 01_M...
## $ work_type_disc         <fct> 04_Private, 04_Private, 04_Private, 05_Self_...
## $ smoking_status_disc    <fct> 04_formerly_smoked, 02_never_smoked, 03_smok...
## $ Residence_type_disc    <fct> 01_Urban, 02_Rural, 01_Urban, 02_Rural, 01_U...

Se limpia el entorno

ls()
rm(list=setdiff(ls(),'df'))

Se crea la caché 2

saveRDS(df,'cache3.rds')

Se carga la caché 2

df <- readRDS(file = 'cache3.rds')

Una vez preparado todo esto, se puede empezar la modelización


5.- MODELIZACIÓN

A la hora de la creación del modelo se va a emplear:

Modelización manual

Se van a crear las funciones de:

  • Matriz confusión

  • Métricas

  • Umbrales

  • ROC

  • AUC

#Matriz de confusión
confusion <- function (real,scoring,umbral){
  conf <- table(real,scoring>=umbral)
  if(ncol(conf)==2)return(conf) else return (NULL)
}

#Metricas
metricas <- function(matriz_conf){
  acierto <- (matriz_conf[1,1]+matriz_conf[2,2])/sum(matriz_conf)*100
  precision <- matriz_conf[2,2]/(matriz_conf[2,2]+matriz_conf[1,2])*100
  cobertura <- matriz_conf[2,2]/(matriz_conf[2,2]+matriz_conf[2,1])*100
  F1<- 2*precision*cobertura/(precision+cobertura)
  salida <- c(acierto,precision,cobertura,F1)
  return(salida)
}

#Umbrales
umbrales<- function(real,scoring){
  umbrales<-data.frame(umbral=rep(0,times=19),acierto=rep(0,times=19),precision=rep(0,times=19),cobertura=rep(0,times=19),F1=rep(0,times=19))
  cont<- 1
  for(cada in seq (0.05,0.95, by=0.05)){
    datos<-metricas(confusion(real,scoring,cada))
    registro <- c(cada,datos)
    umbrales[cont,]<-registro
    cont<- cont +1
  }
  return(umbrales)
}

#ROC
roc <- function(prediction){
  r<-performance(prediction,'tpr','fpr')
  plot(r)
}

#AUC
auc <- function(prediction){
  a<-performance(prediction,'auc')
  return(a@y.values[[1]])
}

Para la hora de la creación de los modelos, es necesario dividir el data frame en dos partes:

  • Train o entrenamiento –> 70% del df

  • Test o comprobación –> 30% del df

df$random <- sample(0:1,size=nrow(df),replace=T,prob=c(0.3,0.7))
train<-filter(df,random==1)
test<-filter(df,random==0)
df$random<-NULL

Se definen las variables tanto independientes como la target

indep <- setdiff(names(df),c('stroke','id'))
target <- 'stroke'

Creación de la formula para emplearla sobre los modelos

formula<-reformulate(indep,target)

Regresión logística

Como se trata de in modelo supervisado y la variable target es dicotómica, es decir, solo toma valores 0 y 1, se realizará la regresión logística en lugar de la lineal

Creación del primer modelo

formula_rl<-formula
rl1<-glm(formula_rl,train,family=binomial(link='logit'))
summary(rl1)
## 
## Call:
## glm(formula = formula_rl, family = binomial(link = "logit"), 
##     data = train)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.0479  -0.3109  -0.1773  -0.0431   3.5736  
## 
## Coefficients:
##                                          Estimate  Std. Error z value
## (Intercept)                             -35.59056  4878.45128  -0.007
## hypertension1                             0.25137     0.22173   1.134
## heart_disease1                            0.28910     0.25808   1.120
## ever_marriedYes                           0.01837     0.32861   0.056
## age_disc02 <= 47                          2.67233     1.06665   2.505
## age_disc03 <= 67                          3.64463     1.02609   3.552
## age_disc04 <= 76                          4.31579     1.03937   4.152
## age_disc05 > 76                           4.98137     1.03189   4.827
## avg_glucose_level_disc02 > 170.95         0.77952     0.20192   3.861
## bmi_disc02 <= 23.7                       14.87753   856.71579   0.017
## bmi_disc03 > 23.7                        15.29775   856.71574   0.018
## gender_disc02_Female                      0.04080     0.18692   0.218
## gender_disc03_Other                     -13.50263 17730.36994  -0.001
## work_type_disc02_children                 0.78338  4862.48247   0.000
## work_type_disc03_Govt_job                13.42681  4802.63741   0.003
## work_type_disc04_Private                 13.41592  4802.63740   0.003
## work_type_disc05_Self_employed           13.28173  4802.63740   0.003
## smoking_status_disc02_never_smoked       -0.06744     0.26680  -0.253
## smoking_status_disc03_smokes              0.27741     0.31021   0.894
## smoking_status_disc04_formerly_smoked     0.07710     0.28509   0.270
## Residence_type_disc02_Rural               0.15657     0.18162   0.862
##                                         Pr(>|z|)    
## (Intercept)                             0.994179    
## hypertension1                           0.256932    
## heart_disease1                          0.262643    
## ever_marriedYes                         0.955428    
## age_disc02 <= 47                        0.012233 *  
## age_disc03 <= 67                        0.000382 ***
## age_disc04 <= 76                      0.00003292 ***
## age_disc05 > 76                       0.00000138 ***
## avg_glucose_level_disc02 > 170.95       0.000113 ***
## bmi_disc02 <= 23.7                      0.986145    
## bmi_disc03 > 23.7                       0.985754    
## gender_disc02_Female                    0.827228    
## gender_disc03_Other                     0.999392    
## work_type_disc02_children               0.999871    
## work_type_disc03_Govt_job               0.997769    
## work_type_disc04_Private                0.997771    
## work_type_disc05_Self_employed          0.997793    
## smoking_status_disc02_never_smoked      0.800454    
## smoking_status_disc03_smokes            0.371185    
## smoking_status_disc04_formerly_smoked   0.786824    
## Residence_type_disc02_Rural             0.388637    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1169.50  on 3425  degrees of freedom
## Residual deviance:  932.14  on 3405  degrees of freedom
## AIC: 974.14
## 
## Number of Fisher Scoring iterations: 19

Conclusiones:

  • Se van a escoger para la creación del modelo, aquellas variables en las cuales al menos uno de sus niveles posea una significación valorada con ***
a_mantener <- c('age_disc','avg_glucose_level_disc')

Creación del segundo modelo

formula_rl<-reformulate(a_mantener,target)
rl2<-glm(formula_rl,train,family=binomial(link='logit'))
summary(rl2)
## 
## Call:
## glm(formula = formula_rl, family = binomial(link = "logit"), 
##     data = train)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -0.8137  -0.2998  -0.1826  -0.0374   3.8123  
## 
## Coefficients:
##                                   Estimate Std. Error z value          Pr(>|z|)
## (Intercept)                        -7.2663     0.9999  -7.267 0.000000000000367
## age_disc02 <= 47                    3.1809     1.0551   3.015           0.00257
## age_disc03 <= 67                    4.1867     1.0107   4.142 0.000034409247856
## age_disc04 <= 76                    4.8737     1.0211   4.773 0.000001813008039
## age_disc05 > 76                     5.4495     1.0144   5.372 0.000000077872588
## avg_glucose_level_disc02 > 170.95   0.8815     0.1950   4.522 0.000006136415885
##                                      
## (Intercept)                       ***
## age_disc02 <= 47                  ** 
## age_disc03 <= 67                  ***
## age_disc04 <= 76                  ***
## age_disc05 > 76                   ***
## avg_glucose_level_disc02 > 170.95 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1169.50  on 3425  degrees of freedom
## Residual deviance:  944.27  on 3420  degrees of freedom
## AIC: 956.27
## 
## Number of Fisher Scoring iterations: 9

Se ha comprobado que las variables seleccionadas tienen una alta significación (***)

Evaluación del modelo

Correlación
pr2_rl <- 1 -(rl2$deviance / rl2$null.deviance)
pr2_rl
## [1] 0.1925915

La correlación que se obtiene es 0.1925915, lo que se traduce en:

  • El modelo es capaz de predecir el 19.25915% de los resultados de la target con dichas variables. Cabe recordar que se han seleccionado las más significativas.

  • Esta correlación se encuentra por debajo del 50%, el límite a partir del cual se considera a un modelo como bueno.

  • Como estas son las variables y los datos que proporciona la empresa, por más que se trabaje en los datos, la correlación no mejorará mucho. Por este motivo, se va a continuar estudiando el modelo.

Comportamiento del modelo en el conjunto de test
rl_predict<-predict(rl2,test,type='response')
plot(rl_predict ~ test$stroke)

El modelo es capaz de separar el scoring de los clientes quese marchan de los que se quedan

Métricas y curvas
#Umbrales
#Para el cálculo de los umbrales se va a maximizar la variable F1=2*(Precision*Cobertura)/(Precision + Cobertura)
umb_rl<- umbrales(test$stroke,rl_predict)
umb_final_rl <- umb_rl[which.max(umb_rl$F1),1]
umb_final_rl
## [1] 0.1
#Matriz de confusion
confusion(test$stroke,rl_predict,umb_final_rl)
##     
## real FALSE TRUE
##    0  3000  286
##    1    76   64
rl_metricas<-filter(umb_rl,umbral==umb_final_rl)
rl_metricas
##   umbral  acierto precision cobertura       F1
## 1    0.1 89.43374  18.28571  45.71429 26.12245
#Evaluación ROC
rl_prediction <- prediction(rl_predict,test$stroke)
roc(rl_prediction)

#Metricas
rl_metricas <- cbind(rl_metricas,AUC=round(auc(rl_prediction),2)*100)
row.names(rl_metricas)<-'rl_metricas'
print(t(rl_metricas))
##           rl_metricas
## umbral        0.10000
## acierto      89.43374
## precision    18.28571
## cobertura    45.71429
## F1           26.12245
## AUC          83.00000

Conclusiones

  • El modelo presenta un AUC = 83 siendo superior al límite de 80, como para considerarlo un buen modelo.

  • Es capaz de diferenciar el Scoring de la variable target

Árbol de decisión

Creación del primer modelo

formula_ar <- formula
ar1<- rpart(formula_ar,train,method='class')
parms=list(split="information")
control=rpart.control(cp=0.00001)

Se observa el resultado

printcp(ar1)
## 
## Classification tree:
## rpart(formula = formula_ar, data = train, method = "class")
## 
## Variables actually used in tree construction:
## character(0)
## 
## Root node error: 140/3426 = 0.040864
## 
## n= 3426 
## 
##   CP nsplit rel error xerror xstd
## 1  0      0         1      0    0
plotcp(ar1)

Los datos y las variables del df no son compatibles con el árbol de decisión, Siendo imposible desarrollar este tipo de modelo

RandomForest

Creación del primer modelo

formula_rf <- formula
rf1 <- randomForest(formula_rf,train,importance=T)
Visualización de la importancia d elas variables
varImpPlot(rf1)#Visualizacion de las variables mas importantes

#Se generan dos variables a la hora de predecir la importancia, tasa de exito y estadistico de Gini. Por eso se va a crear de ambas un unico estadistico
importancia<-importance(rf1)[,3:4]
importancia_norm <- as.data.frame(scale(importancia))
#Se crea una variable única como la suma de las dos variables
importancia_norm <- importancia_norm %>%
  mutate(Variable=rownames(importancia_norm),
         Imp_tot=MeanDecreaseAccuracy+MeanDecreaseGini)%>%
  mutate(Imp_tot=Imp_tot + abs(min(Imp_tot)))%>%
  arrange(desc(Imp_tot))%>%
  select(Variable,Imp_tot,MeanDecreaseAccuracy,MeanDecreaseGini)
importancia_norm
##                  Variable    Imp_tot MeanDecreaseAccuracy MeanDecreaseGini
## 1                age_disc 5.20554839           1.47147446        2.2573877
## 2            ever_married 2.67152714           2.05056351       -0.8557226
## 3     smoking_status_disc 2.40839538          -0.21662651        1.1483357
## 4          work_type_disc 1.36267062          -0.50403371        0.3900181
## 5             gender_disc 1.03103540          -0.07919388       -0.3664570
## 6           heart_disease 0.77442005          -0.22228593       -0.4799803
## 7     Residence_type_disc 0.77269910          -0.36986943       -0.3341177
## 8                bmi_disc 0.49043121          -0.07813149       -0.9081235
## 9  avg_glucose_level_disc 0.05013511          -1.13829595       -0.2882552
## 10           hypertension 0.00000000          -0.91360106       -0.5630852
#Histograma con la importancia de cada variable
ggplot(importancia_norm,aes(reorder(Variable,-Imp_tot),Imp_tot))+
  geom_bar(stat="identity") + theme(axis.text.x=element_text(angle=90,size=12))

El AUC más alto se obtiene sin filtrar las variables, por eso, se van a mantener todas

Comportamiento del modelo en el conjunto de test

rf_predict <- predict(rf1,test,type='prob')[,2]
plot(rf_predict ~ test$stroke)

El modelo es capaz de separar el scoring de los clientes quese marchan de los que se quedan

Métricas y curvas

#Umbrales
#Para el cálculo de los umbrales se va a maximizar la variable F1=2*(Precision*Cobertura)/(Precision + Cobertura)
umb_rf<- umbrales(test$stroke,rf_predict)
umb_final_rf <- umb_rf[which.max(umb_rf$F1),1]
umb_final_rf
## [1] 0.15
#Matriz de confusion
confusion(test$stroke,rf_predict,umb_final_rf)
##     
## real FALSE TRUE
##    0  3257   29
##    1    77   63
rf_metricas<-filter(umb_rf,umbral==umb_final_rf)
rf_metricas
##   umbral  acierto precision cobertura       F1
## 1   0.15 96.90601  68.47826        45 54.31034
#Evaluación ROC
rf_prediction <- prediction(rf_predict,test$stroke)
roc(rf_prediction)

#Metricas
rf_metricas <- cbind(rf_metricas,AUC=round(auc(rf_prediction),2)*100)
row.names(rf_metricas)<-'rf_metricas'
print(t(rf_metricas))
##           rf_metricas
## umbral        0.15000
## acierto      96.90601
## precision    68.47826
## cobertura    45.00000
## F1           54.31034
## AUC          88.00000

Conclusiones

  • El modelo presenta un AUC = 88, siendo superior al límite de 80, como para considerarlo un buen modelo.

  • Es capaz de diferenciar el Scoring de la variable target

Comparativa de los modelos manuales

Se va a generar una tabla, en la cual se van a arecoger las distintas métricas de los 2 modelos creados

comparativa <- rbind(rl_metricas,rf_metricas)
rownames(comparativa)<- c('Regresion Logistica','RandomForest')
comparativa<-t(comparativa)
comparativa%>%knitr::kable()
Regresion Logistica RandomForest
umbral 0.10000 0.15000
acierto 89.43374 96.90601
precision 18.28571 68.47826
cobertura 45.71429 45.00000
F1 26.12245 54.31034
AUC 83.00000 88.00000

Atendiendo al valor del AUC, el modelo que mejor se comporta es el randomForest.

Por este motivo será el seleccionado como mejor modelo dentro de los modelos manuales.

Se guarda el randomForest para dárselo al cliente y lo utilice en futuras ocasiones

saveRDS(rf1,'01_Modelo_final_manual')

Tabla con la probabilidad de sufir un infarto parte de los pacientes que no lo han sufrido

df$SCORING_RF_MANUAL <- predict(rf1,df,type='prob')
df$SCORING_RF_MANUAL <- df$SCORING_RF_MANUAL[,2]
STROKE_RF_MANUAL <- df%>%
  filter(stroke==0)%>%
  select(id,SCORING_RF_MANUAL)%>%
  arrange(desc(SCORING_RF_MANUAL))
head(STROKE_RF_MANUAL)
##      id SCORING_RF_MANUAL
## 1 16066             0.582
## 2  9879             0.492
## 3  7218             0.408
## 4 53697             0.366
## 5 46670             0.354
## 6 44781             0.334
tail(STROKE_RF_MANUAL)
##         id SCORING_RF_MANUAL
## 4695 36901                 0
## 4696 45010                 0
## 4697 22127                 0
## 4698 14180                 0
## 4699 19723                 0
## 4700 37544                 0

Modelización automática

Como se expuso anteriormente, para la realización de la modelización automática se va a emplear H2O

Preparación del entorno

Preparación del cluster de h2o

h2o.init()

Translado de los datos al cluster

df1<-readRDS(file = 'cache3.rds')#Se crea una copia nueva del df
df_h2o<-as.h2o(df1)

Partición del df

split<-h2o.splitFrame(df_h2o)
train_h2o<-split[[1]]
test_h2o<-split[[2]]

Definición de los roles de las variables

y <- 'stroke'
x <- setdiff(names(df_h2o),c('id',y))

Creación de los modelos

Como factor limitante se va a definir la duración en tiempo, concretamente 10 minutos con 3 validaciones cruzadas por modelo

Como en la modelización manual, el ranking se hará en función del AUC

automl_simple <- h2o.automl(
  y=y,
  x=x,
  training_frame = train_h2o,
  validation_frame = test_h2o,
  nfolds=3,
  max_runtime_secs = 600,
  sort_metric='AUC')

Estudio de los modelos generados

automl_simple@leaderboard
##                                     model_id       auc   logloss     aucpr
## 1 GBM_grid__1_AutoML_20210224_114807_model_4 0.8319074 0.1475156 0.1980938
## 2               GBM_5_AutoML_20210224_114807 0.8285628 0.1477594 0.1771652
## 3 GBM_grid__1_AutoML_20210224_114807_model_8 0.8281112 0.1483646 0.1794552
## 4 GBM_grid__1_AutoML_20210224_114807_model_3 0.8270907 0.1479035 0.1879056
## 5 GBM_grid__1_AutoML_20210224_114807_model_6 0.8243547 0.1486819 0.1882354
## 6 GBM_grid__1_AutoML_20210224_114807_model_5 0.8211035 0.1489933 0.1785758
##   mean_per_class_error      rmse        mse
## 1            0.3450292 0.1952802 0.03813435
## 2            0.3268587 0.1960030 0.03841717
## 3            0.3130904 0.1960853 0.03844946
## 4            0.3157515 0.1960401 0.03843171
## 5            0.3475668 0.1961838 0.03848809
## 6            0.3212156 0.1962395 0.03850995
## 
## [27 rows x 7 columns]

Se presenta una tabla donde aparecen los 6 modelos más significativos de un total de 27 modelos que ha generado H2O

#Visualizacion grafica de los modelos
as.data.frame(automl_simple@leaderboard)%>%
  select(model_id,auc)%>%
  ggplot(aes(x=auc,y=reorder(model_id,auc)))+
  geom_point()+geom_label(aes(label=round(auc,3),color=auc),hjust='left')+expand_limits(x=c(0.90,0.820))+theme_bw()

Se visualiza una gráfica donde se ordenan todos los modelos generados en función del AUC. De todos ellos se escogerá el de mayor AUC

Selección del modelo ganador

automl_simple_winner <- automl_simple@leader
#Metricas
automl_simple_winner@model$cross_validation_metrics
## H2OBinomialMetrics: gbm
## ** Reported on cross-validation data. **
## ** 3-fold cross-validation on training data (Metrics computed for combined holdout predictions) **
## 
## MSE:  0.03813435
## RMSE:  0.1952802
## LogLoss:  0.1475156
## Mean Per-Class Error:  0.3450292
## AUC:  0.8319074
## AUCPR:  0.1980938
## Gini:  0.6638148
## R^2:  0.08800596
## 
## Confusion Matrix (vertical: actual; across: predicted) for F1-optimal threshold:
##           0   1    Error       Rate
## 0      3322 199 0.056518  =199/3521
## 1       102  59 0.633540   =102/161
## Totals 3424 258 0.081749  =301/3682
## 
## Maximum Metrics: Maximum metrics at their respective thresholds
##                         metric threshold       value idx
## 1                       max f1  0.138293    0.281623 121
## 2                       max f2  0.049987    0.402109 254
## 3                 max f0point5  0.187395    0.280650  78
## 4                 max accuracy  0.497393    0.956545   0
## 5                max precision  0.497393    1.000000   0
## 6                   max recall  0.008222    1.000000 387
## 7              max specificity  0.497393    1.000000   0
## 8             max absolute_mcc  0.049987    0.261781 254
## 9   max min_per_class_accuracy  0.047562    0.770186 260
## 10 max mean_per_class_accuracy  0.040295    0.774286 279
## 11                     max tns  0.497393 3521.000000   0
## 12                     max fns  0.497393  160.000000   0
## 13                     max fps  0.004906 3521.000000 399
## 14                     max tps  0.008222  161.000000 387
## 15                     max tnr  0.497393    1.000000   0
## 16                     max fnr  0.497393    0.993789   0
## 17                     max fpr  0.004906    1.000000 399
## 18                     max tpr  0.008222    1.000000 387
## 
## Gains/Lift Table: Extract with `h2o.gainsLift(<model>, <data>)` or `h2o.gainsLift(<model>, valid=<T/F>, xval=<T/F>)`
h2o.auc(automl_simple_winner@model$cross_validation_metrics)
## [1] 0.8319074
AUC_winner<-h2o.auc(automl_simple_winner@model$cross_validation_metrics)
#Importancia de las variables
h2o.varimp_plot(automl_simple@leader)

#Se calcula el scoring
SCORING_H2O_STROKE <- as.data.frame(h2o.predict(automl_simple_winner,df_h2o)[3])
df$SCORING_GBM_H2O <- as.numeric(SCORING_H2O_STROKE$p1)
glimpse(df)
## Rows: 4,909
## Columns: 13
## $ id                     <int> 9046, 31112, 60182, 1665, 56669, 53882, 1043...
## $ hypertension           <fct> 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1,...
## $ heart_disease          <fct> 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0,...
## $ ever_married           <fct> Yes, Yes, Yes, Yes, Yes, Yes, No, Yes, Yes, ...
## $ stroke                 <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,...
## $ age_disc               <fct> 03 <= 67, 05 > 76, 03 <= 67, 05 > 76, 05 > 7...
## $ avg_glucose_level_disc <fct> 02 > 170.95, 01 <= 170.95, 02 > 170.95, 02 >...
## $ bmi_disc               <fct> 03 > 23.7, 03 > 23.7, 03 > 23.7, 03 > 23.7, ...
## $ gender_disc            <fct> 01_Male, 01_Male, 02_Female, 02_Female, 01_M...
## $ work_type_disc         <fct> 04_Private, 04_Private, 04_Private, 05_Self_...
## $ smoking_status_disc    <fct> 04_formerly_smoked, 02_never_smoked, 03_smok...
## $ Residence_type_disc    <fct> 01_Urban, 02_Rural, 01_Urban, 02_Rural, 01_U...
## $ SCORING_GBM_H2O        <dbl> 0.15046171, 0.23643685, 0.07776229, 0.262670...

Se guarda el modelo ganador

h2o.saveModel(automl_simple_winner,path='/Users/Alberto/Desktop/Proyecto Data Science/Stroke Prediction')
automl_simple_winner <- h2o.loadModel('/Users/Alberto/Desktop/Proyecto Data Science/Stroke Prediction/GBM_grid_1_AutoML_20210224_114807_model_4')

Tabla con el scoring de abandono de los clientes que no se han ido

STROKE_GBM_H2O <- df%>%
  filter(stroke==0)%>%
  select(id,SCORING_GBM_H2O)%>%
  arrange(desc(SCORING_GBM_H2O))
head(STROKE_GBM_H2O)
##      id SCORING_GBM_H2O
## 1 54353       0.3498449
## 2  1473       0.3400145
## 3 31426       0.3259184
## 4 14000       0.3224669
## 5  8009       0.3224669
## 6 63836       0.3222371
tail(STROKE_GBM_H2O)
##         id SCORING_GBM_H2O
## 4695 56282     0.005695533
## 4696  7868     0.005464337
## 4697  8341     0.005464337
## 4698 42393     0.005464337
## 4699 53489     0.005464337
## 4700 49179     0.005464337

6.-CONCLUSIÓN FINAL

Tanto por la modelización manual como por la automática se han conseguido dos modelos:

Como Random Forest presentan AUC más alto, se va a elegir como modelo ganador

Por lo tanto, al cliente se le entregará:

head(STROKE_RF_MANUAL)%>%knitr::kable()
id SCORING_RF_MANUAL
16066 0.582
9879 0.492
7218 0.408
53697 0.366
46670 0.354
44781 0.334