Capítulo 10 Arboles de decision

Los métodos basados en árboles para regresión y clasificación estratifican o segmentan el espacio predictor en varias regiones. Para hacer una predicción de una observación dada, normalmente utiliza el valor de respuesta promedio de las observaciones de la base de entrenamiento en la región a la que pertenece. En el caso de clasificación se asigna a la categoría mayoritaria dentro del nodo terminal.

10.1 Classification and Regression Tree (CART)

En el caso de árboles de regresión, si \(Y\) es la respuesta y \(X_1\) y \(X_2\) los inputs se parte el espacio \((X_1, X_2)\) en dos regiones, en base a una sola variable (partición horizontal o vertical). Dentro de cada región proponemos como predicción la media muestral de \(Y\).

Se busca elegir la variable y el punto de partición de manera óptima (mejor ajuste global). Es computacionalmente inviable considerar cada posible partición del espacio de atributos en \(J\) regiones. Por lo tanto, toma un enfoque top-down greedy que se conoce como división binaria recursiva. El enfoque es top-down porque comienza en la parte superior del árbol (en cuyo punto todas las observaciones pertenecen a una sola región) y luego divide sucesivamente el espacio predictor; cada división se indica a través de dos nuevas ramas más abajo en el árbol. Es greedy porque en cada paso del proceso de construcción del árbol, la mejor división se hace en ese paso en particular, en lugar de mirar hacia adelante y elegir una división que conducirá a un mejor árbol en algún paso futuro.

El panel izquierdo de la Figura 10.1 muestra un árbol de regresión para predecir el logaritmo del salario (en miles de dólares) de un jugador de béisbol, basado en la cantidad de años que ha jugado en las ligas mayores y la cantidad de hits que hizo en el año anterior. En un nodo interno dado, la etiqueta (de la forma \(X_j < t_k\)) indica la rama izquierda que sale de esa división, y la rama de la derecha corresponde a \(X_j \ge t_k\). Por ejemplo, la división en la parte superior del árbol da como resultado dos ramas grandes. La rama izquierda corresponde a Years < 4,5, y la rama derecha corresponde a Years >= 4,5.19 El árbol tiene dos nodos internos y tres nodos terminales u hojas. El número en cada hoja es la media de la variable de respuesta de las observaciones que caen allí. Por ejemplo, la predicción para el nodo terminal de la izquierda es \(e^{5,107} \times 1.000 = \$165.174\). El panel derecho la Figura 10.1 muestra las regiones en función de Years y Hits.

Arbol de regresión

Figura 10.1: Arbol de regresión

Notar:

  • Cada región tiene su propio modelo.

  • Ciertas variables importan en determinadas regiones y no en otras (Hits).

Dado \(Y\) y \(X\) un vector de \(p\) variables con \(n\) observaciones el algoritmo busca determinar cuál variable usar para la partición y que punto de esa variable usar para la partición. Si \(j\) es la variable de partición y el punto de partición es \(s\), se definen los siguientes semiplanos:

\[\begin{align*} R_1(j,s) = & {X \mid X_j < s} \\ R_2(j,s) = & {X \mid X_j \ge s} \end{align*}\]

Se trata de buscar la variable de partición \(X_j\) y el punto de partición \(s\) que resuelvan (minimizar el \(EMC\) en cada región):

\[\begin{equation} \tag{10.1} \sum_{i: x_i \in R_1(j,s)} (y_i - \hat{y}_{R_1})^2 + \sum_{i: x_i \in R_2(j,s)} (y_i - \hat{y}_{R_2})^2 \end{equation}\]

Donde \(\hat{y}_{R_1}\) y \(\hat{y}_{R_2}\) es el promedio de la respuesta en las regiones \(1\) y \(2\), respectivamente. Para cada variable y punto de partición, la minimización interna se corresponde con la media dentro de cada región.20

¿Cuándo parar de realizar divisiones?

Un árbol demasiado extenso sobreajusta (overfit) los datos. Pero dado que el proceso es secuencial y cada corte no mira lo que puede suceder después, si se detiene el proceso demasiado pronto se puede perder un “gran” corte más abajo. Prunning: ajustar un árbol grande y luego podarlo (prune) usando un criterio de cost-complexity.

Weakest link pruning

Un subárbol \(T \in T_0\) es un árbol que se obtiene colapsando los nodos terminales de otro árbol (cortando ramas).

Cost-complexity del árbol \(T\):

\[\begin{equation} \tag{10.2} C_{\alpha}(T) = \sum_{m=1}^{|T|} n_mQ_m(T) + \alpha[T] \end{equation}\]

