Facultad de Matemáticas · Magíster en Inteligencia Artificial

library(ISLR)
library(MASS)
library(mlbench)
library(car)
library(pROC)
library(corrplot)
library(ggplot2)
library(dplyr)
library(tidyr)
library(knitr)
library(kableExtra)
library(plotly)
library(DT)
library(lmtest)

1 Pregunta 1 — Datos `Smarket`

1.1 a) Análisis descriptivo

Se trabaja con el dataset Smarket (ISLR), que contiene 1250 registros diarios del índice S&P 500 entre 2001 y 2005. La variable de respuesta es \(Y = \texttt{Today}\).

data(Smarket)

1.1.1 Estructura

str(Smarket)

#> 'data.frame':    1250 obs. of  9 variables:
#>  $ Year     : num  2001 2001 2001 2001 2001 ...
#>  $ Lag1     : num  0.381 0.959 1.032 -0.623 0.614 ...
#>  $ Lag2     : num  -0.192 0.381 0.959 1.032 -0.623 ...
#>  $ Lag3     : num  -2.624 -0.192 0.381 0.959 1.032 ...
#>  $ Lag4     : num  -1.055 -2.624 -0.192 0.381 0.959 ...
#>  $ Lag5     : num  5.01 -1.055 -2.624 -0.192 0.381 ...
#>  $ Volume   : num  1.19 1.3 1.41 1.28 1.21 ...
#>  $ Today    : num  0.959 1.032 -0.623 0.614 0.213 ...
#>  $ Direction: Factor w/ 2 levels "Down","Up": 2 2 1 2 2 2 1 2 2 2 ...

1.1.2 Estadísticos descriptivos

Smarket %>%
  select(-Direction) %>%
  pivot_longer(everything(), names_to = "Variable", values_to = "Valor") %>%
  group_by(Variable) %>%
  summarise(
    n       = n(),
    Media   = round(mean(Valor), 4),
    DE      = round(sd(Valor), 4),
    Min     = round(min(Valor), 4),
    Q1      = round(quantile(Valor, 0.25), 4),
    Mediana = round(median(Valor), 4),
    Q3      = round(quantile(Valor, 0.75), 4),
    Max     = round(max(Valor), 4)
  ) %>%
  kable(caption = "Resumen estadístico de las variables numéricas") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Resumen estadístico de las variables numéricas
Variable	n	Media	DE	Min	Q1	Mediana	Q3	Max
Lag1	1250	0.0038	1.1363	-4.9220	-0.6395	0.0390	0.5968	5.7330
Lag2	1250	0.0039	1.1363	-4.9220	-0.6395	0.0390	0.5968	5.7330
Lag3	1250	0.0017	1.1387	-4.9220	-0.6400	0.0385	0.5968	5.7330
Lag4	1250	0.0016	1.1388	-4.9220	-0.6400	0.0385	0.5968	5.7330
Lag5	1250	0.0056	1.1476	-4.9220	-0.6400	0.0385	0.5970	5.7330
Today	1250	0.0031	1.1363	-4.9220	-0.6395	0.0385	0.5968	5.7330
Volume	1250	1.4783	0.3604	0.3561	1.2574	1.4230	1.6417	3.1525
Year	1250	2003.0160	1.4090	2001.0000	2002.0000	2003.0000	2004.0000	2005.0000

1.1.3 Histogramas

Smarket %>%
  select(-Direction) %>%
  pivot_longer(everything(), names_to = "variable", values_to = "valor") %>%
  ggplot(aes(x = valor)) +
  geom_histogram(aes(y = after_stat(density)), bins = 30,
                 fill = "#004B93", alpha = 0.65) +
  geom_density(color = "#5A6C7D", linewidth = 0.8) +
  facet_wrap(~ variable, scales = "free", ncol = 3) +
  labs(x = "Valor", y = "Densidad") +
  theme_minimal(base_size = 11)

1.1.4 Correlaciones

cor_matrix <- cor(Smarket %>% select(-Direction))
corrplot(cor_matrix, method = "color", type = "upper",
         tl.col = "#002D5A", tl.srt = 45, addCoef.col = "#002D5A",
         number.cex = 0.7,
         col = colorRampPalette(c("#8DA0B0", "white", "#004B93"))(200))

