Árboles de decisión

Hasta ahora construimos modelos que ajustan una ecuación: una recta (regresión lineal), una curva en forma de S (regresión logística), o una combinación lineal de variables (LDA, PCA). Todos ellos parten de un supuesto de linealidad, aunque sea de forma implícita.

Hoy vamos a ver un tipo de modelo completamente distinto en su lógica: en lugar de ajustar una ecuación, un árbol de decisión aprende una secuencia de preguntas del tipo “¿esta variable es mayor o menor que tal valor?” para ir separando los datos en grupos cada vez más homogéneos, i.e., el objetivo es agrupar observaciones de tal forma que aquellas dentro de un mismo grupo sean más similares entre sí que con las de otros grupos.

Es un modelo que, como veremos, se puede usar tanto para clasificación (predecir una categoría) como para regresión (predecir un número). Hoy nos enfocaremos en el caso de clasificación.

La idea intuitiva

Antes de tocar datos reales, construyamos un ejemplo de juguete: 45 observaciones, repartidas en 3 clases, descritas por 2 variables numéricas x1 y x2.

set.seed(123)

clase_0 <- data.frame(x1 = rnorm(15, mean = 2, sd = 0.5), x2 = rnorm(15, mean = 1.5, sd = 0.4), clase = "0")
clase_1 <- data.frame(x1 = rnorm(15, mean = 4.5, sd = 0.5), x2 = rnorm(15, mean = 5.2, sd = 0.5), clase = "1")
clase_2 <- data.frame(x1 = rnorm(15, mean = 7, sd = 0.7), x2 = rnorm(15, mean = 3.3, sd = 0.4), clase = "2")

juguete <- rbind(clase_0, clase_1, clase_2)

plot(juguete$x1, juguete$x2,
     col = c("blue", "deeppink", "darkgreen")[as.factor(juguete$clase)],
     pch = 19, xlab = "x1", ylab = "x2",
     main = "3 clases descritas por 2 variables")

legend("topleft", legend = c("Clase 0", "Clase 1", "Clase 2"),
       col = c("blue", "deeppink", "darkgreen"), pch = 19, bty = "n")


Queremos que un modelo aprenda, únicamente a partir de las variables x1 y x2, a qué clase pertenece cada punto nuevo. Un árbol de decisión resuelve esto trazando cortes rectos, uno a la vez, cada uno dividiendo el espacio en dos regiones.

plot(juguete$x1, juguete$x2,
     col = c("blue", "deeppink", "darkgreen")[as.factor(juguete$clase)],
     pch = 19, xlab = "x1", ylab = "x2",
     main = "Cortes que aprendería un árbol de decisión")

abline(v = 3, col = "red", lwd = 2)   # primer corte: x1 <= 3
abline(h = 4, col = "red", lwd = 2)   # segundo corte: x2 <= 4 (solo aplica del lado derecho)

# legend("topleft", legend = c("Clase 0", "Clase 1", "Clase 2"),
#       col = c("blue", "deeppink", "darkgreen"), pch = 19, bty = "n")


Leyendo el gráfico de izquierda a derecha:

  • El primer corte (x1 = 3, línea vertical) ya separa casi perfectamente a la clase 0 (izquierda) del resto (derecha). Este sería el primer corte que el árbol elegiría, porque es el que mejor separa las clases con un solo corte.
  • El segundo corte (x2 = 4, línea horizontal) actúa sobre lo que quedó del lado derecho, y separa la clase 1 (arriba) de la clase 2 (abajo).

Esta secuencia de cortes es exactamente lo que se representa como un árbol: cada corte es un nodo de decisión, y cada región final donde ya no hace falta cortar más es una hoja, que contiene la predicción final para cualquier punto que caiga ahí.

                [x1 < 3?]
               /          \
             sí            no
             /               \
        Clase 0           [x2 < 4?]
                          /          \
                        sí            no
                        /               \
                   Clase 2           Clase 1

¿Cómo decide el árbol dónde cortar? En cada paso, prueba muchas variables y muchos puntos de corte posibles, y elige el corte que deja los dos grupos resultantes lo más “puros” posible (es decir, con la menor mezcla de clases). La métrica más común para medir esa pureza es el índice de Gini, que en términos simples se acerca a 0 cuando un grupo es casi puro (una sola clase domina) y crece mientras más mezcladas están las clases dentro del grupo.

No vamos a entrar en la fórmula exacta del índice de Gini: para este curso nos basta con la intuición de que el árbol busca, en cada paso, el corte que genera grupos lo menos mezclados posible, y que R hace toda esa búsqueda por nosotros.

