Regresión Logística — Enfermedad Cardíaca

Teoría

La regresión logística estima la probabilidad de que ocurra un evento binario (0 o 1) en función de variables predictoras. A diferencia de la regresión lineal, la salida siempre está acotada entre 0 y 1 gracias a la función sigmoide:

\[P(Y=1) = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X_1 + \cdots + \beta_n X_n)}}\]

En este análisis buscamos predecir si un paciente tiene o no enfermedad cardíaca (target = 1 o 0).

---

Instalar paquetes y llamar librerías

# install.packages("tidyverse")
# install.packages("caret")
library(tidyverse)
library(caret)

---

Importar la base de datos

df <- read.csv("heart.csv")

Diccionario de variables:

Variable	Descripción	Tipo
age	Edad del paciente (años)	Numérica
sex	Sexo (1 = masculino, 0 = femenino)	Categórica
cp	Tipo de dolor en el pecho (0–3)	Categórica
trestbps	Presión arterial en reposo (mmHg)	Numérica
chol	Colesterol sérico (mg/dl)	Numérica
fbs	Azúcar en sangre en ayunas > 120 mg/dl (1 = sí)	Categórica
restecg	Resultados del ECG en reposo (0–2)	Categórica
thalach	Frecuencia cardíaca máxima alcanzada	Numérica
exang	Angina inducida por ejercicio (1 = sí)	Categórica
oldpeak	Depresión del ST inducida por ejercicio	Numérica
slope	Pendiente del segmento ST en ejercicio pico (0–2)	Categórica
ca	Nº de vasos principales coloreados por fluoroscopía (0–3)	Categórica
thal	Talasemia (1 = normal, 2 = defecto fijo, 3 = reversible)	Categórica
target	Enfermedad cardíaca (1 = sí, 0 = no)	Objetivo

---

Entender la base de datos

dim(df)

## [1] 1025   14

summary(df)

##       age             sex               cp            trestbps    
##  Min.   :29.00   Min.   :0.0000   Min.   :0.0000   Min.   : 94.0  
##  1st Qu.:48.00   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:120.0  
##  Median :56.00   Median :1.0000   Median :1.0000   Median :130.0  
##  Mean   :54.43   Mean   :0.6956   Mean   :0.9424   Mean   :131.6  
##  3rd Qu.:61.00   3rd Qu.:1.0000   3rd Qu.:2.0000   3rd Qu.:140.0  
##  Max.   :77.00   Max.   :1.0000   Max.   :3.0000   Max.   :200.0  
##       chol          fbs            restecg          thalach     
##  Min.   :126   Min.   :0.0000   Min.   :0.0000   Min.   : 71.0  
##  1st Qu.:211   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:132.0  
##  Median :240   Median :0.0000   Median :1.0000   Median :152.0  
##  Mean   :246   Mean   :0.1493   Mean   :0.5298   Mean   :149.1  
##  3rd Qu.:275   3rd Qu.:0.0000   3rd Qu.:1.0000   3rd Qu.:166.0  
##  Max.   :564   Max.   :1.0000   Max.   :2.0000   Max.   :202.0  
##      exang           oldpeak          slope             ca        
##  Min.   :0.0000   Min.   :0.000   Min.   :0.000   Min.   :0.0000  
##  1st Qu.:0.0000   1st Qu.:0.000   1st Qu.:1.000   1st Qu.:0.0000  
##  Median :0.0000   Median :0.800   Median :1.000   Median :0.0000  
##  Mean   :0.3366   Mean   :1.072   Mean   :1.385   Mean   :0.7541  
##  3rd Qu.:1.0000   3rd Qu.:1.800   3rd Qu.:2.000   3rd Qu.:1.0000  
##  Max.   :1.0000   Max.   :6.200   Max.   :2.000   Max.   :4.0000  
##       thal           target      
##  Min.   :0.000   Min.   :0.0000  
##  1st Qu.:2.000   1st Qu.:0.0000  
##  Median :2.000   Median :1.0000  
##  Mean   :2.324   Mean   :0.5132  
##  3rd Qu.:3.000   3rd Qu.:1.0000  
##  Max.   :3.000   Max.   :1.0000

str(df)

## 'data.frame':    1025 obs. of  14 variables:
##  $ age     : int  52 53 70 61 62 58 58 55 46 54 ...
##  $ sex     : int  1 1 1 1 0 0 1 1 1 1 ...
##  $ cp      : int  0 0 0 0 0 0 0 0 0 0 ...
##  $ trestbps: int  125 140 145 148 138 100 114 160 120 122 ...
##  $ chol    : int  212 203 174 203 294 248 318 289 249 286 ...
##  $ fbs     : int  0 1 0 0 1 0 0 0 0 0 ...
##  $ restecg : int  1 0 1 1 1 0 2 0 0 0 ...
##  $ thalach : int  168 155 125 161 106 122 140 145 144 116 ...
##  $ exang   : int  0 1 1 0 0 0 0 1 0 1 ...
##  $ oldpeak : num  1 3.1 2.6 0 1.9 1 4.4 0.8 0.8 3.2 ...
##  $ slope   : int  2 0 0 2 1 1 0 1 2 1 ...
##  $ ca      : int  2 0 0 1 3 0 3 1 0 2 ...
##  $ thal    : int  3 3 3 3 2 2 1 3 3 2 ...
##  $ target  : int  0 0 0 0 0 1 0 0 0 0 ...

