Во-первых, ваш вопрос некорректен, потому что существует множество алгоритмов для решения Лассо.Самым популярным сейчас является координатный спуск.Вот скелет алгоритма (без критерия остановки).Я использовал numba / jit, потому что циклы for могут быть медленными в питоне.
import numpy as np
from numba import njit
@njit
def ST(x, u):
"Soft thresholding of x at level u"""
return np.sign(x) * np.maximum(np.abs(x) - u, 0.)
@njit
def cd_solver(X, y, alpha, max_iter):
n_samples, n_features = X.shape
beta = np.zeros(n_features)
R = y.copy() # residuals y - X @ beta
lc = (X ** 2).sum(axis=0) # lipschitz constants for coordinate descent
for t in range(max_iter):
for j in range(n_features):
old = beta[j]
beta[j] = ST(old + X[:, j].dot(R) / lc[j], alpha / lc[j])
# keep residuals up to date
if old != beta[j]:
R += (old - beta[j]) * X[:, j]
# I'll leave it up to you to implement a proper stopping criterion
return beta
X = np.random.randn(100, 200)
y = np.random.randn(100)
if not np.isfortran(X):
X = np.asfortranarray(X)
alpha_max = np.max(np.abs(X.T.dot(y)))
cd_solver(X, y, alpha_max / 2., 100)
Вы также можете попробовать использовать проксимальный градиент / ISTA, но по моему опыту это намного медленнее, чем CD.