con \(Q_m(T) = \frac{1}{n_m} \sum_{x_i \in R_m} (y_i - \hat{c}_m)^2\) (impureza) y \(n_m\) cantidad de observaciones en cada partición.

El primer término mide el (mal) ajuste y el segundo la complejidad. Cuando \(\alpha = 0\), entonces el subárbol \(T\) simplemente será igual a \(T_0\), porque en ese caso (10.2) solo mide el error de entrenamiento. Sin embargo, a medida que \(\alpha\) aumenta, hay que pagar un costo por tener un árbol con muchos nodos terminales, por lo que (10.2) tenderá a minimizarse para un subárbol más pequeño.21 Entonces el objetivo es, para un \(\alpha\) dado, encontrar la poda óptima que minimiza \(C_{\alpha}(T)\).

Para cada \(\alpha\) hay un único subárbol \(T_{\alpha}\) que minimiza \(C_{\alpha}(T)\). El mecanismo de búsqueda de \(T_{\alpha}\) (poda óptima dado \(\alpha\)) de weakest link consite en eliminar sucesivamente las ramas que producen el mínimo incremento en \(\sum_{m=1}^{[T]} n_mQ_m(T)\) (impureza). Recordar que un árbol grande aumenta la varianza, por lo tanto, se colapsa la partición menos necesaria. Un árbol más pequeño con menos divisiones (es decir, menos regiones \(R_1,...,R_J\)) tiene menor varianza y es más fácil de interpretar a costa de un pequeño sesgo.

El proceso eventualmente colapsa en el nodo inicial, pero pasa por una sucesión de árboles, desde el más grande, hasta el más chico, por el proceso de weakest link pruning. El árbol óptimo \(T_{\alpha}\) pertenece a esta sucesión.

Classification tree

Un árbol de clasificación es muy similar a un árbol de regresión, excepto que se utiliza para predecir una respuesta cualitativa en lugar de una cuantitativa. Recordar que para un árbol de regresión, la respuesta predicha para una observación esta dada por la respuesta media de las observaciones de entrenamiento que pertenecen al mismo nodo terminal. En contraste, para un árbol de clasificación, predice que cada observación pertenece a la clase que ocurre más comúnmente en las observaciones de entrenamiento en la región a la que pertenece. Se basa en el error de clasificación o índice de Gini (pureza), análogo a \(EMC\) en un árbol de regresión.

10.2 Métodos de Ensamble: Bagging, Random Forest y Boosting

Los métodos de ensamble constituyen una de las principales innovaciones en el campo del aprendizaje estadístico y machine learning. Su fundamento radica en la idea de combinar múltiples modelos débiles —en general, árboles de decisión— con el objetivo de reducir la varianza, el sesgo o ambos simultáneamente, y mejorar la capacidad predictiva del conjunto respecto a cada modelo individual.

10.3 Bagging

Ventajas y desventajas de \(CART\):

  • Forma inteligente de representar no linealidades.

  • Arriba quedan las variables más relevantes entonces es fácil de comunicar. Reproduce proceso decisorio humano.

  • Si la estructura es lineal, \(CART\) no anda bien.

  • Poco robusto, variaciones en los datos modifican el resultado.

Un método de ensemble es un enfoque que combina muchos modelos simples en uno único y potencialmente muy poderoso. Los modelos simples se conocen como modelos de aprendizaje débil, ya que por sí mismos pueden generar predicciones mediocres.

Una posible solución es el bootstrap aggregation que consiste en tomar como predicción el promedio de las predicciones bootstrap.22

\[\begin{equation} \tag{10.3} \hat{f}_{bag} = \frac{1}{B} \sum_{b=1}^{B} \hat{f}^{*b}(x) \end{equation}\]

Esta idea se basa en que la varianza del promedio es menor que la de una predicción sola. Bajo independencia si \(V(x) = \sigma^2\) entonces \(V(\overline{x}) = \frac{\sigma^2}{n}\). Pero existe el problema que si hay un predictor fuerte (siempre va a ser seleccionado primero), los distintos árboles son muy similares entre sí por lo que habrá alta correlación.

10.4 Random Forest

El algoritmo de Random Forest, propuesto por Breiman (2001), constituye una extensión del bagging que introduce un elemento adicional de aleatoriedad en el proceso de crecimiento de los árboles. Además del remuestreo bootstrap, en cada nodo del árbol se selecciona aleatoriamente un subconjunto de variables predictoras de tamaño \(m < p\) (donde \(p\) es el número total de predictores). La división se realiza únicamente entre esas \(m\) variables candidatas, reduciendo así la correlación entre árboles y, por ende, la varianza del ensamble.

