Evaluación de modelos

Fecha de publicación

12 de diciembre de 2024

Objetivo del manual

Familiarizarse con las metricas de evaluación de modelos de aprendizaje estadístico
Ser capaz de evaluar modelos de regresión y clasificación en R

Paquetes a utilizar en este manual:

Código

# instalar/cargar paquetes

sketchy::load_packages(
  c("ggplot2", 
    "viridis", 
    "nnet",
    "caret",
    "glmnet",
    "pROC",
    "nnet",
    "Metrics")
  )

La evaluación del desempeño de un modelo estadístico se basa en comparar las predicciones del modelo con los valores observados para determinar qué tan bien el modelo está capturando la relación entre las variables. Esta comparación se realiza utilizando diversas métricas de desempeño que cuantifican el grado de coincidencia entre las predicciones y las observaciones reales.

La idea es evaluar la capacidad del modelo para hacer predicciones precisas y generalizar bien a datos no vistos. Por lo tanto, evaluar correctamente los modelos permite entender qué tan bien el modelo se ajusta a los datos e idealmente, entender si sus predicciones son útiles para datos no vistos. La naturaleza de la métrica depende de la estructura de los datos predichos (si estos son continuos o categóricos):

flowchart LR
    classDef largeText font-size:18px, padding:15px;

    A(Evaluación de modelos) --> B(Respuesta continua)
    A --> C("Respuesta categórica")
    B --> D(RMSE)
    B --> E(R cuadrado)
    C --> F("Matriz de confusión")
    C --> G("Precisión, sensibilidad, F1")
    C --> H("Pérdida logarítmica")

    style A fill:#40498E66, stroke:#000, stroke-width:2px, color:#FFF, width:180px
    style B fill:#348AA666, stroke:#000, stroke-width:2px, color:#FFF 
    style C fill:#348AA666, stroke:#000, stroke-width:2px, color:#FFF
    style D fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000
    style E fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000
    style F fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000
    style G fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000
    style H fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000

1 Evaluación para variables respuesta continuas

Cuando la variable respuesta $\text{Y}$ es continua el propósito de la evaluación consiste en medir qué tan cerca están las predicciones del modelo ($\hat{Y}$) con respecto a los valores reales observados ($\text{Y}$). Esto implica cuantificar el error entre los valores predichos y los valores observados en una escala que permita comparar la métrica entre modelos.

1.1 Raíz del Error Cuadrático Medio (RMSE)

La raíz del error cuadrático medio (RMSE por sus siglas en inglés) es quizás la medida mas comunmente utilizada para medir la precisión en la predicción de una variable continua. Las unidades de esta métrica son las mismas que la variable de respuesta. La podemos calcular fácilmente:

Código

# Generar datos simulados
set.seed(123)
x <- rnorm(100)
y <- 3 * x + rnorm(100, sd = 0.5)

# Crear un marco de datos
datos <- data.frame(x = x, y = y)

# correr modelo
mod <- lm(y ~ x, data = datos)

# Step 2: Obtener predicciones
y_pred <- predict(mod)

# Step 3: Calcular RMSE
sqrt(mean((datos$y - y_pred)^2))

[1] 0.48048

La raíz del error cuadrático medio penaliza más los errores grandes, ya que los errores se elevan al cuadrado. Esto lo hace útil para detectar errores grandes que podrían ser problemáticos. Sin embargo, es sensible a outliers, ya que los errores grandes tienen un impacto más significativo en el resultado.

1.2 Coeficiente de determinación R²

El R² indica qué proporción de la variabilidad en la variable respuesta es explicada por el modelo. Se calcula así:

Código

r2 <- 1 - sum((datos$y - y_pred)^2) / sum((datos$y - mean(datos$y))^2)
    
# imprimir R^2
r2

[1] 0.96932

Este es el mismo valor que se obtiene al evaluar el modelo con la función summary:

Código

summary(mod)$r.squared

[1] 0.96932

Entre R² sea mas cercano a 1 indica el modelo tiene mejor ajuste. Sin embargo, es importante tener en cuenta que un R² alto no garantiza que el modelo sea útil o que generalice bien a nuevos datos. No es una buena métrica para modelos no lineales, ya que puede no reflejar bien la calidad del ajuste en estos casos.

Ejercicio 1

Los datos mtcars se tomaron de la revista Motor Trend US de 1974, y contienen informacion sobre 10 características del diseño y desempeño de 32 marcas de carros.

Código

# cargar datos
data(mtcars)

Ajusta un modelo de regresión lineal con los datos de mtcars donde la variable respuesta sea mpg (millas por galón) y las variables predictoras sean hp (caballos de fuerza) y wt (peso).

Ajuste 2 modelos más, cada uno con 2 variables predictoras más que el modelo ajustado en el paso 1.
Compare los modelos generados en los puntos 1 y 2 usando la raíz del error cuadrático medio (RMSE) y el R².

2 Evaluación para variables respuesta categóricas

Para demostrar las métricas de evalución para estos modelos utilizaremos un modelo de regresión con los datos del titanic similar al del tutorial de regresión logística.

Código

# cargar datos
data("Titanic")

# dar formato con una observacion por fila
datos_titanic <- as.data.frame(Titanic)
datos_titanic <- datos_titanic[datos_titanic$Freq > 0, ]

datos_tab_titanic <- do.call(rbind, lapply(1:nrow(datos_titanic), function(x) datos_titanic[rep(x, datos_titanic$Freq[x]),]))

# Ajustar el modelo de regresión logística
modelo_log <- glm(Survived ~ Class + Sex, data = datos_tab_titanic, family = binomial)

