Modelos estadisticos

Antes de realizar un modelo hacer un analisis descriptivo del raw data, para comprender de antemano posibles relaciones entre los datos. Por ejemplo, si quiero hacer un modelo de regresión lineal, y al hacer un graf. de dispersión me doy cuenta de que la relación entre los datos no se ve para nada lineal, tengo que hacer ajustes previos o buscar alternativas.

“If you assess hypotheses without examining your data, you risk publishing nonsense” (Wilkinson & Task Force on Statistical Inference, 1999, p. 597).

Paquetes a utilizar:

##### Paquetes

options(repos = c(CRAN = "https://cran.r-project.org"))

#Instalación de paquetes
install.packages("tidyverse") #Para acceder a dplyr, ggplot, etc.
library(tidyverse)

install.packages("texreg") # Permite usar screenreg
library(texreg)  

install.packages("car") # Funciones para análisis de regresión
library(car)  

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

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

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

install.packages("klaR")
library(klaR)

install.packages("rpart")
library(rpart)

install.packages("rpart.plot")
library(rpart.plot)

install.packages("smotefamily") #Sirve para usar la correcion SMOTE
library(smotefamily)

install.packages("pROC") # LIBRERIA PARA CURVAS ROC
library(pROC)

options(warn = -1)
suppressWarnings(install.packages("tidyverse"))
library(tidyverse)
options(warn = 0)

Regresión Lineal

Regresión lineal

La idea de la regresión lineal es explicar el comportamiento de una variable dependiente en función de una o más variables independientes. Esto se hace a través de una recta que pretende generar el menor promedio de residuos posible.

En la RL la variable dependiente siempre es cuantitativa. La variable independiente puede ser cuanti o cuali

  • Independiente cuantitativa: Por cada unidad de aumento en Xn, la variable dependiente aumenta en Bn (unidad de medida).

    Ej: Ingreso = 20.000 + 12000(X) -> El sueldo base es 20.000. Por cada unidad de aumento en X, el ingreso aumenta en $12.000

  • Independiente cuantitativa: Utiliza variables dummie y una categoria de referencia (R la asume automaticamente como la menor, si está en 0 y 1, 0 será la de referencia).

    El coeficiente asociado a la variable cualitativa representa la diferencia promedio en la variable dependiente entre la categoría de referencia y la otra categoría

    Ej. 1: Modelamos el sueldo (variable dependiente) en función del sexo (0:hombre, 1: mujer) La categegoria de referencia será hombre, la dummy mujer (La interpretación se hace en torno a la dummy)

    Sueldo = 20.000 - 2000(D1)

    Para hombres (D1 = 0):Sueldo = 20.000 + 2.000(0) = 20.000

    Para mujeres (D1 = 1):Sueldo = 20.000 - 2.000(1) = 18.000

    Una hombre ganará 20.000 de base. Ser mujer está asociado con un ingreso promedio 2.000 menor al de los hombres.

    Siempre que conviertas una variable categórica en dummies, la categoría de referencia toma el valor 0 en todas las dummies

    Ej. 2: Creemos el mismo modelo pero añadamos el sexo “otro” (0:hombre, 1:mujer, 2:otro)

    Sueldo = 20.000 -2.000(D1) - 5.000(D2)

    Para hombres (D1 = 0, D2=0):Sueldo = 20.000 -2.000(0) - 5.000(0) = 20.000

    Para mujeres (D1 = 1, D2=0):Sueldo = 20.000 -2.000(1) - 5.000(0) = 18.000

    Para otro (D1 = 0, D2=1):Sueldo = 20.000 -2.000(0) - 5.000(1) = 15.000

Interpretación modelos:

Los modelos de regresión lineal tendrán:

  • Un Intercepto y la pendiente para cada variable indepeniente y su p-value. En caso de que haya alguna variable sin significancia, si no afecta el R2 de mi modelo es mejor sacarla, para tener un modelo mas parcemonioso (mas simple).

    Con el intercepto y las pendientes puedo construir la formula y = b0 +b1x para hacer predicciones.

  • R-squared: Nos dice el % de variabilidad de la variable dependiente que se explica en base al modelo. Es ideal que sea sobre 0.8 (La variabilidad de Y se explica en un 80% por este modelo)

  • p-value de anova: No tiene que ver con el test de hipotesis anova. H1: al menos una de las variables independientes tiene un efecto significativo sobre la variable dependiente.

En el siguiente ejemplo crearé un mod. de regresión lineal para explicar el peso a partir de la edad, altura y diabetes (0= no, 1= si)

Recordar que primero debemos revisar y limpiar nuestros datos. Hacer siempre summary para toda la base y table(base$variable) para vairables especificas.

bmi <- readRDS(url("https://github.com/FranvTapia/bases/raw/main/nhanes.rds"))

#modelo
mod_peso <- lm(peso ~ edad + altura + diabetes, data = bmi)
summary(mod_peso)
## 
## Call:
## lm(formula = peso ~ edad + altura + diabetes, data = bmi)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -45.595 -12.139  -2.665   9.037 139.244 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -85.606370   1.509867  -56.70   <2e-16 ***
## edad          0.144367   0.009697   14.89   <2e-16 ***
## altura        0.939364   0.010070   93.28   <2e-16 ***
## diabetesYes  11.272360   0.702579   16.04   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 17.67 on 9625 degrees of freedom
##   (371 observations deleted due to missingness)
## Multiple R-squared:  0.5895, Adjusted R-squared:  0.5894 
## F-statistic:  4607 on 3 and 9625 DF,  p-value: < 2.2e-16

El modelo indica que, manteniendo constantes la edad y la altura, las personas con diabetes pesan en promedio 11,27 kg más que quienes no tienen diabetes; además, por cada centímetro adicional de altura el peso aumenta en 0,94 kg, y por cada año adicional de edad el peso aumenta en 0,14 kg. Al menos una de las variables independientes tiene un efecto significativo sobre la variable dependiente. La variabilidad del peso está explicada en un 58.94% por el modelo.

Para ver los intervalos de confianza confint(modelo)

⚠️ Creo que este modelo tiene problemas de heterocedasticidad y linealidad que no fueron corregidos. Pueden venir dados porque para ciertas alturas y rangos etarios el peso comienza a mantenerse más constante. A pesar de que el R2 no es tan bajo puede que se ajuste mejor a otro tipo de regresiones. Se utilizó como ejemplo unicamente por su simplicidad.

Supuestos de regresión

Si vas a hacer modelos lineales entonces tienen que cumplir con estos supuestos. De no hacerlo se deben tomar medidas correctivas.

1. LINEALIDAD

La relación entre la variable dependiente e independiente sea lineal. Cuando este supuesto se transgrede los coeficientes estimados pueden ser engañosos.