Formalmente, cada árbol \(f^{(b)}(x)\) se construye a partir de un subconjunto de observaciones y un subconjunto aleatorio de predictores. La predicción final es nuevamente un promedio o votación:

\[ \hat{f}_{\text{RF}}(x) = \frac{1}{B} \sum_{b=1}^{B} f^{(b)}(x) \]

Una ventaja crucial del Random Forest es que no requiere poda (pruning), ya que la agregación compensa el sobreajuste de los árboles individuales. Además, el método permite evaluar la importancia de las variables mediante medidas basadas en la reducción promedio de la impureza o en la pérdida de precisión al permutar cada predictor. Desde una perspectiva teórica, el Random Forest combina estabilidad y bajo sesgo.

10.5 Boosting

El boosting persigue el mismo objetivo general —reducir el error de predicción— pero adopta un enfoque secuencial. En lugar de construir múltiples modelos independientes, los árboles se ajustan de manera iterativa, de modo que cada nuevo modelo corrige los errores cometidos por los anteriores. El procedimiento fue formalizado por Freund y Schapire (1997) bajo el algoritmo AdaBoost (Adaptive Boosting), y posteriormente generalizado por Friedman (2001) en el marco de la Gradient Boosting Machine (GBM).

En su forma más general, el boosting minimiza una función de pérdida \(L(y, F(x))\) mediante un proceso aditivo:

\[ F_M(x) = \sum_{m=1}^{M} \gamma_m h_m(x) \]

donde cada \(h_m(x)\) representa un árbol de decisión de poca profundidad (un weak learner), y los pesos \(\gamma_m\) se eligen para minimizar el error residual de la iteración anterior. Este enfoque puede interpretarse como una aproximación de gradiente descendente en el espacio de funciones.

El boosting tiende a reducir el sesgo sin incrementar excesivamente la varianza, lo que lo hace especialmente útil cuando el modelo base es demasiado simple. Sin embargo, su naturaleza secuencial lo vuelve más sensible al sobreajuste, por lo que requiere un control cuidadoso de hiperparámetros como la tasa de aprendizaje (learning rate), la profundidad máxima de los árboles y el número de iteraciones.

Diversas variantes modernas —como XGBoost Chen y Guestrin (2016), LightGBM Ke et al. (2017) y CatBoost Prokhorenkova et al. (2018)— han optimizado el procedimiento original, incorporando regularización explícita, manejo eficiente de datos escasos y paralelización avanzada.

10.5.1 Ada Boost

Adaptive boosting. Dadas las siguientes definiciones:

  • \(y \in -1,1\)

  • Clasificador: \(\hat{y} = G(X)\)

  • Error de predicción = \(\frac{1}{N} \sum_{i=1}^{N}I(y_i \neq G(x_i))\)

Se describe el proceso:

  1. Comienza con con pesos \(w_i = \frac{1}{N}\)

  2. Para \(m = 1,...,M\):

  • Calcula una predicción

  • Calcula el error de predicción agregado

  • Calcula \(\alpha_m = ln[\frac{1 - err_m}{err_m}]\)

  • Actualiza los ponderadores \(w_i\) \(\leftarrow\) \(w_ic_i\)

con \(c_i =\) exp \([\alpha_m \underbrace{I(y_i \neq G(x_i))}_{{0/1}}]\)

  1. Output: \(G(x) =\) sgn \([\sum_{m=1}^{M} \alpha_m G_m(x)]\) (signo del promedio).

Si \(i\) estuvo correctamente predicha, \(c_i = 1\), entonces no hay ajuste. Caso contrario, \(c_i = e^{\alpha_m} = \frac{1 - err_m}{err_m} > 1\). Notar que si siempre se predice la clase mayoritaria la tasa de error nunca puede ser mayor al \(50\%\) y por eso la expresión anterior es mayor a \(1\).

En cada paso el método da más importancia relativa a las observaciones mal predichas. Paso final: promedio ponderado de predicciones en cada paso. La Figura 10.2 muestra el proceso gráficamente.

Ada boost

Figura 10.2: Ada boost

10.6 Aplicacion práctica

El objetivo de esta sección es predecir la probabilidad de default utilizando LightGBM aplicado al dataset Give Me Some Credit (de Kaggle).

10.6.1 Carga de librerías

library(tidyverse)
library(readr)
library(readxl)
library(lightgbm)
library(caret)
library(pROC)
library(data.table)
library(ggplot2)
library(ParBayesianOptimization)
library(shapviz)