Preparando los datos: Titanic

Vamos a trabajar con el dataset Titanic, uno de los más usados para introducir clasificación: para cada pasajero, conocemos varias características y si sobrevivió o no al hundimiento.

# install.packages("titanic")
library(titanic)

data("titanic_train") # Carga en memoria el dataset
titanic <- titanic_train # Copia el dataset en un objeto

dim(titanic)
## [1] 891  12
head(titanic) # Muestra las primeras 6 observaciones
##   PassengerId Survived Pclass
## 1           1        0      3
## 2           2        1      1
## 3           3        1      3
## 4           4        1      1
## 5           5        0      3
## 6           6        0      3
##                                                  Name    Sex Age SibSp Parch
## 1                             Braund, Mr. Owen Harris   male  22     1     0
## 2 Cumings, Mrs. John Bradley (Florence Briggs Thayer) female  38     1     0
## 3                              Heikkinen, Miss. Laina female  26     0     0
## 4        Futrelle, Mrs. Jacques Heath (Lily May Peel) female  35     1     0
## 5                            Allen, Mr. William Henry   male  35     0     0
## 6                                    Moran, Mr. James   male  NA     0     0
##             Ticket    Fare Cabin Embarked
## 1        A/5 21171  7.2500              S
## 2         PC 17599 71.2833   C85        C
## 3 STON/O2. 3101282  7.9250              S
## 4           113803 53.1000  C123        S
## 5           373450  8.0500              S
## 6           330877  8.4583              Q


Nos quedamos únicamente con las variables que vamos a usar como predictores, más la variable respuesta Survived:

library(dplyr)
## Warning: package 'dplyr' was built under R version 4.5.3
## 
## Adjuntando el paquete: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
titanic <- titanic %>%
  dplyr::select(Survived, Pclass, Sex, Age, SibSp, Parch, Fare, Embarked)
  • Survived: 0 = no sobrevivió, 1 = sobrevivió (nuestra variable respuesta).
  • Pclass: clase del boleto (1ª, 2ª o 3ª), como proxy del nivel socioeconómico.
  • Sex, Age: sexo y edad del pasajero.
  • SibSp, Parch: hermanos/cónyuges y padres/hijos a bordo.
  • Fare: precio pagado por el boleto.
  • Embarked: puerto de embarque.

No vamos a detenernos demasiado en el pre-procesamiento, porque el foco de la clase es el modelo. Solo necesitamos dos cosas mínimas antes de poder ajustar un árbol:

1. Convertir a factor las variables categóricas, para que R (y rpart más adelante) las trate como categorías y no como números:

titanic$Survived <- factor(titanic$Survived, levels = c(0, 1), labels = c("Muere", "Vive"))
titanic$Pclass <- factor(titanic$Pclass)
titanic$Sex <- factor(titanic$Sex)
titanic$Embarked <- factor(titanic$Embarked)


Revisemos la estructura de los datos.

str(titanic)
## 'data.frame':    891 obs. of  8 variables:
##  $ Survived: Factor w/ 2 levels "Muere","Vive": 1 2 2 2 1 1 1 1 2 2 ...
##  $ Pclass  : Factor w/ 3 levels "1","2","3": 3 1 3 1 3 3 1 3 3 2 ...
##  $ Sex     : Factor w/ 2 levels "female","male": 2 1 1 1 2 2 2 2 1 1 ...
##  $ Age     : num  22 38 26 35 35 NA 54 2 27 14 ...
##  $ SibSp   : int  1 1 0 1 0 0 0 3 0 1 ...
##  $ Parch   : int  0 0 0 0 0 0 0 1 2 0 ...
##  $ Fare    : num  7.25 71.28 7.92 53.1 8.05 ...
##  $ Embarked: Factor w/ 4 levels "","C","Q","S": 4 2 4 4 4 3 4 4 4 2 ...

2. Eliminar los pocos registros incompletos.

colSums(is.na(titanic))
## Survived   Pclass      Sex      Age    SibSp    Parch     Fare Embarked 
##        0        0        0      177        0        0        0        0


Se observa que Age es la única variable con valores faltantes relevantes, con 177 valores NA sobre un total de 891 observaciones.

Estos datos representan aproximadamente un 20% de la variable Age, pero no afectan a las demás variables.

Eliminar estos registros simplifica el análisis y evita introducir supuestos adicionales, ya que los árboles de decisión son sensibles a la calidad de los datos, por lo que trabajar con datos completos mejora la interpretabilidad.

titanic <- titanic %>% filter(!is.na(Age))

nrow(titanic)
## [1] 714