Ahora podemos predecir la variable respuesta con este modelo:

Código

# predecir valores
pred_vals <- predict(object = modelo_log, newdata = datos_tab_titanic, type = "response")

# binarizar
pred_cat <- ifelse(pred_vals > 0.5, "Yes", "No")

2.1 Matriz de confusión

La matriz de confusión organiza las predicciones y los resultados observados en categorías que indican si el valor predicho es una ocurrencia del evento (positivo vs negativo) y si coincide o no con el valor observado (verdadero vs falso).

Esta matriz la podemos estimar usando la función confusionMatrix del paquete caret:

Código

# hacer la matriz de confusion
mat_conf <-
    confusionMatrix(factor(pred_cat), factor(datos_tab_titanic$Survived))

# imprimir resultado
mat_conf$table

          Reference
Prediction   No  Yes
       No  1364  367
       Yes  126  344

Podemos graficar esta matriz con un gradiente de colores donde los valores esten dado como proporciones, lo cual mas fácil su interpretación:

Código

# convertir a data frame
conf_df <- as.data.frame(mat_conf$table)

# agregar totales por categoria
conf_df$total <-
    sapply(conf_df$Reference, function(x)
        sum(datos_tab_titanic$Survived ==
                x))

# calcular proporciones
conf_df$proportion <- conf_df$Freq / conf_df$total

# graficar
ggplot(conf_df, aes(x = Reference, y = Prediction, fill = proportion)) +
  geom_tile() + 
  coord_equal() + 
  scale_fill_distiller(palette = "Greens", direction = 1) + 
  geom_text(aes(label = round(proportion, 2)), color = "black", size = 3) + 
  labs(x = "Observado", y = "Predicho", fill = "Proporción")

Varias métricas se derivan de los valores contenidos en la matriz de confusion.

2.2 Exactitud (Accuracy)

La exactitud es la proporción de predicciones correctas en relación con el total de predicciones.

La fórmula es:

$\frac{VP + VN}{VP + VN + FP + FN}$

Esta también es estimada por la función confusionMatrix:

Código

mat_conf$overall["Accuracy"]

Accuracy 
 0.77601

2.3 Precisión (Precision) y Sensibilidad (Sensitivity o Recall)

La precisión es la proporción de verdaderos positivos entre todas las predicciones positivas:

$\frac{VP}{VP + FP}$

La sensibilidad es la proporción de verdaderos positivos entre todas las observaciones positivas:

$\frac{VP}{VP + FN}$

Ambas son estimadas por la función confusionMatrix:

Código

mat_conf$byClass[c("Precision", "Sensitivity")]

  Precision Sensitivity 
    0.78798     0.91544

2.4 Índice F1 (o F-Score)

El índice F1 es la media armónica entre la precisión y la sensibilidad:

$2 \times \frac{\text{Precisión} \times \text{Sensibilidad}}{\text{Precisión} + \text{Sensibilidad}}$

Código

precision <- mat_conf$byClass["Precision"]
sensitivity <- mat_conf$byClass["Sensitivity"]


f1_score <- (2 * precision * sensitivity) / (precision + sensitivity)

f1_score

Precision 
  0.84694

Tambien lo podemos extraer del resultado de la función confusionMatrix:

Código

 mat_conf$byClass["F1"]

     F1 
0.84694

2.5 Área Bajo la Curva ROC (AUC-ROC)

El AUC métrica originalmente diseñada para evaluar la precision en la predicción de variables binarias. Esta mide la capacidad del modelo para discriminar entre clases. Se calcula basándose en la curva ROC (Receiver Operating Characteristic). La curva ROC compara la tasa de verdaderos positivos contra la tasa de falsos positivos a diferentes umbrales.La curva ROC es una representación gráfica que evalúa el desempeño de un modelo de clasificación binaria. La curva muestra la relación entre la Tasa de Verdaderos Positivos (True Positive Rate, TPR), también llamada sensibilidad, y la Tasa de Falsos Positivos (False Positive Rate, FPR) para diferentes umbrales de decisión. Los umbrales de decisión son los valores que determinan si una observación se clasifica como positiva o negativa. Podemos ejemplificar su uso con el modelo de regresión logística ajustado mas arriba:

Código

# estimar las probabilidades
pred_vals <- predict(object = modelo_log, newdata = datos_tab_titanic, type = "response")

# Curva ROC
roc_curve <- roc(datos_tab_titanic$Survived, pred_vals, smooth = TRUE)

# Extraer los datos de la curva ROC
roc_data <- data.frame(
  tpr = roc_curve$sensitivities, # True Positive Rate (Sensibilidad)
  fpr = 1 - roc_curve$specificities # False Positive Rate (1 - Especificidad)
)

# Graficar la curva ROC con ggplot2
ggplot(roc_data, aes(x = fpr, y = tpr)) +
  geom_line(color = "tomato") +
  geom_abline(linetype = "dashed", color = "gray") + # Línea diagonal
  labs(
    x = "Tasa de Falsos Positivos (FPR)",
    y = "Tasa de Verdaderos Positivos (TPR)"
  ) +
  scale_x_continuous(limits = c(0, 1), expand = c(0, 0)) +  # Limitar el eje X entre 0 y 1
  scale_y_continuous(limits = c(0, 1), expand = c(0, 0)) +  # Limitar el eje Y entre 0 y 
  theme_classic() +
  theme(panel.grid = element_blank()) +# Remover cuadrícula
  # add polygon showin area under the curve
  geom_polygon(data = roc_data, aes(x = fpr, y = tpr), fill = "tomato", alpha = 0.3) + 