10.6.2 Carga de datos

train_raw <- read_csv("./data/GiveMeSomeCredit/cs-training.csv", show_col_types=FALSE)
test_raw  <- read_csv("./data/GiveMeSomeCredit/cs-test.csv",     show_col_types=FALSE)

10.6.3 Renombramos identificadores

train_raw <- train_raw |> rename(Id = `...1`)
test_raw  <- test_raw  |> rename(Id = `...1`)

10.6.4 Eliminamos registros con edad no razonable

train_raw <- train_raw |> dplyr::filter(age > 18)
test_raw  <- test_raw  |> dplyr::filter(age > 18)

Nota: hasta acá estamos haciendo limpieza estructural (filtros y renombres) aplicada de la misma forma a training y test.

10.6.5 Identificamos variables numéricas

num_vars <- train_raw |>
  dplyr::select(-Id, -SeriousDlqin2yrs) |>
  dplyr::select(where(is.numeric)) |>
  names()

10.6.6 Parámetros de imputación

Calculados solo en training para evitar data leakage (cuando un modelo de machine learning accede a información del futuro o de los datos de test durante el entrenamiento, lo que le permite “hacer trampa” y obtener resultados más precisos).

En caso que sea necesario estandarizar las variables, requerido por algunos algoritmos como redes neuronales para que todas las variables tengan la misma escala, se procede de la misma manera.

medianas <- train_raw |>
  summarise(across(all_of(num_vars), ~ median(.x, na.rm=TRUE)))

p99 <- train_raw |>
  summarise(across(all_of(num_vars), ~ quantile(.x, 0.99, na.rm=TRUE)))

10.6.7 Función de preprocesamiento

Usamos esta función solo como ejemplo de feature engineering sencillo dado que las posibilidades son infinitas: ratios, logs, rezagos, one-hot encoding (convierte una variable categórica en varias columnas binarias (0/1) una por cada categoría posible), etc.

aplicar_preprocesamiento <- function(df, p99, vars){
    df2 <- df
    for(v in vars){
        lim <- p99[[v]][[1]]
        x <- df2[[v]]
        # Winsorización al percentil 99
        x[x > lim] <- lim
        df2[[v]] <- x
    }
    df2
}

10.6.8 Aplicamos preprocesamiento a ambos datasets

train_clean <- aplicar_preprocesamiento(train_raw, p99, num_vars)
test_clean  <- aplicar_preprocesamiento(test_raw,  p99, num_vars)

10.6.9 División en entrenamiento y test (etiquetado)

set.seed(2025)
idx_train   <- createDataPartition(train_clean$SeriousDlqin2yrs, p=0.7, list=FALSE) # divide 70/30
train_model <- train_clean[idx_train, ]
test_labeled <- train_clean[-idx_train, ]

Luego de hacer la partición, se seleccionan las variables para construir las bases de train y test.

X_train <- train_model |> dplyr::select(-Id, -SeriousDlqin2yrs) |> as.matrix()
y_train <- train_model$SeriousDlqin2yrs 

X_test_labeled <- test_labeled |> dplyr::select(-Id, -SeriousDlqin2yrs) |> as.matrix()
y_test_labeled <- test_labeled$SeriousDlqin2yrs

Se convierten los datos a formato dataset de LightGBM.

dtrain <- lgb.Dataset(X_train, label = y_train)

Balanceo por clase:

n_pos <- sum(y_train==1)
n_neg <- sum(y_train==0)
scale_pos_weight <- n_neg / n_pos
scale_pos_weight
## [1] 13.84309
pct_pos <- (n_pos / (n_pos + n_neg)) * 100
pct_pos; 100 - pct_pos
## [1] 6.737143
## [1] 93.26286

Hiperparámetros base:

params_base <- list(
  objective = "binary",                # Define la función de pérdida para clasificación binaria.
  metric = "auc",                      # Establece la métrica de evaluación: área bajo la curva ROC.
  boosting = "gbdt",                   # Define el tipo de boosting a usar: árboles de decisión con gradiente.
  num_leaves = 31,                     # Número máximo de hojas en cada árbol: controla la complejidad del modelo.
  max_depth = -1,                      # Profundidad máxima de cada árbol; -1 significa sin límite explícito.
  min_data_in_leaf = 100,              # Número mínimo de observaciones requeridas en una hoja: evita sobreajuste en hojas pequeñas.
  feature_fraction = 0.8,              # Fracción de características usadas en cada iteración de árbol: para regularización por columnas.
  bagging_fraction = 0.8,              # Fracción de datos usada (submuestra) en cada iteración de árbol: para regularización por filas.
  bagging_freq = 1,                    # Frecuencia de bagging (cada cuántas iteraciones usar muestreo): regula el submuestreo.
  learning_rate = 0.05,                # Tasa de aprendizaje (shrinkage) para actualizar ganancias de cada árbol: elocidad de ajuste.
  scale_pos_weight = scale_pos_weight, # Peso de la clase positiva en la función objetivo: útil para clases desbalanceadas.
  feature_pre_filter = FALSE,          # Indica si se filtrarán automáticamente características irrelevantes antes del entrenamiento.
  verbose = -1                         # Nivel de salida de registro: -1 suprime la mayoría de mensajes de log.
)

