3.1 Árboles de regresión CART

Como se comentó previamente, la construcción del modelo se hace a partir de la muestra de entrenamiento, y consiste en la partición del espacio predictor en \(J\) regiones \(R_1, R_2, \ldots, R_J\), para cada una de las cuales se va a calcular una constante: la media de la variable respuesta \(Y\) para las observaciones de entrenamiento que caen en la región. Estas constantes son las que se van a utilizar para la predicción de nuevas observaciones; para ello solo hay que comprobar cuál es la región que le corresponde.

La cuestión clave es cómo se elige la partición del espacio predictor, para lo que vamos a utilizar como criterio de error la suma de los residuos al cuadrado (RSS, por sus siglas en inglés). Como hemos dicho, vamos a modelizar la respuesta en cada región como una constante, por tanto en la región \(R_j\) nos interesa el \(min_{c_j} \sum_{i\in R_j} (y_i - c_j)^2\), que se alcanza en la media de las respuestas \(y_i\) (de la muestra de entrenamiento) en la región \(R_j\), a la que llamaremos \(\widehat y_{R_j}\). Por tanto, se deben seleccionar las regiones \(R_1, R_2, \ldots, R_J\) que minimicen \[RSS = \sum_{j=1}^{J} \sum_{i\in R_j} (y_i - \widehat y_{R_j})^2\] (Obsérvese el abuso de notación \(i\in R_j\), que significa las observaciones \(i\in N\) que verifican \(x_i \in R_j\)).

Pero este problema es intratable en la práctica, por lo que es necesario simplificarlo. El método CART busca un compromiso entre rendimiento, por una parte, y sencillez e interpretabilidad, por otra, y por ello en lugar de hacer una búsqueda por todas las particiones posibles sigue un proceso iterativo (recursivo) en el que va realizando cortes binarios. En la primera iteración se trabaja con todos los datos:

  • Una variable explicativa \(X_j\) y un punto de corte \(s\) definen dos hiperplanos \(R_1 = \{ X : X_j \le s \}\) y \(R_2 = \{ X : X_j > s \}\).

  • Se seleccionan los valores de \(j\) y \(s\) que minimizen \[ \sum_{i\in R_1} (y_i - \widehat y_{R_1})^2 + \sum_{i\in R_2} (y_i - \widehat y_{R_2})^2\]

A diferencia del problema original, este se soluciona de forma muy rápida. A continuación se repite el proceso en cada una de las dos regiones \(R_1\) y \(R_2\), y así sucesivamente hasta alcanzar un criterio de parada.

Fijémonos en que este método hace dos concesiones importantes: no solo restringe la forma que pueden adoptar las particiones, sino que además sigue un criterio de error codicioso (greedy): en cada iteración busca minimizar el RSS de las dos regiones resultantes, sin preocuparse del error que se va a cometer en iteraciones sucesivas. Y fijémonos también en que este proceso se puede representar en forma de árbol binario (en el sentido de que de cada nodo salen dos ramas, o ninguna cuando se llega al final), de ahí la terminología de hacer crecer el árbol.

¿Y cuándo paramos? Se puede parar cuando se alcance una profundidad máxima, aunque lo más habitual es exigir un número mínimo de observaciones para dividir un nodo.

  • Si el árbol resultante es demasiado grande, va a ser un modelo demasiado complejo, por tanto va a ser difícil de interpretar y, sobre todo, va a provocar un sobreajuste de los datos. Cuando se evalúe el rendimiento utilizando la muestra de validación, los resultados van a ser malos. Dicho de otra manera, tendremos un modelo con poco sesgo pero con mucha varianza y en consecuencia inestable (pequeños cambios en los datos darán lugar a modelos muy distintos). Más adelante veremos que esto justifica la utilización del bagging como técnica para reducir la varianza.

  • Si el árbol es demasiado pequeño, va a tener menos varianza (menos inestable) a costa de más sesgo. Más adelante veremos que esto justifica la utilización del boosting. Los árboles pequeños son más fáciles de interpretar, ya que permiten identificar las variables explicativas que más influyen en la predicción.