Si no se cumple el supuesto de linealidad la corrección es incluir términos cuadráticos para curvar la recta.

Veamoslo mediante un ejemplo. Con set.seed cree un conjunto de datos que NO sigue una tendencia lineal. Si hacemos un grafico de dispersión con geom_point y geom_smooth podemos ver que la linea exponencia podria ajustarse mejor al modelo.

#EJEMPLO RELACION NO LINEAL

# Simular datos con una relación no lineal
set.seed(123)
edad <- seq(10, 90, length.out = 100)
presion_sistolica <- 90 + 0.5 * edad + 0.03 * edad^2 + rnorm(100, 0, 10)  # relación no lineal
datos <- data.frame(edad, presion_sistolica)

# Gráfico comparativo
ggplot(datos, aes(x = edad, y = presion_sistolica)) +
  geom_point(color = "gray40") +
  geom_smooth(method = "lm", formula = y ~ x, se = FALSE, color = "#87CEFA", size = 1.2) +
  geom_smooth(method = "lm", formula = y ~ poly(x, 2), se = FALSE, color = "#EEA2AD", size = 1.2) 

#Al visualizar el raw data de antemano podemos notar que la relación se ajusta mejor a una curva exponencial que a una recta

Podemos verificar de diferentes formas. Primero crearé ambos modelos y los compararé:

# Ajustes de modelos
modelo_lineal <- lm(presion_sistolica ~ edad, data = datos)
modelo_cuadratico <- lm(presion_sistolica ~ edad + I(edad^2), data = datos) #término cuadrático (edad^2)
#ajusté un modelo lineal y otro con el término cuadratico

screenreg(list("Modelo lineal" = modelo_lineal, "Modelo cuadrático" = modelo_cuadratico))
## 
## =============================================
##              Modelo lineal  Modelo cuadrático
## ---------------------------------------------
## (Intercept)   30.67 ***      93.13 ***       
##               (4.29)         (4.26)          
## edad           3.53 ***       0.34           
##               (0.08)         (0.19)          
## edad^2                        0.03 ***       
##                              (0.00)          
## ---------------------------------------------
## R^2            0.95           0.99           
## Adj. R^2       0.95           0.99           
## Num. obs.    100            100              
## =============================================
## *** p < 0.001; ** p < 0.01; * p < 0.05

En el modelo lineal, el coeficiente de x es significativo y parece describir una relación fuerte. Pero al ajustar un modelo cuadrático, vemos que este efecto desaparece, y el término cuadrático toma protagonismo.

No siempre que el coeficiente lineal pierda significancia al incluir el término cuadrático significa que la relación no es lineal. Pero sí es una fuerte señal de que el modelo lineal puede estar mal especificado, a pesar de tener un R² alto.

Ahora revisemos los residuos

residualPlot(modelo_lineal)

residualPlot(modelo_cuadratico)

En los residuos no deberiamos ver ningun patron sistematico

Cuando vemos un patron sistematico (como la curva en el lineal) es un claro indicador de que la relacion no es lineal. En la cuadratica no hay ningun patron, los puntos están dispersos. Eso es lo que deberiamos ver en un buen grafico de residuos.

Por último miremos el error cuadratico medio (MSE) El error cuadrático medio (mean squared error, MSE) es una medida de la calidad del ajuste de un modelo. Se calcula como el promedio de los cuadrados de los residuos. El mejor modelo tendrá el menor MSE.

El primero es el lineal y el segundo el cuadratico:

mse <- function(model) mean(residuals(model)^2) #funcion para calcular el error cuadratico medio
mse(modelo_lineal)  # mayor MSE
## [1] 322.4794
mse(modelo_cuadratico)    # menor MSE
## [1] 81.08254

En este caso, el mejor modelo es el cuadratico.

Todas las pruebas que relizamos (graf dispersión,comparar r2, residualplot, mse) comprueban que nuestros datos no se ajustan a un modelo lineal sino a uno exponencial.

2. HOMOCEDASTICIDAD

Nosotros queremos que nuestro modelo tenga homocedasticidad, no heterocedasticidad.

Cuando hay heterocedasticidad los errores estandar son mas pequeños de lo que deberian, por lo que los intervalos de confianza se agostan, lo que podría resultar en que encuentre significancia estadistcica donde no la hay.

Aumenta la probabilidad de cometer un error de tipo 1: Rechazar H0 cuando es cierta => que mi intervalo NO incluya el 0 cuando SI deberia incluirlo (Recordemos que si mi intervalo incluye el 0 yo debería rechazar H0)

Para ver si hay heterocedasticidad hay que hacer un gráfico de residuos. Si los puntos se van abriendo, hay un problema de heterocedasticidad.

En el siguiente ejemplo se genera un conjunto de datos con heterocedasticidad.Así se ve su gráfico de residuos:

set.seed(123)
n <- 100
x <- runif(n, 0, 100)
e <- rnorm(n, mean = 0, sd = x / 5)  # varianza del error aumenta con x
y <- 2 + 0.5 * x + e

mod_hetero <- lm(y ~ x) #modelo con heterocedasticidad

df_hetero <- data.frame(
  x = x,
  y = y,
  fitted = fitted(mod_hetero),
  resid = resid(mod_hetero)
)

ggplot(df_hetero, aes(x = fitted, y = resid)) +
  geom_point(color = "gray50", size = 2) +
  geom_hline(yintercept = 0, color = "red", linetype = "dashed") +
  labs(title = "Patrón de heterocedasticidad en los residuos",
       x = "Valores ajustados", y = "Residuos") +
  theme_minimal(base_size = 14)  #plot de los residuos

Para corregirlo hay que implementar error estándar robusto, para ensanchar los intervalos de confianza:

coeftest(mod_hetero, vcov. = vcovHC(mod_hetero, type = "HC1"))
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.713120   1.509253  1.1351   0.2591    
## x           0.492107   0.044031 11.1763   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Hacer esta corrección hace una gran diferencia, sobre todo si el intervalo está cerca del 0. Tenemos que hacer la correcion aunque con hetero este muy lejos del 0? En estricto rigor, no. Pero suena más bacan, por lo que se recomienda hacerla de igual manera siempre que haya heterocedasticidad.

3. MULTICOLINEALIDAD

La multicolinealidad se da cuando dos o más variables están explicando lo mismo. Se comprueba con VIF y se corrige eliminando directamente las variables multicolineales del modelo.

Para el sig. ejemplo, en el modelo se incluyeron cuatro variables predictoras: dos altamente correlacionadas entre sí (x1 y x2) y dos no correlacionadas (x3 y x4).