Nota: eliminar filas con NA es la salida más simple, y nos alcanza para esta clase. En un proyecto real valdría la pena evaluar imputar esos valores en lugar de descartarlos, para no perder observaciones valiosas.

División en entrenamiento y prueba

Igual que en las clases anteriores, separamos los datos antes de ajustar cualquier modelo, para poder evaluar después qué tan bien generaliza a observaciones que nunca vio. Usaremos una proporción de 80/20.

Nota: la elección de proporción depende del tamaño del dataset: si los datos son limitados, se prefiere usar 80-20 para aprovechar más información en el entrenamiento; si se dispone de más datos, 70-30 permite una evaluación más sólida.

set.seed(123)

n <- nrow(titanic)

indices_entrenamiento <- sample(1:n, size = 0.8 * n)

train <- titanic[indices_entrenamiento, ]
test <- titanic[-indices_entrenamiento, ]

nrow(train)
## [1] 571
nrow(test)
## [1] 143


Ajustando el primer árbol de decisión

En R, los árboles de decisión no están en una función base, pero rpart (Recursive Partitioning) es la librería estándar para crearlos y viene preinstalada con R.

library(rpart)
library(rpart.plot)

modelo_arbol <- rpart(Survived ~ ., data = train, method = "class")

modelo_arbol
## n= 571 
## 
## node), split, n, loss, yval, (yprob)
##       * denotes terminal node
## 
##  1) root 571 226 Muere (0.60420315 0.39579685)  
##    2) Sex=male 362  73 Muere (0.79834254 0.20165746)  
##      4) Age>=13 329  55 Muere (0.83282675 0.16717325) *
##      5) Age< 13 33  15 Vive (0.45454545 0.54545455)  
##       10) SibSp>=2 15   1 Muere (0.93333333 0.06666667) *
##       11) SibSp< 2 18   1 Vive (0.05555556 0.94444444) *
##    3) Sex=female 209  56 Vive (0.26794258 0.73205742)  
##      6) Pclass=3 86  39 Muere (0.54651163 0.45348837)  
##       12) Fare>=20.8 21   3 Muere (0.85714286 0.14285714) *
##       13) Fare< 20.8 65  29 Vive (0.44615385 0.55384615)  
##         26) Age>=16.5 46  21 Muere (0.54347826 0.45652174)  
##           52) Parch< 0.5 35  14 Muere (0.60000000 0.40000000) *
##           53) Parch>=0.5 11   4 Vive (0.36363636 0.63636364) *
##         27) Age< 16.5 19   4 Vive (0.21052632 0.78947368) *
##      7) Pclass=1,2 123   9 Vive (0.07317073 0.92682927) *


La sintaxis es la misma que ya conocemos de lm() y glm(): respuesta ~ predictores, y data con el data frame. El argumento method = "class" le indica a rpart que se trata de un problema de clasificación (si Survived fuera numérica y quisiéramos predecir un número, usaríamos method = "anova").

La salida en texto es difícil de leer. Visualicemos el árbol con rpart.plot():

rpart.plot(modelo_arbol, type = 2, extra = 104)


Cada caja es un nodo. Para leerlas:

  • La etiqueta superior (Muere o Vive) es la predicción de ese nodo, la clase mayoritaria entre las observaciones que caen ahí.
  • Los dos números de abajo (logrado con extra = 104) son la proporción de cada clase dentro del nodo.
  • El porcentaje inferior indica qué proporción del total de observaciones de entrenamiento cae en ese nodo.
  • Debajo de cada nodo interno aparece la pregunta que usó para dividir (por ejemplo, Sex = male).

El primer corte del árbol es sobre Sex: el sexo del pasajero fue, históricamente, el factor más determinante para la supervivencia (“mujeres y niños primero”). A partir de ahí, el árbol sigue refinando la predicción con Pclass, Age u otras variables, dependiendo de qué grupo se esté analizando.

Esta es una de las grandes ventajas de los árboles frente a modelos como la regresión logística: se leen directamente como un diagrama de flujo, sin necesidad de interpretar coeficientes.

Evaluando el modelo

Prediciendo sobre datos de prueba

# Utilizamos el árbol entrenado para predecir la clase (0 o 1) de cada observación del conjunto de prueba.
pred_clase <- predict(modelo_arbol, newdata = test, type = "class")

head(pred_clase)
##     1     3     9    18    22    27 
## Muere Muere  Vive Muere Muere Muere 
## Levels: Muere Vive


type = "class" le pide a predict() que devuelva directamente la clase predicha (Muere/Vive), en lugar de una probabilidad.

