Exploring New Hyperparameter Dimensions with Laplace Approximated Bayesian Optimization | by Arnaud Capitaine | Jan, 2025

Is it higher than grid search?

Picture by creator from canva

After I discover my mannequin is overfitting, I typically assume, “It’s time to regularize”. However how do I determine which regularization methodology to make use of (L1, L2) and what parameters to decide on? Usually, I carry out hyperparameter optimization by the use of a grid search to pick out the settings. Nevertheless, what occurs if the unbiased variables have totally different scales or various ranges of affect? Can I design a hyperparameter grid with totally different regularization coefficients for every variable? Is any such optimization possible in high-dimensional areas? And are there alternative routes to design regularization? Let’s discover this with a hypothetical instance.

My fictional instance is a binary classification use case with 3 explanatory variables. Every of those variables is categorical and has 7 totally different classes. My reproducible use case is on this pocket book. The operate that generates the dataset is the next:

import numpy as np
import pandas as pd

def get_classification_dataset():
n_samples = 200
cats = ["a", "b", "c", "d", "e", "f"]
X = pd.DataFrame(
knowledge={
"col1": np.random.alternative(cats, dimension=n_samples),
"col2": np.random.alternative(cats, dimension=n_samples),
"col3": np.random.alternative(cats, dimension=n_samples),
}
)
X_preprocessed = pd.get_dummies(X)

theta = np.random.multivariate_normal(
np.zeros(len(cats) * X.form[1]),
np.diag(np.array([1e-1] * len(cats) + [1] * len(cats) + [1e1] * len(cats))),
)

y = pd.Collection(
knowledge=np.random.binomial(1, expit(np.dot(X_preprocessed.to_numpy(), theta))),
index=X_preprocessed.index,
)
return X_preprocessed, y

For info, I intentionally selected 3 totally different values for the theta covariance matrix to showcase the good thing about the Laplace approximated bayesian optimization methodology. If the values had been by some means comparable, the curiosity could be minor.

Together with a easy baseline mannequin that predicts the imply noticed worth on the coaching dataset (used for comparability functions), I opted to design a barely extra complicated mannequin. I made a decision to one-hot encode the three unbiased variables and apply a logistic regression mannequin on high of this primary preprocessing. For regularization, I selected an L2 design and aimed to seek out the optimum regularization coefficient utilizing two strategies: grid search and Laplace approximated bayesian optimization, as you will have anticipated by now. Lastly, I evaluated the mannequin on a take a look at dataset utilizing two metrics (arbitrarily chosen): log loss and AUC ROC.

Earlier than presenting the outcomes, let’s first take a better take a look at the bayesian mannequin and the way we optimize it.

Within the bayesian framework, the parameters are now not mounted constants, however random variables. As a substitute of maximizing the probability to estimate these unknown parameters, we now optimize the posterior distribution of the random parameters, given the noticed knowledge. This requires us to decide on, typically considerably arbitrarily, the design and parameters of the prior. Nevertheless, it is usually doable to deal with the parameters of the prior as random variables themselves — like in Inception, the place the layers of uncertainty maintain stacking on high of one another…

On this examine, I’ve chosen the next mannequin:

I’ve logically chosen a bernouilli mannequin for Y_i | θ, a standard centered prior akin to a L2 regularization for θ | Σ and eventually for Σ_i^{-1}, I selected a Gamma mannequin. I selected to mannequin the precision matrix as an alternative of the covariance matrix as it’s conventional within the literature, like in scikit study consumer information for the Bayesian linear regression [2].

Along with this written mannequin, I assumed the Y_i and Y_j are conditionally (on θ) unbiased in addition to Y_i and Σ.

Probability

In keeping with the mannequin, the probability can consequently be written:

In an effort to optimize, we have to consider almost the entire phrases, except P(Y=y). The phrases within the numerators will be evaluated utilizing the chosen mannequin. Nevertheless, the remaining time period within the denominator can’t. That is the place the Laplace approximation comes into play.

Laplace approximation

In an effort to consider the primary time period of the denominator, we are able to leverage the Laplace approximation. We approximate the distribution of θ | Y, Σ by:

with θ* being the mode of the mode the density distribution of θ | Y, Σ.

Although we have no idea the density operate, we are able to consider the Hessian half because of the next decomposition:

We solely have to know the primary two phrases of the numerator to judge the Hessian which we do.

For these excited by additional clarification, I counsel half 4.4, “The Laplace Approximation”, from Sample Recognition and Machine Studying from Christopher M. Bishop [1]. It helped me lots to grasp the approximation.

Laplace approximated probability

Lastly the Laplace approximated probability to optimize is:

As soon as we approximate the density operate of θ | Y, Σ, we may lastly consider the probability at no matter θ we wish if the approximation was correct in all places. For the sake of simplicity and since the approximation is correct solely near the mode, we consider approximated probability at θ*.

Right here after is a operate that evaluates this loss for a given (scalar) σ²=1/p (along with the given noticed, X and y, and design values, α and β).

import numpy as np
from scipy.stats import gamma

from 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 right here)
# - X: matrix of options
# - 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.form[1]
m_vec = np.array([0] * n_feat)
p_vec = np.array(p * n_feat)