# add polygon coloring lower diagonal 
  geom_polygon(data = data.frame(x = c(0, 1, 1), y = c(0, 0, 1)), aes(x = x, y = y), fill = "tomato", alpha = 0.3) +
  annotate("text", x = 0.5, y = 0.5, label = paste("AUC = ", round(pROC::auc(roc_curve), 2)), size = 6)

Podemos estimar el área bajo la curva con la función auc del paquete pROC:

Código

# estimar AUC
pROC::auc(roc_curve)

Area under the curve: 0.758

Valores mas cercanos a 1 indican un mejor desempeño del modelo en discriminar entre clases tomando en cuenta la incertidumbre relacionada a la escogencia de un umbral de decisión.

En los modelos que predicen más de dos categorías (clasificación multiclase), no se puede aplicar la curva ROC directamente como en el caso binario. Aun así, hay métodos extendidos para calcular una métrica similar a la AUC para clasificación multiclase. Esta métrica se puede extender para modelos que predicen multiples clases. Por ejemplo el metodo de “uno vs el resto” (one vs rest) calcula una curva ROC para cada categoría de manera binaria, considerando una clase como “positiva” y las demás como “negativas”. Luego, se calcula el AUC para cada una de estas curvas ROC y se promedia. Para ejemplificar el uso de este método usaremos una regresión multinomial para predecir la especie de las observaciones del juego de datos iris:

Código

# Ajustar el modelo de regresión multinomial
modelo_multinom <- multinom(Species ~ ., data = iris)

# weights:  18 (10 variable)
initial  value 164.791843 
iter  10 value 16.177348
iter  20 value 7.111438
iter  30 value 6.182999
iter  40 value 5.984028
iter  50 value 5.961278
iter  60 value 5.954900
iter  70 value 5.951851
iter  80 value 5.950343
iter  90 value 5.949904
iter 100 value 5.949867
final  value 5.949867 
stopped after 100 iterations

Código

# Predecir las probabilidades para cada clase
pred_prob <- predict(modelo_multinom, newdata = iris, type = "prob")

# Aplicar multiclass.roc (One-vs-Rest)
roc_multiclass <- multiclass.roc(iris$Species, pred_prob)

# Ver resultados
roc_multiclass


Call:
multiclass.roc.default(response = iris$Species, predictor = pred_prob)

Data: multivariate predictor pred_prob with 3 levels of iris$Species: setosa, versicolor, virginica.
Multi-class area under the curve: 0.999

2.6 Pérdida logarítmica (log-loss)

Mide el rendimiento de un modelo de clasificación al penalizar predicciones erróneas basadas en la probabilidad que asigna el modelo a cada clase.

$-\frac{1}{n} \sum_{i=1}^{n} \left[ y_i \log(p_i) + (1 - y_i) \log(1 - p_i) \right]$

Donde ${y_i}$ es el valor real (0 o 1) y ${p_i}$ es la probabilidad asignada a la clase positiva. Es útil cuando se desea medir la confianza del modelo en sus predicciones y no solo si acierta o falla. Aunque puede ser más complejo de interpretar en comparación con otras métricas.

Código

# Convertir las etiquetas reales a una matriz de indicadores (one-hot encoding)
# note que Species es un factor
y_obs <- model.matrix(~Species-1, data = iris)

# Calcular el log-loss
log_loss_value <- logLoss(y_obs, pred_prob)
log_loss_value

[1] 0.026444

El valor del pérdida logarítmica representa qué tan bien o mal el modelo predice las clases correctas. Un valor más bajo (mas cerano a 0) indica un mejor ajuste, mientras que un valor alto sugiere que el modelo está asignando bajas probabilidades a las clases correctas. El valor de pérdida logarítmica puede tomar cualquier valor real no negativo, y en general se encuentra entre 0 y ∞.

Ejercicio 2

Usaremos los datos de Cortez y Silva (2008) sobre el desempeño académico de estudiantes en la clase de matemáticas. Los datos contienen 395 observaciones y 31 variables. La variable respuesta es G3 (calificación final), y las variables predictoras son las demás variables del conjunto de datos. Pueden econtrar una descripción detallada de los datos aqui. El siguiente código carga los datos y remueve las variables G1 y G2 (calificaciones parciales) para evitar multicolinealidad:

Código

# leer datos
datos_mate <- read.csv("https://raw.githubusercontent.com/maRce10/aprendizaje_estadistico_2024/refs/heads/master/data/student/student-mat.csv", sep = ";")

# remover G1 y G2 (calificaciones parciales)
datos_mate$G1 <- datos_mate$G2 <- NULL

str(datos_mate)