10.6.11 Selección de hiperparámetros (Optimización Bayesiana)

Se implementa una búsqueda basada en optimización bayesiana para maximizar el Área Bajo la Curva ROC (AUC) mediante validación cruzada de 5 folds.

La optimización bayesiana es un método diseñado para encontrar combinaciones óptimas de hiperparámetros cuando evaluar el modelo es costoso —por ejemplo, cuando cada intento requiere entrenar un algoritmo mediante cross-validation. En lugar de probar muchas combinaciones al azar o mediante grillas, esta metodoligía construye un modelo probabilístico (surrogate) de la función objetivo: típicamente un Gaussian Process que aproxima la relación entre los hiperparámetros y el desempeño del modelo (AUC, RMSE, log-loss, etc.). Con esa “función aproximada”, el método usa una función de adquisición (como Upper Confidence Bound o Expected Improvement) que guía la búsqueda hacia regiones prometedoras, equilibrando explotación (probar valores que el surrogate predice como buenos) y exploración (probar valores donde la incertidumbre es alta). El proceso es iterativo: se evalúan algunos puntos iniciales aleatorios, se entrena el surrogate, se elige el siguiente punto según la adquisición, se actualiza el surrogate y así sucesivamente. De este modo, la optimización bayesiana necesita muchas menos evaluaciones reales que una búsqueda tradicional para converger hacia hiperparámetros de alto rendimiento, lo que la vuelve especialmente útil en modelos complejos o costosos de entrenar.

Hiperparámetros base:

params_base_fixed <- list(
  objective = "binary",       # Define la función de pérdida para clasificación binaria.
  metric = "auc",             # Establece la métrica de evaluación: área bajo la curva ROC.
  boosting = "gbdt",          # Define el tipo de boosting a usar: árboles de decisión con gradiente.
  is_unbalance = TRUE,
  feature_pre_filter = FALSE, # Descarta features automáticamente antes de entrenar.
  force_row_wise = TRUE,      # Mejora la estabilidad.    
  verbose = -1                # Nivel de salida de registro: -1 suprime la mayoría de mensajes de log.
)
set.seed(1234)
# Función objetivo (que se busca max o min dependiendo del caso)
bayes_scoring_fun <- function(num_leaves, min_data_in_leaf, feature_fraction, learning_rate) {
  
  pars <- modifyList(params_base_fixed, list(
    num_leaves = as.integer(num_leaves),
    min_data_in_leaf = as.integer(min_data_in_leaf),
    feature_fraction = feature_fraction,
    learning_rate = learning_rate,
    num_threads = 4))

  cv_res <- lgb.cv(
    params = pars,
    data = dtrain,
    nrounds = 1000,
    nfold = 5,
    stratified = TRUE,
    early_stopping_rounds = 50,
    verbose = -1
  )
  
  best_auc <- max(unlist(cv_res$record_evals$valid$auc$eval))
  
  return(list(Score = best_auc, best_iter = cv_res$best_iter))
}

# Espacio de búsqueda
bounds <- list(
  num_leaves       = c(20L, 127L),
  min_data_in_leaf = c(20L, 300L),
  feature_fraction = c(0.4, 1.0),
  learning_rate    = c(0.01, 0.15)
)

# Ejecución de la optimizacion
opt_res <- bayesOpt(
  FUN = bayes_scoring_fun, # Funcion objetivo
  bounds = bounds,         # Espacio de busqueda
  initPoints = 15,         # Cantidad de puntos iniciales aleatorios
  iters.n = 50,            # Cantidad de iteraciones adicionales
  acq = "ucb"              # función de adquisición (Upper Confidence Bound) / 'ei' / 'poi'
)
## 
## Running initial scoring function 15 times in 1 thread(s)...

Hiperparámetros óptimos identificados:

best_pars <- getBestPars(opt_res)
print(best_pars)
## $num_leaves
## [1] 20
## 
## $min_data_in_leaf
## [1] 257
## 
## $feature_fraction
## [1] 0.4968178
## 
## $learning_rate
## [1] 0.01