Construyendo la matriz de confusión

Una matriz de confusión es simplemente una tabla que cruza lo que el modelo predijo contra lo que realmente ocurrió. En R, no necesitamos ninguna función especializada para esto: basta con table().

matriz_confusion <- table(
  Real = test$Survived,
  Predicho = pred_clase
)

matriz_confusion
##        Predicho
## Real    Muere Vive
##   Muere    77    2
##   Vive     23   41


Las filas representan la clase real y las columnas la clase que predijo el modelo. Los cuatro números que aparecen tienen nombre propio:

  • Verdaderos negativos (TN): el modelo predijo “Muere” y el pasajero efectivamente murió.
  • Verdaderos positivos (TP): el modelo predijo “Vive” y el pasajero efectivamente sobrevivió.
  • Falsos positivos (FP): el modelo predijo “Vive” pero el pasajero murió.
  • Falsos negativos (FN): el modelo predijo “Muere” pero el pasajero sobrevivió.

La diagonal principal (TN y TP) son los aciertos del modelo, todo lo que está fuera de la diagonal (FP y FN) son los errores, y cada uno representa un tipo distinto de equivocación.

Visualizando la matriz de confusión

Para exponer esta matriz es mucho más claro mostrarlo como un mapa de calor que leer la tabla cruda. Transformamos la tabla a data frame con as.data.frame() y graficamos con ggplot2.

library(ggplot2)


df_matriz <- as.data.frame(matriz_confusion)

df_matriz$Real <- factor(df_matriz$Real, levels = c("Muere", "Vive"))
df_matriz$Predicho <- factor(df_matriz$Predicho, levels = c("Vive", "Muere"))

ggplot(df_matriz, aes(x = Real, y = Predicho, fill = Freq)) +
  geom_tile(color = "white") +
  geom_text(aes(label = Freq), size = 6, color = "white") +
  scale_fill_gradient(low = "steelblue", high = "darkblue") +
  labs(title = "Matriz de confusión", x = "Valor real", y = "Predicción del modelo") +
  theme_minimal()


Este mapa de calor se lee así:

  • La diagonal (Muere-Muere y Vive-Vive) concentra los aciertos. Mientras más oscuras (más observaciones) estén esas dos celdas en comparación con las otras dos, mejor está funcionando el modelo.
  • Las celdas fuera de la diagonal son los errores. Vale la pena mirarlas por separado, porque no todos los errores son igual de “graves” según el contexto del problema: en este caso, predecir que alguien murió cuando en realidad sobrevivió (falso negativo) no tiene el mismo costo práctico que el error contrario, aunque numéricamente el modelo los trate igual.

En otras aplicaciones esta diferencia es aún más evidente: por ejemplo, en un modelo para detectar una enfermedad suele ser mucho más grave clasificar a un paciente enfermo como sano que generar una falsa alarma en un paciente sano. Por ello, al evaluar un modelo no basta con mirar cuántos errores comete, sino también qué tipo de errores está cometiendo.

Métricas de evaluación

La matriz de confusión es la fuente de toda la información, pero para comparar modelos entre sí conviene resumirla en unos pocos números. Los calculamos indexando directamente la tabla que ya construimos, sin depender de ninguna función adicional.

TN <- matriz_confusion["Muere", "Muere"]
TP <- matriz_confusion["Vive", "Vive"]
FP <- matriz_confusion["Muere", "Vive"]
FN <- matriz_confusion["Vive", "Muere"]

TN; TP; FP; FN
## [1] 77
## [1] 41
## [1] 2
## [1] 23


Exactitud (Accuracy): proporción de predicciones correctas sobre el total.

\[ Accuracy = \frac{TP + TN}{TP + TN + FP + FN} \]

accuracy <- (TP + TN) / (TP + TN + FP + FN)
accuracy
## [1] 0.8251748


Precisión (Precision): de todos los pasajeros que el modelo predijo como sobrevivientes, ¿qué proporción realmente sobrevivió?

\[ Precision = \frac{TP}{TP + FP} \]

precision <- TP / (TP + FP)
precision
## [1] 0.9534884


Sensibilidad (Recall): de todos los pasajeros que realmente sobrevivieron, ¿qué proporción detectó el modelo?

\[ Recall = \frac{TP}{TP + FN} \]

recall <- TP / (TP + FN)
recall
## [1] 0.640625


F1-score: una única métrica que resume precisión y recall a la vez, como su media armónica. Es útil cuando queremos un solo número para comparar modelos, pero sin perder de vista que hay un balance entre ambas métricas.