1.1.5 Boxplots por Direction

Smarket %>%
  pivot_longer(cols = c(Lag1:Lag5, Volume),
               names_to = "variable", values_to = "valor") %>%
  ggplot(aes(x = Direction, y = valor, fill = Direction)) +
  geom_boxplot(alpha = 0.65, outlier.alpha = 0.4) +
  facet_wrap(~ variable, scales = "free_y", ncol = 3) +
  scale_fill_manual(values = c("Down" = "#5A6C7D", "Up" = "#004B93")) +
  theme_minimal(base_size = 11) +
  theme(legend.position = "bottom")

Las variables Lag1–Lag5 no muestran diferencias apreciables entre días de subida y bajada. Lo más notable es la correlación entre Volume y Year (≈ 0.54), lo cual es esperable dado el crecimiento secular del volumen de operaciones.

1.2 b) Selección de modelo lineal para `Today`

Se ajusta primero el modelo completo con todas las variables (excluyendo Direction) y luego se aplica selección stepwise bidireccional con criterio AIC.

modelo_full <- lm(Today ~ . - Direction, data = Smarket)
summary(modelo_full)

#> 
#> Call:
#> lm(formula = Today ~ . - Direction, data = Smarket)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -4.8923 -0.6573  0.0176  0.5904  5.6333 
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -58.170527  54.726633  -1.063    0.288
#> Year          0.029057   0.027365   1.062    0.289
#> Lag1         -0.027620   0.028367  -0.974    0.330
#> Lag2         -0.012707   0.028446  -0.447    0.655
#> Lag3         -0.005425   0.028397  -0.191    0.849
#> Lag4         -0.009786   0.028413  -0.344    0.731
#> Lag5         -0.036272   0.028118  -1.290    0.197
#> Volume       -0.018156   0.107157  -0.169    0.865
#> 
#> Residual standard error: 1.138 on 1242 degrees of freedom
#> Multiple R-squared:  0.003211,   Adjusted R-squared:  -0.002407 
#> F-statistic: 0.5715 on 7 and 1242 DF,  p-value: 0.7795

modelo_step <- step(modelo_full, direction = "both", trace = 0)
summary(modelo_step)

#> 
#> Call:
#> lm(formula = Today ~ 1, data = Smarket)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -4.9251 -0.6426  0.0354  0.5936  5.7299 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept) 0.003138   0.032140   0.098    0.922
#> 
#> Residual standard error: 1.136 on 1249 degrees of freedom

anova(modelo_step) %>%
  kable(digits = 4, caption = "ANOVA del modelo seleccionado") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

ANOVA del modelo seleccionado
	Df	Sum Sq	Mean Sq	F value	Pr(>F)
Residuals	1249	1612.778	1.2913	NA	NA

Se utilizó selección stepwise bidireccional con AIC. El modelo resultante y sus indicadores se muestran arriba.

1.3 c) Regresión logística para `Direction`

Ahora se considera \(Y = \texttt{Direction}\) como variable respuesta y se ajusta un modelo logístico con los rezagos Lag1–Lag5 como predictores.

glm_full <- glm(Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5,
                data = Smarket, family = binomial)
summary(glm_full)

#> 
#> Call:
#> glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5, family = binomial, 
#>     data = Smarket)
#> 
#> Coefficients:
#>              Estimate Std. Error z value Pr(>|z|)
#> (Intercept)  0.074163   0.056674   1.309    0.191
#> Lag1        -0.071325   0.050104  -1.424    0.155
#> Lag2        -0.044136   0.050025  -0.882    0.378
#> Lag3         0.009229   0.049879   0.185    0.853
#> Lag4         0.007211   0.049898   0.145    0.885
#> Lag5         0.009311   0.049490   0.188    0.851
#> 
#> (Dispersion parameter for binomial family taken to be 1)
#> 
#>     Null deviance: 1731.2  on 1249  degrees of freedom
#> Residual deviance: 1728.3  on 1244  degrees of freedom
#> AIC: 1740.3
#> 
#> Number of Fisher Scoring iterations: 3

glm_step <- step(glm_full, direction = "both", trace = 0)
summary(glm_step)