Sin entrar por ahora en métodos combinados (métodos ensemble, tipo bagging o boosting), vamos a explicar cómo encontrar un equilibrio entre sesgo y varianza. Lo que se hace es construir un árbol grande para a continuación empezar a podarlo. Podar un árbol significa colapsar cualquier cantidad de sus nodos internos (no terminales), dando lugar a otro árbol más pequeño al que llamaremos subárbol del árbol original. Sabemos que el árbol completo es el que va a tener menor error si utilizamos la muestra de entrenamiento, pero lo que realmente nos interesa es encontrar el subárbol con un menor error al utilizar la muestra de validación. Lamentablemente, no es una estrategia viable evaluar todos los subárboles: simplemente, hay demasiados. Lo que se hace es, mediante un hiperparámetro (tuning parameter o parámetro de ajuste), controlar el tamaño del árbol, es decir, la complejidad del modelo, seleccionando el subárbol óptimo (para los datos disponibles). Veamos la idea con más detalle.

Dado un subárbol \(T\) con \(R_1, R_2, \ldots, R_t\) nodos terminales, consideramos como medida del error el RSS más una penalización que depende de un hiperparámetro no negativo \(\alpha \ge 0\)

\[\begin{equation} RSS_{\alpha} = \sum_{j=1}^t \sum_{i\in R_j} (y_i - \widehat y_{R_j})^2 + \alpha t \tag{3.1} \end{equation}\]

Para cada valor del parámetro \(\alpha\) existe un único subárbol más pequeño que minimiza este error (obsérvese que aunque hay un continuo de valores distinos de \(\alpha\), solo hay una cantidad finita de subárboles). Evidentemente, cuando \(\alpha = 0\), ese subárbol será el árbol completo, algo que no nos interesa. Pero a medida que se incrementa \(\alpha\) se penalizan los subárboles con muchos nodos terminales, dando lugar a una solución más pequeña (compacta). Encontrarla puede parecer muy costoso computacionalmente, pero lo cierto es que no lo es. El algoritmo consistente en ir colapsando nodos de forma sucesiva, de cada vez el nodo que produzca el menor incremento en el RSS (corregido por un factor que depende del tamaño), da lugar a una sucesión finita de subárboles que contiene, para todo \(\alpha\), la solución.

Para finalizar, solo resta seleccionar un valor de \(\alpha\). Para ello, como se comentó en la Sección 1.3.2, una posible estrategia consiste en dividir la muestra en tres subconjuntos: datos de entrenamiento, de validación y de test. Para cada valor del parámetro de complejidad \(\alpha\) hemos utilizado la muestra de entrenamiento para obtener un árbol (en la jerga, para cada valor del hiperparámetro \(\alpha\) se entrena un modelo). Se emplea la muestra independiente de validación para seleccionar el valor de \(\alpha\) (y por tanto el árbol) con el que nos quedamos. Y por último emplearemos la muestra de test (independiente de las otras dos) para evaluar el rendimiento del árbol seleccionado. No obstante, lo más habitual para seleccionar el valor del hiperparámetro \(\alpha\) es emplear validación cruzada (u otro tipo de remuestreo) en la muestra de entrenamiento en lugar de considerar una muestra adicional de validación.

Hay dos opciones muy utilizadas en la práctica para seleccionar el valor de \(\alpha\): se puede utilizar directamente el valor que minimice el error; o se puede forzar que el modelo sea un poco más sencillo con la regla one-standard-error, que selecciona el árbol más pequeño que esté a una distancia de un error estándar del árbol obtenido mediante la opción anterior.

También es habitual escribir la Ecuación (3.1) reescalando el parámetro de complejidad como \(\tilde \alpha = \alpha / RSS_0\), siendo \(RSS_0 = \sum_{i=1}^{n} (y_i - \bar y)^2\) la variabilidad total (la suma de cuadrados residual del árbol sin divisiones): \[RSS_{\tilde \alpha}=RSS + \tilde \alpha RSS_0 t\]

De esta forma se podría interpretar el hiperparámetro \(\tilde \alpha\) como una penalización en la proporción de variabilidad explicada, ya que dividiendo la expresión anterior por \(RSS_0\) obtendríamos la proporción de variabilidad residual y a partir de ella podríamos definir: \[R^2_{\tilde \alpha}=R^2 - \tilde \alpha t\]