10.6.12 Validación en test etiquetado

Entrenamos modelo candidato con los resultados de Grid Search:

params_grid <- params_base
params_grid$num_leaves       <- best$num_leaves[[1]]
params_grid$min_data_in_leaf <- best$min_data_in_leaf[[1]]
params_grid$feature_fraction <- best$feature_fraction[[1]]
params_grid$learning_rate    <- best$learning_rate[[1]]

best_iter_grid <- best$best_iter[[1]]

candidate_grid <- lgb.train(
  params=params_grid,
  data=dtrain,
  nrounds=best_iter_grid
)

Entrenamos modelo candidato con los resultados de la búsqueda bayesiana:

params_bayes <- modifyList(params_base_fixed, 
    list(num_leaves = best_pars$num_leaves,
    min_data_in_leaf = best_pars$min_data_in_leaf,
    feature_fraction = best_pars$feature_fraction,
    learning_rate = best_pars$learning_rate))

candidate_bayes <- lgb.train(
  params = params_bayes,
  data = dtrain,
  nrounds = 1000
)

10.6.13 Evaluación de los resultados

En el modelo de clasificación aplicado a la base Give Me Some Credit, la matriz de confusión muestra un patrón típico de problemas de la tasa de default con clases muy desbalanceadas (solo ~6,7% de los casos son morosos). El modelo presenta una AUC de 0.87, lo cual indica que es muy bueno para ordenar a los individuos según su probabilidad de incumplimiento: si tomamos un deudor y un no deudor al azar, el modelo asigna un score mayor al deudor el 87% de las veces. Sin embargo, al convertir esos scores en una predicción binaria usando un umbral fijo de 0.5 —que no es ideal en problemas de baja prevalencia— la matriz de confusión revela limitaciones importantes: la sensibilidad (77,6%) muestra que el modelo identifica correctamente a la mayoría de los deudores reales, pero la precisión (Pos Pred Value = 21,8%) es baja, lo que implica que muchos de los casos marcados como default son falsos positivos. A su vez, la accuracy global (0.80) es inferior al clasificador trivial que siempre predice no default (No Information Rate = 0.934), porque la clase positiva es muy rara. Esto ilustra que en contextos de riesgo de crédito, una AUC alta no garantiza buenas predicciones binarias si no se ajusta el umbral de decisión. En estos problemas, métricas como la sensibilidad, la precisión, el F1-score o la balanced accuracy son más informativas que la accuracy tradicional. En práctica, el umbral debe elegirse optimizando un criterio económico (costos relativos de FP y FN) o estadístico (maximizar F1).

# se aplica el modelo entrenado candidate_grid a los datos que separamos para test (X_test_labeled) 
# Se obtienen probabilidades (pred_prob_grid) en el rango [0, 1]
pred_prob_grid <- predict(candidate_grid, X_test_labeled)

#Se genera una clase binaria (requerido para la matriz de confusion) con umbral igual a 0,5.
pred_class_grid <- ifelse(pred_prob_grid >= 0.5, "1", "0")
real_class_grid <- as.character(y_test_labeled)

# Notar que las variables 'character' se convierten a tipo 'factor'
confusionMatrix(
  factor(pred_class_grid, levels=c("0","1")),
  factor(real_class_grid, levels=c("0","1")),
  positive = "1"
)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction     0     1
##          0 33793   671
##          1  8254  2281
##                                           
##                Accuracy : 0.8017          
##                  95% CI : (0.7979, 0.8053)
##     No Information Rate : 0.9344          
##     P-Value [Acc > NIR] : 1               
##                                           
##                   Kappa : 0.2627          
##                                           
##  Mcnemar's Test P-Value : <2e-16          
##                                           
##             Sensitivity : 0.77270         
##             Specificity : 0.80370         
##          Pos Pred Value : 0.21652         
##          Neg Pred Value : 0.98053         
##              Prevalence : 0.06560         
##          Detection Rate : 0.05069         
##    Detection Prevalence : 0.23412         
##       Balanced Accuracy : 0.78820         
##                                           
##        'Positive' Class : 1               
## 
# Se modifica el umbral a 0,6 para comparar resultados.
pred_class_grid2 <- ifelse(pred_prob_grid >= 0.6, "1", "0")

