# sklearn gridsearch lasso regression: find specific number of coefficients

I am using GridSearchCV and Lasso regression in order to fit a dataset composed out of Gaussians. I keep this example similar to this tutorial.

My goal is to find the best solution with a restricted number of non-zero coefficients, e.g. when I know beforehand, the data contains two Gaussians. So far, I used the grid search over the parameter space of number of features (or their spacing) and the width of the features, as well as the alpha parameter. Unfortunately, GridSearchCV does not return the coefficients for each fit, but only for the best one. What is the best way to find the fit which uses exactly two features?

Here my working example

import numpy as np
from matplotlib import pyplot as plt
from sklearn.linear_model import ElasticNet
from sklearn.pipeline import make_pipeline
from sklearn.base import BaseEstimator, TransformerMixin
from sklearn.model_selection import GridSearchCV

class GaussianFeatures(BaseEstimator, TransformerMixin):

def __init__(self, N, width=2.0):
self.N = N
self.width = width

@staticmethod
def _gaussian_basis(x, y, width, axis=None):
arg = (x - y) / width
return np.exp(-0.5 * np.sum(arg ** 2, axis))

def fit(self, X, y=None):
self.centers_ = np.linspace(X.min(), X.max(), self.N)
self.width_ = self.width
return self

def transform(self, X):
return self._gaussian_basis(X[:, :, np.newaxis], self.centers_,
self.width_, axis=1)

# create data
rng = np.random.RandomState(10)
x = 10 * rng.rand(50)
y = np.exp(-np.power((x - 2) / 1.5, 2.) / 2.) + \
np.exp(-np.power((x - 8) / 1.5, 2.) / 2.) + \
0.01 * rng.randn(50)

# initialize pipeline
widths = 2.0
number = 20
gauss_model = make_pipeline(GaussianFeatures(number, widths),
ElasticNet(alpha=0.005, l1_ratio=1, positive=True))

# initialize grid search
grid_width = np.arange(0.5, 5, 0.5)
grid_number = np.arange(1, 50, 1)
grid_alpha = [0.005, 0.004, 0.003]
grid_params = {
'gaussianfeatures__N': grid_number,
'gaussianfeatures__width': grid_width,
'elasticnet__alpha': grid_alpha
}
gs = GridSearchCV(gauss_model, grid_params, scoring='neg_mean_absolute_error')
gs.fit(x[:, np.newaxis], y)

best_widths = gs.best_params_['gaussianfeatures__width']
best_number = gs.best_params_['gaussianfeatures__N']
best_alpha = gs.best_params_['elasticnet__alpha']
best_score = gs.best_score_
print('GridSearch best params: ', gs.best_params_)
print('best score = ', best_score)

# plot best prediction
gauss_model = make_pipeline(GaussianFeatures(best_number, best_widths),
ElasticNet(alpha=best_alpha, l1_ratio=1, positive=True))
gauss_model.fit(x[:, np.newaxis], y)

xfit = np.linspace(0, 10, 1000)
yfit = gauss_model.predict(xfit[:, np.newaxis])

plt.figure()
plt.scatter(x, y)
plt.plot(xfit, yfit)
plt.show(block=False)

plt.figure()
plt.plot(gauss_model.steps[0][1].centers_, gauss_model.steps[1][1].coef_, '.', color='black')
plt.ylabel('theta')
plt.show(block=False)



output:

GridSearch best params:  {'elasticnet__alpha': 0.004, 'gaussianfeatures__N': 9, 'gaussianfeatures__width': 1.0}
best score =  -0.03473800683807379


• maybe SelectKBest with k=2 between GaussianFeatures and Lasso in the pipeline? Commented Mar 4, 2021 at 16:32
• If you want to only consider the models which use exactly two features, couldn't you just fix the N parameter to 2? Or am I misunderstanding you question? Commented Mar 4, 2021 at 16:36
• @Oxbowerce N is a bit misleading here, sorry. The parameter gives the number of features across the interval, and thus the spacing (or resolution) np.linspace(X.min(), X.max(), self.N). I could set it 'N=2' but would have to search the parameter space in terms of the location of the features. It could work, but I am not sure if this is a good approach. Commented Mar 5, 2021 at 10:30
• @BenReiniger thanks! I will look into that. Commented Mar 5, 2021 at 10:30