# computation of theta*
res = decrease(
BayesianLogisticRegression()._loss,
np.array([0] * n_feat),
args=(X, y, m_vec, p_vec),
methodology="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 phrases: the log loss and the regularization time period
loss += baysian_model._loss(theta_star, X, y, m_vec, p_vec)
## third time period: prior distribution over sigma, written p right here
out -= gamma.logpdf(p, a = alpha, scale = 1 / beta)
## fourth time period: Laplace approximated final time period
out += 0.5 * np.linalg.slogdet(H)[1] - 0.5 * n_feat * np.log(2 * np.pi)

return out

In my use case, I’ve chosen to optimize it by the use of Adam optimizer, which code has been taken from this repo.

def adam(
enjoyable,
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.decrease`` appropriate implementation of ADAM -
[http://arxiv.org/pdf/1412.6980.pdf].
Tailored from ``autograd/misc/optimizers.py``.
"""
x = x0
m = np.zeros_like(x)
v = np.zeros_like(x)

for i in vary(startiter, startiter + maxiter):
g = jac(x, *args)

if callback and callback(x):
break

m = (1 - beta1) * g + beta1 * m # first second estimate.
v = (1 - beta2) * (g**2) + beta2 * v # second second 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, enjoyable=enjoyable(x, *args), jac=g, nit=i, nfev=i, success=True)

For this optimization we’d like the spinoff of the earlier loss. We can’t have an analytical kind so I made a decision to make use of a numerical approximation of the spinoff.

As soon as the mannequin is skilled on the coaching dataset, it’s essential to make predictions on the analysis dataset to evaluate its efficiency and examine totally different fashions. Nevertheless, it’s not doable to straight calculate the precise distribution of a brand new level, because the computation is intractable.

It’s doable to approximate the outcomes with:

contemplating:

I selected an uninformative prior over the precision random variable. The naive mannequin performs poorly, with a log lack of 0.60 and an AUC ROC of 0.50. The second mannequin performs higher, with a log lack of 0.44 and an AUC ROC of 0.83, each when hyperoptimized utilizing grid search and bayesian optimization. This means that the logistic regression mannequin, which includes the dependent variables, outperforms the naive mannequin. Nevertheless, there is no such thing as a benefit to utilizing bayesian optimization over grid search, so I’ll proceed with grid seek for now. Thanks for studying.

… However wait, I’m considering. Why are my parameters regularized with the identical coefficient? Shouldn’t my prior depend upon the underlying dependent variables? Maybe the parameters for the primary dependent variable may take larger values, whereas these for the second dependent variable, with its smaller affect, needs to be nearer to zero. Let’s discover these new dimensions.

Up to now now we have thought-about two strategies, the grid search and the bayesian optimization. We will use these similar strategies in larger dimensions.

Contemplating new dimensions may dramatically improve the variety of nodes of my grid. Because of this the bayesian optimization is smart in larger dimensions to get one of the best regularization coefficients. Within the thought-about use case, I’ve supposed there are 3 regularization parameters, one for every unbiased variable. After encoding a single variable, I assumed the generated new variables all shared the identical regularization parameter. Therefore the full regularization parameters of three, even when there are greater than 3 columns as inputs of the logistic regression.

I up to date the earlier loss operate with the next code:

import numpy as np
from scipy.stats import gamma

from 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 right here)
# - X: matrix of options
# - y: vector of observations
# - alpha: prior Gamma distribution alpha parameter over 1/sigma²
# - beta: prior Gamma distribution beta parameter over 1/sigma²
# - X_columns: record of names of X columns
# - col_to_p_id: dictionnary mapping a column identify to a p index
# as a result of many column names can share the identical p worth

n_feat = X.form[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 = decrease(
BayesianLogisticRegression()._loss,
np.array([0] * n_feat),
args=(X, y, m_vec, p_vec),
methodology="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 phrases: the log loss and the regularization time period
loss += baysian_model._loss(theta_star, X, y, m_vec, p_vec)
## third time period: prior distribution over 1/sigma² written p right here
## there may be now a sum as p is now a vector
out -= np.sum(gamma.logpdf(p, a = alpha, scale = 1 / beta))
## fourth time period: Laplace approximated final time period
out += 0.5 * np.linalg.slogdet(H)[1] - 0.5 * n_feat * np.log(2 * np.pi)

return out

With this strategy, the metrics evaluated on the take a look at dataset are the next: 0.39 and 0.88, that are higher than the preliminary mannequin optimized by the use of a grid search and a bayesian strategy with solely a single prior for all of the unbiased variables.

Metrics achieved with the totally different strategies on my use case.

The use case will be reproduced with this pocket book.

I’ve created an instance for instance the usefulness of the approach. Nevertheless, I’ve not been capable of finding an appropriate real-world dataset to completely exhibit its potential. Whereas I used to be working with an precise dataset, I couldn’t derive any important advantages from making use of this system. In case you come throughout one, please let me know — I’d be excited to see a real-world software of this regularization methodology.

In conclusion, utilizing bayesian optimization (with Laplace approximation if wanted) to find out one of the best regularization parameters could also be a great various to conventional hyperparameter tuning strategies. By leveraging probabilistic fashions, bayesian optimization not solely reduces the computational price but in addition enhances the probability of discovering optimum regularization values, particularly in excessive dimension.

  1. Christopher M. Bishop. (2006). Sample Recognition and Machine Studying. Springer.
  2. Bayesian Ridge Regression scikit-learn consumer information: https://scikit-learn.org/1.5/modules/linear_model.html#bayesian-ridge-regression