confusionMatrix(
  factor(pred_class_grid2, levels=c("0","1")),
  factor(real_class_grid, levels=c("0","1")),
  positive = "1"
)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction     0     1
##          0 36490   923
##          1  5557  2029
##                                           
##                Accuracy : 0.856           
##                  95% CI : (0.8527, 0.8592)
##     No Information Rate : 0.9344          
##     P-Value [Acc > NIR] : 1               
##                                           
##                   Kappa : 0.3209          
##                                           
##  Mcnemar's Test P-Value : <2e-16          
##                                           
##             Sensitivity : 0.68733         
##             Specificity : 0.86784         
##          Pos Pred Value : 0.26747         
##          Neg Pred Value : 0.97533         
##              Prevalence : 0.06560         
##          Detection Rate : 0.04509         
##    Detection Prevalence : 0.16858         
##       Balanced Accuracy : 0.77758         
##                                           
##        'Positive' Class : 1               
##

Se aplica el modelo candidato de la búsqueda bayesiana y se calcula la AUC de Grid Search y de la búsquda bayesiana (ahora se utilizan probabilidades):

# se aplica el modelo entrenado candidate_bayes a los datos que separamos para test (X_test_labeled) 
pred_prob_bayes <- predict(candidate_bayes, X_test_labeled)

roc_grid  <- roc(y_test_labeled, pred_prob_grid)
auc_grid  <- auc(roc_grid)
auc_grid
## Area under the curve: 0.8703
roc_bayes <- roc(y_test_labeled, pred_prob_bayes)
auc_bayes <- auc(roc_bayes)
auc_bayes
## Area under the curve: 0.8709

Se utilizan los resultados anteriores para armar una base de datos (en formato ‘long’) para graficar.

# Curva ROC del modelo Grid
df_grid <- data.frame(
  fpr = 1 - roc_grid$specificities,
  tpr = roc_grid$sensitivities,
  model = "Grid Search"
)

# Curva ROC del modelo Bayesiano
df_bayes <- data.frame(
  fpr = 1 - roc_bayes$specificities,
  tpr = roc_bayes$sensitivities,
  model = "Bayesian Search"
)

# Unimos en formato largo
roc_long <- rbind(df_grid, df_bayes)

Curva ROC:

10.6.14 Entrenamiento final

Dado que ahora se busca entrenar un modelo final para luego aplicar a la predicción final (objetivo del ejercicio) se utilizan TODOS los datos disponibles para entrenar con los mejores hiperparámetros identificados.

X_full <- train_clean |> dplyr::select(-Id, -SeriousDlqin2yrs) |> as.matrix()
y_full <- train_clean$SeriousDlqin2yrs
dtrain_full <- lgb.Dataset(X_full, label = y_full)
final_model_grid <- lgb.train(
  params=params_grid,
  data=dtrain_full,
  nrounds=best_iter_grid
)

final_model_bayes <- lgb.train(
  params=params_bayes,
  data=dtrain_full,
  nrounds=1000
)

10.6.15 Análisis de resultados

En LightGBM, la importancia de variables se resume mediante tres métricas complementarias: Gain, Cover y Frequency. Gain es la más informativa y la que suele usarse como referencia principal: mide cuánto contribuye una variable a reducir la función de pérdida a lo largo de todo el modelo. En otras palabras, captura el “poder explicativo” real de cada predictor y su peso en la construcción del score final. Cover indica cuántas observaciones, en promedio, pasan por las divisiones donde interviene una variable; refleja qué tan “amplio” es el segmento del dataset donde esa variable influye. Una variable puede tener bajo Gain pero alto Cover si aparece en splits tempranos que afectan a muchos registros, pero sin mejorar demasiado la pérdida. Finalmente, Frequency simplemente cuenta cuántas veces una variable aparece en un split; es útil para detectar variables muy utilizadas, aunque su aporte marginal pueda ser pequeño. En la práctica, Gain es la métrica dominante, porque resume la mejora real del modelo, mientras que Cover y Frequency se interpretan como medidas auxiliares para entender el rol estructural de cada variable dentro de los árboles.