\[ F1 = \frac{2 \times Precision \times Recall}{Precision + Recall} \]

f1 <- 2 * (precision * recall) / (precision + recall)
f1
## [1] 0.7663551


Juntamos todo en una sola tabla, para que sea fácil de leer y de comparar más adelante contra otros modelos:

data.frame(
  Accuracy = accuracy,
  Precision = precision,
  Recall = recall,
  F1 = f1
)
##    Accuracy Precision   Recall        F1
## 1 0.8251748 0.9534884 0.640625 0.7663551


En conjunto, estas métricas muestran que el modelo tiene un desempeño razonablemente bueno, pero también evidencian sus limitaciones. La accuracy indica que aproximadamente el 82.5% de las predicciones fueron correctas; sin embargo, esta medida por sí sola no describe completamente el comportamiento del modelo. La precision elevada (95.3%) nos dice que, cuando el árbol predice que un pasajero sobrevivió, suele acertar. En cambio, el recall (64.1%) revela que el modelo deja sin identificar una proporción importante de los pasajeros que realmente sobrevivieron. El F1-score resume este equilibrio entre precisión y sensibilidad, mostrando que aún existe margen de mejora.

Aunque la accuracy es una métrica fácil de interpretar, solo indica el porcentaje total de aciertos y no distingue entre los diferentes tipos de errores. En cambio, el F1-score obliga al modelo a encontrar un equilibrio entre identificar correctamente los casos positivos (recall) y evitar falsas alarmas (precision). Un modelo puede tener una accuracy elevada simplemente porque la mayoría de las observaciones pertenecen a una misma clase, mientras que el F1 disminuirá si el modelo ignora muchos casos de la clase de interés o genera demasiadas predicciones incorrectas. Por ello, cuando el objetivo es evaluar qué tan bien se identifica una clase específica, el F1-score suele ser una medida más representativa que la accuracy por sí sola.

Con esto ya tenemos un primer árbol ajustado y evaluado de punta a punta. En la siguiente parte de la clase vamos a ver qué pasa si dejamos crecer el árbol sin restricciones (overfitting), cómo controlarlo, y por qué combinar muchos árboles (Random Forest) suele superar a un árbol único.

Controlando la complejidad del árbol

Si dejamos que un árbol crezca sin ninguna restricción, seguirá dividiendo los datos una y otra vez, hasta que cada hoja quede casi perfectamente pura (en el extremo, con una sola observación). Esto suena bien, pero es exactamente el mismo problema de sobreajuste que ya vimos con la regresión: un árbol así memoriza las particularidades del conjunto de entrenamiento en lugar de aprender un patrón que generalice.

Vamos a comprobarlo directamente, comparando el desempeño en train contra el desempeño en test para tres árboles con distinto nivel de restricción.

Un árbol sin restricciones

modelo_profundo <- rpart(
  Survived ~ ., data = train, method = "class",
  control = rpart.control(cp = 0, minsplit = 2, maxdepth = 30)
)


  • cp = 0: el parámetro de complejidad (complexity parameter) controla qué tan grande debe ser la mejora en pureza para que valga la pena hacer un corte. En cp = 0, aceptamos absolutamente cualquier corte, por mínima que sea la mejora.
  • minsplit = 2: permite intentar dividir un nodo aunque tenga solamente 2 observaciones.
  • maxdepth = 30: prácticamente sin límite de profundidad.

En conjunto, estos tres argumentos le dan al árbol total libertad para crecer.

Un árbol demasiado simple

modelo_simple <- rpart(
  Survived ~ ., data = train, method = "class",
  control = rpart.control(maxdepth = 1)
)


maxdepth = 1 fuerza al árbol a hacer un único corte y detenerse ahí, sin importar qué tan mezclados sigan los grupos resultantes.

Nota: a este tipo de árboles se les conoce como decision stumps, árboles de decisión con profundidad 1 (maxdepth = 1), i.e., con un único nodo de decisión, un solo corte y que produce únicamente dos nodos hoja.

Comparando train vs. test

Definimos una pequeña función para no repetir el mismo cálculo tres veces: simplemente compara las predicciones del modelo contra la clase real y devuelve la proporción de aciertos (accuracy).

calcular_accuracy <- function(modelo, datos) {
  predicciones <- predict(modelo, newdata = datos, type = "class")
  mean(predicciones == datos$Survived)
}