#> 
#> Call:
#> glm(formula = Direction ~ 1, family = binomial, data = Smarket)
#> 
#> Coefficients:
#>             Estimate Std. Error z value Pr(>|z|)
#> (Intercept)  0.07363    0.05661   1.301    0.193
#> 
#> (Dispersion parameter for binomial family taken to be 1)
#> 
#>     Null deviance: 1731.2  on 1249  degrees of freedom
#> Residual deviance: 1731.2  on 1249  degrees of freedom
#> AIC: 1733.2
#> 
#> Number of Fisher Scoring iterations: 3

glm_best <- if (length(coef(glm_step)) <= 1) glm_full else glm_step
cat("Variables en el modelo final:", length(coef(glm_best)) - 1, "\n")

#> Variables en el modelo final: 5

Se aplicó stepwise con AIC. Dado que los retornos pasados tienen escaso poder predictivo sobre la dirección futura (resultado consistente con la hipótesis de mercados eficientes), el procedimiento puede reducir el modelo al intercepto. En tal caso, se reportan los odds ratios del modelo completo para facilitar la interpretación.

1.4 d) Odds Ratios

ci_raw <- confint(glm_best)
ci_matrix <- matrix(ci_raw, ncol = 2)
rownames(ci_matrix) <- names(coef(glm_best))

data.frame(
  Variable    = names(coef(glm_best)),
  Log_odds    = round(coef(glm_best), 4),
  OR          = round(exp(coef(glm_best)), 4),
  IC_2.5      = round(exp(ci_matrix[, 1]), 4),
  IC_97.5     = round(exp(ci_matrix[, 2]), 4)
) %>%
  kable(col.names = c("Variable", "Coef (log-odds)", "OR", "IC 2.5%", "IC 97.5%"),
        caption = "Odds Ratios con intervalos de confianza al 95%") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Odds Ratios con intervalos de confianza al 95%
	Variable	Coef (log-odds)	OR	IC 2.5%	IC 97.5%
(Intercept)	(Intercept)	0.0742	1.0770	0.9638	1.2036
Lag1	Lag1	-0.0713	0.9312	0.8436	1.0270
Lag2	Lag2	-0.0441	0.9568	0.8671	1.0553
Lag3	Lag3	0.0092	1.0093	0.9152	1.1131
Lag4	Lag4	0.0072	1.0072	0.9133	1.1109
Lag5	Lag5	0.0093	1.0094	0.9160	1.1124

Un OR mayor a 1 indica que un aumento unitario en el predictor incrementa las probabilidades de que el mercado suba. Cuando el intervalo de confianza incluye el 1, el efecto no es significativo al 5%. En este caso, ninguno de los Lags muestra significancia, lo cual es coherente con la dificultad de predecir retornos bursátiles a partir de información pasada.

2 Pregunta 2 — Datos `Boston`

data(Boston)

El dataset Boston (librería MASS) contiene información sobre 506 vecindarios del área de Boston. Se busca predecir medv, el valor mediano de las viviendas (en miles de USD).

2.1 a) Modelo de regresión lineal múltiple

lm_boston <- lm(medv ~ ., data = Boston)
summary(lm_boston)

#> 
#> Call:
#> lm(formula = medv ~ ., data = Boston)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -15.595  -2.730  -0.518   1.777  26.199 
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  3.646e+01  5.103e+00   7.144 3.28e-12 ***
#> crim        -1.080e-01  3.286e-02  -3.287 0.001087 ** 
#> zn           4.642e-02  1.373e-02   3.382 0.000778 ***
#> indus        2.056e-02  6.150e-02   0.334 0.738288    
#> chas         2.687e+00  8.616e-01   3.118 0.001925 ** 
#> nox         -1.777e+01  3.820e+00  -4.651 4.25e-06 ***
#> rm           3.810e+00  4.179e-01   9.116  < 2e-16 ***
#> age          6.922e-04  1.321e-02   0.052 0.958229    
#> dis         -1.476e+00  1.995e-01  -7.398 6.01e-13 ***
#> rad          3.060e-01  6.635e-02   4.613 5.07e-06 ***
#> tax         -1.233e-02  3.760e-03  -3.280 0.001112 ** 
#> ptratio     -9.527e-01  1.308e-01  -7.283 1.31e-12 ***
#> black        9.312e-03  2.686e-03   3.467 0.000573 ***
#> lstat       -5.248e-01  5.072e-02 -10.347  < 2e-16 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 4.745 on 492 degrees of freedom
#> Multiple R-squared:  0.7406, Adjusted R-squared:  0.7338 
#> F-statistic: 108.1 on 13 and 492 DF,  p-value: < 2.2e-16