set.seed(123)
n <- 100
# Variables correlacionadas
x1 <- rnorm(n)
x2 <- x1 + rnorm(n, sd = 0.05)  # alta colinealidad con x1
# Variables no correlacionadas
x3 <- rnorm(n)
x4 <- rnorm(n)
# Variable dependiente
y <- 3 + 1.5 * x1 + 2 * x2 + 0.5 * x3 + 0.7 * x4 + rnorm(n, sd = 1)


mod_full <- lm(y ~ x1 + x2 + x3 + x4)

# Diagnóstico de colinealidad
vif(mod_full)
##         x1         x2         x3         x4 
## 357.538469 357.318411   1.020348   1.006429

VALOR DE VIF ES ACEPTABLE < 10

El diagnóstico de colinealidad mediante el VIF muestra que x1 y x2 presentan valores muy elevados, lo cual indica que la varianza de sus coeficientes está inflada artificialmente debido a la redundancia entre ellos. Esto genera estimaciones inestables,con errores estándar grandes y posible pérdida de significancia.

La solución es sacar una variable del modelo. Esto queda a elección propia, debemos elegir la que más nos haga sentido.

mod_corregido <- lm(y ~ x1 + x3 + x4) #sacamos x2
vif(mod_corregido)
##       x1       x3       x4 
## 1.019554 1.019625 1.004563

Ya no hay problemas de multicolinealidad.

** Existen otros supuestos como la normalidad de los residuos o la independencia de los errores. Que se usan CUANDO NO ME ACUERDO BUSCAR

Train/Test/C.Matrix

Siempre al crear modelos para predicción hay que separarlos datos entre train(~70%), para entrenar el modelo, y test (~30%), para comparar datos reales con las predicciones del modelo, con una matriz de confusión.

index <- createDataPartition(datos$diagnosis, p = 0.7, list = F)

train <- datos[index,] #train queda con el 70%

test <- datos[-index,] #test queda con el 30%

Si no separamos podemos tener problemas de overfit (Que el modelo se ajuste demasiado a los datos de entrenamiento y no sea eficiente al clasificar datos externos)

  1. Separaremos los datos en train y test.
  2. Creamos el modelo usando train. En el method va el tipo de modelo. Ej:“nb” para naive bayes.
  3. Realizamos la predicción en la base de testeo
  4. Comparamos con una matriz de confusión. Podemos especificar cual es el positivo. En este caso queremos que los positivos sean los tumores malignos (positive= “M”)

🔵 Usamos la matriz de confusión para ver el rendimiento de nuestro modelo, pero esta no es la única métrica. También es importante la curva ROC/AUC, sobre todo cuando hay desbalance de clases y el acuraccy no me entrega mucha información.

El output incluirá:

  • Matriz de confusión: Compara predicciones con los datos reales de test.

  • Accuracy y su CI: Proporción de aciertos del modelo.

  • No Information Rate: Es la precisión que se obtendría si siempre se predijera la clase mayoritaria, sin usar un modelo. Es, en la práctica, simplemente la proporción de la clase mayoritaria en tus datos. Por ejemplo, si la prev de una enfermedad donde evaluamos benigno y maligno es un 60%, el NIR sería 0.6, porque esta seria la clase mayoritaria. Con dos categorias el varlor min del NIR es 0.5 (Si la prev < 50% entonces NIR= prev - 1), con tres sería 0.33, porque siempre considera la clase mayoritaria.

  • P-Value [Acc > NIR]: H1:Acc > NIR La idea es que la accuracy de nuestro modelo sea mayor al NIR, que viene dado simplemente por lo random. Nuestro modelo debería ser mejor que lo random. Si nir > accuracy mi modelo no está aportando ninguna valor real, predice igual o peor que simplemente elegir la clase mayoritaria

    Hay que tener ojo con el NIR cuando tengo un fuerte desbalance de clases, es decir, que una clase aparece mucho más que la otra (ej: 95% de los casos son benignos, 5% son malignos), ya que en estos casos, nuestro NIR será alto.

🔴 Ante desbalance de clases se recomienda usar SMOTE, una técnica que genera nuevos ejemplos sintéticos de la clase minoritaria para balancear el conjunto de entrenamiento.

  • Kappa: Es una métrica que evalúa el grado de acuerdo entre el modelo y la realidad, ajustando por el acuerdo que podría darse por azar (calcula el azar y lo compara con la accuracy). A diferencia del accuracy, que solo mide el porcentaje total de aciertos, Kappa calcula cuánto del acuerdo es realmente significativo, considerando la distribución de clases.

    Mientras más cercano a 1, mejor.

  • Mcnemar’s Test P-Value: H0: Variables NO discrepan/ H1: Variables discrepan. Nos interesa que no discrepen, NO rechazar H0 (pvalue>0.05)

  • Sensibilidad, especificidad, vpp, vpn

  • Prevalencia : (Enfermos/Total)

  • Detection rate: (VP/total)

  • Detection Prevalence: (+/total)

  • Balanced Accuracy: Promedio entre sensibilidad y especificidad

🔴 Modelos como Random Forest, SVM, KNN, Redes Neuronales tienen HIPERPARAMETROS. Son configuraciones externas al modelo que tú defines al entrenarlo y que afectan cómo aprende.

Para encontrar los mejores hiperparametros puedes llevar a cabo un Grid Search (con tuneLength y tuneGrid en Caret). Si no los especificas R usa una combinación de hiperparámetros por defecto. Para modelos mejor optimizados se recomienda hacer Grid Search.

En los ejemplos encontraremos accuracy altas porque son bases diseñadas para hacer modelos. No son realistas. Si en la realidad encontramos accuracy sobre 95% es un poco sospechosos

Naive Bayes

Naive Bayes es un método de clasificación para variables dependientes categoricas. Es naive porque asume que las variables predictoras son independientes entre sí, por esto se recomienda hacer test de correlación antes de realizar los modelos, y eliminar variables altamente correlacionadas.

Su clasificación se basa en el teorema de bayes:

Dado un conjunto de clases C1, C2, …, Ck, el objetivo de Naive Bayes es predecir a qué clase pertenece un nuevo dato X. Para ello, sigue los siguientes pasos:

  1. Calcula la probabilidad a priori de cada clase, es decir, la proporción de veces que aparece cada clase en los datos de entrenamiento.

  2. Calcula la probabilidad condicional de cada atributo del dato dado cada clase, asumiendo independencia entre atributos.

  3. Aplica el teorema de Bayes para obtener la probabilidad de que el dato pertenezca a cada clase posible.

  4. Compara las probabilidades obtenidas y asigna al dato la clase que tenga la probabilidad más alta.

De esta manera, Naive Bayes permite clasificar nuevos datos eligiendo la clase más probable según la evidencia observada.

