Predicción de enfermedad cardiovasculares cardiovasculares para la competición de Drivendata: https://www.drivendata.org/competitions/54/machine-learning-with-a-heart/.
Se trata de un problema de clasificación binaria, para intentar predecir la presencia o no de enfermedad cardiovascular cardiovascular en un paciente según una serie de variables clínicas.
El dataset proviende de la Cleveland Heart Disease Database disponible en el UCI Machine Learning Repository
Para el modelado utilizaremos la librería H2O.AI, con la función de AutoML.
Bibliografía: Aha, D., and Dennis Kibler. “Instance-based prediction of heart-disease presence with the Cleveland database.” University of California 3.1 (1988): 3-2.
El dataset contiene los siguientes archivos:
El orden de pacientes es el mismo para los dos archivos del conjunto de entrenamiento, por lo que podemos simplemente unir la variable target al dataframe inicial.
Uniremos los conjuntos de entrenamiento y validación para hacer el análisis inicial y transformación de variables, que luego separaremos en la fase de modelización.
Del sitio web de la competición encontramos una descripcción de las variables independientes y de la target. Apuntamos también las acciones que tenemos que llevar a cabo en cada variable.
| Variable | Descripción | Acción |
|---|---|---|
| patient_id (type:factor) | Identificación del paciente | No usar |
| slope_of_peak_exercise_st_segment (type: int) | Inclinación del segmento ST del EKG durante el ejercicio. | Pasar a factor |
| thal (type:categorical) | resultado del test de stress con talio. Los niveles son: “normal”,“fixed_defect”,“reversible_defect” | No precisa |
| “resting_blood_pressure (type: int)” | Tensión arterial sistólica en reposo | No precisa |
| “chest_pain_type (type: int)” | Tipo de dolor torácico, no tiene un texto asociado, entendemos que 0 es la ausencia de dolor, y el resto de valores algún tipo de escala | Pasar a factor |
| “num_major_vessels (type: int)” | Número de vasos principales coloreados por fluoroscopia | Pasar a factor |
| “fasting_blood_sugar_gt_120_mg_per_dl(type: binary)” | Glucemia en ayunas mayor de 120mg/dl, variable dicotómica (0/1), aunque se importa como numérica | Pasar a factor |
| “resting_ekg_results” | Resultado del EKG en reposo. Tres niveles numéricos sin descripción, suponemos que 0 representa la ausencia de alteraciones | Pasar a factor |
| “serum_cholesterol_mg_per_dl (type: int)” | Nivel de colesterol en sangre expresado en mg/dl. | Mantener |
| “oldpeak_eq_st_depression” (type: float) | Descenso del segmento ST en el EKG durante el ejercicio en relación con el EKG en reposo | Mantener |
| “sex (type: binary)” | Sexo, variable dicotómica (0/1) | Pasar a factor |
| “age (type: int)” | oldpeak_eq_st_depression, variable continua | Valorar discretización |
| “max_heart_rate_achieved (type: int)” | Frecuencia cardiaca máxima alcanzada en la prueba de esfuerzo | Mantener |
| “exercise_induced_angina (type: binary)” | Angina producida por el ejercicio, variable dicotómica (0/1) | Pasar a factor |
| TARGET | Presencia o no de enfermedad cardiovascular cardiovascular, originalmente llamada “heart_disease_present” en el listado “train_labels.csv” | Pasar a factor |
Transformamos las variables según este primer análisis:
# pasamos las variables numéricas a factor
a_factores <- c("slope_of_peak_exercise_st_segment","chest_pain_type","num_major_vessels","resting_ekg_results","fasting_blood_sugar_gt_120_mg_per_dl","sex","age","exercise_induced_angina")
df[a_factores] <- lapply(df[a_factores],as.factor)
Hacemos un primer vistazo con skimr:
| variable | missing | complete | n | n_unique | top_counts | ordered |
|---|---|---|---|---|---|---|
| age | 0 | 270 | 270 | 41 | 54: 16, 58: 15, 51: 12, 57: 12 | FALSE |
| chest_pain_type | 0 | 270 | 270 | 4 | 4: 129, 3: 79, 2: 42, 1: 20 | FALSE |
| exercise_induced_angina | 0 | 270 | 270 | 2 | 0: 181, 1: 89, NA: 0 | FALSE |
| fasting_blood_sugar_gt_120_mg_per_dl | 0 | 270 | 270 | 2 | 0: 230, 1: 40, NA: 0 | FALSE |
| num_major_vessels | 0 | 270 | 270 | 4 | 0: 160, 1: 58, 2: 33, 3: 19 | FALSE |
| patient_id | 0 | 270 | 270 | 270 | 02c: 1, 08u: 1, 0g1: 1, 0n5: 1 | FALSE |
| resting_ekg_results | 0 | 270 | 270 | 3 | 2: 137, 0: 131, 1: 2, NA: 0 | FALSE |
| sex | 0 | 270 | 270 | 2 | 1: 183, 0: 87, NA: 0 | FALSE |
| slope_of_peak_exercise_st_segment | 0 | 270 | 270 | 3 | 1: 130, 2: 122, 3: 18, NA: 0 | FALSE |
| thal | 0 | 270 | 270 | 3 | nor: 152, rev: 104, fix: 14, NA: 0 | FALSE |
| variable | missing | complete | n | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|---|
| max_heart_rate_achieved | 0 | 270 | 270 | 149.68 | 23.17 | 71 | 133 | 153.5 | 166 | 202 | ▁▁▂▃▆▇▅▁ |
| resting_blood_pressure | 0 | 270 | 270 | 131.34 | 17.86 | 94 | 120 | 130 | 140 | 200 | ▂▇▇▆▃▁▁▁ |
| serum_cholesterol_mg_per_dl | 0 | 270 | 270 | 249.66 | 51.69 | 126 | 213 | 245 | 280 | 564 | ▂▇▇▃▁▁▁▁ |
| variable | missing | complete | n | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|---|
| oldpeak_eq_st_depression | 0 | 270 | 270 | 1.05 | 1.15 | 0 | 0 | 0.8 | 1.6 | 6.2 | ▇▃▂▁▁▁▁▁ |
| TARGET | 0 | 270 | 270 | 0.3 | 0.46 | 0 | 0 | 0 | 1 | 1 | ▇▁▁▁▁▁▁▃ |
No observamos variables con valores nulos. La distribución de las variables numéricas parecen tener distribuciones que podríamos considerar dentro de la normalidad para el caso que estamos estudiando.
Encontramos en la variable de colesterol en sangre valores mayores de 400, que si bien son muy altos, pueden darse en la realidad, por lo que mantenemos estos registros.
Utilizamos corrplot para hacer un gráfico de correlación entre variables
corrplot(cor((sapply(df,as.integer))),method="number")
## Warning in lapply(X = X, FUN = FUN, ...): NAs introducidos por coerción
Observamos que la correlación entre variables es débil, por lo que no haremos ninguna transformación.
table(df$TARGET)
##
## 0 1
## 190 80
ggplot (data=df,aes(x=TARGET,fill=TARGET)) + geom_bar() + guides(fill=FALSE)
La penetración de la variable target es del 44%, por lo que no existe desbalanceo de clases y no será necesario ninguna actuación en este sentido.
Transformamos la variable TARGET a factor
df <- df %>%
mutate (TARGET = as.factor(as.numeric(TARGET)))
La fase de preselección de variables independientes en este caso la dejaremos, ya que el conjunto de datos presenta un número muy reducido de las mismas y todas son conocidas como factores relacionados con la patología cardiovascular, tanto variables personales como los hallazgos electrocardiográficos.
Discretizamos la variable “resting_blood_pressure” según la clasificación de tensión arterial.
# discretización de "resting_blood_pressure"
df <- df %>% mutate(resting_blood_pressure_DISC = as.factor(case_when(
between(resting_blood_pressure,-Inf,120) ~ '01_MENOS_DE_120',
between(resting_blood_pressure,120,130) ~ '02_DE_120_A_130',
between(resting_blood_pressure,130,140) ~ '03_DE_130_A_140',
between(resting_blood_pressure,140,160) ~ '04_DE_140_A_160',
between(resting_blood_pressure,160,180) ~ '04_DE_160_A_180',
between(resting_blood_pressure,180,Inf) ~ '05_MAYOR_DE_180',
TRUE ~ '00_ERROR')))%>%
select(-resting_blood_pressure)
Para la modelización utilizaremos el algoritmo AutoML de H2O, que crea modelos sobre un conjunto de datos, con diversos hiperparámetros y nos devuelve un ranking con los mejores modelos creados.
Preparación de datos para la modelización con H2O
# quitamos los datos de test del conjunto de entrenamiento
df <- df %>%
filter(CONJUNTO=="train")
Cargamos el cluster H2O y le pasamos el dataframe para que lo pueda utilizar.
Creamos el listado de variables independientes.
quitar <- c("patient_id","TARGET","CONJUNTO")
independientes <- setdiff(names(df),quitar)
Iniciamos H2O AutoML, seleccionamos una validación cruzada de 5 partes. Dejamos que cree modelos durante 5 minutos y que los puntúe según la AUC obtenida.
automl <- h2o.automl(
y = "TARGET",
x = independientes,
training_frame = dfh,
nfolds = 5,
max_runtime_secs = 300,
sort_metric = 'AUC'
)
##
|
| | 0%
|
| | 1%
|
|= | 1%
|
|= | 2%
|
|== | 3%
|
|=== | 4%
|
|=== | 5%
|
|===== | 7%
|
|===== | 8%
|
|========= | 13%
|
|========== | 16%
|
|============ | 18%
|
|============== | 21%
|
|=============== | 24%
|
|=================== | 29%
|
|==================== | 31%
|
|===================== | 32%
|
|====================== | 34%
|
|======================= | 35%
|
|======================= | 36%
|
|======================== | 36%
|
|======================== | 37%
|
|======================== | 38%
|
|========================= | 38%
|
|======================================= | 61%
|
|======================================== | 61%
|
|================================================== | 76%
|
|================================================== | 77%
|
|================================================== | 78%
|
|=================================================== | 78%
|
|=================================================== | 79%
|
|==================================================== | 79%
|
|==================================================== | 80%
|
|==================================================== | 81%
|
|===================================================== | 81%
|
|===================================================== | 82%
|
|====================================================== | 82%
|
|====================================================== | 83%
|
|====================================================== | 84%
|
|======================================================= | 84%
|
|======================================================= | 85%
|
|======================================================== | 86%
|
|======================================================== | 87%
|
|========================================================= | 87%
|
|========================================================= | 88%
|
|========================================================== | 89%
|
|========================================================== | 90%
|
|=========================================================== | 90%
|
|=========================================================== | 91%
|
|============================================================ | 92%
|
|============================================================== | 96%
|
|=================================================================| 100%
Visualizamos la AUC de los 10 mejores modelos
as.data.frame(automl@leaderboard[1:10,]) %>%
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.936,0.947)) +
theme_bw()
Extraemos el mejor modelo del objeto automl
ganador <- automl@leader
La evaluación del modelo tenemos que hacerla sobre la validación cruzada que hace H2O, ya que el dataset original solo incluye las etiquetas de los datos para el conjunto de entrenamiento
Visualizamos las métricas del modelo:
#Sacamos solo las metricas de validación cruzada
ganador@model$cross_validation_metrics
## H2OBinomialMetrics: xgboost
## ** Reported on cross-validation data. **
## ** 5-fold cross-validation on training data (Metrics computed for combined holdout predictions) **
##
## MSE: 0.1364117
## RMSE: 0.3693396
## LogLoss: 0.4317517
## Mean Per-Class Error: 0.1875
## AUC: 0.88425
## pr_auc: 0.8565171
## Gini: 0.7685
##
## Confusion Matrix (vertical: actual; across: predicted) for F1-optimal threshold:
## 0 1 Error Rate
## 0 80 20 0.200000 =20/100
## 1 14 66 0.175000 =14/80
## Totals 94 86 0.188889 =34/180
##
## Maximum Metrics: Maximum metrics at their respective thresholds
## metric threshold value idx
## 1 max f1 0.418090 0.795181 84
## 2 max f2 0.286911 0.856481 110
## 3 max f0point5 0.603685 0.828488 64
## 4 max accuracy 0.603685 0.822222 64
## 5 max precision 0.950208 1.000000 0
## 6 max recall 0.085134 1.000000 164
## 7 max specificity 0.950208 1.000000 0
## 8 max absolute_mcc 0.603685 0.641889 64
## 9 max min_per_class_accuracy 0.459205 0.800000 80
## 10 max mean_per_class_accuracy 0.510898 0.817500 74
##
## Gains/Lift Table: Extract with `h2o.gainsLift(<model>, <data>)` or `h2o.gainsLift(<model>, valid=<T/F>, xval=<T/F>)`
#Sacamos concretamente el AUC
h2o.auc(ganador@model$cross_validation_metrics)
## [1] 0.88425
Y podemos ver las 10 variables que están resultando más importantes
kable(as.data.frame(h2o.varimp(automl@leader)[1:10,]))
| variable | relative_importance | scaled_importance | percentage |
|---|---|---|---|
| chest_pain_type.4 | 182.05032 | 1.0000000 | 0.1951036 |
| thal.reversible_defect | 144.59659 | 0.7942671 | 0.1549644 |
| num_major_vessels.0 | 117.90881 | 0.6476715 | 0.1263631 |
| oldpeak_eq_st_depression | 117.66848 | 0.6463514 | 0.1261055 |
| thal.normal | 106.26064 | 0.5836883 | 0.1138797 |
| exercise_induced_angina.0 | 82.08660 | 0.4509006 | 0.0879723 |
| max_heart_rate_achieved | 40.32474 | 0.2215032 | 0.0432161 |
| serum_cholesterol_mg_per_dl | 27.14366 | 0.1490997 | 0.0290899 |
| slope_of_peak_exercise_st_segment.2 | 23.82738 | 0.1308835 | 0.0255358 |
| exercise_induced_angina.1 | 20.12483 | 0.1105454 | 0.0215678 |
Encontramos que variables como el EKG en reposo, la tensión arterial, la edad, el número de vasos y la glucemia en ayunas son las que el modelo considera más importantes a la hora de predecir la enfermedad cardiovascular, resultados concordantes con los factores de riesgo cardiovascular.