anova(lm_boston) %>%
  kable(digits = 4, caption = "Tabla ANOVA del modelo completo") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Tabla ANOVA del modelo completo
	Df	Sum Sq	Mean Sq	F value	Pr(>F)
crim	1	6440.7831	6440.7831	286.0300	0.0000
zn	1	3554.3362	3554.3362	157.8452	0.0000
indus	1	2551.2364	2551.2364	113.2984	0.0000
chas	1	1529.8479	1529.8479	67.9393	0.0000
nox	1	76.2476	76.2476	3.3861	0.0664
rm	1	10938.1166	10938.1166	485.7530	0.0000
age	1	90.2679	90.2679	4.0087	0.0458
dis	1	1779.5011	1779.5011	79.0262	0.0000
rad	1	34.1343	34.1343	1.5159	0.2188
tax	1	329.5541	329.5541	14.6352	0.0001
ptratio	1	1309.3093	1309.3093	58.1454	0.0000
black	1	593.3376	593.3376	26.3496	0.0000
lstat	1	2410.8387	2410.8387	107.0634	0.0000
Residuals	492	11078.7846	22.5179	NA	NA

El estadístico \(F\) global contrasta \(H_0\): todos los coeficientes son cero simultáneamente. Con un p-valor extremadamente bajo, rechazamos \(H_0\) — al menos un predictor tiene relación lineal significativa con medv. A nivel individual, variables como rm, lstat, dis y ptratio presentan p-valores muy por debajo de 0.05.

2.2 b) Selección stepwise con AIC

lm_step <- step(lm_boston, direction = "both", trace = 0)
summary(lm_step)

#> 
#> Call:
#> lm(formula = medv ~ crim + zn + chas + nox + rm + dis + rad + 
#>     tax + ptratio + black + lstat, data = Boston)
#> 
#> Residuals:
#>      Min       1Q   Median       3Q      Max 
#> -15.5984  -2.7386  -0.5046   1.7273  26.2373 
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  36.341145   5.067492   7.171 2.73e-12 ***
#> crim         -0.108413   0.032779  -3.307 0.001010 ** 
#> zn            0.045845   0.013523   3.390 0.000754 ***
#> chas          2.718716   0.854240   3.183 0.001551 ** 
#> nox         -17.376023   3.535243  -4.915 1.21e-06 ***
#> rm            3.801579   0.406316   9.356  < 2e-16 ***
#> dis          -1.492711   0.185731  -8.037 6.84e-15 ***
#> rad           0.299608   0.063402   4.726 3.00e-06 ***
#> tax          -0.011778   0.003372  -3.493 0.000521 ***
#> ptratio      -0.946525   0.129066  -7.334 9.24e-13 ***
#> black         0.009291   0.002674   3.475 0.000557 ***
#> lstat        -0.522553   0.047424 -11.019  < 2e-16 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 4.736 on 494 degrees of freedom
#> Multiple R-squared:  0.7406, Adjusted R-squared:  0.7348 
#> F-statistic: 128.2 on 11 and 494 DF,  p-value: < 2.2e-16

coef_tbl <- summary(lm_step)$coefficients %>%
  as.data.frame() %>%
  mutate(Variable = rownames(.)) %>%
  select(Variable, everything())
rownames(coef_tbl) <- NULL

coef_tbl %>%
  kable(digits = 4, caption = "Coeficientes del modelo seleccionado (stepwise AIC)") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Coeficientes del modelo seleccionado (stepwise AIC)