head(df)

##   age sex cp trestbps chol fbs restecg thalach exang oldpeak slope ca thal
## 1  52   1  0      125  212   0       1     168     0     1.0     2  2    3
## 2  53   1  0      140  203   1       0     155     1     3.1     0  0    3
## 3  70   1  0      145  174   0       1     125     1     2.6     0  0    3
## 4  61   1  0      148  203   0       1     161     0     0.0     2  1    3
## 5  62   0  0      138  294   1       1     106     0     1.9     1  3    2
## 6  58   0  0      100  248   0       0     122     0     1.0     1  0    2
##   target
## 1      0
## 2      0
## 3      0
## 4      0
## 5      0
## 6      1

---

Limpieza y preparación

Valores inválidos

thal = 0 y ca = 4 están fuera del rango clínico válido (errores de registro). Los convertimos a NA y eliminamos esas filas.

df$thal[df$thal == 0] <- NA
df$ca[df$ca == 4]     <- NA

cat("Filas antes de limpiar: ", nrow(df), "\n")

## Filas antes de limpiar:  1025

df <- na.omit(df)
cat("Filas después de limpiar:", nrow(df), "\n")

## Filas después de limpiar: 1000

Conversión de variables categóricas a factores

df$target  <- as.factor(df$target)
df$sex     <- as.factor(df$sex)
df$cp      <- as.factor(df$cp)
df$fbs     <- as.factor(df$fbs)
df$restecg <- as.factor(df$restecg)
df$exang   <- as.factor(df$exang)
df$slope   <- as.factor(df$slope)
df$ca      <- as.factor(df$ca)
df$thal    <- as.factor(df$thal)

str(df)

## 'data.frame':    1000 obs. of  14 variables:
##  $ age     : int  52 53 70 61 62 58 58 55 46 54 ...
##  $ sex     : Factor w/ 2 levels "0","1": 2 2 2 2 1 1 2 2 2 2 ...
##  $ cp      : Factor w/ 4 levels "0","1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
##  $ trestbps: int  125 140 145 148 138 100 114 160 120 122 ...
##  $ chol    : int  212 203 174 203 294 248 318 289 249 286 ...
##  $ fbs     : Factor w/ 2 levels "0","1": 1 2 1 1 2 1 1 1 1 1 ...
##  $ restecg : Factor w/ 3 levels "0","1","2": 2 1 2 2 2 1 3 1 1 1 ...
##  $ thalach : int  168 155 125 161 106 122 140 145 144 116 ...
##  $ exang   : Factor w/ 2 levels "0","1": 1 2 2 1 1 1 1 2 1 2 ...
##  $ oldpeak : num  1 3.1 2.6 0 1.9 1 4.4 0.8 0.8 3.2 ...
##  $ slope   : Factor w/ 3 levels "0","1","2": 3 1 1 3 2 2 1 2 3 2 ...
##  $ ca      : Factor w/ 4 levels "0","1","2","3": 3 1 1 2 4 1 4 2 1 3 ...
##  $ thal    : Factor w/ 3 levels "1","2","3": 3 3 3 3 2 2 1 3 3 2 ...
##  $ target  : Factor w/ 2 levels "0","1": 1 1 1 1 1 2 1 1 1 1 ...
##  - attr(*, "na.action")= 'omit' Named int [1:25] 15 53 84 129 209 243 291 320 330 341 ...
##   ..- attr(*, "names")= chr [1:25] "15" "53" "84" "129" ...

Distribución de la variable objetivo

table(df$target)

## 
##   0   1 
## 492 508

prop.table(table(df$target))

## 
##     0     1 
## 0.492 0.508

---

Exploración visual

etiquetas <- c("0" = "Sin enfermedad", "1" = "Con enfermedad")

# Distribución del target
ggplot(df, aes(x = target, fill = target)) +
  geom_bar() +
  scale_x_discrete(labels = etiquetas) +
  scale_fill_manual(values = c("0" = "#74a9d8", "1" = "#e05c5c")) +
  labs(title = "Distribución de la variable objetivo", x = "", y = "Pacientes") +
  theme_minimal() + theme(legend.position = "none")