data.frame(
  modelo = c("Simple (maxdepth = 1)", "Por defecto", "Profundo (sin restricciones)"),
  accuracy_train = c(
    calcular_accuracy(modelo_simple, train),
    calcular_accuracy(modelo_arbol, train),
    calcular_accuracy(modelo_profundo, train)
  ),
  accuracy_test = c(
    calcular_accuracy(modelo_simple, test),
    calcular_accuracy(modelo_arbol, test),
    calcular_accuracy(modelo_profundo, test)
  )
)
##                         modelo accuracy_train accuracy_test
## 1        Simple (maxdepth = 1)      0.7740806     0.8041958
## 2                  Por defecto      0.8406305     0.8251748
## 3 Profundo (sin restricciones)      0.9877408     0.7832168


Interpretamos:

  • Modelo simple: accuracy baja tanto en train como en test, y bastante parecida entre ambos. Es un caso de underfitting: el modelo es tan restringido que ni siquiera logra aprender bien el patrón dentro de sus propios datos de entrenamiento.
  • Modelo profundo: accuracy en train muy cercana a 1 (casi memoriza el dataset completo), pero accuracy en test claramente menor, incluso por debajo del árbol por defecto. Es la firma clásica de overfitting: aprendió ruido específico del conjunto de entrenamiento que no generaliza.
  • Modelo por defecto: un punto intermedio, sin restricciones manuales tan extremas, y por eso generaliza mejor que el árbol profundo.

rpart() por defecto configura cp distinto de 0 (usualmente 0.01), lo cual explica por qué el “modelo por defecto” no sufre el mismo sobreajuste extremo que el árbol sin restricciones. rpart está, en cierta forma, haciendo ya un poco de este trabajo por nosotros.

Visualizando el efecto de la complejidad

rpart trae una función para visualizar cómo cambia el error (estimado por validación cruzada interna) según el valor de cp usado:

plotcp(modelo_profundo)


Este gráfico muestra cómo cambia el error del árbol a medida que permitimos realizar más cortes. El eje vertical representa un error relativo, es decir, una medida del error expresada en comparación con el árbol más simple posible (un árbol sin divisiones). Valores más pequeños indican un mejor desempeño. Lo que buscamos es el punto donde agregar más cortes apenas reduce el error: a partir de ahí el árbol se vuelve más complejo, pero la mejora obtenida es muy pequeña, aumentando el riesgo de sobreajuste. Esta es la misma idea del “codo” que vimos al seleccionar el número de componentes principales en PCA.

La conclusión práctica se motiva de la misma forma que en regresión, debe haber un balance entre qué tan flexible dejamos que sea el modelo y qué tan bien generaliza, y ese balance no lo elegimos mirando solo el desempeño en entrenamiento.

Al principio, cada nueva división ayuda al árbol a aprender patrones importantes y el error disminuye. Sin embargo, llega un momento en el que seguir agregando cortes ya no aporta información útil; el árbol comienza a memorizar particularidades del conjunto de entrenamiento en lugar de aprender patrones generales. Cuando esto ocurre, el modelo pierde capacidad para predecir correctamente datos nuevos y el error de validación vuelve a aumentar. Por esta razón, no siempre el árbol más grande es el mejor: buscamos un árbol lo suficientemente complejo para capturar los patrones de los datos, pero no tanto como para sobreajustarse.

Hasta un cp cercano a 0.003–0.004 el error disminuye. A partir de ese punto el error deja de mejorar y empieza a aumentar gradualmente. Esa es la razón por la que se suele podar el árbol en esa zona, en lugar de quedarse con el árbol más grande (137 nodos terminales), que ya está aprendiendo demasiado del conjunto de entrenamiento.

Las limitaciones de un solo árbol

Los árboles de decisión tienen una debilidad importante: son inestables. Un cambio relativamente pequeño en los datos de entrenamiento puede hacer que el árbol elija un corte distinto desde el primer nivel, y como cada corte posterior depende del anterior, todo el árbol puede terminar siendo muy diferente.

Podemos verlo directamente: ajustemos el árbol usando dos muestras de entrenamiento ligeramente distintas (obtenidas re-muestreando train con reemplazo) y comparemos.

set.seed(1)
muestra_a <- train[sample(1:nrow(train), replace = TRUE), ]
arbol_a <- rpart(Survived ~ ., data = muestra_a, method = "class")

set.seed(2)
muestra_b <- train[sample(1:nrow(train), replace = TRUE), ]
arbol_b <- rpart(Survived ~ ., data = muestra_b, method = "class")

rpart.plot(arbol_a, main = "Árbol A")

rpart.plot(arbol_b, main = "Árbol B")


Aunque ambas muestras salieron de los mismos datos (train), los dos árboles tienen distintos cortes, distinta profundidad, e incluso el orden de las variables usadas cambia. Esta sensibilidad a los datos específicos es justamente lo que hace que un solo árbol tienda a sobreajustar.