Variable	Estimate	Std. Error	t value	Pr(>\|t\|)
(Intercept)	36.3411	5.0675	7.1714	0.0000
crim	-0.1084	0.0328	-3.3074	0.0010
zn	0.0458	0.0135	3.3902	0.0008
chas	2.7187	0.8542	3.1826	0.0016
nox	-17.3760	3.5352	-4.9151	0.0000
rm	3.8016	0.4063	9.3562	0.0000
dis	-1.4927	0.1857	-8.0370	0.0000
rad	0.2996	0.0634	4.7255	0.0000
tax	-0.0118	0.0034	-3.4925	0.0005
ptratio	-0.9465	0.1291	-7.3337	0.0000
black	0.0093	0.0027	3.4746	0.0006
lstat	-0.5226	0.0474	-11.0187	0.0000

Interpretación de dos coeficientes seleccionados:

rm: cada habitación adicional en promedio se asocia con un aumento de aproximadamente 3.8 miles de USD en el valor mediano, manteniendo el resto constante. Es el predictor con mayor efecto positivo.

lstat: un punto porcentual más de población de bajos ingresos reduce medv en cerca de 0.52 miles de USD. Esto refleja la relación inversa entre estatus socioeconómico y precio de vivienda.

2.3 c) Supuestos del modelo

2.3.1 Multicolinealidad (VIF)

vif_vals <- vif(lm_step)
vif_df <- data.frame(Variable = names(vif_vals), VIF = round(vif_vals, 2))
rownames(vif_df) <- NULL

vif_df %>%
  arrange(desc(VIF)) %>%
  ggplot(aes(x = reorder(Variable, VIF), y = VIF)) +
  geom_col(fill = "#004B93", alpha = 0.75, width = 0.65) +
  geom_hline(yintercept = 5, linetype = "dashed", color = "#5A6C7D") +
  geom_hline(yintercept = 10, linetype = "dotted", color = "#333") +
  coord_flip() +
  labs(x = NULL, y = "VIF",
       caption = "Línea discontinua: VIF = 5 | Línea punteada: VIF = 10") +
  theme_minimal(base_size = 11)

VIF > 5 sugiere multicolinealidad moderada y VIF > 10 multicolinealidad severa. Las variables rad y tax suelen mostrar valores elevados por su fuerte asociación mutua; sin embargo, el modelo puede mantenerse si el objetivo es predicción más que inferencia causal.

2.3.2 Homocedasticidad

par(mfrow = c(1, 2))
plot(lm_step, which = 1, col = "#004B93", pch = 16, cex = 0.6)
plot(lm_step, which = 3, col = "#004B93", pch = 16, cex = 0.6)

par(mfrow = c(1, 1))

bp_test <- lmtest::bptest(lm_step)
cat("Test de Breusch-Pagan\n")

#> Test de Breusch-Pagan

cat("  BP =", round(bp_test$statistic, 4), "\n")

#>   BP = 59.9073

cat("  p  =", format.pval(bp_test$p.value, digits = 4), "\n")

#>   p  = 9.647e-09

El p-valor del test de Breusch-Pagan es inferior a 0.05, por lo que se rechaza la hipótesis nula de homocedasticidad. Esto indica que la varianza de los residuos no es constante, lo cual podría afectar la inferencia sobre los coeficientes.

2.4 d) Observaciones influyentes

2.4.1 Influence Plot

car::influencePlot(lm_step, main = "Gráfico de influencia")

2.4.2 Distancia de Cook

cook_d <- cooks.distance(lm_step)
n <- nrow(Boston)
umbral <- 4 / n

data.frame(Obs = 1:n, Cook = cook_d) %>%
  ggplot(aes(x = Obs, y = Cook)) +
  geom_segment(aes(xend = Obs, yend = 0),
               color = ifelse(cook_d > umbral, "#004B93", "#CBD5E0"), linewidth = 0.4) +
  geom_point(color = ifelse(cook_d > umbral, "#004B93", "#CBD5E0"), size = 0.8) +
  geom_hline(yintercept = umbral, linetype = "dashed", color = "#5A6C7D") +
  labs(x = "Observación", y = "Distancia de Cook",
       caption = paste0("Umbral 4/n = ", round(umbral, 4))) +
  theme_minimal(base_size = 11)