'data.frame':   395 obs. of  31 variables:
 $ school    : chr  "GP" "GP" "GP" "GP" ...
 $ sex       : chr  "F" "F" "F" "F" ...
 $ age       : int  18 17 15 15 16 16 16 17 15 15 ...
 $ address   : chr  "U" "U" "U" "U" ...
 $ famsize   : chr  "GT3" "GT3" "LE3" "GT3" ...
 $ Pstatus   : chr  "A" "T" "T" "T" ...
 $ Medu      : int  4 1 1 4 3 4 2 4 3 3 ...
 $ Fedu      : int  4 1 1 2 3 3 2 4 2 4 ...
 $ Mjob      : chr  "at_home" "at_home" "at_home" "health" ...
 $ Fjob      : chr  "teacher" "other" "other" "services" ...
 $ reason    : chr  "course" "course" "other" "home" ...
 $ guardian  : chr  "mother" "father" "mother" "mother" ...
 $ traveltime: int  2 1 1 1 1 1 1 2 1 1 ...
 $ studytime : int  2 2 2 3 2 2 2 2 2 2 ...
 $ failures  : int  0 0 3 0 0 0 0 0 0 0 ...
 $ schoolsup : chr  "yes" "no" "yes" "no" ...
 $ famsup    : chr  "no" "yes" "no" "yes" ...
 $ paid      : chr  "no" "no" "yes" "yes" ...
 $ activities: chr  "no" "no" "no" "yes" ...
 $ nursery   : chr  "yes" "no" "yes" "yes" ...
 $ higher    : chr  "yes" "yes" "yes" "yes" ...
 $ internet  : chr  "no" "yes" "yes" "yes" ...
 $ romantic  : chr  "no" "no" "no" "yes" ...
 $ famrel    : int  4 5 4 3 4 5 4 4 4 5 ...
 $ freetime  : int  3 3 3 2 3 4 4 1 2 5 ...
 $ goout     : int  4 3 2 2 2 2 4 4 2 1 ...
 $ Dalc      : int  1 1 2 1 1 1 1 1 1 1 ...
 $ Walc      : int  1 1 3 1 2 2 1 1 1 1 ...
 $ health    : int  3 3 3 5 5 5 3 1 1 5 ...
 $ absences  : int  6 4 10 2 4 10 0 6 0 0 ...
 $ G3        : int  6 6 10 15 10 15 11 6 19 15 ...

Ajuste un modelo de regresión lineal con los datos de datos_mate donde la variable respuesta sea G3. Utilice su conocimiento del sistema (como estudiante) para escoger las variables predictoras que usted espere tengan un mayor efecto sobre el desempeño académico. Calcule el R² para este modelo.

Ajuste un modelo de regresión multinomial análogo al del punto anterior (e.g. mismos predictores y misma respuesta).

Calcule la matriz de confusión, la exactitud, el área bajo la curva y pérdida logarítmica para el modelo del punto anterior.

El siguiente código binariza la variable respuesta G3. Si la calificación es menor a 9, se considera que el estudiante reprobó (Reprobado), de lo contrario, se considera que aprobó (Aprobado):

Código

datos_mate$G3_bin <- ifelse(datos_mate$G3 < 9, 0, 1)

Utilizando la variable binaria G3_bin como respuesta, ajuste un modelo de regresión logística con las mismas variables predictoras que el modelo ajustado en el punto 2.

5.Calcule la matriz de confusión, la exactitud, precisión, indice F1, area bajo la curva y pérdida logarítmica para el modelo del punto anterior.

Compare los modelos generados en los puntos 2 y 4. ¿Cual modelo parece tener un mejor desempeño? ¿Por qué?
Explore la sensibilidad de los modelos multinomial (punto 2) y logístico (punto 4). Para esto puede simplemente imprimir en la consola el resultado de la función confusionMatrix. ¿Cómo difieren los valores y la forma en que estos se estructuran entre los 2 modelos?

3 Referencias

Cortez, P., & Silva, A.M. (2008). Using data mining to predict secondary school student performance. (enlace)

Información de la sesión

R version 4.3.2 (2023-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 22.04.2 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0 
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

time zone: America/Costa_Rica
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] Metrics_0.1.4     nnet_7.3-17       pROC_1.18.5       caret_6.0-94     
 [5] lattice_0.20-45   glmnet_4.1-8      Matrix_1.6-5      viridis_0.6.5    
 [9] viridisLite_0.4.2 ggplot2_3.5.1     knitr_1.48       

loaded via a namespace (and not attached):
 [1] tidyselect_1.2.1     timeDate_4032.109    dplyr_1.1.4         
 [4] farver_2.1.2         fastmap_1.2.0        digest_0.6.37       
 [7] rpart_4.1.16         timechange_0.2.0     lifecycle_1.0.4     
[10] survival_3.2-13      magrittr_2.0.3       compiler_4.3.2      
[13] rlang_1.1.4          tools_4.3.2          utf8_1.2.4          
[16] yaml_2.3.10          data.table_1.14.10   labeling_0.4.3      
[19] htmlwidgets_1.6.4    RColorBrewer_1.1-3   plyr_1.8.9          
[22] withr_3.0.1          purrr_1.0.2          grid_4.3.2          
[25] stats4_4.3.2         fansi_1.0.6          e1071_1.7-16        
[28] colorspace_2.1-1     future_1.34.0        globals_0.16.3      
[31] scales_1.3.0         iterators_1.0.14     MASS_7.3-55         
[34] cli_3.6.3            rmarkdown_2.28       crayon_1.5.3        
[37] generics_0.1.3       remotes_2.5.0        rstudioapi_0.16.0   
[40] future.apply_1.11.2  reshape2_1.4.4       proxy_0.4-27        
[43] stringr_1.5.1        splines_4.3.2        parallel_4.3.2      
[46] vctrs_0.6.5          hardhat_1.3.0        jsonlite_1.8.9      
[49] listenv_0.9.1        packrat_0.9.2        foreach_1.5.2       
[52] gower_1.0.1          recipes_1.0.9        glue_1.8.0          
[55] parallelly_1.38.0    codetools_0.2-18     xaringanExtra_0.8.0 
[58] lubridate_1.9.3      stringi_1.8.4        shape_1.4.6         
[61] gtable_0.3.5         munsell_0.5.1        tibble_3.2.1        
[64] pillar_1.9.0         htmltools_0.5.8.1    ipred_0.9-14        
[67] lava_1.7.3           R6_2.5.1             evaluate_1.0.0      
[70] sketchy_1.0.3        class_7.3-20         Rcpp_1.0.13         
[73] gridExtra_2.3        nlme_3.1-155         prodlim_2023.08.28  
[76] xfun_0.47            pkgconfig_2.0.3      ModelMetrics_1.2.2.2