# Edad según diagnóstico
ggplot(df, aes(x = target, y = age, fill = target)) +
  geom_boxplot() +
  scale_x_discrete(labels = etiquetas) +
  scale_fill_manual(values = c("0" = "#74a9d8", "1" = "#e05c5c")) +
  labs(title = "Edad según diagnóstico", x = "", y = "Edad (años)") +
  theme_minimal() + theme(legend.position = "none")

# Frecuencia cardíaca máxima según diagnóstico
ggplot(df, aes(x = target, y = thalach, fill = target)) +
  geom_boxplot() +
  scale_x_discrete(labels = etiquetas) +
  scale_fill_manual(values = c("0" = "#74a9d8", "1" = "#e05c5c")) +
  labs(title = "Frecuencia cardíaca máxima según diagnóstico",
       x = "", y = "Freq. cardíaca máx. (bpm)") +
  theme_minimal() + theme(legend.position = "none")

# Tipo de dolor de pecho según diagnóstico
ggplot(df, aes(x = cp, fill = target)) +
  geom_bar(position = "fill") +
  scale_fill_manual(values = c("0" = "#74a9d8", "1" = "#e05c5c"),
                    labels = etiquetas) +
  labs(title = "Proporción de diagnóstico por tipo de dolor",
       x = "Tipo de dolor en el pecho (0–3)", y = "Proporción", fill = "") +
  theme_minimal()

---

División en entrenamiento y prueba

Usamos 70% para entrenar el modelo y 30% para evaluarlo sobre datos que no vio.

set.seed(42)
idx   <- sample(1:nrow(df), size = 0.7 * nrow(df))
train <- df[idx, ]
test  <- df[-idx, ]

cat("Train:", nrow(train), "filas\n")

## Train: 700 filas

cat("Test: ", nrow(test),  "filas\n")

## Test:  300 filas

---

Modelo completo

Ajustamos un primer modelo con todas las variables para ver cuáles son significativas.

modelo_completo <- glm(target ~ ., data = train, family = binomial)
summary(modelo_completo)

## 
## Call:
## glm(formula = target ~ ., family = binomial, data = train)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  3.943258   2.004327   1.967 0.049140 *  
## age          0.010232   0.017034   0.601 0.548059    
## sex1        -2.074022   0.385453  -5.381 7.42e-08 ***
## cp1          0.730420   0.377098   1.937 0.052751 .  
## cp2          1.914424   0.363680   5.264 1.41e-07 ***
## cp3          2.761132   0.499559   5.527 3.25e-08 ***
## trestbps    -0.021819   0.008100  -2.694 0.007068 ** 
## chol        -0.003331   0.002852  -1.168 0.242862    
## fbs1         0.406982   0.406241   1.002 0.316429    
## restecg1     0.150416   0.276272   0.544 0.586132    
## restecg2    -0.849470   2.288134  -0.371 0.710451    
## thalach      0.016745   0.008158   2.053 0.040109 *  
## exang1      -1.071410   0.312255  -3.431 0.000601 ***
## oldpeak     -0.491236   0.171170  -2.870 0.004106 ** 
## slope1      -0.860702   0.556722  -1.546 0.122101    
## slope2       0.892960   0.589836   1.514 0.130048    
## ca1         -2.370383   0.348406  -6.804 1.02e-11 ***
## ca2         -3.680048   0.531857  -6.919 4.54e-12 ***
## ca3         -2.437445   0.737874  -3.303 0.000955 ***
## thal2       -0.478688   0.580266  -0.825 0.409402    
## thal3       -1.900330   0.570886  -3.329 0.000872 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 970.38  on 699  degrees of freedom
## Residual deviance: 395.26  on 679  degrees of freedom
## AIC: 437.26
## 
## Number of Fisher Scoring iterations: 6

---

Selección de variables (stepwise)

Usamos step() con dirección bidireccional para quedarnos solo con las variables que mejoran el modelo (criterio AIC).

modelo_step <- step(modelo_completo, direction = "both", trace = 0)
summary(modelo_step)

## 
## Call:
## glm(formula = target ~ sex + cp + trestbps + thalach + exang + 
##     oldpeak + slope + ca + thal, family = binomial, data = train)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  4.150988   1.616938   2.567 0.010253 *  
## sex1        -1.971277   0.369515  -5.335 9.57e-08 ***
## cp1          0.742618   0.374395   1.984 0.047310 *  
## cp2          1.991191   0.358751   5.550 2.85e-08 ***
## cp3          2.852750   0.493898   5.776 7.65e-09 ***
## trestbps    -0.020774   0.007505  -2.768 0.005637 ** 
## thalach      0.014098   0.007459   1.890 0.058757 .  
## exang1      -1.049898   0.308871  -3.399 0.000676 ***
## oldpeak     -0.538231   0.167603  -3.211 0.001321 ** 
## slope1      -0.903916   0.546633  -1.654 0.098207 .  
## slope2       0.852821   0.580680   1.469 0.141925    
## ca1         -2.332281   0.337419  -6.912 4.77e-12 ***
## ca2         -3.534259   0.507011  -6.971 3.15e-12 ***
## ca3         -2.388648   0.699304  -3.416 0.000636 ***
## thal2       -0.608255   0.566956  -1.073 0.283341    
## thal3       -2.050888   0.556518  -3.685 0.000229 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 970.38  on 699  degrees of freedom
## Residual deviance: 398.31  on 684  degrees of freedom
## AIC: 430.31
## 
## Number of Fisher Scoring iterations: 6