top5 <- sort(cook_d, decreasing = TRUE)[1:5]
data.frame(
  Observación = as.integer(names(top5)),
  Cook_D = round(top5, 5)
) %>%
  kable(caption = "Las 5 observaciones con mayor distancia de Cook") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Las 5 observaciones con mayor distancia de Cook
	Observación	Cook_D
369	369	0.16121
373	373	0.10858
365	365	0.07896
366	366	0.07400
370	370	0.06142

Las observaciones listadas arriba combinan residuos grandes y/o alto leverage, lo que les otorga una influencia desproporcionada sobre las estimaciones. Conviene verificar si corresponden a registros atípicos genuinos o a errores de medición.

3 Pregunta 3 — Datos `PimaIndiansDiabetes`

data(PimaIndiansDiabetes)

Este dataset (mlbench) contiene 768 registros de mujeres de la comunidad Pima con variables clínicas. La variable respuesta diabetes indica presencia (pos) o ausencia (neg) de la enfermedad.

3.1 a) Regresión logística completa

glm_pima <- glm(diabetes ~ ., data = PimaIndiansDiabetes, family = binomial)
summary(glm_pima)

#> 
#> Call:
#> glm(formula = diabetes ~ ., family = binomial, data = PimaIndiansDiabetes)
#> 
#> Coefficients:
#>               Estimate Std. Error z value Pr(>|z|)    
#> (Intercept) -8.4046964  0.7166359 -11.728  < 2e-16 ***
#> pregnant     0.1231823  0.0320776   3.840 0.000123 ***
#> glucose      0.0351637  0.0037087   9.481  < 2e-16 ***
#> pressure    -0.0132955  0.0052336  -2.540 0.011072 *  
#> triceps      0.0006190  0.0068994   0.090 0.928515    
#> insulin     -0.0011917  0.0009012  -1.322 0.186065    
#> mass         0.0897010  0.0150876   5.945 2.76e-09 ***
#> pedigree     0.9451797  0.2991475   3.160 0.001580 ** 
#> age          0.0148690  0.0093348   1.593 0.111192    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> (Dispersion parameter for binomial family taken to be 1)
#> 
#>     Null deviance: 993.48  on 767  degrees of freedom
#> Residual deviance: 723.45  on 759  degrees of freedom
#> AIC: 741.45
#> 
#> Number of Fisher Scoring iterations: 5

coef_pima <- summary(glm_pima)$coefficients %>%
  as.data.frame() %>%
  mutate(Variable = rownames(.), OR = exp(Estimate)) %>%
  select(Variable, Estimate, OR, `Std. Error`, `z value`, `Pr(>|z|)`)
rownames(coef_pima) <- NULL

coef_pima %>%
  kable(digits = 4,
        col.names = c("Variable", "Log-odds", "OR", "Err. Std.", "z", "p-valor"),
        caption = "Coeficientes del modelo logístico") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Coeficientes del modelo logístico
Variable	Log-odds	OR	Err. Std.	z	p-valor
(Intercept)	-8.4047	0.0002	0.7166	-11.7280	0.0000
pregnant	0.1232	1.1311	0.0321	3.8401	0.0001
glucose	0.0352	1.0358	0.0037	9.4814	0.0000
pressure	-0.0133	0.9868	0.0052	-2.5404	0.0111
triceps	0.0006	1.0006	0.0069	0.0897	0.9285
insulin	-0.0012	0.9988	0.0009	-1.3223	0.1861
mass	0.0897	1.0938	0.0151	5.9453	0.0000
pedigree	0.9452	2.5733	0.2991	3.1596	0.0016
age	0.0149	1.0150	0.0093	1.5929	0.1112

glucose: cada unidad adicional de glucosa incrementa el log-odds en 0.0352, lo que equivale a multiplicar los odds por 1.0358. Aunque el efecto unitario parece pequeño, el rango de glucosa es amplio (44–199), por lo que el impacto acumulado es considerable.

mass (IMC): OR = 1.0938. Un punto más de IMC aumenta los odds de diabetes en un 9.38%.

Respecto a la significancia: se compara cada p-valor contra \(\alpha = 0.05\).

3.2 b) Bondad de ajuste (Devianza)

dev_nula  <- glm_pima$null.deviance
dev_resid <- glm_pima$deviance
pct_exp   <- (1 - dev_resid / dev_nula) * 100