--- title: Evaluación de modelos editor_options: chunk_output_type: console --- ```{r,echo=FALSE,message=FALSE} options("digits"=5) options("digits.secs"=3) ``` ```{r, echo = FALSE, message=FALSE} library(knitr) library(ggplot2) library(viridis) library(glmnet) library(caret) library(pROC) library(nnet) library(Metrics) theme_set(theme_classic(base_size = 20)) ``` ::: {.alert .alert-info} # Objetivo del manual {.unnumbered .unlisted} - Familiarizarse con las metricas de evaluación de modelos de aprendizaje estadístico - Ser capaz de evaluar modelos de regresión y clasificación en R ::: Paquetes a utilizar en este manual: ```{r} # instalar/cargar paquetes sketchy::load_packages( c("ggplot2", "viridis", "nnet", "caret", "glmnet", "pROC", "nnet", "Metrics") ) ``` ::: grid ::: g-col-8 La evaluación del desempeño de un modelo estadístico se basa en comparar las predicciones del modelo con los valores observados para determinar qué tan bien el modelo está capturando la relación entre las variables. Esta comparación se realiza utilizando diversas métricas de desempeño que cuantifican el grado de coincidencia entre las predicciones y las observaciones reales. ::: ::: g-col-4 ![](images/dart_target.png){fig-align="right"} ::: ::: La idea es evaluar la capacidad del modelo para hacer predicciones precisas y generalizar bien a datos no vistos. Por lo tanto, evaluar correctamente los modelos permite entender qué tan bien el modelo se ajusta a los datos e idealmente, entender si sus predicciones son útiles para datos no vistos. La naturaleza de la métrica depende de la estructura de los datos predichos (si estos son continuos o categóricos): ```{mermaid} flowchart LR classDef largeText font-size:18px, padding:15px; A(Evaluación de modelos) --> B(Respuesta continua) A --> C("Respuesta categórica") B --> D(RMSE) B --> E(R cuadrado) C --> F("Matriz de confusión") C --> G("Precisión, sensibilidad, F1") C --> H("Pérdida logarítmica") style A fill:#40498E66, stroke:#000, stroke-width:2px, color:#FFF, width:180px style B fill:#348AA666, stroke:#000, stroke-width:2px, color:#FFF style C fill:#348AA666, stroke:#000, stroke-width:2px, color:#FFF style D fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000 style E fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000 style F fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000 style G fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000 style H fill:#49C1AD66, stroke:#000, stroke-width:2px, color:#000 ``` # Evaluación para variables respuesta continuas Cuando la variable respuesta $\text{Y}$ es continua el propósito de la evaluación consiste en medir qué tan cerca están las predicciones del modelo ($\hat{Y}$) con respecto a los valores reales observados ($\text{Y}$). Esto implica cuantificar el error entre los valores predichos y los valores observados en una escala que permita comparar la métrica entre modelos. ## Raíz del Error Cuadrático Medio (RMSE) La raíz del error cuadrático medio (RMSE por sus siglas en inglés) es quizás la medida mas comunmente utilizada para medir la precisión en la predicción de una variable continua. Las unidades de esta métrica son las mismas que la variable de respuesta. La podemos calcular fácilmente: ```{r} # Generar datos simulados set.seed(123) x <- rnorm(100) y <- 3 * x + rnorm(100, sd = 0.5) # Crear un marco de datos datos <- data.frame(x = x, y = y) # correr modelo mod <- lm(y ~ x, data = datos) # Step 2: Obtener predicciones y_pred <- predict(mod) # Step 3: Calcular RMSE sqrt(mean((datos$y - y_pred)^2)) ``` La raíz del error cuadrático medio penaliza más los errores grandes, ya que los errores se elevan al cuadrado. Esto lo hace útil para detectar errores grandes que podrían ser problemáticos. Sin embargo, es sensible a outliers, ya que los errores grandes tienen un impacto más significativo en el resultado. ## Coeficiente de determinación R^2^ El R^2^ indica qué proporción de la variabilidad en la variable respuesta es explicada por el modelo. Se calcula así: ```{r} r2 <- 1 - sum((datos$y - y_pred)^2) / sum((datos$y - mean(datos$y))^2) # imprimir R^2 r2 ``` Este es el mismo valor que se obtiene al evaluar el modelo con la función `summary`: ```{r} summary(mod)$r.squared ``` Entre R^2^ sea mas cercano a 1 indica el modelo tiene mejor ajuste. Sin embargo, es importante tener en cuenta que un R^2^ alto no garantiza que el modelo sea útil o que generalice bien a nuevos datos. No es una buena métrica para modelos no lineales, ya que puede no reflejar bien la calidad del ajuste en estos casos. ::: {.alert .alert-info} ### Ejercicio 1 {.unnumbered} Los datos `mtcars` se tomaron de la revista Motor Trend US de 1974, y contienen informacion sobre 10 características del diseño y desempeño de 32 marcas de carros. ```{r} # cargar datos data(mtcars) ``` 1. Ajusta un modelo de regresión lineal con los datos de `mtcars` donde la variable respuesta sea `mpg` (millas por galón) y las variables predictoras sean `hp` (caballos de fuerza) y `wt` (peso). ```{r solucion del ejercicio 1.1} #| eval: false #| echo: false # cargar datos data(mtcars) # ajustar modelo mod_mtcars <- lm(mpg ~ hp + wt, data = mtcars) # imprimir resumen summary(mod_mtcars) ``` 2. Ajuste 2 modelos más, cada uno con 2 variables predictoras más que el modelo ajustado en el paso 1. 3. Compare los modelos generados en los puntos 1 y 2 usando la raíz del error cuadrático medio (RMSE) y el R^2^. ```{r solucion del ejercicio 1.2} #| eval: false #| echo: false # obtener predicciones pred_mtcars <- predict(mod_mtcars) # calcular RMSE rmse_mtcars <- sqrt(mean((mtcars$mpg - pred_mtcars)^2)) # calcular R^2^ r2_mtcars <- summary(mod_mtcars)$r.squared ``` ::: # Evaluación para variables respuesta categóricas Para demostrar las métricas de evalución para estos modelos utilizaremos un modelo de regresión con los datos del titanic similar al del tutorial de regresión logística. ```{r} # cargar datos data("Titanic") # dar formato con una observacion por fila datos_titanic <- as.data.frame(Titanic) datos_titanic <- datos_titanic[datos_titanic$Freq > 0, ] datos_tab_titanic <- do.call(rbind, lapply(1:nrow(datos_titanic), function(x) datos_titanic[rep(x, datos_titanic$Freq[x]),])) # Ajustar el modelo de regresión logística modelo_log <- glm(Survived ~ Class + Sex, data = datos_tab_titanic, family = binomial) ``` Ahora podemos predecir la variable respuesta con este modelo: ```{r} # predecir valores pred_vals <- predict(object = modelo_log, newdata = datos_tab_titanic, type = "response") # binarizar pred_cat <- ifelse(pred_vals > 0.5, "Yes", "No") ``` ## Matriz de confusión La matriz de confusión organiza las predicciones y los resultados observados en categorías que indican si el valor predicho es una ocurrencia del evento (positivo vs negativo) y si coincide o no con el valor observado (verdadero vs falso). ```{r} #| echo: false outcomes <- c("Verdadero Positivos\n(VP)", "Falso Positivos\n(FP)", "Falso Negativos\n(FN)", "Verdadero Negativos\n(VN)") # Crear los datos para la matriz de confusión confusion_matrix <- data.frame( Predicted = rep(c("Positivo", "Negativo"), each = 2), Actual = rep(c("Positivo", "Negativo"), 2), Count = c(50, 10, 5, 35), # Ejemplo de números Outcome = outcomes ) cols <- mako(2, begin = 0.4, end = 0.9, alpha = 0.7)[c(1, 2, 2, 1)] names(cols) <- outcomes # Crear el gráfico con ggplot ggplot(confusion_matrix, aes(x = Actual, y = Predicted, fill = Outcome)) + geom_tile(color = "white", linewidth = 0.7) + geom_text(aes(label = Outcome), color = "black", size = 5) + scale_fill_manual(values = cols, guide = "none") + labs( x = "Valores observados", y = "Valores predichos" ) + theme_minimal() + theme( axis.title = element_text(size = 14), axis.text = element_text(size = 12), plot.title = element_text(hjust = 0.5, size = 16), panel.grid = element_blank() ) ``` Esta matriz la podemos estimar usando la función `confusionMatrix` del paquete `caret`: ```{r} # hacer la matriz de confusion mat_conf <- confusionMatrix(factor(pred_cat), factor(datos_tab_titanic$Survived)) # imprimir resultado mat_conf$table ``` Podemos graficar esta matriz con un gradiente de colores donde los valores esten dado como proporciones, lo cual mas fácil su interpretación: ```{r} # convertir a data frame conf_df <- as.data.frame(mat_conf$table) # agregar totales por categoria conf_df$total <- sapply(conf_df$Reference, function(x) sum(datos_tab_titanic$Survived == x)) # calcular proporciones conf_df$proportion <- conf_df$Freq / conf_df$total # graficar ggplot(conf_df, aes(x = Reference, y = Prediction, fill = proportion)) + geom_tile() + coord_equal() + scale_fill_distiller(palette = "Greens", direction = 1) + geom_text(aes(label = round(proportion, 2)), color = "black", size = 3) + labs(x = "Observado", y = "Predicho", fill = "Proporción") ``` Varias métricas se derivan de los valores contenidos en la matriz de confusion. ## Exactitud (Accuracy) La exactitud es la proporción de predicciones correctas en relación con el total de predicciones. La fórmula es: <center> $\frac{VP + VN}{VP + VN + FP + FN}$ </center> Esta también es estimada por la función `confusionMatrix`: ```{r} mat_conf$overall["Accuracy"] ``` ## Precisión (Precision) y Sensibilidad (Sensitivity o Recall) La precisión es la proporción de verdaderos positivos entre todas las predicciones positivas: <center> $\frac{VP}{VP + FP}$ </center> La sensibilidad es la proporción de verdaderos positivos entre todas las observaciones positivas: <center> $\frac{VP}{VP + FN}$ </center> Ambas son estimadas por la función `confusionMatrix`: ```{r} mat_conf$byClass[c("Precision", "Sensitivity")] ``` ## Índice F1 (o F-Score) El índice F1 es la media armónica entre la precisión y la sensibilidad: <center> $2 \times \frac{\text{Precisión} \times \text{Sensibilidad}}{\text{Precisión} + \text{Sensibilidad}}$ </center> ```{r} precision <- mat_conf$byClass["Precision"] sensitivity <- mat_conf$byClass["Sensitivity"] f1_score <- (2 * precision * sensitivity) / (precision + sensitivity) f1_score ``` Tambien lo podemos extraer del resultado de la función `confusionMatrix`: ```{r} mat_conf$byClass["F1"] ``` ## Área Bajo la Curva ROC (AUC-ROC) El AUC métrica originalmente diseñada para evaluar la precision en la predicción de variables binarias. Esta mide la capacidad del modelo para discriminar entre clases. Se calcula basándose en la curva ROC (Receiver Operating Characteristic). La curva ROC compara la tasa de verdaderos positivos contra la tasa de falsos positivos a diferentes umbrales.La curva ROC es una representación gráfica que evalúa el desempeño de un modelo de clasificación binaria. La curva muestra la relación entre la Tasa de Verdaderos Positivos (True Positive Rate, TPR), también llamada sensibilidad, y la Tasa de Falsos Positivos (False Positive Rate, FPR) para diferentes umbrales de decisión. Los umbrales de decisión son los valores que determinan si una observación se clasifica como positiva o negativa. Podemos ejemplificar su uso con el modelo de regresión logística ajustado mas arriba: ```{r} #| message: false # estimar las probabilidades pred_vals <- predict(object = modelo_log, newdata = datos_tab_titanic, type = "response") # Curva ROC roc_curve <- roc(datos_tab_titanic$Survived, pred_vals, smooth = TRUE) # Extraer los datos de la curva ROC roc_data <- data.frame( tpr = roc_curve$sensitivities, # True Positive Rate (Sensibilidad) fpr = 1 - roc_curve$specificities # False Positive Rate (1 - Especificidad) ) # Graficar la curva ROC con ggplot2 ggplot(roc_data, aes(x = fpr, y = tpr)) + geom_line(color = "tomato") + geom_abline(linetype = "dashed", color = "gray") + # Línea diagonal labs( x = "Tasa de Falsos Positivos (FPR)", y = "Tasa de Verdaderos Positivos (TPR)" ) + scale_x_continuous(limits = c(0, 1), expand = c(0, 0)) + # Limitar el eje X entre 0 y 1 scale_y_continuous(limits = c(0, 1), expand = c(0, 0)) + # Limitar el eje Y entre 0 y theme_classic() + theme(panel.grid = element_blank()) +# Remover cuadrícula # add polygon showin area under the curve geom_polygon(data = roc_data, aes(x = fpr, y = tpr), fill = "tomato", alpha = 0.3) + # add polygon coloring lower diagonal geom_polygon(data = data.frame(x = c(0, 1, 1), y = c(0, 0, 1)), aes(x = x, y = y), fill = "tomato", alpha = 0.3) + annotate("text", x = 0.5, y = 0.5, label = paste("AUC = ", round(pROC::auc(roc_curve), 2)), size = 6) ``` Podemos estimar el área bajo la curva con la función `auc` del paquete `pROC`: ```{r} # estimar AUC pROC::auc(roc_curve) ``` Valores mas cercanos a 1 indican un mejor desempeño del modelo en discriminar entre clases tomando en cuenta la incertidumbre relacionada a la escogencia de un umbral de decisión. En los modelos que predicen más de dos categorías (clasificación multiclase), no se puede aplicar la curva ROC directamente como en el caso binario. Aun así, hay métodos extendidos para calcular una métrica similar a la AUC para clasificación multiclase. Esta métrica se puede extender para modelos que predicen multiples clases. Por ejemplo el metodo de "uno vs el resto" (one vs rest) calcula una curva ROC para cada categoría de manera binaria, considerando una clase como "positiva" y las demás como "negativas". Luego, se calcula el AUC para cada una de estas curvas ROC y se promedia. Para ejemplificar el uso de este método usaremos una regresión multinomial para predecir la especie de las observaciones del juego de datos iris: ```{r} # Ajustar el modelo de regresión multinomial modelo_multinom <- multinom(Species ~ ., data = iris) # Predecir las probabilidades para cada clase pred_prob <- predict(modelo_multinom, newdata = iris, type = "prob") # Aplicar multiclass.roc (One-vs-Rest) roc_multiclass <- multiclass.roc(iris$Species, pred_prob) # Ver resultados roc_multiclass ``` ## Pérdida logarítmica (log-loss) Mide el rendimiento de un modelo de clasificación al penalizar predicciones erróneas basadas en la probabilidad que asigna el modelo a cada clase. <center> $-\frac{1}{n} \sum_{i=1}^{n} \left[ y_i \log(p_i) + (1 - y_i) \log(1 - p_i) \right]$ </center> Donde ${y_i}$ es el valor real (0 o 1) y ${p_i}$ es la probabilidad asignada a la clase positiva. Es útil cuando se desea medir la confianza del modelo en sus predicciones y no solo si acierta o falla. Aunque puede ser más complejo de interpretar en comparación con otras métricas. ```{r} # Convertir las etiquetas reales a una matriz de indicadores (one-hot encoding) # note que Species es un factor y_obs <- model.matrix(~Species-1, data = iris) # Calcular el log-loss log_loss_value <- logLoss(y_obs, pred_prob) log_loss_value ``` El valor del pérdida logarítmica representa qué tan bien o mal el modelo predice las clases correctas. Un valor más bajo (mas cerano a 0) indica un mejor ajuste, mientras que un valor alto sugiere que el modelo está asignando bajas probabilidades a las clases correctas. El valor de pérdida logarítmica puede tomar cualquier valor real no negativo, y en general se encuentra entre 0 y ∞. ::: {.alert .alert-info} ## Ejercicio 2 {.unnumbered} Usaremos los datos de Cortez y Silva (2008) sobre el desempeño académico de estudiantes en la clase de matemáticas. Los datos contienen 395 observaciones y 31 variables. La variable respuesta es `G3` (calificación final), y las variables predictoras son las demás variables del conjunto de datos. Pueden econtrar una [descripción detallada de los datos aqui](https://archive.ics.uci.edu/dataset/320/student+performance). El siguiente código carga los datos y remueve las variables `G1` y `G2` (calificaciones parciales) para evitar multicolinealidad: ```{r} # leer datos datos_mate <- read.csv("https://raw.githubusercontent.com/maRce10/aprendizaje_estadistico_2024/refs/heads/master/data/student/student-mat.csv", sep = ";") # remover G1 y G2 (calificaciones parciales) datos_mate$G1 <- datos_mate$G2 <- NULL str(datos_mate) ``` 1. Ajuste un modelo de regresión lineal con los datos de `datos_mate` donde la variable respuesta sea `G3`. Utilice su conocimiento del sistema (como estudiante) para escoger las variables predictoras que usted espere tengan un mayor efecto sobre el desempeño académico. Calcule el R^2^ para este modelo. ```{r solucion del ejercicio 1} #| eval: false #| echo: false # ajustar modelo mod_mate <- lm(G3 ~ school + sex + age + address + famsize + Pstatus, data = datos_mate) # get summary summary(mod_mate) ``` 2. Ajuste un modelo de regresión multinomial análogo al del punto anterior (e.g. mismos predictores y misma respuesta). ```{r} #| eval: false #| echo: false mod_mult <- multinom(G3 ~ school + sex + age + address + famsize + Pstatus, data = datos_mate) ``` 3. Calcule la matriz de confusión, la exactitud, el área bajo la curva y pérdida logarítmica para el modelo del punto anterior. ```{r} #| eval: false #| echo: false # predecir probabilidades pred_prob <- predict(mod_mult, newdata = datos_mate, type = "probs") # calcular log-loss log_loss_value <- logLoss(datos_mate$G3, pred_prob) log_loss_value # convertir respuesta a factor datos_mate$G3_f <- as.factor(datos_mate$G3) # Convertir las etiquetas reales a una matriz de indicadores (one-hot encoding) y_obs <- model.matrix(~G3_f-1, data = datos_mate) # Calcular el log-loss log_loss_value <- logLoss(y_obs, pred_prob) log_loss_value # hacer la matriz de confusion # predecir categorias pred_cat <- predict(mod_mult, newdata = datos_mate) mat_conf <- confusionMatrix(factor(pred_cat), factor(datos_mate$G3)) mat_conf$table mat_conf$overall["Accuracy"] mat_conf$overall # area bajo la curva roc_multiclass <- multiclass.roc(datos_mate$G3, pred_prob) roc_multiclass ``` El siguiente código binariza la variable respuesta `G3`. Si la calificación es menor a 9, se considera que el estudiante reprobó (`Reprobado`), de lo contrario, se considera que aprobó (`Aprobado`): ```{r} datos_mate$G3_bin <- ifelse(datos_mate$G3 < 9, 0, 1) ``` 4. Utilizando la variable binaria `G3_bin` como respuesta, ajuste un modelo de regresión logística con las mismas variables predictoras que el modelo ajustado en el punto 2. ```{r} #| eval: false #| echo: false # ajustar modelo mod_log <- glm(G3_bin ~ school + sex + age + address + famsize + Pstatus, data = datos_mate, family = binomial) ``` 5.Calcule la matriz de confusión, la exactitud, precisión, indice F1, area bajo la curva y pérdida logarítmica para el modelo del punto anterior. ```{r} #| eval: false #| echo: false # predecir probabilidades pred_prob <- predict(mod_log, newdata = datos_mate, type = "response") # calcular log-loss log_loss_value <- logLoss(datos_mate$G3_bin, pred_prob) log_loss_value # hacer la matriz de confusion # binarizar prediccion pred_cat <- ifelse(pred_prob > 0.5, 1, 0) mat_conf <- confusionMatrix(factor(pred_cat), factor(datos_mate$G3_bin)) mat_conf$table mat_conf$overall["Accuracy"] mat_conf$overall # area bajo la curva roc_curve <- roc(datos_mate$G3_bin, pred_prob, smooth = TRUE) roc_curve ``` 6. Compare los modelos generados en los puntos 2 y 4. ¿Cual modelo parece tener un mejor desempeño? ¿Por qué? 7. Explore la sensibilidad de los modelos multinomial (punto 2) y logístico (punto 4). Para esto puede simplemente imprimir en la consola el resultado de la función `confusionMatrix`. ¿Cómo difieren los valores y la forma en que estos se estructuran entre los 2 modelos? ::: # Referencias Cortez, P., & Silva, A.M. (2008). Using data mining to predict secondary school student performance. ([enlace](https://www.semanticscholar.org/paper/Using-data-mining-to-predict-secondary-school-Cortez-Silva/61d468d5254730bbecf822c6b60d7d6595d9889c)) ------------------------------------------------------------------------ # Información de la sesión {.unnumbered .unlisted} ```{r session info, echo=F} sessionInfo() ```

Objetivo del manual

1 Evaluación para variables respuesta continuas

1.1 Raíz del Error Cuadrático Medio (RMSE)

1.2 Coeficiente de determinación R2

Ejercicio 1

2 Evaluación para variables respuesta categóricas

2.1 Matriz de confusión

2.2 Exactitud (Accuracy)

2.3 Precisión (Precision) y Sensibilidad (Sensitivity o Recall)

2.4 Índice F1 (o F-Score)

2.5 Área Bajo la Curva ROC (AUC-ROC)

2.6 Pérdida logarítmica (log-loss)

Ejercicio 2

3 Referencias

Información de la sesión

1.2 Coeficiente de determinación R²