Regresión Logística
Lic. Miguel Angel
2024-03-10
Ejemplo de Regresión Logística
Un estudio considera que existe relación entre el hecho de que un estudiante asista a clases de repaso de lectura (sí = 1, no = 0), la nota que obtiene en un examen de lectura estándar (realizado antes de iniciar las clases de repaso) y el sexo (hombre = 1, mujer = 0). Se quiere generar un modelo en el que a partir de las variables puntuación del examen y sexo, prediga la probabilidad de que el estudiante tenga que asistir a clases de repaso.
datos <- data.frame(sexo = c("hombre", "hombre", "mujer", "mujer", "mujer", "hombre",
"mujer", "hombre", "mujer", "mujer", "hombre", "hombre",
"hombre", "hombre", "mujer", "mujer", "hombre", "mujer",
"hombre", "mujer", "hombre", "mujer", "mujer", "hombre",
"hombre", "mujer", "mujer", "mujer", "hombre", "hombre",
"hombre", "mujer", "hombre", "mujer", "hombre", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "mujer",
"hombre", "mujer", "mujer", "mujer", "mujer", "hombre",
"mujer", "hombre", "mujer", "hombre", "mujer", "mujer",
"hombre", "hombre", "hombre", "hombre", "hombre", "hombre",
"hombre", "hombre", "hombre", "mujer", "hombre", "hombre",
"hombre", "hombre", "mujer", "hombre", "mujer", "hombre",
"hombre", "hombre", "mujer", "hombre", "mujer", "mujer",
"hombre", "mujer", "mujer", "mujer", "hombre", "hombre",
"hombre", "hombre", "hombre", "mujer", "mujer", "mujer",
"mujer", "hombre", "mujer", "mujer", "mujer", "mujer",
"mujer", "mujer", "mujer", "mujer","mujer", "mujer",
"hombre", "mujer", "hombre", "hombre", "mujer", "mujer",
"mujer", "hombre","mujer", "hombre", "mujer", "mujer",
"mujer", "hombre", "mujer", "hombre", "mujer", "hombre",
"mujer", "hombre", "mujer", "mujer", "mujer", "mujer",
"mujer", "mujer", "mujer", "mujer", "hombre", "mujer",
"hombre", "hombre", "hombre", "hombre", "hombre", "hombre",
"hombre", "mujer","mujer", "mujer", "hombre", "hombre",
"mujer", "mujer", "hombre", "mujer", "hombre", "hombre",
"hombre", "mujer", "mujer", "mujer", "mujer", "hombre",
"hombre", "mujer", "hombre", "hombre", "mujer", "hombre",
"hombre", "hombre", "hombre", "mujer", "hombre", "hombre",
"mujer", "mujer", "hombre", "hombre", "hombre", "hombre",
"hombre", "mujer", "mujer", "mujer", "mujer", "hombre",
"hombre", "hombre", "mujer", "hombre", "mujer", "hombre",
"hombre", "hombre", "mujer"),
examen_lectura = c(91, 77.5, 52.5, 54, 53.5, 62, 59, 51.5,
61.5, 56.5, 47.5, 75, 47.5, 53.5, 50, 50,
49, 59, 60, 60, 60.5, 50, 101, 60, 60,
83.5, 61, 75, 84, 56.5, 56.5, 45, 60.5,
77.5, 62.5, 70, 69, 62, 107.5, 54.5, 92.5,
94.5, 65, 80, 45, 45, 66, 66, 57.5, 42.5,
60, 64, 65, 47.5, 57.5, 55, 55, 76.5,
51.5, 59.5, 59.5, 59.5, 55, 70, 66.5,
84.5, 57.5, 125, 70.5, 79, 56, 75, 57.5,
56, 67.5, 114.5, 70, 67, 60.5, 95, 65.5,
85, 55, 63.5, 61.5, 60, 52.5, 65, 87.5,
62.5, 66.5, 67, 117.5, 47.5, 67.5, 67.5,
77, 73.5, 73.5, 68.5, 55, 92, 55, 55, 60,
120.5, 56, 84.5, 60, 85, 93, 60, 65, 58.5,
85, 67, 67.5, 65, 60, 47.5, 79, 80, 57.5,
64.5, 65, 60, 85, 60, 58, 61.5, 60, 65,
93.5, 52.5, 42.5, 75, 48.5, 64, 66, 82.5,
52.5, 45.5, 57.5, 65, 46, 75, 100, 77.5,
51.5, 62.5, 44.5, 51, 56, 58.5, 69, 65,
60, 65, 65, 40, 55, 52.5, 54.5, 74, 55,
60.5, 50, 48, 51, 55, 93.5, 61, 52.5,
57.5, 60, 71, 65, 60, 55, 60, 77, 52.5,
95, 50, 47.5, 50, 47, 71, 65),
clases_repaso = c("0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "0", "0", "0", "0", "0", "0",
"0", "0", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1", "1", "1", "1",
"1", "1", "1", "1", "1"))
datos$sexo <- as.factor(datos$sexo)
datos$clases_repaso <- as.factor(datos$clases_repaso)library(ggplot2)
tabla <- table(datos$clases_repaso, datos$sexo,
dnn = c("clases de repaso","sexo"))
addmargins(tabla)## sexo
## clases de repaso hombre mujer Sum
## 0 57 73 130
## 1 36 23 59
## Sum 93 96 189tabla_frecuencias <- prop.table(tabla)*100
addmargins(tabla_frecuencias)## sexo
## clases de repaso hombre mujer Sum
## 0 30.15873 38.62434 68.78307
## 1 19.04762 12.16931 31.21693
## Sum 49.20635 50.79365 100.00000ggplot(data = datos, aes(x = clases_repaso, y = examen_lectura, colour = sexo)) +
geom_boxplot() +
theme_bw() +
theme(legend.position = "bottom")
El número de estudiantes en la muestra es semejante para ambos sexos (96, 93). Parece ser mayor el porcentaje de hombres que necesitan clases de repaso (19.04762, 12.16931). El promedio de las notas de lectura de los estudiantes que requieren de clases particulares es menor que el de los que no requieren clases. En vista de estos datos, tiene sentido considerar el sexo y la nota como posibles predictores.
Generar el modelo de regresión logística
modelo_glm <- glm(clases_repaso ~ examen_lectura + sexo, data = datos,
family = "binomial")
summary(modelo_glm)##
## Call:
## glm(formula = clases_repaso ~ examen_lectura + sexo, family = "binomial",
## data = datos)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.18365 0.78559 1.507 0.1319
## examen_lectura -0.02617 0.01223 -2.139 0.0324 *
## sexomujer -0.64749 0.32484 -1.993 0.0462 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 234.67 on 188 degrees of freedom
## Residual deviance: 224.64 on 186 degrees of freedom
## AIC: 230.64
##
## Number of Fisher Scoring iterations: 4Acorde al modelo, el logaritmo de odds de que un estudiante necesite clases de repaso esta negativamente relacionado con la puntuación obtenida en el examen de lectura (coeficiente parcial = -0.02617), siento significativa esta relación (p-value = 0.0324). También existe una relación significativa positiva entre el logaritmo de odds de necesitar clases de repaso y el género del estudiante (p-value = 0.0462), siendo, para un mismo resultado en el examen de lectura, mayor si el estudiante es hombre. En concreto los odds de que un hombre requiera clases de repaso son e^(0.64749)=1.910739 mayores que los de las mujeres. (Esto se puede ver gráficamente representando el modelo para hombres y mujeres).
Además del valor estimado de los coeficientes parciales de correlación calculados por el modelo, es conveniente generar sus correspondientes intervalos de confianza. En el caso de regresión logística, estos intervalos suelen calcularse basados en profile likelihood (en R es el método por defecto si se tiene instalado el paquete MASS).
confint(modelo_glm)## Waiting for profiling to be done...## 2.5 % 97.5 %
## (Intercept) -0.2973940 2.799267537
## examen_lectura -0.0518029 -0.003536647
## sexomujer -1.2931216 -0.015658932# En caso de querer los intervalos basados en el error estándar.
confint.default(modelo_glm)## 2.5 % 97.5 %
## (Intercept) -0.35608462 2.723378445
## examen_lectura -0.05015001 -0.002192983
## sexomujer -1.28416605 -0.010807125Representación gráfica del modelo
Al tratarse de un modelo con 2 predictores, no se puede obtener una representación en 2D en la que se incluyan ambos predictores a la vez. Sí es posible representar la curva del modelo logístico cuando se mantiene constante uno de los dos predictores. Por ejemplo, al representar las predicciones del modelo diferenciando entre hombres y mujeres (fijando el valor del predictor sexo) se aprecia que la curva de los hombres siempre está por encima. Esto se debe a que, como indica el coeficiente parcial de correlación del predictor sexo, para una misma nota en el examen de lectura, el logaritmo de ODDs de necesitar clases de repaso es 0.64749 veces mayor en hombres.
library(ggplot2)
# Para graficar los valores en ggplot junto con la curva, la variable respuesta
# tiene que ser numérica en lugar de factor.
datos$clases_repaso <- as.numeric(as.character(datos$clases_repaso))
# Se crea un dataframe que contenga la probabilidad de que se necesiten clases
# de repaso dada una determinada nota en el examen de lectura y siendo hombre.
# Vector con nuevos valores interpolados en el rango de observaciones.
nuevos_valores_examen <- seq(from = min(datos$examen_lectura),
to = max(datos$examen_lectura), by = 0.5)
sexo <- as.factor(rep(x = "hombre", length(nuevos_valores_examen)))
# Predicciones de los nuevos puntos según el modelo. type = "response" devuelve
# las predicciones en forma de probabilidad en lugar de en log_ODDs.
predicciones <- predict(object = modelo_glm,
newdata=data.frame(examen_lectura=nuevos_valores_examen,
sexo = sexo),
type = "response")
# Se crea un data frame con los nuevos puntos y sus predicciones para graficar
# la curva.
datos_curva_hombre <- data.frame(examen_lectura = nuevos_valores_examen,
sexo = sexo,
clases_repaso = predicciones)
# Mismo proceso para mujeres (sexo = 0).
nuevos_valores_examen <- seq(from = min(datos$examen_lectura),
to = max(datos$examen_lectura), by = 0.5)
sexo <- as.factor(rep("mujer", length(nuevos_valores_examen)))
predicciones <- predict(object = modelo_glm,
newdata=data.frame(examen_lectura=nuevos_valores_examen,
sexo = sexo),
type = "response")
datos_curva_mujer <- data.frame(examen_lectura = nuevos_valores_examen,
sexo = sexo, clases_repaso = predicciones)
# Se unifican los dos dataframe.
datos_curva <- rbind(datos_curva_hombre, datos_curva_mujer)
ggplot(data = datos, aes(x = examen_lectura, y = as.numeric(clases_repaso),
color = sexo)) +
geom_point() +
geom_line(data = datos_curva, aes(y = clases_repaso)) +
geom_line(data = datos_curva, aes(y = clases_repaso)) +
theme_bw() +
labs(title = "P. clases repaso en función de nota lectura y sexo",
y = "P(clase de repaso)") +
theme(plot.title = element_text(size = 10))
La imagen muestra un gráfico de dispersión con líneas de tendencia que representan la probabilidad de que los estudiantes tomen clases de repaso en función de su nota de lectura y diferenciado por sexo.
Líneas de Tendencia: Las líneas sólidas en el gráfico representan la tendencia general en los datos para cada sexo. La línea roja se refiere a los hombres y la azul a las mujeres. Estas líneas indican cómo cambia la probabilidad de tomar clases de repaso con la variación de las notas del examen de lectura.
Interpretación de las Tendencias:
A medida que las notas del examen de lectura aumentan (moviéndose de izquierda a derecha en el eje x), la probabilidad de asistir a clases de repaso disminuye para ambos sexos. Esto sugiere que los estudiantes con mejores notas de lectura sienten menos necesidad de asistencia adicional.
La línea de tendencia de las mujeres (azul) parece estar consistentemente por debajo de la de los hombres (rojo) en todo el rango de notas. Esto podría interpretarse como que, a igualdad de notas, las mujeres tienen una menor probabilidad de tomar clases de repaso en comparación con los hombres. Esto puede deberse a varios factores, como la confianza en la materia, las prácticas de estudio, o diferencias en la percepción de la necesidad de ayuda adicional.
Evaluación del modelo
Likelihood ratio:
# Diferencia de residuos
dif_residuos <- modelo_glm$null.deviance - modelo_glm$deviance
# Grados libertad
df <- modelo_glm$df.null - modelo_glm$df.residual
# p-value
p_value <- pchisq(q = dif_residuos,df = df, lower.tail = FALSE)
paste("Diferencia de residuos:", round(dif_residuos, 4))## [1] "Diferencia de residuos: 10.0334"paste("Grados de libertad:", df)## [1] "Grados de libertad: 2"paste("p-value:", round(p_value, 4))## [1] "p-value: 0.0066"El modelo en conjunto sí es significativo y, acorde a los p-values mostrados en el summary(), también es significativa la contribución al modelo de ambos predictores.
Comparación de las predicciones con las observaciones
Se va a emplear un threshold de 0.5. Si la probabilidad predicha de asistir a clases de repaso es superior a 0.5 se asigna al nivel 1 (sí asiste), si es menor se asigna al nivel 0 (no clases de repaso).
predicciones <- ifelse(test = modelo_glm$fitted.values > 0.5, yes = 1, no = 0)
matriz_confusion <- table(modelo_glm$model$clases_repaso, predicciones,
dnn = c("observaciones", "predicciones"))
matriz_confusion## predicciones
## observaciones 0 1
## 0 129 1
## 1 56 3El modelo es capaz de clasificar correctamente
\[ \frac{129+3}{129+3+56+1} = 0.698 \, (69.8\%) \]
de las observaciones de entrenamiento. Si se analiza en detalle cómo se distribuye el error, se aprecia que el modelo solo ha sido capaz de identificar correctamente a 3 de los 59 alumnos que realmente asisten a clases de repaso. El porcentaje de falsos negativos es muy alto. Seleccionar otro threshold puede mejorar la exactitud del modelo.
#install.packages("vcd")
library(vcd)## Warning: package 'vcd' was built under R version 4.3.3## Loading required package: gridpredicciones <- ifelse(test = modelo_glm$fitted.values > 0.45, yes = 1, no = 0)
matriz_confusion <- table(modelo_glm$model$clases_repaso, predicciones,
dnn = c("observaciones", "predicciones"))
matriz_confusion## predicciones
## observaciones 0 1
## 0 122 8
## 1 45 14El modelo logístico creado para predecir la probabilidad de que un alumno tenga que asistir a clases de repaso a partir de la nota obtenida en un examen de lectura y el sexo del alumno es en conjunto significativo acorde al Likelihood ratio (p-value = 0.0066). El p-value de ambos predictores es significativo (examen_lectura = 0.0324, sexo1 = 0.0462). El ratio de error obtenido empleando las observaciones con las que se ha entrenado el modelo muestra un porcentaje de falsos negativos muy alto.
\[ \text{logit}(\text{clases de repaso}) = 0.53616 - 0.02617 \times \text{examen\_lectura} + 0.64749 \times \text{sexo} \]
\[ P(\text{clases de repaso}) = \frac{e^{0.53616 - 0.02617 \times \text{examen\_lectura} + 0.64749 \times \text{sexo}}}{1 + e^{0.53616 - 0.02617 \times \text{examen\_lectura} + 0.64749 \times \text{sexo}}} \]
No hay que olvidar que los errores calculados son de entrenamiento, por lo que no son generalizables a nuevas observaciones. Para obtener una estimación más realista hay que calcular el error de test.