Modelo

Creemos un modelo de predicción de cancer de mama.

Primero limpiamos la base. Naive Bayes asume independencia entre las variables, por lo que para trabajar con este modelo es eliminar de ante mano las variables correlacionadas, esto lo podemos ver con una matriz de correlación. En este caso, deje que R eligiera las variables a eliminar, pero en un proyecto más serio esto se debería hacer mas cautelosamente. Estas son las variables que eliminé:

datos <- read.csv(url("https://github.com/FranvTapia/bases/raw/main/cancer.csv"))
#limpieza
datos <- datos[, !names(datos) %in% c( "X", "id")] #x no tenia nada, id no lo usaremos
names(datos) <- gsub(" ", "_", names(datos)) #hay algunos que no tienen _

#eliminar datos correlacionados
# Calcular la matriz de correlación
cor_matrix <- cor(datos[, !names(datos) %in% "diagnosis"])
#tambien se puede hacer asi pero markdown no me deja cor_matrix <- cor(select(datos, -diagnosis))

# Convertimos a matriz superior para evitar duplicados
high_cor_pairs <- which(abs(cor_matrix) > 0.9 & abs(cor_matrix) < 1, arr.ind = TRUE)
high_cor_pairs <- as.data.frame(high_cor_pairs)

# Añadimos nombres de variables
high_cor_pairs$var1 <- rownames(cor_matrix)[high_cor_pairs$row]
high_cor_pairs$var2 <- colnames(cor_matrix)[high_cor_pairs$col]

# Eliminamos duplicados (e.g., A vs B y B vs A)
high_cor_pairs <- high_cor_pairs[high_cor_pairs$row < high_cor_pairs$col, ]

#eliminar todas las var 2
vars_to_remove <- unique(high_cor_pairs$var2)

vars_to_remove
##  [1] "perimeter_mean"       "area_mean"            "concave.points_mean" 
##  [4] "perimeter_se"         "area_se"              "radius_worst"        
##  [7] "texture_worst"        "perimeter_worst"      "area_worst"          
## [10] "concave.points_worst"
# Creamos nueva base sin esas variables
datos <- datos[, !names(datos) %in% vars_to_remove]

Ahora trabajaré en el modelo, partiré por la partición de datos, y terminaré con la matriz de confusión.

#necesitamos caret,
#separar datos
index <- createDataPartition(datos$diagnosis, p = 0.7, list = F)
train <- datos[index,] #train queda con el 70%
test <- datos[-index,] #test queda con el 30%

#crear modelo
modelo_nb <- train(diagnosis ~ ., data = train, method = "nb")

# PREDECIR EN LA BASE DE TESTEO
p_nb <- predict(modelo_nb, test) #el predic siempre se hace en la base de testeo
# MATRIZ DE CONFUSIÓN
test$diagnosis <- factor(test$diagnosis, levels = levels(p_nb)) #aqui tuve que decirle que  factorizara
#aquí puedo especificar cual va a ser mi positivo. Puedo dejar los malignos coomo positivos
confusionMatrix(p_nb, test$diagnosis, positive= "M")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   B   M
##          B 101  12
##          M   6  51
##                                          
##                Accuracy : 0.8941         
##                  95% CI : (0.8378, 0.936)
##     No Information Rate : 0.6294         
##     P-Value [Acc > NIR] : 5.13e-15       
##                                          
##                   Kappa : 0.7685         
##                                          
##  Mcnemar's Test P-Value : 0.2386         
##                                          
##             Sensitivity : 0.8095         
##             Specificity : 0.9439         
##          Pos Pred Value : 0.8947         
##          Neg Pred Value : 0.8938         
##              Prevalence : 0.3706         
##          Detection Rate : 0.3000         
##    Detection Prevalence : 0.3353         
##       Balanced Accuracy : 0.8767         
##                                          
##        'Positive' Class : M              
##

Nuestro modelo se equivocó 18 veces, tiene una proporción de acierto del 89,4%. P-Value [Acc > NIR], Kappa, Mcnemar’s Test son favorables.

El 80.95% de los casos positivos reales fueron detectados (VP). El 94.39% de los casos negativos reales fueron correctamente clasificados (VN). De todas las predicciones positivas, el 89.47% realmente eran positivas. De todas las predicciones negativas, el 89.38% realmente eran negativas.

Arboles de decisión

Utiliza las probabilidades en varias etapas, utilizando el teorema de bayes Consiste en la representación de las ramificaciones de decisiones de manera fácil de la trayectoria histórica para cada una de las alternativas posibles.

Los arboles de decisión pueden ser categoricos o de regresión (dependiente cuanti). R elige las particiones (para las cuanti) y las variables que usa para el modelo.

Es por esto que Los árboles grandes tienen poco sentido intuitivo y las predicciones tienen cierto aire de “cajas negras”, no sabemos exactamente por qué el arbol hizo las particiones donde las hizo, ni por qué escogió las variables que escogió.

Sigamos con el ejemplo del cáncer de mama.

Primero veamos las divisiones y variables que conforman el árbol usando rpart.plot. Recordar que estas divisiones y variables las escoge R.

#los datos ya están limpios
#separar datos
index <- createDataPartition(datos$diagnosis, p = 0.7, list = F)
train <- datos[index,] #train queda con el 70%
test <- datos[-index,] #test queda con el 30%

#crear modelo
modelo_dt <- rpart(diagnosis ~ ., data = train, method = "class") 

#dibujar modelo

rpart.plot(modelo_dt)

Los porcentajes represental total de los datos que cae en ese nodo. Si un nodo tiene 31%, significa que el 31% de todos los casos del dataset caen en ese nodo.

Las proporciones son la proporción de los casos en ese nodo que pertenecen a la clase predicha. Si un nodo muestra M, 0.95, eso significa que el 95% de los casos que llegaron a ese nodo son malignos, por eso el árbol predice “M” ahí.

El nodo de los que tienen area_worst < 885 abarca un 68% de los datos, mientras que el de los que no, un 32%. De este 32%, un 95% de los casos son malignos, lo que indica una alta probabilidad de malignidad cuando area_worst es elevado. Por otro lado, del 68% restante (con area_worst < 885), se genera una nueva división según concavity_mean: si es menor a 0.1, se agrupa el 12% de los datos, de los cuales el 69% son benignos; este grupo se subdivide nuevamente: si radius_mean < 15, se clasifica como benigno (60% benignos en ese nodo, con un 3% de representación total); si no, el modelo evalúa smoothness_worst: si es menor a 0.14, se clasifica como benigno (4% del total); si es mayor, se predice maligno (6% del total). En cambio, si concavity_mean ≥ 0.1, se clasifica directamente como maligno, con una proporción del 63% y una representación del 10% del total.

