Regresión Logística - Heart

Teoría

La regresión logística es un modelo estadístico de clasificación binaria, que estima la probabilidad de que ocurra un evento (valor 1) frente a que no ocurra (valor 0), en función de variables independientes. # Instalar paquetes y llamar librerías

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

## Cargando paquete requerido: ggplot2

## Cargando paquete requerido: lattice

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

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.2.0     ✔ readr     2.1.6
## ✔ forcats   1.0.1     ✔ stringr   1.6.0
## ✔ lubridate 1.9.5     ✔ tibble    3.3.1
## ✔ purrr     1.2.1     ✔ tidyr     1.3.2

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ✖ purrr::lift()   masks caret::lift()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library (dplyr)

Crear la base de datos

#file.choose()
df <- read.csv("C:\\Users\\lucia\\Downloads\\concentración lit\\heart.csv")

Entender la base de datos

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 ...

df <- df %>% select(-slope, -fbs)
df <- na.omit(df)
df$cp <- as.factor(df$cp)
df$sex <- as.factor(df$sex)
df$restecg <- as.factor(df$restecg)
df$exang <- as.factor(df$exang)
df$ca <- as.factor(df$ca)
df$thal <- as.factor(df$thal)
df$target <- as.factor(df$target)

Crear el modelo

modelo <- glm(target ~ ., data=df, family=binomial)
summary(modelo)

## 
## Call:
## glm(formula = target ~ ., family = binomial, data = df)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.725058   1.717509  -0.422 0.672910    
## age          0.020809   0.013851   1.502 0.132996    
## sex1        -1.684486   0.296195  -5.687 1.29e-08 ***
## cp1          0.983753   0.301704   3.261 0.001112 ** 
## cp2          1.902427   0.270775   7.026 2.13e-12 ***
## cp3          2.102385   0.368328   5.708 1.14e-08 ***
## trestbps    -0.021882   0.006284  -3.482 0.000498 ***
## chol        -0.005749   0.002230  -2.578 0.009947 ** 
## restecg1     0.478813   0.213759   2.240 0.025093 *  
## restecg2    -0.702944   1.599559  -0.439 0.660327    
## thalach      0.028230   0.006228   4.532 5.83e-06 ***
## exang1      -0.790824   0.244171  -3.239 0.001200 ** 
## oldpeak     -0.617427   0.121957  -5.063 4.13e-07 ***
## ca1         -2.009666   0.264895  -7.587 3.28e-14 ***
## ca2         -3.002481   0.397583  -7.552 4.29e-14 ***
## ca3         -1.962875   0.483527  -4.059 4.92e-05 ***
## ca4          0.971714   0.935472   1.039 0.298925    
## thal1        2.260062   1.129027   2.002 0.045309 *  
## thal2        2.112730   1.068003   1.978 0.047905 *  
## thal3        0.580400   1.069570   0.543 0.587372    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1420.24  on 1024  degrees of freedom
## Residual deviance:  635.44  on 1005  degrees of freedom
## AIC: 675.44
## 
## Number of Fisher Scoring iterations: 6

Probar el modelo

# Probar el modelo (incluyendo age)

prueba <- data.frame(cp=as.factor(c(1,1)),
  age     = 52,
  sex     = factor(1, levels = levels(df$sex)),
  cp      = factor(0, levels = levels(df$cp)),
  trestbps= 125,
  chol    = 212,
  restecg = factor(1, levels = levels(df$restecg)),
  thalach = 168,
  exang   = factor(0, levels = levels(df$exang)),
  oldpeak = 1,
  ca      = factor(2, levels = levels(df$ca)),
  thal    = factor(3, levels = levels(df$thal))
)

probabilidad <- predict(modelo, newdata = prueba, type = "response")
clase <- ifelse(probabilidad >= 0.5, 1, 0)

cbind(prueba, Probabilidad_Target1 = probabilidad, Prediccion = clase)

##   cp age sex cp.1 trestbps chol restecg thalach exang oldpeak ca thal
## 1  1  52   1    0      125  212       1     168     0       1  2    3
## 2  1  52   1    0      125  212       1     168     0       1  2    3
##   Probabilidad_Target1 Prediccion
## 1            0.1075632          0
## 2            0.1075632          0