---

Evaluación del modelo

Matriz de confusión

pred_prob  <- predict(modelo_step, newdata = test, type = "response")
pred_clase <- factor(ifelse(pred_prob >= 0.5, 1, 0), levels = c(0, 1))

confusionMatrix(pred_clase, test$target, positive = "1")

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 118  20
##          1  26 136
##                                           
##                Accuracy : 0.8467          
##                  95% CI : (0.8008, 0.8855)
##     No Information Rate : 0.52            
##     P-Value [Acc > NIR] : <2e-16          
##                                           
##                   Kappa : 0.6923          
##                                           
##  Mcnemar's Test P-Value : 0.461           
##                                           
##             Sensitivity : 0.8718          
##             Specificity : 0.8194          
##          Pos Pred Value : 0.8395          
##          Neg Pred Value : 0.8551          
##              Prevalence : 0.5200          
##          Detection Rate : 0.4533          
##    Detection Prevalence : 0.5400          
##       Balanced Accuracy : 0.8456          
##                                           
##        'Positive' Class : 1               
## 

---

Probar el modelo con casos nuevos

Creamos dos perfiles de pacientes para predecir su probabilidad de enfermedad cardíaca:

• Paciente A: Hombre de 62 años, dolor típico, alta frecuencia cardíaca, sin angina.
• Paciente B: Mujer de 45 años, sin dolor, baja frecuencia cardíaca, con angina.

nuevos <- data.frame(
  age      = c(62, 45),
  sex      = factor(c(1, 0),   levels = levels(df$sex)),
  cp       = factor(c(3, 0),   levels = levels(df$cp)),
  trestbps = c(130, 150),
  chol     = c(254, 280),
  fbs      = factor(c(0, 1),   levels = levels(df$fbs)),
  restecg  = factor(c(1, 0),   levels = levels(df$restecg)),
  thalach  = c(170, 120),
  exang    = factor(c(0, 1),   levels = levels(df$exang)),
  oldpeak  = c(0.5, 2.5),
  slope    = factor(c(2, 0),   levels = levels(df$slope)),
  ca       = factor(c(0, 2),   levels = levels(df$ca)),
  thal     = factor(c(2, 3),   levels = levels(df$thal))
)

probabilidad <- predict(modelo_step, newdata = nuevos, type = "response")

resultado <- cbind(
  Paciente = c("A (hombre, 62)", "B (mujer, 45)"),
  nuevos[, c("age", "sex", "cp", "thalach", "exang")],
  Prob_Enfermedad = round(probabilidad, 3)
)
print(resultado)

##         Paciente age sex cp thalach exang Prob_Enfermedad
## 1 A (hombre, 62)  62   1  3     170     0           0.991
## 2  B (mujer, 45)  45   0  0     120     1           0.005

---

Conclusiones

El modelo final (tras selección stepwise) retuvo 9 de 13 variables — descartó age, chol, fbs y restecg por no mejorar el AIC. Sobre el conjunto de prueba (30% de los datos, n = 300):

• Accuracy: 84.7% — el modelo clasifica correctamente 8 de cada 10 pacientes.
• Sensibilidad: 87.2% — detecta correctamente el 87% de los pacientes con enfermedad (falsos negativos bajos).
• Especificidad: 81.9% — identifica correctamente el 82% de los pacientes sin enfermedad.

Variables con mayor impacto (coeficiente más alto en valor absoluto):

• ca (vasos coloreados): El predictor más fuerte. Tener 1, 2 o 3 vasos afectados reduce la log-odds en −2.3 a −3.5 puntos respecto a 0 vasos.
• cp (tipo de dolor): El tipo 3 suma +2.85 a la log-odds — el más asociado a diagnóstico positivo.
• thal3 (defecto reversible): Coeficiente de −2.05, fuerte predictor negativo respecto a la referencia.
• exang1 (angina por ejercicio): Reduce la log-odds en −1.05 al estar presente.

En los casos de prueba, el Paciente A (hombre, 62 años, sin angina, 0 vasos) obtuvo una probabilidad de enfermedad del 99.1%, mientras que el Paciente B (mujer, 45 años, con angina, 2 vasos) obtuvo apenas un 0.5%.