Navigate back to the homepage

ML Part 3: Gradient Boosting

Cristobal Aguirre
June 26th, 2019 · 5 min read
  • Part 1: Árboles de decisión
  • Part 2: Random Forests
  • Repositorio con el código

Introducción

Finalmente hemos llegado a Gradient Boosting Decision Trees (GBDT)—el final en esta serie que introduce el modelo de machine learning que por 5 años ya (una eternidad a la velocidad que avanzan estas cosas) ha estado liderando las competencias de data science y potenciando varias de las aplicaciones que usas día a día. Su éxito se debe a una combinación única de calidad predictiva, velocidad de entrenamiento (requiriendo menor intensidad computacional), y capacidad de recibir cualquier función objetivo (mientras sea diferenciable) lo que le permite resolver problemas de regresión o clasificación sin muchos ajustes necesarios.

Esta técnica se me hizo bastante más complicada de entender en profundidad, y por eso no voy a tratar de cubrirla desde 0 en este post—hay recursos mucho mejores que creo que hacen un mejor trabajo de lo que podría hacer acá y son necesarios para entender bien la teoría. Recomiendo partir con este post que construye el algoritmo de manera intuitiva con un ejemplo simple, y luego esta serie de 3 posts que también toman un approach gradual y ahondan un poco más en la matemática detrás.

En vez, voy a tomar un enfoque un poco distinto que asume un conocimiento básico, y es lo que me ayudó a mi a cimentar la teoría: construir el algoritmo en Python desde 0. Ver en código el proceso en vivo me dejó muy claro cada paso y los efectos de distintas funciones de pérdida, estimadores base, etc.

Repaso rápido de la teória

La técnica consiste básicamente en estimar una serie de estimadores débiles de manera secuencial, con cada estimador corrigiendo los errores del estimador anterior. Esta imágen que tomé de uno de los posts recomendados más arriba es una buena forma de visualizar el proceso:

Ejemplo de gradient descent con golfista

En términos matemáticos, en cada paso mm se va a estimar un estimador fm(x)f_m(x) sobre los residuos entre la variable objetivo y la estimación del paso anterior fm1(x)f_{m-1}(x). Al final del proceso, el estimador final Fm(x)F_m(x) consistirá en un estimador base f0f_0---que puede ser una serie de 0’s, el promedio, la mediana, etc.---y una sumatoria de mm estimadores débiles fm(x)f_m(x), ajustados por un parámetro de regularización η\eta que afecta la contribución de cada estimador al estimador final (conocido en inglés como shrinkage o learning rate) para controlar el overfitting.

FM(x)=f0+ηmMfm(x)F_M(x) = f_0 + \eta \sum_m^M f_m(x)

Donde el algoritmo secuencial implica que cada estimador durante la etapa mm toma la siguiente forma:

Fm(x)=Fm1(x)+ηfm(x)F_m(x) = F_{m-1}(x) + \eta f_m(x)

Uno de los principales atractivos de ésta técnica es que es te permite usar distintas funciones de pérdida, lo que le da un montón de flexibilidad; sin embargo, yo creo que la forma más fácil de entenderlo es tomando el caso de la regresión (variable dependiente contínua), tomando como función de pérdida el error cuadrático medio---esta función tiene la gracia de que su gradiente equivale a los residuos comunes yf(x)y - f(x).

En este caso, el algoritmo consiste en:

  1. Partir de un estimador base F0F_{0} (por ejemplo el promedio de la variable objetivo)
  2. Calcular el gradiente (residuos entre la variable objetivo y este estimador): yF0(x)y - F_0(x)
  3. Estimar un estimador débil a estos residuos f1(x)=yF0(x)f_1(x) = y - F_0(x)
  4. Actualizar el estimador agregado con esta nueva estimación F1(x)=F0+ηf1(x)F_1(x) = F_0 + \eta f_1(x)
  5. Repetir los pasos 2-4 mm veces. El siguiente paso calcularía un estimador débil f2(x)f_2(x) sobre los residuos yF1(x)y - F_1(x), para obtener F2(x)=F1+ηf2(x)F_2(x) = F_1 + \eta f_2(x), y así sucesivamente.

Show me the code

Este proceso es mucho más fácil de apreciar mirando una implementación simple en código. Sin embargo, mirar los códigos base de implementaciones populares como las de Scikit-learn o XGBoost son de poca ayuda porque hacen muchos ajustes de optimización y versatilidad que mejoran la calidad del estimador pero implican demasiadas líneas de código, lo que hace difícil apreciar el algoritmo de fondo.

Lo que sigue es una implementación que quita todo el fluff del proceso y deja los puros huesos—sorprendentemente (para mi) logrando sin embargo capacidad predictiva bastante competitiva. Pueden ver un notebook con algunas comparaciones en el mismo repositorio, o directamente acá.

Empecemos con la función de pérdida. Este ejemplo simple va a soportar tres tipos: error cuadrático, error absoluto, y la función logística (para clasificación binaria). Todos los ejemplos que vi haciendo research para este post se concentran en regresión, porque es tanto más fácil interpretar el algoritmo con la función de error cuadrática como hice más arriba, pero quiero que este post demuestre también como el algoritmo es extendible a otros problemas, dándole la versatilidad que lo hace famoso.