data.frame(
  Indicador = c("Devianza nula (solo intercepto)", "Devianza residual (modelo completo)", "% Devianza explicada"),
  Valor = c(round(dev_nula, 2), round(dev_resid, 2), paste0(round(pct_exp, 2), "%"))
) %>%
  kable(caption = "Comparación de devianzas") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Comparación de devianzas
Indicador	Valor
Devianza nula (solo intercepto)	993.48
Devianza residual (modelo completo)	723.45
% Devianza explicada	27.18%

El modelo explica aproximadamente un 27.2% de la devianza total. Este valor puede parecer modesto, pero es habitual en estudios biomédicos donde la variabilidad entre individuos es alta y factores no capturados (genéticos, ambientales) juegan un rol importante. La reducción es estadísticamente significativa, lo cual se confirma con el test de razón de verosimilitudes implícito en la diferencia de devianzas.

3.3 c) Curva ROC y punto de corte óptimo (Youden)

prob_pima <- predict(glm_pima, type = "response")
roc_pima  <- roc(PimaIndiansDiabetes$diabetes, prob_pima, levels = c("neg", "pos"))

coords_opt <- coords(roc_pima, x = "best", best.method = "youden",
                     ret = c("threshold", "sensitivity", "specificity"))

roc_df <- data.frame(
  FPR = 1 - roc_pima$specificities,
  TPR = roc_pima$sensitivities
)

ggplot(roc_df, aes(x = FPR, y = TPR)) +
  geom_line(color = "#004B93", linewidth = 1) +
  geom_abline(slope = 1, intercept = 0, linetype = "dashed", color = "#CBD5E0") +
  geom_point(aes(x = 1 - coords_opt$specificity, y = coords_opt$sensitivity),
             color = "#5A6C7D", size = 3.5, shape = 18) +
  annotate("text", x = 0.55, y = 0.3,
           label = paste0("AUC = ", round(auc(roc_pima), 4)),
           size = 4.5, color = "#004B93") +
  annotate("text", x = 1 - coords_opt$specificity + 0.12,
           y = coords_opt$sensitivity - 0.04,
           label = paste0("Youden (", round(coords_opt$threshold, 3), ")"),
           size = 3.2, color = "#5A6C7D") +
  labs(x = "1 - Especificidad", y = "Sensibilidad") +
  theme_minimal(base_size = 11)

cat("AUC:", round(auc(roc_pima), 4), "\n")

#> AUC: 0.8394

cat("Punto de corte óptimo (Youden):", round(coords_opt$threshold, 4), "\n")

#> Punto de corte óptimo (Youden): 0.3537

cat("  Sensibilidad:", round(coords_opt$sensitivity, 4), "\n")

#>   Sensibilidad: 0.7388

cat("  Especificidad:", round(coords_opt$specificity, 4), "\n")

#>   Especificidad: 0.784

El AUC de 0.839 indica una capacidad discriminativa aceptable. El índice de Youden (\(J = \text{Sens} + \text{Esp} - 1\)) identifica el umbral que balancea mejor la detección de positivos y negativos. En este caso, el punto de corte óptimo resulta ser 0.354, distinto del 0.5 habitual, lo que refleja el desbalance entre clases y el trade-off inherente a la detección temprana.

3.4 d) Matriz de confusión y métricas

umbral_opt <- coords_opt$threshold
pred_class <- factor(ifelse(prob_pima >= umbral_opt, "pos", "neg"), levels = c("neg", "pos"))
cm <- table(Predicho = pred_class, Real = PimaIndiansDiabetes$diabetes)

as.data.frame(cm) %>%
  ggplot(aes(x = Real, y = Predicho, fill = Freq)) +
  geom_tile(color = "white", linewidth = 1.5) +
  geom_text(aes(label = Freq), size = 7, fontface = "bold", color = "white") +
  scale_fill_gradient(low = "#8DA0B0", high = "#004B93") +
  labs(title = paste0("Matriz de confusión (umbral = ", round(umbral_opt, 4), ")")) +
  theme_minimal(base_size = 11) +
  theme(legend.position = "none",
        axis.text = element_text(size = 11))