importance_tbl <- lgb.importance(final_model_grid)
importance_tbl |> head(15)
##                                  Feature        Gain      Cover  Frequency
##                                   <char>       <num>      <num>      <num>
##  1: RevolvingUtilizationOfUnsecuredLines 0.418344921 0.18261028 0.15544872
##  2:              NumberOfTimes90DaysLate 0.159013264 0.11306596 0.03974359
##  3: NumberOfTime30-59DaysPastDueNotWorse 0.150011588 0.10402867 0.05512821
##  4: NumberOfTime60-89DaysPastDueNotWorse 0.072762270 0.09127634 0.03621795
##  5:                                  age 0.049887113 0.11360970 0.13685897
##  6:                            DebtRatio 0.049664173 0.11068942 0.20865385
##  7:                        MonthlyIncome 0.039347900 0.11195220 0.16602564
##  8:      NumberOfOpenCreditLinesAndLoans 0.030908778 0.09007919 0.10929487
##  9:         NumberRealEstateLoansOrLines 0.023280020 0.05987491 0.05416667
## 10:                   NumberOfDependents 0.006779972 0.02281334 0.03846154
lgb.plot.importance(importance_tbl, top_n = 15, measure = "Gain") +
  ggplot2::labs(title = "Feature Importance — Modelo LightGBM",
                subtitle = "Medida: Gain") +
  ggplot2::theme_minimal(base_size = 14)

## NULL

10.6.16 Análisis de interpretabilidad (SHAP)

Los valores SHAP (SHapley Additive exPlanations) son una forma de entender qué aportan las variables dentro de un modelo complejo, como LightGBM. La idea viene de la teoría de juegos cooperativos: imaginar que cada variable es un “jugador” y que la predicción final del modelo es el resultado de un “juego” donde todos colaboran. Los valores SHAP calculan cuál es la contribución marginal promedio de cada variable a través de todas las posibles combinaciones con el resto de las variables. Es decir, miden cuánto cambia la predicción cuando se agrega una variable a un conjunto cualquiera de otras variables. Esto permite explicar predicciones individuales (por ejemplo, por qué para un caso puntual el modelo predijo un valor más alto o más bajo) y al mismo tiempo construir resúmenes globales que muestran las variables más influyentes de manera consistente

shap_obj <- shapviz(final_model_bayes, X_pred = X_test_labeled)

sv_importance(shap_obj, kind = "bar", show_numbers = TRUE, max_display = 15) +
  labs(title = "Importancia Global (Mean |SHAP|)") +
  theme_minimal()

10.6.17 Aplicación final para predicción

Se aplica el modelo final entrenado en la base de datos de entrenamiento completa utilizando los mejores hiperparámetros identificados por cada método.

test_proc <- test_clean

X <- test_proc |>
  dplyr::select(all_of(num_vars)) |> as.matrix()

prob_default_grid  <- predict(final_model_grid, X)
prob_default_bayes <- predict(final_model_bayes, X)

resultados_default <- tibble(
  Id                 = test_proc$Id,
  prob_default_grid  = prob_default_grid,
  prob_default_bayes = prob_default_bayes,
)

head(resultados_default, 10)
## # A tibble: 10 × 3
##       Id prob_default_grid prob_default_bayes
##    <dbl>             <dbl>              <dbl>
##  1     1            0.473              0.484 
##  2     2            0.429              0.388 
##  3     3            0.144              0.179 
##  4     4            0.494              0.515 
##  5     5            0.638              0.596 
##  6     6            0.265              0.305 
##  7     7            0.279              0.321 
##  8     8            0.339              0.319 
##  9     9            0.0317             0.0300
## 10    10            0.924              0.921

Bibliografia

Breiman, Leo. 2001. «Random Forests». Machine Learning 45 (1): 5-32.
Chen, Tianqi, y Carlos Guestrin. 2016. «XGBoost: A Scalable Tree Boosting System». En Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, 785-94.
Freund, Yoav, y Robert E. Schapire. 1997. «A Decision-Theoretic Generalization of On-Line Learning and an Application to Boosting». Journal of Computer and System Sciences 55 (1): 119-39.
Friedman, Jerome H. 2001. «Greedy Function Approximation: A Gradient Boosting Machine». Annals of Statistics 29 (5): 1189-1232.
Ke, Guolin, Qi Meng, Thomas Finley, Taifeng Wang, Wei Chen, Weidong Ma, Qiwei Ye, y Tie-Yan Liu. 2017. «LightGBM: A Highly Efficient Gradient Boosting Decision Tree». En Advances in Neural Information Processing Systems. Vol. 30.
Prokhorenkova, Liudmila, Gleb Gusev, Aleksandr Vorobev, Anna Veronika Dorogush, y Andrey Gulin. 2018. «CatBoost: Unbiased Boosting with Categorical Features». En Advances in Neural Information Processing Systems. Vol. 31.

  1. Al estar arriba, Years es la variable más importante para explicar el salario.↩︎

  2. Recordar que si se quiere predecir una variable aleatoria \(Y\) con una constante \(m\) el mejor predictor es su esperanza, es decir, \(m = E(Y)\).↩︎

  3. Idea similar a Lasso.↩︎

  4. Recordar muestras con remplazo Sección 6.8.2.↩︎