Definiremos el gradiente de cada función de pérdida de manera simple, con las siguientes funciones:

1def mse_gradient(y, f):
2 """
3 Error cuadrático, cuyo gradiente
4 es equivalente a los residuos
5 """
6 return y - f
7
8
9def lad_gradient(y, f):
10 """
11 Error absoluto, cuyo gradiente
12 es 1.0 si y - pred > 0.0 y sino -1.0
13 """
14 return 2.0 * (y - f > 0) - 1.0
15
16
17def logistic_gradient(y, f):
18 """
19 Función logistica para clasificación, cuyo gradiente es
20 y - expit(x). expit(x) es la funcion logistica: 1/(1+exp(-x))
21 """
22 return y - expit(f)

Con las distintas funciones para los pseudo-residuos listas, podemos pasar a implementar el algoritmo. Para el cuerpo principal, tendremos:

1# Voy a ignorar los imports en los segmentos de codigo mas abajo
2from sklearn.base import clone, BaseEstimator
3from sklearn.tree import DecisionTreeRegressor
4from sklearn.dummy import DummyRegressor
5
6import numpy as np
7from scipy.special import expit
8
9# BaseEstimator provee algunos metodos de utilidad no usados aca
10class GBoost(BaseEstimator):
11 def __init__(
12 self,
13 n_estimators=100, # Tambien referido como M arriba
14 learning_rate=0.1,
15 loss='ls',
16 tol=0.5 # Solo relevante para clasificacion
17 ):
18
19 # Inicializar argumentos
20 self.n_estimators = n_estimators
21 self.learning_rate = learning_rate
22 self.tol = tol
23 # Definir estimadores inicial y base
24 self.first_estimator = DummyRegressor(
25 strategy='constant',
26 constant=0
27 )
28 self.base_estimator = DecisionTreeRegressor(max_depth=1)
29 # Definir función de pérdida
30 if loss == 'lad':
31 self.gradient = lad_gradient
32 elif loss == 'ls':
33 self.gradient = mse_gradient
34 elif loss == 'logistic':
35 self.gradient = logistic_gradient
36 else:
37 raise AttributeError(f'Unknown loss function: {loss}')

Como estimador inicial usaremos simplemente un vector de 0’s que tomará el tamaño correspondiente a la variable objetivo. Esta es la primera gran diferencia con las implementaciones reales, que usan uno especifico para distintas funciones de pérdida (promedio para el error cuadrático, mediana para error absoluto, etc.). Luego, como estimador base se define un árbol de decisión simple que podrá realizar solamente un split---ciertamente muy restrictivo pero encarnando el espíritu de tomar estimadores débiles como base.

Luego vienen las implementaciones del algoritmo, repartidas entre una función para estimar: fit, y otra para predecir: predict, que como veremos ahora están ligadas por el proceso de estimación secuencial. Voy a empezar describiendo el algoritmo para el caso de regresión, y luego describiré los ajustes necesarios para realizar clasificación.

1def fit(self, X, y):
2 # Inicializar lista para estimadores debiles
3 self._estimators = [self.first_estimator]
4 # Paso 0
5 F = self.first_estimator.fit(X, y)
6 # Empezar proceso secuencial
7 for m in range(self.n_estimators):
8 # Paso 1: predecir usando todos los estimadores
9 # debiles entrenados hasta este momento
10 F = self.predict(X)
11 # Paso 2: Calcular pseudo-residuos del paso anterior
12 residuals = self.gradient(y, F)
13 # Paso 3: Estimar nuevo estimador débil a estos residuos
14 f = clone(self.base_estimator).fit(X, residuals)
15 # Paso 4: Guardar este estimador para
16 # incluirlo en la prediccion del paso 1
17 self._estimators.append(f)
18 return self
19
20def predict(self, X):
21 # Tomar el estimador base
22 f_0 = self._estimators[0].predict(X)
23 # Y sumarle las predicciones ajustadas de cada estimador debil
24 boosting = np.sum(
25 [ self.learning_rate * f.predict(X)
26 for f in self._estimators[1:] ],
27 axis=0
28 )
29 return f_0 + boosting

Toda la magia está realmente en las líneas resaltadas. El paso 1 consiste en llamar el método predict, que va a tomar todos los estimadores almacenados en la lista _estimators para calcular la predicción, es decir, equivale a Fm=F0+ηF1+ηF2+...F_m = F_0 + \eta F_1 + \eta F_2 + .... Luego se toma esta predicción y se llama al método gradient definido previamente para calcular los residuos en el paso 2, pero ahora debería quedar más claro que en realidad no tienen por qué ser necesariamente los residuos, sino que cualquier función que describa una divergencia entre la variable objetivo y la estimación. El paso 3 estima el nuevo estimador, pero usando los residuos como variable objetivo (la función auxiliar clone de scikit-learn solamente sirve para crear una copia).

Extendiendo a clasificación