Así, el árbol revela rutas claras asociadas al riesgo de malignidad, destacando especialmente a los casos con alto area_worst o altos valores combinados de concavity_mean y smoothness_worst.

Ahora veamos la matriz de confusón

# PREDECIR EN LA BASE DE TESTEO
p_dt <- predict(modelo_dt, test, type = "class") #hay que poner class en los dt
# MATRIZ DE CONFUSIÓN
test$diagnosis <- factor(test$diagnosis, levels = levels(p_dt)) #aqui tuve que decirle que  factorizara
#aquí puedo especificar cual va a ser mi positivo. Puedo dejar los malignos coomo positivos
confusionMatrix(p_dt, test$diagnosis, positive= "M") #señalar la clase positiva
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  B  M
##          B 98  9
##          M  9 54
##                                          
##                Accuracy : 0.8941         
##                  95% CI : (0.8378, 0.936)
##     No Information Rate : 0.6294         
##     P-Value [Acc > NIR] : 5.13e-15       
##                                          
##                   Kappa : 0.773          
##                                          
##  Mcnemar's Test P-Value : 1              
##                                          
##             Sensitivity : 0.8571         
##             Specificity : 0.9159         
##          Pos Pred Value : 0.8571         
##          Neg Pred Value : 0.9159         
##              Prevalence : 0.3706         
##          Detection Rate : 0.3176         
##    Detection Prevalence : 0.3706         
##       Balanced Accuracy : 0.8865         
##                                          
##        'Positive' Class : M              
##

Random Forest

Random Forest combina varios árboles de decisión para generar una predicción más robusta y precisa. Los árboles que forman parte del “bosque” se construyen de manera aleatoria, tanto en los datos de entrenamiento como en las variables utilizadas en cada división. Por esto puede que los resultados varien en diferentes simulaciones.

#Recordar que los datos deben estar separados en train y test. No lo hago aca pk ya lo hice arriba
modelo_rf <- train(diagnosis~., train, method = "rf") #modelo
p_rf <- predict(modelo_rf, test) #predicción
test$diagnosis <- factor(test$diagnosis, levels = levels(p_rf)) #aqui tuve que decirle que  factorizara
confusionMatrix(p_rf,test$diagnosis)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   B   M
##          B 104   7
##          M   3  56
##                                           
##                Accuracy : 0.9412          
##                  95% CI : (0.8945, 0.9714)
##     No Information Rate : 0.6294          
##     P-Value [Acc > NIR] : <2e-16          
##                                           
##                   Kappa : 0.8722          
##                                           
##  Mcnemar's Test P-Value : 0.3428          
##                                           
##             Sensitivity : 0.9720          
##             Specificity : 0.8889          
##          Pos Pred Value : 0.9369          
##          Neg Pred Value : 0.9492          
##              Prevalence : 0.6294          
##          Detection Rate : 0.6118          
##    Detection Prevalence : 0.6529          
##       Balanced Accuracy : 0.9304          
##                                           
##        'Positive' Class : B               
##

Support Vector Machine

Sirve para modelos de regresión y clasificación. Su objetivo principal es encontrar un hiperplano que separe de mejor manera las diferentes clases de datos. Los support vector son los primeros datos de cada categoria en tocar el margen del hiperplano.

SVM siempre busca una línea recta (o un plano, en más dimensiones) para separar los datos, pero a veces los datos están mezclados de una forma que no se pueden separar con una línea recta en el espacio original, como en en el ejemplo anterior. Es aquí donde entra el concepto de Kernel.

El kernel en SVM transforma los datos a otro espacio (a veces con muchas dimensiones) donde sí se puede trazar una línea recta que los separe bien. Esa línea recta en el nuevo espacio se puede ver como una curva en el espacio original. Es como si tomaras una hoja de papel con puntos mezclados y la doblaras o estiraras (con el kernel), para que luego sí puedas trazar una línea recta que los separe.

Existen diferentes tipos de Kernel. Los principales para SVM son:

  • Kernel Lineal -> C

  • Kernel Polinómico -> C y D

  • Kernel Gaussiano -> C y gamma

Tienen diferentes hiperparametros (señalados con las flechas) que pueden ser ajustados con Grid Search. Si no se ajustan R usará hiperparametros por defecto.

Kernel Lineal

Es el más simple. Utiliza la técnica de máximos de distancia de hiperplanos.
Se puede ajusta la “Holgura” (Slack Variable) con el hiperparametro C , modificando el tamaño del margen que nos permite clasificar.

  • C grande: el modelo quiere evitar errores a toda costa. Permite márgenes más estrechos. Puede llevar a overfitting, porque se ajusta demasiado a los datos de entrenamiento.

  • C pequeño: el modelo tolera más errores.Permite márgenes más amplios, aunque eso implique clasificar mal algunos puntos. Ayuda a generalizar mejor, reduciendo el riesgo de overfitting.

🔹 Todos los kernel se pueden ajustar por C

El siguiente ejemplo muestra un modelo de SVM lineal para cáncer de mama (No se ajustó con GridSearch, se dejó el C por defecto)

#Modelo lineal sin ajustar con Grid Search
#Recordar que los datos deben estar separados en train y test. No lo hago aca pk ya lo hice arriba

svm1 <- train(diagnosis ~., data = train, method = "svmLinear")

p_svm1 <- predict(svm1, test)
test$diagnosis <- factor(test$diagnosis, levels = levels(p_svm1))
confusionMatrix(p_svm1, test$diagnosis, positive= "M")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   B   M
##          B 105   4
##          M   2  59
##                                           
##                Accuracy : 0.9647          
##                  95% CI : (0.9248, 0.9869)
##     No Information Rate : 0.6294          
##     P-Value [Acc > NIR] : <2e-16          
##                                           
##                   Kappa : 0.9238          
##                                           
##  Mcnemar's Test P-Value : 0.6831          
##                                           
##             Sensitivity : 0.9365          
##             Specificity : 0.9813          
##          Pos Pred Value : 0.9672          
##          Neg Pred Value : 0.9633          
##              Prevalence : 0.3706          
##          Detection Rate : 0.3471          
##    Detection Prevalence : 0.3588          
##       Balanced Accuracy : 0.9589          
##                                           
##        'Positive' Class : M               
##
Kernel Polinómico

Al igual que el lineal, se puede ajustar por C pero también por D, correspondiente al degree del polinomio. D = 2 crea una transformación cuadrática, D = 3 cúbica, etc. (Las lineales tienen un D=1)

Modelo polinomial cáncer de mama:

#Modelo poli sin ajustar con Grid Search
#Recordar que los datos deben estar separados en train y test. No lo hago aca pk ya lo hice arriba
svm2 <- train(diagnosis ~ ., data = train, method = "svmPoly") #kernel polinomial

p_svm2 <- predict(svm2, test)
test$diagnosis <- factor(test$diagnosis, levels = levels(p_svm2))
confusionMatrix(p_svm2, test$diagnosis, positive= "M")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   B   M
##          B 106   6
##          M   1  57
##                                          
##                Accuracy : 0.9588         
##                  95% CI : (0.917, 0.9833)
##     No Information Rate : 0.6294         
##     P-Value [Acc > NIR] : <2e-16         
##                                          
##                   Kappa : 0.9103         
##                                          
##  Mcnemar's Test P-Value : 0.1306         
##                                          
##             Sensitivity : 0.9048         
##             Specificity : 0.9907         
##          Pos Pred Value : 0.9828         
##          Neg Pred Value : 0.9464         
##              Prevalence : 0.3706         
##          Detection Rate : 0.3353         
##    Detection Prevalence : 0.3412         
##       Balanced Accuracy : 0.9477         
##                                          
##        'Positive' Class : M              
##
Kernel Gaussiano

Usa Gamma, que determina cuánta influencia tiene cada punto de entrenamiento al construir la frontera de decisión.

  • Gamma muy pequeño: Cada punto tiene una influencia muy amplia, como si su efecto se “difundiera” mucho.El modelo genera una frontera muy suave.

  • Gamma muy grande: Cada punto influye solo en su zona muy cercana, y el modelo trata de ajustarse fuertemente a cada punto. Esto puede causar overfitting

A medida que aumenta gamma, la frontera de decisión se vuelve más precisa pero también más arriesgada, porque puede sobreajustarse. Para la imagen anterior, el mejor modelo sería el gamma = 0.01

En cuanto al C, cumple la misma función que en los otros kernel: controla qué tan duro es el modelo con los errores mediante los margenes, pero su efecto está influenciado por gamma, por lo que su impacto en el modelo depende de cómo ambos hiperparámetros interactúan.

Para el primero, la frontera es muy suave: el modelo tolera errores para mantener un margen amplio. Para el último: el modelo intenta clasificar todo correctamente, generando overfitting.

#Modelo gaussiano sin ajustar con Grid Search
#Recordar que los datos deben estar separados en train y test. No lo hago aca pk ya lo hice arriba

svm3 <- train(diagnosis ~ ., data = train, method = "svmRadial") #kernel gaussiano

#GAUSSIANO
p_svm3 <- predict(svm3, test)
test$diagnosis <- factor(test$diagnosis, levels = levels(p_svm3))
confusionMatrix(p_svm3, test$diagnosis, positive= "M")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   B   M
##          B 105   5
##          M   2  58
##                                          
##                Accuracy : 0.9588         
##                  95% CI : (0.917, 0.9833)
##     No Information Rate : 0.6294         
##     P-Value [Acc > NIR] : <2e-16         
##                                          
##                   Kappa : 0.9109         
##                                          
##  Mcnemar's Test P-Value : 0.4497         
##                                          
##             Sensitivity : 0.9206         
##             Specificity : 0.9813         
##          Pos Pred Value : 0.9667         
##          Neg Pred Value : 0.9545         
##              Prevalence : 0.3706         
##          Detection Rate : 0.3412         
##    Detection Prevalence : 0.3529         
##       Balanced Accuracy : 0.9510         
##                                          
##        'Positive' Class : M              
##

Regresión Logistica

Modelo de clasificación (Aprendizaje supervisado).Se utiliza cuando la variable dependiente es binaria.

ODDS

Antes de adentrarnos en el modelo, es necesario entender el concepto de odds: el ratio entre los casos favorables y los casos no favorables

Para comparar las odds entre dos grupos usamos calculamos el Odds Ratio (OR), que puede ir desde 0 a infinito pusitivo.

  • OR > 1 Factor de riesgo (Para interpretar veo cuanto le sobra a 1)

  • OR = 1 No hay efecto

  • OR < 1 Factor protector (Para interpretar veo cuanto le falta para llegar a 1)

Se interpreta con porcentajes

Ej 1: Se investiga el efecto del consumo de alcohol en la probabilidad de tener un accidente (Comparamos dos grupos: Consumidores y no consumidores de alcohol) Se obtiene un OR= 2.8.

Ya que OR >1 veo cuanto le sobra al 1 (2.8-1 = 1.8) -> Las odds de tener un accidente son 180% mayores para aquellos que consumen alcohol en comparación a los que no consumen alcohol.

Ej 2: Se investiga el efecto de la dieta en la probabilidad de tener cancer estomacal, comparando a un grupo vegetariano con uno omnivoro. OR = 0.9

Ya que OR<1 veo cuánto le falta para llegar a 1 (1-0.9 = 0.1) -> Las odds de tener cancer estomacal son un 10% menores para los vegetarianos.

El logaritmo de las Odds nos permite generar simetría.

Ej: En una semana (7 días) llueve solo 1.

Odds de que llueva 1/6=0.17 -> log(0,17)=−1,79

Odds de que no llueva 6/1=6 -> log(6)=1,79

Esta caracteristica es deseable en estadistica. Por eso los modelos de regresión logit se calculan con el log(odds) y no directamente con odds.

Modelo de regresión Logit

El objetivo de la regresión logística es determinar P(Y =1|X1…Xk) suponiendo que 1 es presencia y 0 ausencia.

Siempre tiene forma de s, ya que la linea busca capturar las tendencias de dos grupos diferentes.

El modelo se define como:

Técnicamente, los coeficientes beta están midiendo el cambio en las log(odds). Por esto para interpreatarlos deberemos sacarles la exponencial (exp(coeficientes)).

Al interpretar

Predictores cualis: Ej: OR de tener diabetes: 0.8 al comparar entre hacer actividad física (1= ref) o no hacer actividad física (0=dummy)

Las odds de tener diabetes disminuyen en un 20% al hacer actividad física

Predictor cuanti:Ej: OR diabetes: 1.2 en base a los mg/dL de glucosa en sangre

Por cada mg/dL adicional de glucosa, las odds de tener diabetes aumentan en un 20%

Los modelos logit pertenecen a la familia de los Modelos Lineales Generalizados (GLM), una extensión de los modelos lineales que permite abordar situaciones donde la variable dependiente no sigue una distribución normal. En un GLM, se especifica una distribución para la variable respuesta (como binomial, Poisson, entre otras) y se utiliza una función de enlace que transforma la relación no lineal en una forma linealizable. En el caso de la regresión logística, se asume una distribución binomial y se emplea la función de enlace logit.