¿Y si, en lugar de confiar en un único árbol, entrenamos muchos árboles (cada uno con una muestra distinta) y promediamos sus predicciones?

Esa es la idea central de Bagging (Bootstrap Aggregating): generar muchas muestras con reemplazo a partir del dataset original (igual que hicimos arriba con muestra_a y muestra_b), ajustar un árbol distinto en cada una, y combinar sus predicciones (por votación, en clasificación). Como cada árbol individual ve una versión ligeramente distinta de los datos, sus errores tienden a ser distintos entre sí, y al promediar muchos árboles, esos errores individuales tienden a equilibrarse.

Random Forest lleva esta idea un paso más allá: además de usar una muestra distinta para cada árbol, en cada corte del árbol solo le permite considerar un subconjunto aleatorio de las variables disponibles (no todas). Esto fuerza a que los árboles del “bosque” sean menos parecidos entre sí (si no lo hiciéramos, muchos árboles elegirían el mismo primer corte, como Sex, y terminarían siendo muy similares entre ellos), lo cual mejora todavía más el beneficio de promediarlos.

En resumen, un Random Forest es, literalmente, un conjunto de muchos árboles como los que ya sabemos construir, combinados con dos fuentes de aleatoriedad (muestras con reemplazo + subconjuntos de variables) para hacerlos menos parecidos entre sí, y cuya predicción final es un promedio (o votación) de todos ellos.

Ajustando un Random Forest en R

R no incluye Random Forest en sus paquetes base, pero randomForest es la implementación estándar y es liviana de instalar.

# install.packages("randomForest")
library(randomForest)

set.seed(123)

modelo_rf <- randomForest(
  Survived ~ ., data = train,
  ntree = 500,
  mtry = 3,
  importance = TRUE
)

modelo_rf
## 
## Call:
##  randomForest(formula = Survived ~ ., data = train, ntree = 500,      mtry = 3, importance = TRUE) 
##                Type of random forest: classification
##                      Number of trees: 500
## No. of variables tried at each split: 3
## 
##         OOB estimate of  error rate: 18.74%
## Confusion matrix:
##       Muere Vive class.error
## Muere   306   39   0.1130435
## Vive     68  158   0.3008850


Los argumentos que sí vale la pena controlar:

  • ntree: el número de árboles a entrenar. Más árboles casi nunca perjudica el desempeño (solo el tiempo de cómputo); 500 suele ser un valor razonable por defecto.
  • mtry: cuántas variables se sortean como candidatas en cada corte. Es el parámetro que introduce la segunda fuente de aleatoriedad que mencionamos antes. Para clasificación, un valor de referencia habitual es la raíz cuadrada del número total de predictores.
  • importance = TRUE: le pide al modelo que además calcule qué tan importante fue cada variable (lo veremos en un momento).

Al imprimir el modelo aparece algo llamado OOB estimate of error rate (out-of-bag). Cada árbol del bosque se entrena con una muestra con reemplazo del ~63% de las observaciones de train; el ~37% restante (las que ese árbol nunca vio) se llaman observaciones “out-of-bag” para ese árbol. Random Forest aprovecha esto: para cada observación, promedia las predicciones únicamente de los árboles que no la usaron para entrenar, obteniendo así una estimación del error de generalización sin necesidad de un conjunto de prueba separado. Aun así, nosotros seguiremos evaluando también sobre nuestro propio test, para poder comparar de forma consistente contra el árbol individual.

Nota: estos porcentajes son una consecuencia matemática del muestreo con reemplazo cuando se toman n observaciones de un conjunto que también tiene tamaño n.

Evaluando el Random Forest

Repetimos exactamente el mismo procedimiento que usamos para el árbol individual: predecir sobre test, construir la matriz de confusión con table(), y calcular las métricas a mano.

pred_rf <- predict(modelo_rf, newdata = test, type = "class")

matriz_confusion_rf <- table(
  Real = test$Survived,
  Predicho = pred_rf
)

matriz_confusion_rf
##        Predicho
## Real    Muere Vive
##   Muere    78    1
##   Vive     19   45


TN_rf <- matriz_confusion_rf["Muere", "Muere"]
TP_rf <- matriz_confusion_rf["Vive", "Vive"]
FP_rf <- matriz_confusion_rf["Muere", "Vive"]
FN_rf <- matriz_confusion_rf["Vive", "Muere"]

