Code
library(tidyverse)
library(bayesQR)
library(broom)
library(corrplot)
library(yardstick)
library(pROC)
#library(ggrepel)
# library(purrr)
library(plotly)Unidad 3(RLog): Determinación de Sistemas de Aprendizaje automático
Regresión Logística
Algoritmos de Machine Learning
Regresión Lineal
Regresión Logística
Árboles de Decisión
k Nearest Neighbor
SVM
Naive Bayes
Clustering jerárquico
K-Means
PCA
Redes Neuronales
Aprendizaje profundo
¿Qué es la regresión logística?
Esta introducción ha sido extraída del texto Aurora Gónzalez Vidal (2022).
La regresión logística es el conjunto de modelos estadísticos utilizados cuando se desea conocer la relación entre:
Una variable dependiente cualitativa, dicotómica (regresión logística binaria o binomial) o con más de dos categorías (regresión logística multinomial).
Una o más variables explicativas, llamadas covariables, ya sean cualitativas o cuantitativas.
Las covariables cualitativas deben ser dicotómicas, tomando valor 0 para su ausencia y 1 para su presencia. Si la covariable tuviera más de dos categorías debemos realizar una transformación de la misma en varias covariables cualitativas dicotómicas ficticias (variables dummy). Al hacer esta transformación cada categoría de la variable entraría en el modelo de forma individual.
Los modelos de regresión logística tienen tres finalidades:
Cuantificar la importancia de la relación existente entre cada una de las covariables y la variable dependiente.
Clarificar la existencia de interacción y confusión entre covariables respecto a la variable dependiente (es decir, los odds ratio para cada covariable).
Clasificar individuos dentro de las categorías (presente/ausente) de la variable dependiente.
Por tanto, el objetivo de la regresión logística no es, como en regresión lineal, predecir el valor de la variable \(Y\) a partir de una o varias variables predictoras (\(X_s\)), sino que queremos predecir la probabilidad de que ocurra Y conocidos los valores de las variables \(X_s\). La ecuación general es de la forma:
\[\begin{equation} P(Y) = \frac{1}{1 + e^{-(b_0 + b_1 X_1 + b_2 X_2 + \dots + b_n X_n)}}\end{equation}\] En su forma más sencilla, cuando tenemos sólo una variable predictora \(X_1\), la ecuación de la regresión logística (simple) viene dada por:
\[\begin{equation} P(Y) = \frac{1}{1 + e^{-(b_0 + b_1 X_1)}}\end{equation}\]
Los valores posibles de estas ecuaciones varían entre 0 y 1. Un valor cercano a 0 significa que es muy improbable que Y haya ocurrido, y un valor cercano a 1 significa que es muy probable que tuviese lugar.
Como en la regresión lineal cada variable predictora de la ecuación logística tiene su propio coeficiente. Los valores de los parámetros se estiman utilizando el método de máxima verosimilitud que selecciona los coeficientes que hacen más probable que los valores observados ocurran.
Fávero, Belfiore, and de Freitas Souza (2023)
Ejemplo Regresión logística simple (binaria)
El dataset de ejemplo es sobre una compañía de servicios financieros. Los datos se encuentran en la librería bayesQR.
# Cargamos los datos
data(Churn) #bayesQRCode
glimpse(Churn)Rows: 400
Columns: 5
$ churn <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ gender <int> 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 1, 0, 1, …
$ Social_Class_Score <dbl> -0.11567485, -0.31342471, -1.32543870, 1.96651484, …
$ lor <dbl> -1.08922097, 1.18298297, -0.84615637, 0.08694165, -…
$ recency <dbl> -0.72132150, 3.63443539, -0.42758226, -0.53567170, …Hay 400 filas y cada fila representa un cliente. Los datos que vamos a utilizar como predictores son el tiempo desde la primera compra (lor - length of relationship) y el tiempo desde la última compra (recency - recency of activity). La variable de respuesta será churn, es decir si el cliente ha abandonado la compañía. Los datos de tiempo están normalizados por z-score, con lo que no representan unidades de tiempo y por supuesto, hay valores negativos.
Code
Churn %>%
cor() %>%
corrplot(method = "number")
Modelo lineal simple: churn vs recency
La variable de respuesta churn en función de recency:
mdl_churn_vs_recency_lm <- lm(churn ~ recency, Churn)
summary(mdl_churn_vs_recency_lm)
Call:
lm(formula = churn ~ recency, data = Churn)
Residuals:
Min 1Q Median 3Q Max
-0.8689 -0.4636 -0.1329 0.5132 0.5647
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 0.49078 0.02503 19.610 < 2e-16 ***
recency 0.06378 0.02252 2.832 0.00486 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 0.4963 on 398 degrees of freedom
Multiple R-squared: 0.01975, Adjusted R-squared: 0.01729
F-statistic: 8.02 on 1 and 398 DF, p-value: 0.00486Code
ggplot(Churn, aes(recency, churn, color = churn)) +
geom_point(alpha = 0.5) +
geom_smooth(method = "lm", se = FALSE) +
theme_bw()
Code
ggplot(Churn, aes(recency, churn, color = churn)) +
geom_jitter(height = 0.1, alpha = 0.5) +
geom_smooth(method = "lm", se = FALSE) +
theme_bw()
Obviamente, la variable de respuesta no puede ser un número decimal. Lo que nos da el modelo de regresión lineal, realmente no son las predicciones, sino la probabilidad de que el cliente abandone o no la compañía, en función del tiempo transcurrido desde la última actividad. A mayor tiempo transcurrido, más probable que el cliente haya abandonado.
Si sustituimos geom_smooth por abline y modificamos los límites de los ejes, vemos la gráfica desde fuera, para ver la línea de tendencia completa:
Code
coeffs <- coefficients(mdl_churn_vs_recency_lm)
intercept <- coeffs[1]
slope <- coeffs[2]
ggplot(Churn, aes(recency, churn, color = churn)) +
geom_point(alpha = 0.1) +
geom_abline(intercept = intercept, slope = slope) +
xlim(-10, 10) +
ylim(-0.2, 1.2) +
theme_bw()
Modelo lineal generalizado: glm()
glm(churn ~ recency, Churn, family = gaussian)
Call: glm(formula = churn ~ recency, family = gaussian, data = Churn)
Coefficients:
(Intercept) recency
0.49078 0.06378
Degrees of Freedom: 399 Total (i.e. Null); 398 Residual
Null Deviance: 100
Residual Deviance: 98.02 AIC: 578.7El término family indica la distribución de los residuos del modelo.
Para un modelo de regresión logística binaria, usaremos binomial
mdl_churn_vs_recency_glm_bin <- glm(churn ~ recency, Churn, family = binomial)
mdl_churn_vs_recency_glm_bin
Call: glm(formula = churn ~ recency, family = binomial, data = Churn)
Coefficients:
(Intercept) recency
-0.03502 0.26921
Degrees of Freedom: 399 Total (i.e. Null); 398 Residual
Null Deviance: 554.5
Residual Deviance: 546.4 AIC: 550.4Predicciones como probabilidad
Es importante especificar en las predicciones el atributo type = “response”
explanatory_data <- data.frame(recency = -10:10)
predicted <- explanatory_data %>%
mutate(churn_probability = predict(mdl_churn_vs_recency_glm_bin,
newdata = ., type = "response"),
churn_log_odds = predict(mdl_churn_vs_recency_glm_bin,
newdata = .))Code
head(predicted) recency churn_probability churn_log_odds
1 -10 0.06138909 -2.727169
2 -9 0.07885883 -2.457954
3 -8 0.10076629 -2.188739
4 -7 0.12791464 -1.919524
5 -6 0.16106717 -1.650309
6 -5 0.20083332 -1.381094Code
plt <- ggplot(Churn) +
geom_point(data = Churn, aes(recency, churn, color = churn), alpha = 0.1) +
geom_line(data = predicted, aes(recency, churn_probability), color = "blue") +
geom_point(data = predicted, aes(recency, churn_probability), color = "blue") +
theme_bw()
plt
De probabilidades a predicción de churn: el más probable
Se escoge arbitrariamente un valor límite a partir del cual consideraremos churn o no churn. De momento y por simplificar, diremos que es 0.5, si la probabilidad de churn > 0.5 churn_predicted = 1
Code
predicted <- predicted %>%
mutate(churn_predicted = round(churn_probability))
plt <- plt +
geom_point(data = predicted, aes(recency, churn_predicted), color = "brown", alpha = 0.5)
plt
Ratios de probabilidad y relación con los coeficientes
\[odds\_ratio = \frac{p}{1-p}\] \[\log\left(\frac{1 - p}{p}\right) = \beta_0 + \beta_1X\]
Exponenciando ambos lados para eliminar el logaritmo:
\[\frac{p}{1 - p} = e^{\beta_0 + \beta_1X}\]
Resolver para p:
\[p = \frac{e^{\beta_0 + \beta_1X}}{1 + e^{\beta_0 + \beta_1X}}\quad \text{o de forma equivalente} \quad p = \frac{1}{1 + e^{-(\beta_0 + \beta_1X)}}\]
Code
predicted <- predicted %>%
mutate(churn_odds_ratio = churn_probability/(1-churn_probability))
predicted recency churn_probability churn_log_odds churn_predicted churn_odds_ratio
1 -10 0.06138909 -2.72716897 0 0.06540419
2 -9 0.07885883 -2.45795401 0 0.08560993
3 -8 0.10076629 -2.18873906 0 0.11205796
4 -7 0.12791464 -1.91952411 0 0.14667675
5 -6 0.16106717 -1.65030916 0 0.19199054
6 -5 0.20083332 -1.38109421 0 0.25130342
7 -4 0.24752070 -1.11187926 0 0.32894022
8 -3 0.30097395 -0.84266430 0 0.43056185
9 -2 0.36044129 -0.57344935 0 0.56357810
10 -1 0.42452268 -0.30423440 0 0.73768793
11 0 0.49124603 -0.03501945 0 0.96558664
12 1 0.55828273 0.23419550 1 1.26389156
13 2 0.62326047 0.50341045 1 1.65435376
14 3 0.68408855 0.77262541 1 2.16544397
15 4 0.73920495 1.04184036 1 2.83442858
16 5 0.78768969 1.31105531 1 3.71008694
17 6 0.82924279 1.58027026 1 4.85626810
18 7 0.86406665 1.84948521 1 6.35654641
19 8 0.89270749 2.11870016 1 8.32031542
20 9 0.91590112 2.38791512 1 10.89076428
21 10 0.93444909 2.65713007 1 14.25531855Cuantificando el ajuste del modelo: matriz de confusión
TP: True Positive. La predicción es correcta y el resultado es 1.
TN: True Negative. La predicción es correcta y el resultado es 0.
FP: False Positive. Predicción fallida y resultado real = 0 vs predicción = 1.
FN: False Negative. Predicción fallida y resultado real = 1 vs predicción = 0.
Si el modelo es multinomial o multiclase, la matriz de confusión se hace más compleja.
# Matriz de confusión manual:
real_churn <- Churn$churn
predicted_churn <- round(fitted(mdl_churn_vs_recency_glm_bin))
outcomes <- table(predicted_churn, real_churn)confusion <- conf_mat(outcomes)
autoplot(confusion)
Métricas
summary(confusion, event_level = "second") # el atributo es para considerar positivos# A tibble: 13 × 3
.metric .estimator .estimate
<chr> <chr> <dbl>
1 accuracy binary 0.575
2 kap binary 0.15
3 sens binary 0.445
4 spec binary 0.705
5 ppv binary 0.601
6 npv binary 0.560
7 mcc binary 0.155
8 j_index binary 0.150
9 bal_accuracy binary 0.575
10 detection_prevalence binary 0.37
11 precision binary 0.601
12 recall binary 0.445
13 f_meas binary 0.511Exactitud (Accuracy): Es una medida general de cuántas predicciones fueron correctas. Un modelo es generalmente mejor cuando la exactitud es alta. \(\text{Accuracy} = \frac{TP + TN}{TP + TN + FP + FN}\)
Sensibilidad (Sens): También conocida como tasa de verdaderos positivos o recall. Indica la proporción de casos positivos reales que fueron correctamente identificados por el modelo. Una sensibilidad alta es preferible, especialmente en situaciones donde no detectar los casos positivos (como enfermedades) es crítico. \(\text{Recall} = \frac{TP}{TP + FN}\)
Especificidad (Spec): También conocida como tasa de verdaderos negativos. Muestra la proporción de casos negativos reales que fueron correctamente identificados. Al igual que con la sensibilidad, una especificidad alta es generalmente mejor, particularmente en situaciones donde los falsos positivos son problemáticos. \(\text{Specificity} = \frac{TN}{TN + FP}\)
ROC Área bajo la curva (AUC)
ROC: La curva ROC es un gráfico que muestra el rendimiento de un modelo de clasificación en todos los umbrales de clasificación posibles. Esta curva traza dos parámetros: la Tasa de Verdaderos Positivos (Sensibilidad o Recall) en el eje Y y la Tasa de Falsos Positivos (1 - Especificidad) en el eje X.
Cómo se usa: Al cambiar el umbral de decisión (cutoff) para la clasificación (el punto en el que decidimos clasificar un resultado como positivo o negativo), la Sensibilidad y la Especificidad del modelo cambiarán, creando un trazado en el gráfico. Esta curva proporciona una medida clara de la compensación entre Sensibilidad y Especificidad en diferentes umbrales, y es especialmente útil cuando las clases están desequilibradas.
Qué se necesita para crearla: necesitamos el vector de la variable de respuesta y un vector de probabilidades.
AUC: es una métrica que resume la curva ROC en un solo valor, calculando el área bajo la curva ROC. El valor de AUC varía entre 0 y 1. Un modelo con un AUC de 0.5 no tiene capacidad discriminativa y es equivalente a una decisión aleatoria. Un modelo perfecto tiene un AUC de 1.0, donde identifica correctamente todos los verdaderos positivos y verdaderos negativos. Por lo tanto, cuanto más alto sea el valor del AUC, mejor será el modelo en la diferenciación entre las clases positivas y negativas.
ROC <- roc(Churn$churn, fitted(mdl_churn_vs_recency_glm_bin))Setting levels: control = 0, case = 1Setting direction: controls < casesplot(ROC, col = "blue")
Code
auc(ROC)Area under the curve: 0.5945Ejercicio 1: Regresión logística simple
Repite todos los pasos usando como variable explicativa lor. ¿Qué modelo obtiene mejores resultados en términos de exactitud y AUC?
Ejercicio 2: Regresión logística múltiple
Cuando usamos dos o más variables predictoras, las visualizaciones se complican. La función ya no es una sigmoide con un par de coeficientes, sino que al igual que ocurre con la regresión lineal múltiple, se convierten en gráficas de superficie en el espacio o fórmulas más complejas que no tiene sentido representar gráficamente. Además pueden aparecer interacciones, lo cual complica aún más el modelo añadiendo más coeficientes.
Crea un modelo usando las dos variables explicativas: lor y recency. Repite todos los pasos excepto las visualizaciones. Usa la fórmula churn ~ recency+lor
Ejercicio 3: Regresión logística múltiple
Repite el ejercicio 2, usa la fórmula con interacciones
References
Aurora Gónzalez Vidal, Antonio Maurandi López &. 2022. Análisis de Datos y Métodos Estadísticos Con r. Ediciones de la Universidad de Murcia. https://doi.org/https://doi.org/10.6018/editum.2967.
Fávero, Luiz Paulo, Patrícia Belfiore, and Rafael de Freitas Souza. 2023. Chapter 15 - Binary and Multinomial Logistic Regression Models. Edited by Luiz Paulo Fávero, Patrícia Belfiore, and Rafael de Freitas Souza. Academic Press. https://doi.org/https://doi.org/10.1016/B978-0-12-824271-1.00008-1.