Veamos un ejemplo en R. Voy a crear una base de datos sobre diabetes (0=no, 1= si) con edad, imc, presion y glucosa. Luego crearé un modelo logit binomial.

set.seed(123)

n <- 1000
edad <- rnorm(n, mean = 50, sd = 10)
imc <- rnorm(n, mean = 27, sd = 4)
presion <- rnorm(n, mean = 130, sd = 15)
glucosa <- rnorm(n, mean = 100, sd = 25)

# Modelo más realista (coeficientes más suaves)
logit_p <- -12 + 0.03*edad + 0.07*imc + 0.02*presion + 0.05*glucosa
p <- 1 / (1 + exp(-logit_p))
diabetes <- rbinom(n, size = 1, prob = p)

datos_diabetes <- data.frame(
  diabetes = factor(diabetes, levels = c(0, 1)),
  edad = edad,
  imc = imc,
  presion = presion,
  glucosa = glucosa
)

#AJUSTAR MODELO
modelo_diab <- glm(diabetes ~ ., data = datos_diabetes, family = binomial)
coef <- coef(modelo_diab)
exp(coef)
##  (Intercept)         edad          imc      presion      glucosa 
## 5.486063e-06 1.036758e+00 1.074902e+00 1.017178e+00 1.051182e+00

Por cada año adicional de edad, las odds de tener diabetes aumentan en un 3.67%.

Por cada unidad adicional de IMC, las odds de tener diabetes aumentan en un 7.49%.

Por cada mmHg adicional de presión, las odds aumentan en un 1.71%.

Por cada mg/dL adicional de glucosa, las odds aumentan en un 5.11%.

Si queremos probar la capacidad previctiva de nuestro modelo con una matriz de confusión, deberíamos separar nuestros datos en train/test, y luego al probar en test tenemos que usar type = “response” para transformar el logit a probabilidad. Nosotros debemos asignar un punto de corte (podemos determinarlo con la curva ROC) sobre el cual se clasificará como diabetes o no (ej: (p_diab > 0.5, 1, 0) Si la prob es sobre 0.5 se clasifica como diabetico)

#separar datos

index_d <- createDataPartition(datos_diabetes$diabetes, p = 0.7, list = F)
train_d <- datos_diabetes[index_d,] #train queda con el 70%
test_d <- datos_diabetes[-index_d,] #test queda con el 30%

#AJUSTAR MODELO
modelo_diab2 <- glm(diabetes ~ ., data = train_d, family = binomial)

#Predicciones

# Paso 1: Asegúrate de predecir probabilidades
prob_diab <- predict(modelo_diab2, test_d, type = "response") # transforma el logit a probabilidad mediante la función logística inversa.

# Paso 2: Convertir a clases usando un punto de corte (ej. 0.5) Hay metodos para elegir el mjr punto de corte. Si las clases están desbalanceadas puede que 0.5 no lo sea
pred_clase <- ifelse(prob_diab > 0.5, 1, 0)

# Paso 3: Convertir a factor con mismos niveles que test$diabetes
pred_clase <- factor(pred_clase, levels = c(0, 1))
test_d$diabetes <- factor(test_d$diabetes, levels = c(0, 1))

# Paso 4: Matriz de confusión
confusionMatrix(pred_clase, test_d$diabetes, positive = "1")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 193  41
##          1  17  48
##                                           
##                Accuracy : 0.806           
##                  95% CI : (0.7566, 0.8493)
##     No Information Rate : 0.7023          
##     P-Value [Acc > NIR] : 3.078e-05       
##                                           
##                   Kappa : 0.497           
##                                           
##  Mcnemar's Test P-Value : 0.002527        
##                                           
##             Sensitivity : 0.5393          
##             Specificity : 0.9190          
##          Pos Pred Value : 0.7385          
##          Neg Pred Value : 0.8248          
##              Prevalence : 0.2977          
##          Detection Rate : 0.1605          
##    Detection Prevalence : 0.2174          
##       Balanced Accuracy : 0.7292          
##                                           
##        'Positive' Class : 1               
##

Otra forma de probar los modelos es con la Curva ROC y el AUC, se recomienda siempre usar este método de evaluación

Curva ROC/AUC

La curva ROC es un gráfico que representa visualmente el rendimiento de un modelo de clasificación en todos sus posibles puntos (thresholds). Grafica el VPR (Sensibilidad) x FPR (Especificidad - 1). El AUC es el area debajo de la curva ROC

Al medio siempre habrá una linea recta, que representa un modelo sin poder de discriminación (AUC = 0.5).

Una curva ROC perfecta tendría Sen = 1 y Esp = 1 (FP = 0), y su AUC= 1. Mientras más se acerque nuestro modelo a esa curva perfecta.

Podemos usar el AUC y curva roc para hacer comparaciones entre modelos. Sobre todo cuando hay desbalances de clases, aunque se recomienda hacerlo siempre que se hagan modelos predictivos de clasificación. Mientras más AUC, mejor.

Otra funcionalidad es poder determinar el punto de corte ideal (threshold) de probabilidad para la clasificación, ya que dependiendo de este tendré más o menos sensibilidad y especificiada.

En general, en coords, la función de R que nos permite obtener el threshold, especificamos “best” para buscar el threshold que maximiza la suma de sensibilidad y especificidad (es decir, el punto más equilibrado entre ambos). Esto es útil cuando ambos errores (falsos positivos y falsos negativos) son igual de graves.

Hay casos especificos en los que puede que necesitemos más sensibilidad (cuando no nos podemos arriesgar a falsos negativos, por ejemplo, para diagnosticar enfermedades) o más especificidad (cuando no nos podemos arriesgar a un falso positivo, por ejemplo, para la detección de fraudes bancarios). Ante esto también podemos utilizar coords a modo de que favorezca lo que necesitamos.

  • Threshold bajo, más casos superan ese valor → clasificas más positivos, incluso algunos que no lo eran realmente → favorece sensibilidad, pero baja especificidad.

  • Threshold alto, menos casos superan ese valor → clasificas menos positivos, pero con mayor seguridad → aumenta especificidad, pero baja sensibilidad.

Veamoslo en el ejemplo de predicción de diabetes con regresión logit. En el apartado anterior creamos un modelo al que le asignamos (prob_diab > 0.5). Veamos su AUC

#Curva roc
curva_rocd <- roc(test_d$diabetes, prob_diab)

plot(curva_rocd)

auc(curva_rocd) #Buen AUC
## Area under the curve: 0.8162

Es un buen AUC, aún así debemos mirar el threshold

Asignaremos “best” para buscar el punto más equilibrado entre sensibilidad y especificidad

coords(curva_rocd, "best")
##   threshold specificity sensitivity
## 1 0.3424608   0.7666667   0.7303371

