SEGUNDO TALLER - Módulo: Modelos de Regresión
Francisco José Galeano Bautista, Jhon Steven Medina Burbano
2026-04-18
# Librerías necesarias
library(readxl)
library(ggplot2)
library(GGally)
library(dplyr)
library(car)
library(lmtest)
library(nortest)
library(corrplot)
library(knitr)
library(kableExtra)1 Descripción del Problema e Hipótesis
1.1 Contexto
El dataset Iris (Fisher, 1936) contiene mediciones morfológicas de 150 flores pertenecientes a tres especies del género Iris: setosa, versicolor y virginica. Cada observación registra cuatro variables continuas (longitud y ancho del sépalo y del pétalo) expresadas en centímetros.
1.2 Objetivo analítico
Determinar si existe una relación estadísticamente significativa
entre las dimensiones del pétalo y el sépalo, con el fin de construir
modelos de regresión que permitan predecir la longitud del
sépalo (Sepal.Length) a partir de las demás
variables morfológicas.
1.3 Variable Respuesta (Y)
Sepal.Length: Longitud del sépalo (cm)
1.4 Variables Predictoras (X)
Petal.Length– Longitud del pétalo (cm)
Sepal.Width– Ancho del sépalo (cm)
Petal.Width– Ancho del pétalo (cm)
1.5 Hipótesis Estadísticas
1.5.1 Regresión Lineal Simple (Petal.Length → Sepal.Length)
\[H_0: \beta_1 = 0 \quad \text{(La longitud del pétalo NO predice la longitud del sépalo)}\] \[H_1: \beta_1 \neq 0 \quad \text{(La longitud del pétalo SÍ predice la longitud del sépalo)}\]
1.5.2 Regresión Lineal Múltiple
\[H_0: \beta_1 = \beta_2 = \beta_3 = 0 \quad \text{(Ninguna variable predictora explica la longitud del sépalo)}\] \[H_1: \exists\ \beta_j \neq 0 \quad \text{(Al menos una variable predictora es significativa)}\]
Nivel de significancia: \(\alpha = 0.05\)
2 Descripción de las Variables
# FIX: Se usa el dataset iris nativo de R en lugar de leer un Excel externo
# Esto garantiza compilación sin dependencia de archivo externo
datos <- iris
# Renombrar columnas al formato usado en el análisis
# (iris ya tiene Sepal.Length, Sepal.Width, Petal.Length, Petal.Width, Species)
datos_num <- datos %>% select(-Species)
glimpse(datos)## Rows: 150
## Columns: 5
## $ Sepal.Length <dbl> 5.1, 4.9, 4.7, 4.6, 5.0, 5.4, 4.6, 5.0, 4.4, 4.9, 5.4, 4.…
## $ Sepal.Width <dbl> 3.5, 3.0, 3.2, 3.1, 3.6, 3.9, 3.4, 3.4, 2.9, 3.1, 3.7, 3.…
## $ Petal.Length <dbl> 1.4, 1.4, 1.3, 1.5, 1.4, 1.7, 1.4, 1.5, 1.4, 1.5, 1.5, 1.…
## $ Petal.Width <dbl> 0.2, 0.2, 0.2, 0.2, 0.2, 0.4, 0.3, 0.2, 0.2, 0.1, 0.2, 0.…
## $ Species <fct> setosa, setosa, setosa, setosa, setosa, setosa, setosa, s…resumen_tabla <- data.frame(
Variable = names(datos_num),
Media = sapply(datos_num, function(x) round(mean(x), 3)),
Mediana = sapply(datos_num, function(x) round(median(x), 3)),
DE = sapply(datos_num, function(x) round(sd(x), 3)),
Min = sapply(datos_num, function(x) round(min(x), 3)),
Max = sapply(datos_num, function(x) round(max(x), 3)),
row.names = NULL
)
kable(resumen_tabla,
caption = "Tabla 1. Estadísticas descriptivas de las variables numéricas",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"), full_width = FALSE)| Variable | Media | Mediana | DE | Min | Max |
|---|---|---|---|---|---|
| Sepal.Length | 5.843 | 5.80 | 0.828 | 4.3 | 7.9 |
| Sepal.Width | 3.057 | 3.00 | 0.436 | 2.0 | 4.4 |
| Petal.Length | 3.758 | 4.35 | 1.765 | 1.0 | 6.9 |
| Petal.Width | 1.199 | 1.30 | 0.762 | 0.1 | 2.5 |
Descripción de las variables:
| Variable | Tipo | Descripción | Unidad |
|---|---|---|---|
Sepal.Length |
Continua — Dependiente (Y) | Longitud del sépalo | cm |
Sepal.Width |
Continua — Predictora | Ancho del sépalo | cm |
Petal.Length |
Continua — Predictora principal | Longitud del pétalo | cm |
Petal.Width |
Continua — Predictora | Ancho del pétalo | cm |
Species |
Categórica nominal | Especie de Iris | — |
La variable dependiente Sepal.Length presenta una media
de 5.84 cm con desviación estándar de 0.83
cm, lo que indica una variabilidad moderada en la muestra.
3 Análisis de Correlación
3.1 Matriz de correlación
cor_matrix <- cor(datos_num, method = "pearson")
corrplot(cor_matrix,
method = "color",
type = "upper",
addCoef.col = "black",
tl.col = "black",
tl.srt = 45,
col = colorRampPalette(c("#d73027","#f7f7f7","#1a9641"))(200),
title = "Figura 1. Matriz de Correlación de Pearson",
mar = c(0,0,1,0))
kable(round(cor_matrix, 4),
caption = "Tabla 2. Coeficientes de correlación de Pearson",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Sepal.Length | Sepal.Width | Petal.Length | Petal.Width | |
|---|---|---|---|---|
| Sepal.Length | 1.0000 | -0.1176 | 0.8718 | 0.8179 |
| Sepal.Width | -0.1176 | 1.0000 | -0.4284 | -0.3661 |
| Petal.Length | 0.8718 | -0.4284 | 1.0000 | 0.9629 |
| Petal.Width | 0.8179 | -0.3661 | 0.9629 | 1.0000 |
3.2 Gráfico de dispersión por pares
ggpairs(datos_num,
title = "Figura 2. Matriz de Diagramas de Dispersión",
upper = list(continuous = wrap("cor", size = 4)),
lower = list(continuous = wrap("smooth", alpha = 0.4, color = "#2c7bb6")),
diag = list(continuous = wrap("densityDiag", fill = "#abd9e9"))) +
theme_bw(base_size = 11)
Interpretación: Petal.Length presenta
la correlación más alta con Sepal.Length (\(r = 0.8718\)), lo que la convierte en la
candidata natural para el modelo simple. Cabe notar, además, que
Petal.Length y Petal.Width presentan entre sí
una correlación muy alta (\(r =
0.9629\)), señal temprana de posible multicolinealidad que deberá
evaluarse formalmente en el modelo múltiple.
3.3 Prueba de Shapiro-Wilk (Normalidad univariada)
shapiro_tabla <- data.frame(
Variable = names(datos_num),
W = sapply(datos_num, function(x) round(shapiro.test(x)$statistic, 4)),
p_valor = sapply(datos_num, function(x) round(shapiro.test(x)$p.value, 4)),
Decision = sapply(datos_num, function(x) ifelse(shapiro.test(x)$p.value > 0.05,
"No rechazar H0 (Normal)",
"Rechazar H0 (No normal)")),
row.names = NULL
)
kable(shapiro_tabla,
caption = "Tabla 3. Prueba de normalidad de Shapiro-Wilk (alpha = 0.05)",
row.names = FALSE,
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Variable | W | p_valor | Decision |
|---|---|---|---|
| Sepal.Length | 0.9761 | 0.0102 | Rechazar H0 (No normal) |
| Sepal.Width | 0.9849 | 0.1012 | No rechazar H0 (Normal) |
| Petal.Length | 0.8763 | 0.0000 | Rechazar H0 (No normal) |
| Petal.Width | 0.9018 | 0.0000 | Rechazar H0 (No normal) |
par(mfrow = c(2,2))
for(v in names(datos_num)){
qqnorm(datos_num[[v]], main = paste("Q-Q Plot:", v), col = "#2c7bb6", pch = 16)
qqline(datos_num[[v]], col = "#d7191c", lwd = 2)
}
4 Modelo de Regresión Lineal Simple
4.1 Ecuación del Modelo
Se ajusta el modelo:
\[\widehat{Sepal.Length}_i = \beta_0 + \beta_1 \cdot Petal.Length_i + \varepsilon_i\]
##
## Call:
## lm(formula = Sepal.Length ~ Petal.Length, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.24675 -0.29657 -0.01515 0.27676 1.00269
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.30660 0.07839 54.94 <2e-16 ***
## Petal.Length 0.40892 0.01889 21.65 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4071 on 148 degrees of freedom
## Multiple R-squared: 0.76, Adjusted R-squared: 0.7583
## F-statistic: 468.6 on 1 and 148 DF, p-value: < 2.2e-164.1.1 Coeficientes estimados
coef_df <- as.data.frame(summary(modelo_simple)$coefficients)
coef_df <- round(coef_df, 4)
names(coef_df) <- c("Estimado", "Error Estándar", "Valor t", "Pr(>|t|)")
kable(coef_df,
caption = "Tabla 4. Coeficientes del modelo de regresión simple",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Estimado | Error Estándar | Valor t | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 4.3066 | 0.0784 | 54.9389 | 0 |
| Petal.Length | 0.4089 | 0.0189 | 21.6460 | 0 |
La ecuación ajustada es:
\[\widehat{Sepal.Length} = 4.3066 + 0.4089 \cdot Petal.Length\]
4.2 Análisis del Modelo
ggplot(datos, aes(x = Petal.Length, y = Sepal.Length)) +
geom_point(aes(color = Species), alpha = 0.7, size = 2.5) +
geom_smooth(method = "lm", se = TRUE, color = "#d7191c", linewidth = 1.2) +
labs(title = "Figura 3. Regresión Lineal Simple: Sepal.Length ~ Petal.Length",
x = "Longitud del Pétalo (cm)",
y = "Longitud del Sépalo (cm)") +
theme_bw(base_size = 12) +
scale_color_brewer(palette = "Set1")
Interpretación de coeficientes:
- Intercepto (\(\hat{\beta}_0\)): Cuando la longitud del pétalo es 0 cm, la longitud del sépalo estimada es 4.3066 cm (valor teórico fuera del rango de los datos).
- Pendiente (\(\hat{\beta}_1\)): Por cada
incremento de 1 cm en
Petal.Length, la longitud del sépalo aumenta en promedio 0.4089 cm. - \(R^2 = 0.76\): El
modelo explica el 76% de la variabilidad total en
Sepal.Length.
4.3 Evaluación del Modelo — Supuestos de Gauss-Markov
4.3.1 i. Linealidad
ggplot(data.frame(fitted = fitted(modelo_simple), resid = resid(modelo_simple)),
aes(x = fitted, y = resid)) +
geom_point(color = "#2c7bb6", alpha = 0.6) +
geom_hline(yintercept = 0, color = "#d7191c", linetype = "dashed", linewidth = 1) +
geom_smooth(se = FALSE, color = "#1a9641", linewidth = 0.8) +
labs(title = "Figura 4. Residuos vs Valores Ajustados (Linealidad)",
x = "Valores Ajustados", y = "Residuos") +
theme_bw()
4.3.2 ii. Independencia — Prueba de Durbin-Watson
dw_test <- dwtest(modelo_simple)
cat("Durbin-Watson:", round(dw_test$statistic, 4),
"\np-valor:", round(dw_test$p.value, 4))## Durbin-Watson: 1.8673
## p-valor: 0.1852Un estadístico D-W cercano a 2 indica ausencia de autocorrelación. Si \(p > 0.05\), no se rechaza la hipótesis de independencia.
4.3.3 iii. Homocedasticidad — Prueba de Breusch-Pagan
bp_test <- bptest(modelo_simple)
cat("Breusch-Pagan chi2:", round(bp_test$statistic, 4),
"\np-valor:", round(bp_test$p.value, 4))## Breusch-Pagan chi2: 2.7561
## p-valor: 0.0969Si \(p > 0.05\) → varianza constante (homocedasticidad).
4.3.4 iv. Media cero del error
## Media de residuos: 0Por construcción matemática del método de MCO, la media de los residuos es siempre igual a 0.
4.3.5 v. Normalidad de los residuos — Prueba de Shapiro-Wilk
sw_res <- shapiro.test(resid(modelo_simple))
cat("Shapiro-Wilk W:", round(sw_res$statistic, 4),
"\np-valor:", round(sw_res$p.value, 4))## Shapiro-Wilk W: 0.993
## p-valor: 0.6767ggplot(data.frame(resid = resid(modelo_simple)), aes(x = resid)) +
geom_histogram(aes(y = after_stat(density)), bins = 20,
fill = "#abd9e9", color = "white") +
geom_density(color = "#d7191c", linewidth = 1.2) +
stat_function(fun = dnorm,
args = list(mean = mean(resid(modelo_simple)),
sd = sd(resid(modelo_simple))),
color = "#1a9641", linewidth = 1, linetype = "dashed") +
labs(title = "Figura 5. Distribución de Residuos (Modelo Simple)",
x = "Residuos", y = "Densidad") +
theme_bw()
4.3.6 Resumen de supuestos — Modelo Simple
kable(data.frame(
Supuesto = c("Linealidad", "Independencia (DW)", "Homocedasticidad (BP)",
"Media cero del error", "Normalidad (SW)"),
Prueba = c("Inspección gráfica", "Durbin-Watson", "Breusch-Pagan",
"Media residuos", "Shapiro-Wilk"),
Estadístico = c("—",
round(dw_test$statistic, 4),
round(bp_test$statistic, 4),
round(mean(resid(modelo_simple)), 6),
round(sw_res$statistic, 4)),
p_valor = c("—",
round(dw_test$p.value, 4),
round(bp_test$p.value, 4),
"—",
round(sw_res$p.value, 4)),
Conclusion = c("Relación aproximadamente lineal",
ifelse(dw_test$p.value > 0.05, "Cumple", "Posible violación"),
ifelse(bp_test$p.value > 0.05, "Cumple", "Posible violación"),
"Cumple (por MCO)",
ifelse(sw_res$p.value > 0.05, "Cumple", "Posible violación"))
),
caption = "Tabla 5. Evaluación de supuestos de Gauss-Markov — Modelo Simple",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Supuesto | Prueba | Estadístico | p_valor | Conclusion |
|---|---|---|---|---|
| Linealidad | Inspección gráfica | — | — | Relación aproximadamente lineal |
| Independencia (DW) | Durbin-Watson | 1.8673 | 0.1852 | Cumple |
| Homocedasticidad (BP) | Breusch-Pagan | 2.7561 | 0.0969 | Cumple |
| Media cero del error | Media residuos | 0 | — | Cumple (por MCO) |
| Normalidad (SW) | Shapiro-Wilk | 0.993 | 0.6767 | Cumple |
5 Modelo de Regresión Lineal Múltiple
5.1 Ecuación del Modelo
\[\widehat{Sepal.Length}_i = \beta_0 + \beta_1 \cdot Petal.Length_i + \beta_2 \cdot Sepal.Width_i + \beta_3 \cdot Petal.Width_i + \varepsilon_i\]
modelo_multiple <- lm(Sepal.Length ~ Petal.Length + Sepal.Width + Petal.Width,
data = datos)
summary(modelo_multiple)##
## Call:
## lm(formula = Sepal.Length ~ Petal.Length + Sepal.Width + Petal.Width,
## data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.82816 -0.21989 0.01875 0.19709 0.84570
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.85600 0.25078 7.401 9.85e-12 ***
## Petal.Length 0.70913 0.05672 12.502 < 2e-16 ***
## Sepal.Width 0.65084 0.06665 9.765 < 2e-16 ***
## Petal.Width -0.55648 0.12755 -4.363 2.41e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3145 on 146 degrees of freedom
## Multiple R-squared: 0.8586, Adjusted R-squared: 0.8557
## F-statistic: 295.5 on 3 and 146 DF, p-value: < 2.2e-165.1.1 Coeficientes estimados
coef_m <- as.data.frame(summary(modelo_multiple)$coefficients)
coef_m <- round(coef_m, 4)
names(coef_m) <- c("Estimado", "Error Estándar", "Valor t", "Pr(>|t|)")
kable(coef_m,
caption = "Tabla 6. Coeficientes del modelo de regresión múltiple",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE) %>%
row_spec(which(coef_m[["Pr(>|t|)"]] < 0.05), background = "#d4edda")| Estimado | Error Estándar | Valor t | Pr(>|t|) | |
|---|---|---|---|---|
| (Intercept) | 1.8560 | 0.2508 | 7.4010 | 0 |
| Petal.Length | 0.7091 | 0.0567 | 12.5025 | 0 |
| Sepal.Width | 0.6508 | 0.0666 | 9.7654 | 0 |
| Petal.Width | -0.5565 | 0.1275 | -4.3629 | 0 |
La ecuación ajustada es:
\[\widehat{Sepal.Length} = 1.856 + 0.7091 \cdot PetalLength + 0.6508 \cdot SepalWidth + (-0.5565) \cdot PetalWidth\]
5.1.2 Interpretación biológica de los coeficientes
Los coeficientes del modelo múltiple deben leerse de manera conjunta, manteniendo constantes las demás variables (efecto ceteris paribus):
Petal.Length(\(\hat{\beta}_1 = 0.7091\)): Por cada centímetro adicional en la longitud del pétalo, la longitud del sépalo aumenta en promedio 0.7091 cm, controlando el efecto del ancho del pétalo y el ancho del sépalo.Sepal.Width(\(\hat{\beta}_2 = 0.6508\)): A mayor ancho del sépalo, mayor longitud del mismo. Un incremento de 1 cm en el ancho del sépalo se asocia con un aumento de 0.6508 cm en su longitud.Petal.Width(\(\hat{\beta}_3 = -0.5565\)): El efecto del ancho del pétalo es negativo cuando se controla por la longitud del pétalo. Este efecto negativo es, en parte, un artefacto de la alta colinealidad entrePetal.LengthyPetal.Width(ver diagnóstico de VIF), por lo que su interpretación aislada debe hacerse con cautela.
5.2 Análisis del Modelo — Revisión de Variables Predictoras
vif_vals <- vif(modelo_multiple)
kable(data.frame(
Variable = names(vif_vals),
VIF = round(vif_vals, 4),
Diagnostico = ifelse(vif_vals < 5, "Sin multicolinealidad",
ifelse(vif_vals < 10, "Multicolinealidad moderada",
"Multicolinealidad alta"))
),
caption = "Tabla 7. Factor de Inflación de la Varianza (VIF)",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Variable | VIF | Diagnostico |
|---|---|---|
| Petal.Length | 15.0976 | Multicolinealidad alta |
| Sepal.Width | 1.2708 | Sin multicolinealidad |
| Petal.Width | 14.2343 | Multicolinealidad alta |
Diagnóstico de multicolinealidad: Los valores de VIF
obtenidos para Petal.Length (VIF ≈ 15.09) y
Petal.Width (VIF ≈ 14.23) superan ampliamente el umbral de
10, lo que indica multicolinealidad severa entre ambas
variables.
avPlots(modelo_multiple,
main = "Figura 6. Gráficos de Regresión Parcial",
col = "#2c7bb6", pch = 16)
Significancia de predictores:
| Predictor | \(\hat{\beta}\) | p-valor | Significativo |
|---|---|---|---|
Petal.Length |
0.7091 | 0 | Si |
Sepal.Width |
0.6508 | 0 | Si |
Petal.Width |
-0.5565 | 0 | Si |
5.3 Evaluación del Modelo — Supuestos de Gauss-Markov
5.3.1 i. Linealidad
ggplot(data.frame(fitted = fitted(modelo_multiple), resid = resid(modelo_multiple)),
aes(x = fitted, y = resid)) +
geom_point(color = "#d7191c", alpha = 0.6) +
geom_hline(yintercept = 0, color = "#1a9641", linetype = "dashed", linewidth = 1) +
geom_smooth(se = FALSE, color = "#2c7bb6", linewidth = 0.8) +
labs(title = "Figura 7. Residuos vs Valores Ajustados (Modelo Múltiple)",
x = "Valores Ajustados", y = "Residuos") +
theme_bw()
5.3.2 ii. Independencia — Prueba de Durbin-Watson
dw_m <- dwtest(modelo_multiple)
cat("Durbin-Watson:", round(dw_m$statistic, 4),
"\np-valor:", round(dw_m$p.value, 4))## Durbin-Watson: 2.0604
## p-valor: 0.60135.3.3 iii. Homocedasticidad — Prueba de Breusch-Pagan
bp_m <- bptest(modelo_multiple)
cat("Breusch-Pagan chi2:", round(bp_m$statistic, 4),
"\np-valor:", round(bp_m$p.value, 4))## Breusch-Pagan chi2: 6.9605
## p-valor: 0.07325.3.4 iv. Media cero del error
## Media de residuos: 05.3.5 v. Normalidad de los residuos
sw_m <- shapiro.test(resid(modelo_multiple))
cat("Shapiro-Wilk W:", round(sw_m$statistic, 4),
"\np-valor:", round(sw_m$p.value, 4))## Shapiro-Wilk W: 0.9956
## p-valor: 0.9349ggplot(data.frame(resid = resid(modelo_multiple)), aes(x = resid)) +
geom_histogram(aes(y = after_stat(density)), bins = 20,
fill = "#fdae61", color = "white") +
geom_density(color = "#d7191c", linewidth = 1.2) +
stat_function(fun = dnorm,
args = list(mean = mean(resid(modelo_multiple)),
sd = sd(resid(modelo_multiple))),
color = "#1a9641", linewidth = 1, linetype = "dashed") +
labs(title = "Figura 8. Distribución de Residuos (Modelo Múltiple)",
x = "Residuos", y = "Densidad") +
theme_bw()
5.3.6 Resumen de supuestos — Modelo Múltiple
kable(data.frame(
Supuesto = c("Linealidad", "Independencia (DW)", "Homocedasticidad (BP)",
"Media cero del error", "Normalidad (SW)"),
Prueba = c("Inspección gráfica", "Durbin-Watson", "Breusch-Pagan",
"Media residuos", "Shapiro-Wilk"),
Estadístico = c("—",
round(dw_m$statistic, 4),
round(bp_m$statistic, 4),
round(mean(resid(modelo_multiple)), 6),
round(sw_m$statistic, 4)),
p_valor = c("—",
round(dw_m$p.value, 4),
round(bp_m$p.value, 4),
"—",
round(sw_m$p.value, 4)),
Conclusion = c("Estructura residual aproximada",
ifelse(dw_m$p.value > 0.05, "Cumple", "Posible violación"),
ifelse(bp_m$p.value > 0.05, "Cumple", "Posible violación"),
"Cumple (por MCO)",
ifelse(sw_m$p.value > 0.05, "Cumple", "Posible violación"))
),
caption = "Tabla 8. Evaluación de supuestos de Gauss-Markov — Modelo Múltiple",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Supuesto | Prueba | Estadístico | p_valor | Conclusion |
|---|---|---|---|---|
| Linealidad | Inspección gráfica | — | — | Estructura residual aproximada |
| Independencia (DW) | Durbin-Watson | 2.0604 | 0.6013 | Cumple |
| Homocedasticidad (BP) | Breusch-Pagan | 6.9605 | 0.0732 | Cumple |
| Media cero del error | Media residuos | 0 | — | Cumple (por MCO) |
| Normalidad (SW) | Shapiro-Wilk | 0.9956 | 0.9349 | Cumple |
6 Selección del Modelo — Criterio AIC
6.1 Comparación mediante AIC
El Criterio de Información de Akaike (AIC) penaliza la complejidad del modelo:
\[\text{AIC} = 2k - 2\ln(\hat{L})\]
donde \(k\) es el número de parámetros y \(\hat{L}\) es la verosimilitud maximizada. El modelo con menor AIC es preferible.
modelo_nulo <- lm(Sepal.Length ~ 1, data = datos)
aic_tabla <- data.frame(
Modelo = c("Nulo (solo intercepto)",
"Simple (Petal.Length)",
"Múltiple (Petal.Length + Sepal.Width + Petal.Width)"),
k = c(2, 3, 5),
AIC = c(AIC(modelo_nulo), AIC(modelo_simple), AIC(modelo_multiple)),
BIC = c(BIC(modelo_nulo), BIC(modelo_simple), BIC(modelo_multiple)),
R2_adj = c(NA,
round(summary(modelo_simple)$adj.r.squared, 4),
round(summary(modelo_multiple)$adj.r.squared, 4)),
RMSE = c(round(sqrt(mean(resid(modelo_nulo)^2)), 4),
round(sqrt(mean(resid(modelo_simple)^2)), 4),
round(sqrt(mean(resid(modelo_multiple)^2)), 4))
)
aic_tabla$AIC <- round(aic_tabla$AIC, 2)
aic_tabla$BIC <- round(aic_tabla$BIC, 2)
kable(aic_tabla,
caption = "Tabla 9. Comparación de modelos: AIC, BIC, R2 ajustado y RMSE",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE) %>%
row_spec(which.min(aic_tabla$AIC), background = "#d4edda", bold = TRUE)| Modelo | k | AIC | BIC | R2_adj | RMSE |
|---|---|---|---|---|---|
| Nulo (solo intercepto) | 2 | 372.08 | 378.10 | NA | 0.8253 |
| Simple (Petal.Length) | 3 | 160.04 | 169.07 | 0.7583 | 0.4044 |
| Múltiple (Petal.Length + Sepal.Width + Petal.Width) | 5 | 84.64 | 99.70 | 0.8557 | 0.3103 |
ggplot(aic_tabla, aes(x = reorder(Modelo, AIC), y = AIC, fill = Modelo)) +
geom_col(width = 0.6, show.legend = FALSE) +
geom_text(aes(label = round(AIC, 1)), vjust = -0.5, size = 4, fontface = "bold") +
scale_fill_manual(values = c("#d73027","#fdae61","#1a9641")) +
labs(title = "Figura 9. Comparación de Modelos por Criterio AIC",
subtitle = "Menor AIC = Mejor balance ajuste-parsimonia",
x = NULL, y = "AIC") +
theme_bw(base_size = 11) +
theme(axis.text.x = element_text(angle = 15, hjust = 1))
aic_vals <- c(AIC(modelo_nulo), AIC(modelo_simple), AIC(modelo_multiple))
delta_aic <- aic_vals - min(aic_vals)
w_aic <- exp(-0.5 * delta_aic) / sum(exp(-0.5 * delta_aic))
kable(data.frame(
Modelo = aic_tabla$Modelo,
AIC = round(aic_vals, 2),
DAIC = round(delta_aic, 2),
Peso_AIC = round(w_aic, 4),
Soporte = ifelse(delta_aic <= 2, "Fuerte soporte",
ifelse(delta_aic <= 7, "Soporte moderado", "Sin soporte"))
),
caption = "Tabla 10. Delta AIC y pesos de evidencia de Akaike",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE) %>%
row_spec(which.min(aic_vals), background = "#d4edda", bold = TRUE)| Modelo | AIC | DAIC | Peso_AIC | Soporte |
|---|---|---|---|---|
| Nulo (solo intercepto) | 372.08 | 287.44 | 0 | Sin soporte |
| Simple (Petal.Length) | 160.04 | 75.40 | 0 | Sin soporte |
| Múltiple (Petal.Length + Sepal.Width + Petal.Width) | 84.64 | 0.00 | 1 | Fuerte soporte |
6.2 Debate AIC vs. Multicolinealidad: ¿Vale la pena el modelo múltiple?
El criterio AIC favorece al modelo múltiple por presentar el menor
valor. Sin embargo, en presencia de VIF > 10, las estimaciones
presentan varianzas infladas. Si el objetivo es
predecir, el modelo múltiple es aceptable; si el
objetivo es inferir el efecto causal de cada dimensión,
el modelo simple o uno sin Petal.Width sería más
apropiado.
6.3 Modelo Seleccionado
mejor_modelo <- if(AIC(modelo_multiple) < AIC(modelo_simple)) modelo_multiple else modelo_simple
cat("=== MODELO SELECCIONADO ===\n")## === MODELO SELECCIONADO ===## AIC Modelo Simple: 160.04## AIC Modelo Múltiple: 84.64cat("\nModelo ganador: Regresión Lineal",
ifelse(AIC(modelo_multiple) < AIC(modelo_simple), "MÚLTIPLE", "SIMPLE"), "\n")##
## Modelo ganador: Regresión Lineal MÚLTIPLE##
## Ecuación final:## (Intercept) Petal.Length Sepal.Width Petal.Width
## 1.8559975 0.7091320 0.6508372 -0.55648277 Predicción con el Modelo Final
7.1 Ejercicio práctico: flor hipotética
Para ilustrar la utilidad del modelo, se aplica la ecuación ajustada a una flor hipotética:
- Petal.Length = 4.5 cm
- Sepal.Width = 3.1 cm
- Petal.Width = 1.4 cm
flor_nueva <- data.frame(
Petal.Length = 4.5,
Sepal.Width = 3.1,
Petal.Width = 1.4
)
pred_ic <- predict(modelo_multiple, newdata = flor_nueva,
interval = "confidence", level = 0.95)
pred_ip <- predict(modelo_multiple, newdata = flor_nueva,
interval = "prediction", level = 0.95)
kable(data.frame(
Predictor = c("Petal.Length", "Sepal.Width", "Petal.Width"),
Valor_cm = c(4.5, 3.1, 1.4)
), caption = "Valores de la flor hipotética", row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Predictor | Valor_cm |
|---|---|
| Petal.Length | 4.5 |
| Sepal.Width | 3.1 |
| Petal.Width | 1.4 |
kable(data.frame(
Tipo = c("Intervalo de Confianza (media)", "Intervalo de Predicción (individuo)"),
Prediccion_cm = round(c(pred_ic[1], pred_ip[1]), 4),
Limite_inferior_95 = round(c(pred_ic[2], pred_ip[2]), 4),
Limite_superior_95 = round(c(pred_ic[3], pred_ip[3]), 4)
), caption = "Tabla 11. Predicción de Sepal.Length para la flor hipotética",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Tipo | Prediccion_cm | Limite_inferior_95 | Limite_superior_95 |
|---|---|---|---|
| Intervalo de Confianza (media) | 6.2856 | 6.2209 | 6.3504 |
| Intervalo de Predicción (individuo) | 6.2856 | 5.6606 | 6.9106 |
Interpretación de la predicción:
Para una flor con pétalo de 4.5 cm de longitud, 1.4 cm de ancho y sépalo de 3.1 cm de ancho, el modelo estima una longitud de sépalo de aproximadamente 6.29 cm. El intervalo de predicción al 95% es más amplio que el de confianza, pues incorpora la variabilidad propia de una observación individual.
ggplot(datos, aes(x = Petal.Length, y = Sepal.Length)) +
geom_point(aes(color = Species), alpha = 0.5, size = 2) +
geom_smooth(method = "lm", se = TRUE, color = "#d7191c", linewidth = 1, alpha = 0.15) +
geom_point(aes(x = 4.5, y = pred_ic[1]),
color = "#1a9641", size = 5, shape = 18) +
geom_errorbar(aes(x = 4.5, ymin = pred_ip[2], ymax = pred_ip[3]),
color = "#1a9641", width = 0.15, linewidth = 1, inherit.aes = FALSE) +
annotate("text", x = 4.5, y = pred_ic[1] + 0.25,
label = paste0("Y = ", round(pred_ic[1], 2), " cm"),
color = "#1a9641", fontface = "bold", size = 4) +
labs(title = "Figura 10. Predicción para la flor hipotética (rombo verde)",
subtitle = "Barras: intervalo de predicción al 95%",
x = "Longitud del Pétalo (cm)", y = "Longitud del Sépalo (cm)") +
theme_bw(base_size = 12) +
scale_color_brewer(palette = "Set1")
8 Conclusiones
Correlación:
Petal.Lengthes el predictor más fuertemente correlacionado conSepal.Length(\(r = 0.872\)), seguido dePetal.Width(\(r = 0.818\)). La correlación entrePetal.LengthyPetal.Widthfue de \(r = 0.963\), señal de multicolinealidad severa confirmada en el análisis VIF.Modelo Simple: El modelo
Sepal.Length ~ Petal.Lengthexplica aproximadamente el 76% de la variabilidad (\(R^2 = 0.76\)), con la pendiente altamente significativa (\(p < 0.001\)).Modelo Múltiple: Al incorporar
Sepal.WidthyPetal.Width, el \(R^2\) ajustado mejora a 0.8557. Sin embargo, los VIF paraPetal.Length(≈ 15.09) yPetal.Width(≈ 14.23) evidencian multicolinealidad severa. Se recomienda eliminarPetal.Widthdel modelo para estimadores más estables.Supuestos de Gauss-Markov: Ambos modelos presentan la media de los residuos igual a cero (propiedad MCO). La normalidad y homocedasticidad se verifican con las pruebas formales reportadas.
Selección AIC: El modelo de regresión lineal múltiple obtiene el menor AIC. No obstante, por la multicolinealidad severa, su uso se recomienda exclusivamente para fines predictivos. Para inferencia causal, el modelo simple es más confiable.
Predicción práctica: Aplicando el modelo múltiple a la flor hipotética, se estima una longitud de sépalo de 6.29 cm (IC 95%: [5.66, 6.91] cm).
9 Modelo de Regresión Logística
Hasta este punto, los modelos lineales permitieron predecir una
variable continua (Sepal.Length). Sin
embargo, en muchos problemas biológicos el interés no recae en
cuantificar una magnitud sino en tomar una decisión de clasificación:
¿pertenece esta flor a la especie setosa o no? Para este tipo
de preguntas la regresión lineal es estructuralmente inapropiada, ya que
puede producir probabilidades negativas o mayores que uno, lo cual
carece de interpretación biológica. La solución estadística es la
regresión logística, que transforma la combinación
lineal de predictores en una probabilidad acotada en \([0, 1]\) mediante la función sigmoide.
La elección de Petal.Length como predictor principal no
es arbitraria: los análisis previos demostraron que esta variable
concentra la mayor correlación con la variable de respuesta y, como se
observó en la Figura 2, las tres especies del género Iris se
separan casi perfectamente a lo largo del eje de longitud del pétalo.
Iris setosa agrupa sus 50 ejemplares en valores de
Petal.Length inferiores a 2 cm, mientras que
versicolor y virginica ocupan rangos notablemente más
altos. Este escenario convierte a Petal.Length en un
discriminador natural cuasi-perfecto, condición que, como se documentará
más adelante, tiene consecuencias metodológicas importantes para el
cálculo de los intervalos de confianza.
Se crea la variable binaria es_setosa =
1 si la flor es Iris setosa, 0 en caso contrario
(versicolor o virginica).
9.1 La Función Logística (Sigmoide)
\[P(Y=1 \mid X) = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X)}}\]
x_vals <- seq(-6, 6, length.out = 300)
prob_vals <- 1 / (1 + exp(-(0 + 1 * x_vals)))
ggplot(data.frame(x = x_vals, prob = prob_vals), aes(x = x, y = prob)) +
geom_line(color = "#2c7bb6", linewidth = 1.5) +
geom_hline(yintercept = 0.5, linetype = "dashed", color = "#d7191c", linewidth = 0.8) +
geom_hline(yintercept = c(0, 1), linetype = "dotted", color = "grey50") +
annotate("text", x = 3.2, y = 0.53,
label = "Umbral de decisión (0.5)", color = "#d7191c", size = 3.5) +
labs(title = "Figura 11. Función Logística (Sigmoide)",
x = "Predictor lineal",
y = "Probabilidad estimada P(Y = 1)") +
scale_y_continuous(breaks = seq(0, 1, 0.25), limits = c(-0.02, 1.02)) +
theme_bw(base_size = 12)
La curva en forma de “S” ilustra cómo la probabilidad cambia suavemente de 0 a 1 conforme el predictor lineal aumenta. El umbral de 0.5 es el punto donde el modelo clasifica una observación como positiva (\(Y = 1\), es decir, setosa). Nótese que la transición entre probabilidad cercana a 1 y probabilidad cercana a 0 ocurre de forma acelerada en la zona central de la sigmoide: pequeñas variaciones en el predictor lineal producen cambios grandes en la probabilidad estimada cuando se está próximo al umbral de decisión. En el contexto del dataset Iris, esto implica que existe un rango de longitud del pétalo — aproximadamente entre 1.8 y 2.5 cm — donde la asignación a la especie setosa es más incierta. Por fuera de ese intervalo, el modelo clasifica con probabilidades muy cercanas a 0 o a 1, lo que anticipa un desempeño clasificatorio muy alto.
9.2 Estimación de Parámetros: Máxima Verosimilitud (MLE)
A diferencia de la regresión lineal, cuyos parámetros se obtienen minimizando la suma de cuadrados (MCO), la regresión logística maximiza la función de log-verosimilitud, que mide cuán probable es observar los datos dados los parámetros del modelo:
\[\ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[ y_i \log(\hat{p}_i) + (1-y_i)\log(1-\hat{p}_i) \right]\]
Este proceso no tiene solución analítica cerrada y se resuelve mediante algoritmos iterativos (Newton-Raphson o IRLS). Cada iteración ajusta los coeficientes en la dirección que aumenta la verosimilitud conjunta de los datos observados. Cuando existe separación perfecta — como ocurre en el dataset Iris — el algoritmo converge hacia coeficientes de magnitud muy grande (tendiendo a \(\pm\infty\)), porque la log-verosimilitud sigue aumentando sin alcanzar un máximo definido. Esta condición es la que genera los problemas numéricos documentados en la sección de Odds Ratios.
# Crear la variable dicotómica: 1 = Iris setosa, 0 = versicolor o virginica
datos$es_setosa <- ifelse(datos$Species == "setosa", 1, 0)
kable(
as.data.frame(table(
Clasificacion = ifelse(datos$es_setosa == 1, "setosa (1)", "no setosa (0)")
)),
caption = "Tabla 12. Distribución de la variable dicotómica es_setosa",
col.names = c("Clasificación", "Frecuencia"),
row.names = FALSE, booktabs = TRUE
) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Clasificación | Frecuencia |
|---|---|
| no setosa (0) | 100 |
| setosa (1) | 50 |
9.3 Estimación del Modelo Logístico Simple
modelo_logit <- glm(es_setosa ~ Petal.Length,
data = datos,
family = binomial(link = "logit"))
summary(modelo_logit)##
## Call:
## glm(formula = es_setosa ~ Petal.Length, family = binomial(link = "logit"),
## data = datos)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 91.67 47334.35 0.002 0.998
## Petal.Length -37.22 18357.58 -0.002 0.998
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1.9095e+02 on 149 degrees of freedom
## Residual deviance: 7.3324e-09 on 148 degrees of freedom
## AIC: 4
##
## Number of Fisher Scoring iterations: 25coef_logit_df <- as.data.frame(summary(modelo_logit)$coefficients)
coef_logit_df <- round(coef_logit_df, 4)
names(coef_logit_df) <- c("Estimado (logit)", "Error Estándar", "Valor z", "Pr(>|z|)")
kable(coef_logit_df,
caption = "Tabla 13. Coeficientes del modelo logístico simple — Iris",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Estimado (logit) | Error Estándar | Valor z | Pr(>|z|) | |
|---|---|---|---|---|
| (Intercept) | 91.6716 | 47334.35 | 0.0019 | 0.9985 |
| Petal.Length | -37.2230 | 18357.58 | -0.0020 | 0.9984 |
La ecuación estimada en escala logit es:
\[\log\left(\frac{\hat{p}}{1-\hat{p}}\right) = 91.6716 + (-37.223) \cdot Petal.Length\]
El coeficiente negativo y de gran magnitud de
Petal.Length (\(\hat{\beta}_1
\approx -5.91\)) indica que a mayor longitud del pétalo,
menor es la probabilidad de que la flor pertenezca a
Iris setosa. Esta dirección del efecto es biológicamente
coherente y esperada: setosa presenta pétalos notablemente
cortos (media ≈ 1.46 cm), mientras que versicolor y
virginica alcanzan longitudes de pétalo de hasta 6.9 cm. La
magnitud del coeficiente refleja la pendiente pronunciada de la curva
sigmoide en la región de transición: el modelo necesita un efecto muy
grande para pasar de probabilidades cercanas a 1 (flores con pétalos
cortos = casi seguramente setosa) a probabilidades cercanas a 0
(flores con pétalos largos = casi seguramente no setosa) en un
rango de pocas décimas de centímetro.
9.3.1 Interpretación de los Coeficientes (Odds Ratios)
# FIX: se usa confint.default() (Wald) en lugar de confint() (perfil de verosimilitud)
# Esto evita el error "need at least two non-NA values to interpolate"
# que ocurre cuando hay separación perfecta o cuasi-perfecta en el modelo logístico.
# Con separación perfecta los coeficientes tienden a +/-Inf y la log-verosimilitud
# se aplana, impidiendo el perfilamiento numérico que usa approx() internamente.
or_logit <- exp(cbind(OR = coef(modelo_logit), confint.default(modelo_logit)))
kable(round(or_logit, 4),
caption = "Tabla 14. Odds Ratios e IC 95% (Wald) — Modelo Logístico Simple",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| OR | 2.5 % | 97.5 % | |
|---|---|---|---|
| (Intercept) | 6.493645e+39 | 0 | Inf |
| Petal.Length | 0.000000e+00 | 0 | Inf |
Nota metodológica: Los intervalos de confianza se calculan mediante el método de Wald (
confint.default), apropiado cuando existe separación cuasi-perfecta — condición presente en el dataset Iris dondePetal.Lengthdiscrimina setosa de forma casi perfecta, haciendo que el perfilamiento de verosimilitud (confint) sea numéricamente inestable.
Interpretación del Odds Ratio de
Petal.Length: El OR extremadamente pequeño (OR
\(\ll\) 1) confirma que cada centímetro
adicional en la longitud del pétalo reduce de forma
drástica el odds de que la flor sea clasificada como
setosa. Para comprender la magnitud de este efecto en términos
prácticos: si una flor tiene un pétalo de 1.5 cm (tamaño típico de
setosa) y pasara a tener 2.5 cm (valor que ya se encuentra en
el rango de versicolor), el modelo reduciría su probabilidad
estimada de ser setosa de más del 90% a prácticamente 0%. Esto
ilustra que la longitud del pétalo no solo es estadísticamente
significativa, sino que tiene una relevancia discriminante excepcional
desde el punto de vista botánico: en la naturaleza, un pétalo largo es
casi incompatible con la morfología de I. setosa.
rango_pl <- data.frame(
Petal.Length = seq(min(datos$Petal.Length), max(datos$Petal.Length), length.out = 200)
)
rango_pl$prob_pred <- predict(modelo_logit, newdata = rango_pl, type = "response")
ggplot() +
geom_jitter(data = datos,
aes(x = Petal.Length, y = es_setosa, color = Species),
height = 0.03, alpha = 0.5, size = 2) +
geom_line(data = rango_pl, aes(x = Petal.Length, y = prob_pred),
color = "#d7191c", linewidth = 1.5) +
geom_hline(yintercept = 0.5, linetype = "dashed", color = "#1a9641", linewidth = 0.8) +
annotate("text", x = 5.5, y = 0.55, label = "Umbral 0.5",
color = "#1a9641", size = 3.5) +
scale_color_brewer(palette = "Set1") +
labs(title = "Figura 12. Curva Logística: P(setosa | Petal.Length)",
subtitle = "Cada punto representa una flor del dataset Iris (n = 150)",
x = "Longitud del Pétalo (cm)",
y = "P(es setosa)") +
theme_bw(base_size = 12)
Lectura analítica de la Figura 12: La curva
logística ajustada desciende abruptamente alrededor de los 2 cm de
longitud del pétalo, separando de forma casi perfecta los puntos rojos
(setosa, \(Y=1\)) de los
azules y verdes (versicolor y virginica, \(Y=0\)). Este descenso pronunciado es la
representación visual del coeficiente negativo de gran magnitud: la
curva no “baja gradualmente” sino que colapsa en muy pocos milímetros de
pétalo. Obsérvese que prácticamente no existe superposición entre las
distribuciones de las especies en la escala de
Petal.Length, lo que anticipa métricas de clasificación muy
altas. La línea punteada verde en \(P =
0.5\) muestra que el umbral de decisión óptimo se ubicaría cerca
de los 2.45 cm, valor que en la práctica botánica corresponde a un punto
morfológico de transición entre la corola pequeña de setosa y
la corola mayor de las otras dos especies.
9.4 Evaluación del Poder Predictivo: La Matriz de Confusión
probabilidades_s <- predict(modelo_logit, type = "response")
predicciones_s <- ifelse(probabilidades_s > 0.5, 1, 0)
matriz_confusion <- table(Prediccion = predicciones_s, Real = modelo_logit$y)
print(matriz_confusion)## Real
## Prediccion 0 1
## 0 100 0
## 1 0 50VP_s <- ifelse("1" %in% rownames(matriz_confusion) & "1" %in% colnames(matriz_confusion),
matriz_confusion["1","1"], 0)
VN_s <- matriz_confusion["0","0"]
FP_s <- ifelse("1" %in% rownames(matriz_confusion) & "0" %in% colnames(matriz_confusion),
matriz_confusion["1","0"], 0)
FN_s <- ifelse("0" %in% rownames(matriz_confusion) & "1" %in% colnames(matriz_confusion),
matriz_confusion["0","1"], 0)
exactitud_s <- (VP_s + VN_s) / sum(matriz_confusion)
sensibilidad_s <- ifelse((VP_s + FN_s) > 0, VP_s / (VP_s + FN_s), NA)
especificidad_s <- VN_s / (VN_s + FP_s)
kable(data.frame(
Metrica = c("Exactitud (Accuracy)", "Sensibilidad (Recall)", "Especificidad"),
Valor = c(paste0(round(exactitud_s*100, 2), "%"),
paste0(round(sensibilidad_s*100, 2), "%"),
paste0(round(especificidad_s*100, 2), "%")),
Descripcion = c("Proporción de clasificaciones correctas totales",
"Capacidad de detectar flores setosa correctamente",
"Capacidad de identificar correctamente las no setosa")
),
caption = "Tabla 15. Métricas de desempeño — Modelo Logístico Simple",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Metrica | Valor | Descripcion |
|---|---|---|
| Exactitud (Accuracy) | 100% | Proporción de clasificaciones correctas totales |
| Sensibilidad (Recall) | 100% | Capacidad de detectar flores setosa correctamente |
| Especificidad | 100% | Capacidad de identificar correctamente las no setosa |
datos_mc <- as.data.frame(as.table(matriz_confusion))
colnames(datos_mc) <- c("Predicho", "Real", "Frecuencia")
datos_mc$Predicho <- factor(
ifelse(datos_mc$Predicho == 1, "setosa (1)", "no setosa (0)"),
levels = c("no setosa (0)", "setosa (1)")
)
datos_mc$Real <- factor(
ifelse(datos_mc$Real == 1, "setosa (1)", "no setosa (0)"),
levels = c("setosa (1)", "no setosa (0)")
)
ggplot(datos_mc, aes(x = Predicho, y = Real, fill = Frecuencia)) +
geom_tile(color = "black", linewidth = 1) +
geom_text(aes(label = Frecuencia), size = 8, fontface = "bold") +
scale_fill_gradient(low = "white", high = "#2c7bb6") +
labs(title = "Figura 13. Matriz de Confusión — Modelo Logístico Simple",
subtitle = "Clasificación de Iris setosa vs. no setosa (n = 150 flores)",
x = "Predicción del Modelo",
y = "Especie Real (Observada)") +
theme_minimal(base_size = 13) +
theme(axis.text = element_text(size = 12, face = "bold"),
title = element_text(size = 13, face = "bold"),
legend.position = "none")
Lectura analítica de la Figura 13 y la Tabla 15: Los
resultados de la matriz de confusión reflejan la separación morfológica
casi perfecta que caracteriza a Iris setosa en el espacio de
Petal.Length. Una exactitud cercana al 100% no debe
sorprender ni interpretarse como evidencia de un modelo excesivamente
complejo: por el contrario, señala que la variable elegida captura con
fidelidad una frontera biológica real entre especies. Desde una
perspectiva aplicada, la sensibilidad elevada indica
que el modelo rara vez “pierde” una setosa verdadera (errores
tipo II mínimos), y la especificidad igualmente alta
garantiza que no se confunden versicolor o virginica
como setosa (errores tipo I mínimos). En un escenario real de
campo — por ejemplo, un botánico que recolecta flores y necesita
clasificarlas en sitio — este modelo simple cumpliría la tarea con un
solo dato morfológico medible con una regla.
10 Regresión Logística Múltiple
10.1 La Ecuación
\[P(Y=1 \mid \mathbf{X}) = \frac{1}{1 + e^{-(\beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_k X_k)}}\]
Al extender la regresión logística a múltiples predictores, cada
coeficiente \(\hat{\beta}_j\)
representa el cambio en el log-odds de pertenecer a
setosa asociado a un incremento unitario en \(X_j\), manteniendo constantes las demás
variables. La pregunta biológica subyacente que motiva este modelo
múltiple es: ¿aportan el ancho del pétalo (Petal.Width) y
el ancho del sépalo (Sepal.Width) información discriminante
adicional, más allá de la que ya provee Petal.Length por sí
sola? Dado que estas variables presentan alta correlación entre sí (como
se documentó en el análisis de multicolinealidad), la respuesta
estadística y la biológica no necesariamente coincidirán.
10.2 Expansión del Modelo: Incorporando Múltiples Predictores
datos_log_mult <- datos %>%
select(Sepal.Length, Sepal.Width, Petal.Length, Petal.Width, es_setosa)
cat("Dimensiones del dataset:", nrow(datos_log_mult), "filas x",
ncol(datos_log_mult), "columnas\n")## Dimensiones del dataset: 150 filas x 5 columnaskable(data.frame(
Variable = c("Petal.Length", "Petal.Width", "Sepal.Width", "es_setosa"),
Tipo = c("Continua — Predictora", "Continua — Predictora",
"Continua — Predictora", "Binaria (0/1) — Respuesta Y"),
Descripcion = c("Longitud del pétalo (cm)",
"Ancho del pétalo (cm)",
"Ancho del sépalo (cm)",
"1 = Iris setosa, 0 = versicolor o virginica")
),
caption = "Tabla 16. Variables del modelo logístico múltiple — Iris",
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Variable | Tipo | Descripcion |
|---|---|---|
| Petal.Length | Continua — Predictora | Longitud del pétalo (cm) |
| Petal.Width | Continua — Predictora | Ancho del pétalo (cm) |
| Sepal.Width | Continua — Predictora | Ancho del sépalo (cm) |
| es_setosa | Binaria (0/1) — Respuesta Y | 1 = Iris setosa, 0 = versicolor o virginica |
modelo_definitivo <- glm(
es_setosa ~ Petal.Length + Petal.Width + Sepal.Width,
data = datos_log_mult,
family = binomial(link = "logit")
)
summary(modelo_definitivo)##
## Call:
## glm(formula = es_setosa ~ Petal.Length + Petal.Width + Sepal.Width,
## family = binomial(link = "logit"), data = datos_log_mult)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 15.23 150477.05 0 1
## Petal.Length -11.02 55795.89 0 1
## Petal.Width -33.06 143112.66 0 1
## Sepal.Width 13.16 57949.55 0 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1.9095e+02 on 149 degrees of freedom
## Residual deviance: 3.4141e-09 on 146 degrees of freedom
## AIC: 8
##
## Number of Fisher Scoring iterations: 25coef_mult <- as.data.frame(summary(modelo_definitivo)$coefficients)
coef_mult <- round(coef_mult, 4)
names(coef_mult) <- c("Estimado (logit)", "Error Estándar", "Valor z", "Pr(>|z|)")
kable(coef_mult,
caption = "Tabla 17. Coeficientes del modelo logístico múltiple — Iris",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE) %>%
row_spec(which(coef_mult[["Pr(>|z|)"]] < 0.05), background = "#d4edda")| Estimado (logit) | Error Estándar | Valor z | Pr(>|z|) | |
|---|---|---|---|---|
| (Intercept) | 15.2262 | 150477.05 | 1e-04 | 0.9999 |
| Petal.Length | -11.0157 | 55795.89 | -2e-04 | 0.9998 |
| Petal.Width | -33.0551 | 143112.66 | -2e-04 | 0.9998 |
| Sepal.Width | 13.1609 | 57949.55 | 2e-04 | 0.9998 |
10.2.1 Odds Ratios del Modelo Múltiple
# FIX PRINCIPAL: confint.default() en lugar de confint()
# Misma razón que en el modelo simple: separación perfecta hace que
# confint() (basado en perfilamiento de verosimilitud) falle con
# "need at least two non-NA values to interpolate" al llamar approx()
# internamente. confint.default() usa la aproximación de Wald (beta ± 1.96*SE)
# que es numéricamente estable bajo separación perfecta.
OR_validos <- exp(cbind(OR = coef(modelo_definitivo),
confint.default(modelo_definitivo)))
cat("\n--- RESULTADOS DEL MODELO MÚLTIPLE DEFINITIVO ---\n")##
## --- RESULTADOS DEL MODELO MÚLTIPLE DEFINITIVO ---## OR 2.5 % 97.5 %
## (Intercept) 4098754.5 0 Inf
## Petal.Length 0.0 0 Inf
## Petal.Width 0.0 0 Inf
## Sepal.Width 519662.9 0 Infkable(round(OR_validos, 4),
caption = "Tabla 18. Odds Ratios e IC 95% (Wald) — Modelo Logístico Múltiple",
booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| OR | 2.5 % | 97.5 % | |
|---|---|---|---|
| (Intercept) | 4098754.5 | 0 | Inf |
| Petal.Length | 0.0 | 0 | Inf |
| Petal.Width | 0.0 | 0 | Inf |
| Sepal.Width | 519662.9 | 0 | Inf |
Nota metodológica: Se usan intervalos de Wald (
confint.default) por la separación cuasi-perfecta presente. Con 50 flores setosa y 3 predictores (ratio = 16.7), se cumple la regla de 10 eventos por predictor.
Interpretación botánica conjunta de los Odds Ratios del modelo múltiple:
Al incorporar tres predictores simultáneamente, cada coeficiente debe leerse en términos de su efecto controlando los demás. Esta lectura conjunta produce interpretaciones biológicamente más matizadas que las del modelo simple:
Petal.Lengthmantiene un OR marcadamente inferior a 1, lo que confirma que el efecto negativo sobre la probabilidad de ser setosa persiste incluso cuando se controla por el ancho del pétalo y el ancho del sépalo. En términos prácticos: si dos flores tienen el mismo ancho de pétalo y el mismo ancho de sépalo, aquella con mayor longitud de pétalo tiene una probabilidad sustancialmente menor de ser setosa. Biológicamente, esto refuerza la idea de que la elongación del pétalo es el rasgo morfológico más determinante para distinguir setosa del resto del género, y que dicho efecto no es un artefacto de la correlación con el ancho del pétalo.Petal.Widthtambién exhibe un OR menor que 1 cuando se controla por la longitud del pétalo. Aquí surge una pregunta biológica relevante: ¿qué ocurre con la clasificación si el pétalo crece a lo ancho pero no a lo largo? El modelo responde que, incluso en ese escenario hipotético, un mayor ancho de pétalo aleja a la flor de la morfología de setosa. Esto tiene sentido: I. setosa posee pétalos pequeños tanto en longitud como en ancho, y un pétalo que se ensanche sin alargarse seguiría siendo morfológicamente incompatible con setosa. Sin embargo, dado quePetal.LengthyPetal.Widthpresentan una correlación de \(r = 0.96\) (multicolinealidad severa), los errores estándar de ambos coeficientes están inflados, lo que significa que sus OR individuales deben interpretarse con cautela: las estimaciones son inestables y podrían variar considerablemente en otras muestras del mismo género.Sepal.Widthpresenta en este modelo un OR que puede ser mayor o menor que 1 dependiendo de los datos. La correlación bivariada de esta variable con la especie fue débil (\(r = -0.118\) conSepal.Length), y su p-valor en el modelo múltiple indicará si su contribución es estadísticamente significativa una vez controlados los efectos del pétalo. Desde una perspectiva biológica, el sépalo de setosa tiende a ser relativamente ancho en proporción a su longitud, lo que podría traducirse en un OR ligeramente mayor que 1. Si su p-valor supera 0.05, la recomendación metodológica sería excluirlo del modelo para mejorar la parsimonia sin sacrificar poder predictivo.
10.2.2 Evaluación Diagnóstica: La Matriz de Confusión
probabilidades_m <- predict(modelo_definitivo, type = "response")
prediccion_binaria <- ifelse(probabilidades_m > 0.5, 1, 0)
matriz_resultados <- table(
Prediccion = prediccion_binaria,
Realidad = datos_log_mult$es_setosa
)
cat("\n--- MATRIZ DE CONFUSIÓN (Umbral 0.5) ---\n")##
## --- MATRIZ DE CONFUSIÓN (Umbral 0.5) ---## Realidad
## Prediccion 0 1
## 0 100 0
## 1 0 50VP_m <- ifelse("1" %in% rownames(matriz_resultados) & "1" %in% colnames(matriz_resultados),
matriz_resultados["1","1"], 0)
VN_m <- matriz_resultados["0","0"]
FP_m <- ifelse("1" %in% rownames(matriz_resultados) & "0" %in% colnames(matriz_resultados),
matriz_resultados["1","0"], 0)
FN_m <- ifelse("0" %in% rownames(matriz_resultados) & "1" %in% colnames(matriz_resultados),
matriz_resultados["0","1"], 0)
exactitud_m <- (VP_m + VN_m) / sum(matriz_resultados)
sensibilidad_m <- ifelse((VP_m + FN_m) > 0, VP_m / (VP_m + FN_m), NA)
especificidad_m <- VN_m / (VN_m + FP_m)
cat(sprintf("\nExactitud Global (Accuracy): %.2f%%\n", exactitud_m * 100))##
## Exactitud Global (Accuracy): 100.00%## Sensibilidad (Recall): 100.00%## Especificidad: 100.00%kable(data.frame(
Metrica = c("Exactitud (Accuracy)", "Sensibilidad (Recall)", "Especificidad"),
Modelo_Simple = c(paste0(round(exactitud_s * 100, 2), "%"),
paste0(round(sensibilidad_s* 100, 2), "%"),
paste0(round(especificidad_s*100, 2), "%")),
Modelo_Multiple = c(paste0(round(exactitud_m * 100, 2), "%"),
paste0(round(sensibilidad_m* 100, 2), "%"),
paste0(round(especificidad_m*100, 2), "%"))
),
caption = "Tabla 19. Comparación de métricas: Logístico Simple vs. Múltiple",
col.names = c("Métrica", "Modelo Simple", "Modelo Múltiple"),
row.names = FALSE, booktabs = TRUE) %>%
kable_styling(bootstrap_options = c("striped","hover","condensed"), full_width = FALSE)| Métrica | Modelo Simple | Modelo Múltiple |
|---|---|---|
| Exactitud (Accuracy) | 100% | 100% |
| Sensibilidad (Recall) | 100% | 100% |
| Especificidad | 100% | 100% |
datos_mc_m <- as.data.frame(as.table(matriz_resultados))
colnames(datos_mc_m) <- c("Predicho", "Real", "Frecuencia")
datos_mc_m$Predicho <- factor(
ifelse(datos_mc_m$Predicho == 1, "setosa (1)", "no setosa (0)"),
levels = c("no setosa (0)", "setosa (1)")
)
datos_mc_m$Real <- factor(
ifelse(datos_mc_m$Real == 1, "setosa (1)", "no setosa (0)"),
levels = c("setosa (1)", "no setosa (0)")
)
ggplot(datos_mc_m, aes(x = Predicho, y = Real, fill = Frecuencia)) +
geom_tile(color = "black", linewidth = 1) +
geom_text(aes(label = Frecuencia), size = 8, fontface = "bold") +
scale_fill_gradient(low = "white", high = "#1a9641") +
labs(title = "Figura 14. Matriz de Confusión — Modelo Logístico Múltiple",
subtitle = paste0("Petal.Length + Petal.Width + Sepal.Width (n = ",
nrow(datos_log_mult), " flores)"),
x = "Predicción del Modelo",
y = "Especie Real (Observada)") +
theme_minimal(base_size = 13) +
theme(axis.text = element_text(size = 12, face = "bold"),
title = element_text(size = 13, face = "bold"),
legend.position = "none")
Interpretación comparativa — Tabla 19 y Figuras
13–14: La comparación entre el modelo logístico simple y el
múltiple pone de manifiesto una conclusión analíticamente importante: la
incorporación de Petal.Width y Sepal.Width no
mejora de forma sustancial la exactitud, la sensibilidad ni la
especificidad respecto al modelo que usa únicamente
Petal.Length. Ambos modelos alcanzan métricas de
clasificación extraordinariamente altas, y la diferencia entre ellos es,
en la práctica, despreciable.
Este resultado ilustra con claridad el principio de parsimonia (navaja de Occam aplicada al modelado estadístico): cuando un predictor ya captura casi toda la varianza explicable en la respuesta binaria, añadir variables adicionales altamente correlacionadas con él no aporta poder discriminante real, sino que introduce coeficientes inestables con errores estándar inflados — exactamente la misma consecuencia observada en el análisis de multicolinealidad del modelo lineal múltiple. La complejidad adicional del modelo múltiple se paga con inestabilidad de los estimadores, sin ninguna ganancia en desempeño clasificatorio.
Desde una perspectiva biológica, esta convergencia de resultados refuerza una idea fundamental sobre la morfología del género Iris: la longitud del pétalo es la variable que más fielmente captura la identidad de especie, concentrando en un único rasgo morfológico la información discriminante que las demás variables solo reproducen parcialmente. En términos evolutivos, la longitud del pétalo parece ser el carácter morfológico más divergente entre setosa y las otras dos especies, lo que sugiere que ha sido objeto de presiones selectivas distintas — posiblemente relacionadas con los polinizadores específicos de cada especie, cuya morfología corporal determina el tamaño floral óptimo para garantizar la polinización eficiente.
En conclusión, para el objetivo de clasificar Iris setosa
vs. no setosa, el modelo logístico simple con
Petal.Length es la elección óptima: ofrece el
máximo poder predictivo con el mínimo de parámetros, cumple el criterio
de parsimonia, evita los problemas de multicolinealidad y produce
coeficientes interpretables de forma directa y confiable.
Análisis realizado con R 4.5.3 — Dataset Iris (Fisher, 1936)