TP <- cm["pos", "pos"]; TN <- cm["neg", "neg"]
FP <- cm["pos", "neg"]; FN <- cm["neg", "pos"]

sens <- TP / (TP + FN)
espec <- TN / (TN + FP)
acc <- (TP + TN) / sum(cm)

data.frame(
  Métrica = c("Sensibilidad", "Especificidad", "Accuracy"),
  Valor   = round(c(sens, espec, acc), 4)
) %>%
  kable(caption = "Métricas de clasificación") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)

Métricas de clasificación
Métrica	Valor
Sensibilidad	0.7388
Especificidad	0.7840
Accuracy	0.7682

En un problema de detección de enfermedades, generalmente conviene priorizar la sensibilidad por sobre la especificidad. Un falso negativo implica que un paciente diabético no recibe diagnóstico ni tratamiento oportuno, con riesgo de complicaciones como retinopatía, neuropatía o enfermedad cardiovascular. Un falso positivo, en cambio, se resuelve con exámenes confirmatorios adicionales. Por esta razón, un punto de corte menor a 0.5 (como el obtenido por Youden) puede ser clínicamente preferible.

EPG4001 — Aprendizaje Supervisado · Prof. Jorge Luis Bazán
Pontificia Universidad Católica de Chile · Facultad de Matemáticas
Santiago, 14 de June de 2026

Laboratorio 1 — Aprendizaje Supervisado

EPG4001 · Pontificia Universidad Católica de Chile

Armando Arredondo Valle & Julynho Merlin Toledo Rios

14 de June de 2026

1 Pregunta 1 — Datos `Smarket`

1.1 a) Análisis descriptivo

1.1.1 Estructura

1.1.2 Estadísticos descriptivos

1.1.3 Histogramas

1.1.4 Correlaciones

1.1.5 Boxplots por Direction

1.2 b) Selección de modelo lineal para `Today`

1.3 c) Regresión logística para `Direction`

1.4 d) Odds Ratios

2 Pregunta 2 — Datos `Boston`

2.1 a) Modelo de regresión lineal múltiple

2.2 b) Selección stepwise con AIC

2.3 c) Supuestos del modelo

2.3.1 Multicolinealidad (VIF)

2.3.2 Homocedasticidad

2.4 d) Observaciones influyentes

2.4.1 Influence Plot

2.4.2 Distancia de Cook

3 Pregunta 3 — Datos `PimaIndiansDiabetes`

3.1 a) Regresión logística completa

3.2 b) Bondad de ajuste (Devianza)

3.3 c) Curva ROC y punto de corte óptimo (Youden)

3.4 d) Matriz de confusión y métricas

Laboratorio 1 — Aprendizaje Supervisado

EPG4001 · Pontificia Universidad Católica de Chile

Armando Arredondo Valle & Julynho Merlin Toledo Rios

14 de June de 2026

1 Pregunta 1 — Datos Smarket

1.1 a) Análisis descriptivo

1.1.1 Estructura

1.1.2 Estadísticos descriptivos

1.1.3 Histogramas

1.1.4 Correlaciones

1.1.5 Boxplots por Direction

1.2 b) Selección de modelo lineal para Today

1.3 c) Regresión logística para Direction

1.4 d) Odds Ratios

2 Pregunta 2 — Datos Boston

2.1 a) Modelo de regresión lineal múltiple

2.2 b) Selección stepwise con AIC

2.3 c) Supuestos del modelo

2.3.1 Multicolinealidad (VIF)

2.3.2 Homocedasticidad

2.4 d) Observaciones influyentes

2.4.1 Influence Plot

2.4.2 Distancia de Cook

3 Pregunta 3 — Datos PimaIndiansDiabetes

3.1 a) Regresión logística completa

3.2 b) Bondad de ajuste (Devianza)

3.3 c) Curva ROC y punto de corte óptimo (Youden)

3.4 d) Matriz de confusión y métricas

1 Pregunta 1 — Datos `Smarket`

1.2 b) Selección de modelo lineal para `Today`

1.3 c) Regresión logística para `Direction`

2 Pregunta 2 — Datos `Boston`

3 Pregunta 3 — Datos `PimaIndiansDiabetes`