Cuando noto que mi modelo está sobreajustadoa menudo pienso, “Es hora de regularizar”. Pero, ¿cómo decido qué método de regularización usar (L1, L2) y qué parámetros elegir? Normalmente, realizo una optimización de hiperparámetros mediante una búsqueda en cuadrícula para seleccionar la configuración. Sin embargo, ¿qué sucede si las variables independientes tienen diferentes escalas o distintos niveles de influencia? ¿Puedo diseñar una cuadrícula de hiperparámetros con diferentes coeficientes de regularización para cada variable? ¿Es factible este tipo de optimización en espacios de alta dimensión? ¿Y existen formas alternativas de diseñar la regularización? Exploremos esto con un ejemplo hipotético.
Mi ejemplo ficticio es un caso de uso de clasificación binaria con 3 variables explicativas. Cada una de estas variables es categórica y tiene 7 categorías diferentes. Mi caso de uso reproducible está en este computadora portátil. La función que genera el conjunto de datos es la siguiente:
import numpy as np
import pandas as pddef get_classification_dataset():
n_samples = 200
cats = ["a", "b", "c", "d", "e", "f"]
X = pd.DataFrame(
data={
"col1": np.random.choice(cats, size=n_samples),
"col2": np.random.choice(cats, size=n_samples),
"col3": np.random.choice(cats, size=n_samples),
}
)
X_preprocessed = pd.get_dummies(X)
theta = np.random.multivariate_normal(
np.zeros(len(cats) * X.shape[1]),
np.diag(np.array([1e-1] * len(cats) + [1] * len(cats) + [1e1] * len(cats))),
)
y = pd.Series(
data=np.random.binomial(1, expit(np.dot(X_preprocessed.to_numpy(), theta))),
index=X_preprocessed.index,
)
return X_preprocessed, y
A modo de información, elegí deliberadamente 3 valores diferentes para la matriz de covarianza theta para mostrar el beneficio del método de optimización bayesiano aproximado de Laplace. Si los valores fueran de alguna manera similares, el interés sería menor.
Junto con un modelo de referencia simple que predice el valor medio observado En el conjunto de datos de entrenamiento (utilizado con fines de comparación), opté por diseñar un modelo un poco más complejo. Decidí codificar en caliente las tres variables independientes y aplicar un modelo de regresión logística además de este preprocesamiento básico. Para la regularización, elegí un diseño L2 y busqué encontrar el coeficiente de regularización óptimo utilizando dos técnicas: búsqueda de cuadrícula y Laplace aproxima la optimización bayesianacomo ya habrás anticipado. Finalmente, evalué el modelo en un conjunto de datos de prueba utilizando dos métricas (seleccionadas arbitrariamente): pérdida de registro y AUC ROC.
Antes de presentar los resultados, primero echemos un vistazo más de cerca al modelo bayesiano y cómo lo optimizamos.
En el marco bayesiano, los parámetros ya no son constantes fijas, sino variables aleatorias. En lugar de maximizar la probabilidad de estimar estos parámetros desconocidos, ahora optimizamos la distribución posterior de los parámetros aleatorios, dados los datos observados. Esto nos obliga a elegir, a menudo de forma un tanto arbitraria, el diseño y los parámetros del anterior. Sin embargo, también es posible tratar los parámetros del anterior como variables aleatorias en sí mismas, como en Comienzodonde las capas de incertidumbre se siguen acumulando unas sobre otras…
En este estudio, he elegido el siguiente modelo:
Lógicamente he elegido un modelo de bernouilli para Y_i | θ, una normal centrada anterior correspondiente a una regularización L2 para θ | Σ y finalmente para Σ_i^{-1}, elegí un modelo Gamma. Elegí modelar la matriz de precisión en lugar de la matriz de covarianza como es tradicional en la literatura, como en la guía del usuario de scikit learn para la regresión lineal bayesiana. [2].
Además de este modelo escrito, supuse que Y_i e Y_j son condicionalmente (en θ) independientes, así como Y_i y Σ.
Probabilidad
Según el modelo, la probabilidad se puede escribir:
Para optimizar, necesitamos evaluar casi todos los términos, con la excepción de P(Y=y). Los términos de los numeradores se pueden evaluar utilizando el modelo elegido. Sin embargo, el término restante en el denominador no puede. Aquí es donde el Aproximación de Laplace entra en juego.
Aproximación de Laplace
Para evaluar el primer término del denominador, podemos aprovechar la aproximación de Laplace. Aproximamos la distribución de θ | Y, Σ por:
siendo θ* la moda de la distribución de densidad de θ | Y, Σ.
Aunque no conocemos la función de densidad, podemos evaluar la parte hessiana gracias a la siguiente descomposición:
Sólo necesitamos conocer los dos primeros términos del numerador para evaluar el hessiano, lo cual hacemos.
Para aquellos interesados en obtener más explicaciones, recomiendo la parte 4.4, “La aproximación de Laplace”, de Reconocimiento de patrones y aprendizaje automático de Christopher M. Bishop. [1]. Me ayudó mucho a entender la aproximación.
Probabilidad aproximada de Laplace
Finalmente, la probabilidad aproximada de optimización de Laplace es:
Una vez que aproximamos la función de densidad de θ | Y, Σ, finalmente podríamos evaluar la probabilidad en cualquier θ que queramos si la aproximación fuera precisa en todas partes. En aras de la simplicidad y debido a que la aproximación es precisa sólo cerca de la moda, evaluamos la probabilidad aproximada en θ*.
A continuación se muestra una función que evalúa esta pérdida para un dado (escalar) σ²=1/p (además de los valores observados dados, X e y, y de diseño, α y β).
import numpy as np
from scipy.stats import gammafrom module.bayesian_model import BayesianLogisticRegression
def loss(p, X, y, alpha, beta):
# computation of the loss for given values:
# - 1/sigma² (named p for precision here)
# - X: matrix of features
# - y: vector of observations
# - alpha: prior Gamma distribution alpha parameter over 1/sigma²
# - beta: prior Gamma distribution beta parameter over 1/sigma²
n_feat = X.shape[1]
m_vec = np.array([0] * n_feat)
p_vec = np.array(p * n_feat)
# computation of theta*
res = minimize(
BayesianLogisticRegression()._loss,
np.array([0] * n_feat),
args=(X, y, m_vec, p_vec),
method="BFGS",
jac=BayesianLogisticRegression()._jac,
)
theta_star = res.x
# computation the Hessian for the Laplace approximation
H = BayesianLogisticRegression()._hess(theta_star, X, y, m_vec, p_vec)
# loss
loss = 0
## first two terms: the log loss and the regularization term
loss += baysian_model._loss(theta_star, X, y, m_vec, p_vec)
## third term: prior distribution over sigma, written p here
out -= gamma.logpdf(p, a = alpha, scale = 1 / beta)
## fourth term: Laplace approximated last term
out += 0.5 * np.linalg.slogdet(H)[1] - 0.5 * n_feat * np.log(2 * np.pi)
return out
En mi caso de uso, he optado por optimizarlo mediante Optimizador de Adán, qué código se ha tomado de este repositorio.
def adam(
fun,
x0,
jac,
args=(),
learning_rate=0.001,
beta1=0.9,
beta2=0.999,
eps=1e-8,
startiter=0,
maxiter=1000,
callback=None,
**kwargs
):
"""``scipy.optimize.minimize`` compatible implementation of ADAM -
[http://arxiv.org/pdf/1412.6980.pdf].
Adapted from ``autograd/misc/optimizers.py``.
"""
x = x0
m = np.zeros_like(x)
v = np.zeros_like(x)for i in range(startiter, startiter + maxiter):
g = jac(x, *args)
if callback and callback(x):
break
m = (1 - beta1) * g + beta1 * m # first moment estimate.
v = (1 - beta2) * (g**2) + beta2 * v # second moment estimate.
mhat = m / (1 - beta1**(i + 1)) # bias correction.
vhat = v / (1 - beta2**(i + 1))
x = x - learning_rate * mhat / (np.sqrt(vhat) + eps)
i += 1
return OptimizeResult(x=x, fun=fun(x, *args), jac=g, nit=i, nfev=i, success=True)
Para esta optimización necesitamos la derivada de la pérdida anterior. No podemos tener una forma analítica así que decidí usar una aproximación numérica de la derivada.
Una vez que el modelo se entrena en el conjunto de datos de entrenamiento, es necesario hacer predicciones en el conjunto de datos de evaluación para evaluar su desempeño y comparar diferentes modelos. Sin embargo, no es posible calcular directamente la distribución real de un nuevo punto, ya que el cálculo es intratable.
Es posible aproximar los resultados con:
en vista de:
Elegí una variable previa no informativa en lugar de la variable aleatoria de precisión. El modelo ingenuo tiene un rendimiento deficiente, con una pérdida logarítmica de 0,60 y un AUC ROC de 0,50. El segundo modelo funciona mejor, con una pérdida logarítmica de 0,44 y un AUC ROC de 0,83, ambos cuando se hiperoptimiza mediante búsqueda de cuadrícula y optimización bayesiana. Esto indica que el modelo de regresión logística, que incorpora las variables dependientes, supera al modelo ingenuo. Sin embargo, no hay ninguna ventaja en utilizar la optimización bayesiana sobre la búsqueda en cuadrícula, por lo que continuaré con la búsqueda en cuadrícula por ahora. Gracias por leer.
… Pero espera, estoy pensando. ¿Por qué mis parámetros se regularizan con el mismo coeficiente? ¿No debería mi anterior depender de las variables dependientes subyacentes? Quizás los parámetros de la primera variable dependiente podrían tomar valores más altos, mientras que los de la segunda variable dependiente, con su menor influencia, deberían estar más cerca de cero. Exploremos estas nuevas dimensiones.
Hasta ahora hemos considerado dos técnicas, la búsqueda en grilla y la optimización bayesiana. Podemos utilizar estas mismas técnicas en dimensiones superiores.
Considerar nuevas dimensiones podría aumentar drásticamente la cantidad de nodos de mi red. Es por eso que la optimización bayesiana tiene sentido en dimensiones superiores para obtener los mejores coeficientes de regularización. En el caso de uso considerado, supuse que hay 3 parámetros de regularización, uno para cada variable independiente. Después de codificar una sola variable, supuse que todas las nuevas variables generadas compartían el mismo parámetro de regularización. De ahí los parámetros de regularización totales de 3, incluso si hay más de 3 columnas como entradas de la regresión logística.
Actualicé la función de pérdida anterior con el siguiente código:
import numpy as np
from scipy.stats import gammafrom module.bayesian_model import BayesianLogisticRegression
def loss(p, X, y, alpha, beta, X_columns, col_to_p_id):
# computation of the loss for given values:
# - 1/sigma² vector (named p for precision here)
# - X: matrix of features
# - y: vector of observations
# - alpha: prior Gamma distribution alpha parameter over 1/sigma²
# - beta: prior Gamma distribution beta parameter over 1/sigma²
# - X_columns: list of names of X columns
# - col_to_p_id: dictionnary mapping a column name to a p index
# because many column names can share the same p value
n_feat = X.shape[1]
m_vec = np.array([0] * n_feat)
p_list = []
for col in X_columns:
p_list.append(p[col_to_p_id[col]])
p_vec = np.array(p_list)
# computation of theta*
res = minimize(
BayesianLogisticRegression()._loss,
np.array([0] * n_feat),
args=(X, y, m_vec, p_vec),
method="BFGS",
jac=BayesianLogisticRegression()._jac,
)
theta_star = res.x
# computation the Hessian for the Laplace approximation
H = BayesianLogisticRegression()._hess(theta_star, X, y, m_vec, p_vec)
# loss
loss = 0
## first two terms: the log loss and the regularization term
loss += baysian_model._loss(theta_star, X, y, m_vec, p_vec)
## third term: prior distribution over 1/sigma² written p here
## there is now a sum as p is now a vector
out -= np.sum(gamma.logpdf(p, a = alpha, scale = 1 / beta))
## fourth term: Laplace approximated last term
out += 0.5 * np.linalg.slogdet(H)[1] - 0.5 * n_feat * np.log(2 * np.pi)
return out
Con este enfoque, las métricas evaluadas en el conjunto de datos de prueba son las siguientes: 0,39 y 0,88, que son mejores que el modelo inicial optimizado mediante una búsqueda en grilla y un enfoque bayesiano con un solo previo para todas las variables independientes.
El caso de uso se puede reproducir con este computadora portátil.
He creado un ejemplo para ilustrar la utilidad de la técnica. Sin embargo, no he podido encontrar un conjunto de datos del mundo real adecuado para demostrar plenamente su potencial. Mientras trabajaba con un conjunto de datos real, no pude obtener ningún beneficio significativo al aplicar esta técnica. Si encuentra uno, hágamelo saber. Me encantaría ver una aplicación en el mundo real de este método de regularización.
En conclusión, utilizar la optimización bayesiana (con la aproximación de Laplace si es necesario) para determinar los mejores parámetros de regularización puede ser una buena alternativa a los métodos tradicionales de ajuste de hiperparámetros. Al aprovechar los modelos probabilísticos, la optimización bayesiana no solo reduce el costo computacional sino que también mejora la probabilidad de encontrar valores de regularización óptimos, especialmente en dimensiones altas.
- Cristóbal M. Obispo. (2006). Reconocimiento de patrones y aprendizaje automático. Saltador.
- Guía del usuario de aprendizaje de ciencia ficción de regresión de crestas bayesianas: https://scikit-learn.org/1.5/modules/linear_model.html#bayesian-ridge-regression