accuracy_rf <- (TP_rf + TN_rf) / (TP_rf + TN_rf + FP_rf + FN_rf)
precision_rf <- TP_rf / (TP_rf + FP_rf)
recall_rf <- TP_rf / (TP_rf + FN_rf)
f1_rf <- 2 * (precision_rf * recall_rf) / (precision_rf + recall_rf)


Comparando árbol individual vs. Random Forest

data.frame(
  modelo = c("Árbol de decisión", "Random Forest"),
  Accuracy = c(accuracy, accuracy_rf),
  Precision = c(precision, precision_rf),
  Recall = c(recall, recall_rf),
  F1 = c(f1, f1_rf)
)
##              modelo  Accuracy Precision   Recall        F1
## 1 Árbol de decisión 0.8251748 0.9534884 0.640625 0.7663551
## 2     Random Forest 0.8601399 0.9782609 0.703125 0.8181818


Observamos que Random Forest supera al árbol individual en todas las métricas. La razón detrás de este fenómeno es que un solo árbol es una herramienta inestable y propensa a sobreajustar, mientras que el bosque promedia muchos árboles entrenados con variaciones distintas de los datos, lo cual reduce esa inestabilidad.

La contraparte de esta mejora: perdimos la posibilidad de graficar “el árbol” con rpart.plot(), porque ya no existe un único árbol, sino 500. Este es otro ejemplo del mismo tipo de trade-off que vimos con PCA: ganamos desempeño y estabilidad, a costa de interpretabilidad directa.

Importancia de variables

Aunque ya no podemos visualizar un árbol específico, ganamos que Random Forest sí nos puede decir qué variables aportaron más a las predicciones del bosque en su conjunto, gracias al argumento importance = TRUE que usamos al ajustar el modelo.

importance(modelo_rf)
##              Muere      Vive MeanDecreaseAccuracy MeanDecreaseGini
## Pclass   16.064510 29.780521            34.023065        25.766917
## Sex      52.230164 74.998293            78.503638        66.369472
## Age      19.587268 20.745682            28.965659        59.063550
## SibSp    17.061106 -1.019895            15.272184        12.278517
## Parch     3.920499  3.226074             5.536004         9.049881
## Fare     15.716131 16.689069            24.613800        59.050298
## Embarked  1.116075  0.818380             1.504953         7.212165


Esta tabla trae dos columnas relevantes:

  • MeanDecreaseAccuracy: cuánto empeora la accuracy del modelo, en promedio, si “revolvemos” aleatoriamente los valores de esa variable (rompiendo su relación real con Survived). Mientras más cae la accuracy al desordenar una variable, más importante era esa variable para las predicciones.
  • MeanDecreaseGini: cuánto contribuye esa variable, en promedio y sumando sobre todos los árboles del bosque, a mejorar la pureza (índice de Gini) de los nodos donde fue usada para cortar.

Podemos visualizar ambas de una vez con varImpPlot():

varImpPlot(modelo_rf, main = "Importancia de variables")


Sex aparece como la variable más importante en ambos criterios, seguida por Fare o Pclass (relacionadas con la clase socioeconómica) y Age. Esto es consistente con lo que históricamente se sabe del hundimiento del Titanic, y también con lo que vimos como primer corte en nuestro árbol individual al inicio de la clase.

El árbol individual nos daba un diagrama de flujo fácil de leer pero inestable; el Random Forest nos da un modelo más robusto y preciso, pero a cambio perdemos esa lectura directa, y varImpPlot() es la principal herramienta que nos queda para entender, de forma agregada, qué variables está usando el modelo.

Cierre

En la práctica, encontrar el mejor modelo requiere ajustar y evaluar sistemáticamente los parámetros que controlan su complejidad y comportamiento, como la profundidad del árbol, el tamaño mínimo de muestras en los nodos o los umbrales de clasificación.

Este proceso combina técnicas como validación cruzada y búsqueda de hiperparámetros, midiendo el impacto de cada ajuste en métricas clave como precisión, recall o F1, según el objetivo del problema.

Dominar el análisis de datos requiere práctica y claridad sobre lo que buscamos descubrir, entendiendo cómo nuestras decisiones influyen en los resultados. R se convierte en una herramienta poderosa en este proceso, permitiendo explorar, visualizar y analizar datos de manera ágil y precisa. Con sus capacidades para manipular información, generar gráficos intuitivos y aplicar modelos estadísticos, R facilita identificar patrones y construir modelos confiables. Con práctica, aprender a aprovechar R nos permite transformar datos en conocimiento sólido y tomar decisiones fundamentadas.

¡Cada dato es una oportunidad para aprender algo nuevo! ʕ•́ᴥ•̀ʔっ♡