El punto donde se maximiza la sensibilidad y especificidad es 0.34. Veamos qué pasa con la matriz de confusión al elegir este nuevo punto de corte

# Paso 1: Asegúrate de predecir probabilidades
prob_diab <- predict(modelo_diab2, test_d, type = "response") # transforma el logit a probabilidad mediante la función logística inversa.

# Paso 2:Cambiar punto de corte
pred_clase <- ifelse(prob_diab > 0.34, 1, 0)

# Convertir a factor con mismos niveles que test$diabetes
pred_clase <- factor(pred_clase, levels = c(0, 1))
test_d$diabetes <- factor(test_d$diabetes, levels = c(0, 1))

# Paso 3: Matriz de confusión
confusionMatrix(pred_clase, test_d$diabetes, positive = "1")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 160  24
##          1  50  65
##                                           
##                Accuracy : 0.7525          
##                  95% CI : (0.6996, 0.8004)
##     No Information Rate : 0.7023          
##     P-Value [Acc > NIR] : 0.031687        
##                                           
##                   Kappa : 0.454           
##                                           
##  Mcnemar's Test P-Value : 0.003659        
##                                           
##             Sensitivity : 0.7303          
##             Specificity : 0.7619          
##          Pos Pred Value : 0.5652          
##          Neg Pred Value : 0.8696          
##              Prevalence : 0.2977          
##          Detection Rate : 0.2174          
##    Detection Prevalence : 0.3846          
##       Balanced Accuracy : 0.7461          
##                                           
##        'Positive' Class : 1               
##

Bajó el accuracy, pero mejoró la sensibilidad

Ajustes/Correcciones

SMOTE

Ante desbalance de clases se recomienda usar SMOTE, una técnica que genera nuevos ejemplos sintéticos de la clase minoritaria para balancear el conjunto de entrenamiento.

SMOTE no genera datos aleatorios, sino que hace interpolaciones realistas entre ejemplos reales de la misma clase minoritaria. No conviene aplicar SMOTE antes de hacer el train/test split → se debe aplicar dentro del trainControl()

Para ejemplificarlo simule un conjunto de datos de una enfermedad rara,con una prevalencia del 0.05, por lo que hay un desbalance de clases:

set.seed(123)

# Generamos 1000 observaciones
n <- 1000
data <- data.frame(
  edad = rnorm(n, 45, 15),
  sexo = rbinom(n, 1, 0.5),  # 0 = mujer, 1 = hombre
  presion = rnorm(n, 120, 15),
  colesterol = rnorm(n, 200, 30),
  enfermedadX = sample(c(0, 1), n, replace = TRUE, prob = c(0.95, 0.05))  # 1 = tiene la enfermedad
)

#antes de smote
ggplot(data, aes(x = enfermedadX)) +
  geom_bar(fill = "steelblue") +
  labs(title = "Distribución antes de SMOTE",
       x = "Tiene la enfermedad", y = "Frecuencia")

Cuando usamos SMOTE este desbalance se corrige, y se generan nuevos datos para la clase minoritaria

#APLICAR SMOTE
# Dividir datos en variables predictoras (X) y la variable objetivo (y)
X <- data[, -which(names(data) == "enfermedadX")]
y <- data$enfermedadX

# Aplicar SMOTE: duplica clase minoritaria 19 veces (ajustable)
smote_data <- SMOTE(X, y, dup_size = 19)$data

# Renombrar la clase como enfermedadX
smote_data$enfermedadX <- as.factor(smote_data$class)
smote_data$class <- NULL

#Despues de smote
ggplot(smote_data, aes(x = enfermedadX)) +
  geom_bar(fill = "darkgreen") +
  labs(title = "Después de SMOTE", x = "Tiene la enfermedad", y = "Frecuencia")

Problemas de escala

Este problema se da en situaciones en que las unidades de medida de las variables explicativas son muy pequeñas (por ejemplo, milímetros o fracciones decimales), lo que provoca que, al interpretar los coeficientes del modelo logístico como efectos por “una unidad de cambio”, los resultados parezcan exageradamente grandes o poco realistas.

Esto no significa que el modelo esté mal, sino que la escala de las variables distorsiona la interpretación de las pendientes.

Para ejemplificar usemos el modelo de cancer de mama, donde la diferencia entre el tamaño de los tumores está en unidades de medida muy pequeñas.

Sin corregir las escalas veremos valores extraños, muy extremos:

modelo_log <- glm(factor(diagnosis) ~ .-smoothness_mean, data = datos, family = binomial)

coeficientes <- coef(modelo_log)

exp(coeficientes)
##             (Intercept)             radius_mean            texture_mean 
##            9.272327e-29            3.887911e+00            1.653684e+00 
##        compactness_mean          concavity_mean           symmetry_mean 
##            1.142697e-25            1.207148e+33            1.796275e-15 
##  fractal_dimension_mean               radius_se              texture_se 
##            2.214863e-06            5.250296e+09            7.128854e-01 
##           smoothness_se          compactness_se            concavity_se 
##            4.736145e+71            3.551773e+69            7.314839e-49 
##       concave.points_se             symmetry_se    fractal_dimension_se 
##           6.411961e+188           6.001345e-119            0.000000e+00 
##        smoothness_worst       compactness_worst         concavity_worst 
##            7.641034e+19            1.657350e-11            7.121752e+04 
##          symmetry_worst fractal_dimension_worst 
##            1.095300e+22           1.504683e+119

Ahora, al aplicar la corrección, estos valores no serán tan extremos:

# Estandarizar todas las columnas numéricas excepto el diagnóstico
datos_escalados <- datos
datos_escalados[, -1] <- scale(datos[, -1])  # Asumiendo que diagnosis está en la 1ra columna

#modelo
modelo_log2 <- glm(factor(diagnosis) ~ .-smoothness_mean, data = datos_escalados, family = binomial)

coeficientes <- coef(modelo_log2)

exp(coeficientes)
##             (Intercept)             radius_mean            texture_mean 
##              0.92325881            119.72611986              8.70105104 
##        compactness_mean          concavity_mean           symmetry_mean 
##              0.04816600            433.78210406              0.39423867 
##  fractal_dimension_mean               radius_se              texture_se 
##              0.91217066            496.05601951              0.82969514 
##           smoothness_se          compactness_se            concavity_se 
##              1.64137098             17.60037191              0.03523465 
##       concave.points_se             symmetry_se    fractal_dimension_se 
##             14.62153528              0.10537439              0.00237280 
##        smoothness_worst       compactness_worst         concavity_worst 
##              2.84433449              0.02012892             10.28886421 
##          symmetry_worst fractal_dimension_worst 
##             23.09564454            142.06781369