Este estimador, tal cual está definido, puede realizar clasificación en lugar de regresión sin ninguna modificación necesaria---el algoritmo es el mismo. Sin embargo, la predicción va a ser imposible de interpretar, ya que debe ser convertida nuevamente por la función logística en una probabilidad, y luego en la clase correspondiente. Esto se puede lograr agregando los siguientes métodos a la clase:

1def _proba_to_class(self, sample):
2 """
3 Método auxiliar para convertir
4 probabilidad en clase, dada una tolerancia
5 """
6 return int(sample > self.tol)
7
8def predict_class(self, X):
9 """ Convierte prediccion en clase binaria {0, 1} """
10 # Convertir prediccion en probabilidad
11 predicted_probas = expit(self.predict(X))
12 # Convertir probabilidad en clase
13 return np.array(
14 [self._proba_to_class(sample) for sample in predicted_probas]
15 )

Con estos dos métodos finales terminamos la implementación, que con toda la funcionalidad básica y desempeño comparable a las de full scale, pero sin los controles ni optimizaciones, alcanza menos de 50 líneas.

1from sklearn.base import clone, BaseEstimator
2from sklearn.tree import DecisionTreeRegressor
3from sklearn.dummy import DummyRegressor
4
5import numpy as np
6from scipy.special import expit
7
8class GBoost(BaseEstimator):
9 def __init__(self, n_estimators=100, learning_rate=0.1, loss='ls', tol=0.5):
10 self.n_estimators = n_estimators
11 self.learning_rate = learning_rate
12 self.tol = tol # Relevante solo para clasificación
13 self.first_estimator = DummyRegressor(strategy='constant', constant=0)
14 self.base_estimator = DecisionTreeRegressor(max_depth=1)
15 if loss == 'lad':
16 self.gradient = lad_gradient
17 elif loss == 'ls':
18 self.gradient = mse_gradient
19 elif loss == 'logistic':
20 self.gradient = logistic_gradient
21 else:
22 raise AttributeError('Unknown loss function: {}'.format(loss))
23
24 def fit(self, X, y):
25 self._estimators = [self.first_estimator]
26 F = self.first_estimator.fit(X, y)
27 for m in range(self.n_estimators):
28 F = self.predict(X)
29 residuals = self.gradient(y, F)
30 f = clone(self.base_estimator).fit(X, residuals)
31 self._estimators.append(f)
32 return self
33
34 def predict(self, X):
35 f_0 = self._estimators[0].predict(X)
36 boosting = np.sum([self.learning_rate * f.predict(X) for f in self._estimators[1:]], axis=0)
37 return f_0 + boosting
38
39 def _proba_to_class(self, sample):
40 """ método auxiliar para convertir probabilidad en clase, dada una tolerancia """
41 return int(sample > self.tol)
42
43 def predict_class(self, X):
44 """ Convierte prediccion en clase binaria {0, 1} """
45 predicted_probas = expit(self.predict(X))
46 return np.array([self._proba_to_class(sample) for sample in predicted_probas])

Conclusiones

Espero que este proceso haya esclarecido un poco el algoritmo detrás de estos estimadores, y te de más confianza a la hora de usar una implementación real como XBoost, LightGBM o CatBoost. Cada una de estas aplica diferentes técnicas para mejorar la calidad predictiva, velocidad de entrenamiento, tratamiento de variables categóricas, etc. que le dan su merecido reconocimiento; sin embargo, en el fondo el algoritmo es el mismo y con una buena comprensión de este y los mecanismos que los subyacen se hace más abordable el proceso de descubrir sus funcionalidades y características.

En términos de estrategias generales a la hora de estimar, me gustaría escribir otro post dedicado a estas implementaciones y como usarlas con problemas de la vida real. En términos generales, esta técnica es mucho más susceptible a overfitting que random forests, por lo que es más importante tener cuidado—una técnica común es early stopping, que consiste en entregar al estimador una función de pérdida adicional junto con un set de validación que no es usado para entrenar (o usar out-of-bag como vimos en el post anterior en la serie sobre random forests). Con estos, se calcula un test score o test error que describe la capacidad de generalización del modelo durante cada etapa mm de entrenamiento, y al momento que este empieza a empeorar significa que el modelo esta overfitting y detiene el proceso. Esta técnica tiene 0 costo de implementar, y en cambio no solo garantiza la mejor generalización posible con los parámetros dados, sino que además disminuye potencialmente los tiempos de estimación y predicción, por lo que es altamente recomendada.

Ejemplo de overfitting

Otra consideración importante es que a menor learning rate, menor es el riesgo de overfitting también, sin embargo se van a requerir más estimadores para converger a un optimo.

More articles from Cristobal Aguirre

Scraping Twitter

Un programa sencillo escrito en Python para descarga información histórica de Twitter

June 4th, 2019 · 5 min read

ML Part 2: Random Forests

El segundo artículo en esta serie introductoria de Machine Learning, donde veremos en detalle el algoritmo conocido como Random Forests.

May 29th, 2019 · 5 min read
© 2019 Cristobal Aguirre
Link to $https://twitter.com/jcaguirre89Link to $https://github.com/jcaguirre89