5、LASSO模型选择：交叉验证-AIC-BIC

Posted 2023-03-21

参考技术A 5、LASSO模型选择：交叉验证-AIC-BIC

import time

import numpy as np

import matplotlib.pyplot as plt

from sklearn.linear_model import LassoCV, LassoLarsCV, LassoLarsIC

from sklearn import datasets

plt.rcParams['font.sans-serif'] = ['SimHei']

plt.rcParams['axes.unicode_minus'] = False

# 这是为了在执行np.log10时避免被零除

EPSILON = 1e-4

X, y = datasets.load_diabetes(return_X_y=True)

rng = np.random.RandomState(42)

X = np.c_[X, rng.randn(X.shape[0], 14)]  # 增加一些不好的特性

# 按照Lars的方法对数据进行规范化，以便进行比较

X /= np.sqrt(np.sum(X ** 2, axis=0))

# LassoLarsIC: 基于BIC/AIC准则的最小角回归

model_bic = LassoLarsIC(criterion='bic')

t1 = time.time()

model_bic.fit(X, y)

t_bic = time.time() - t1

alpha_bic_ = model_bic.alpha_

model_aic = LassoLarsIC(criterion='aic')

model_aic.fit(X, y)

alpha_aic_ = model_aic.alpha_

def plot_ic_criterion(model, name, color):

    criterion_ = model.criterion_

    plt.semilogx(model.alphas_ + EPSILON, criterion_, '--', color=color,

                linewidth=3, label='%s criterion' % name)

    plt.axvline(model.alpha_ + EPSILON, color=color, linewidth=3,

                label='alpha: %s estimate' % name)

    plt.xlabel(r'$\alpha$')

    plt.ylabel('criterion')

plt.figure()

plot_ic_criterion(model_aic, 'AIC', 'b')

plot_ic_criterion(model_bic, 'BIC', 'r')

plt.legend()

plt.title('信息-模型选择的标准 (训练时间: %.3fs)'

          % t_bic)

# LassoCV: 坐标下降

# 计算路径

print("Computing regularization path using the coordinate descent lasso...")

t1 = time.time()

model = LassoCV(cv=20).fit(X, y)

t_lasso_cv = time.time() - t1

# 显示结果

plt.figure()

ymin, ymax = 2300, 3800

plt.semilogx(model.alphas_ + EPSILON, model.mse_path_, ':')

plt.plot(model.alphas_ + EPSILON, model.mse_path_.mean(axis=-1), 'k',

        label='Average across the folds', linewidth=2)

plt.axvline(model.alpha_ + EPSILON, linestyle='--', color='k',

            label='alpha: CV estimate')

plt.legend()

plt.xlabel(r'$\alpha$')

plt.ylabel('Mean square error')

plt.title('每个折叠上的均方误差：坐标下降'

          '(训练时间 : %.2fs)' % t_lasso_cv)

plt.axis('tight')

plt.ylim(ymin, ymax)

# LassoLarsCV:最小角回归

# 计算路径

print("Computing regularization path using the Lars lasso...")

t1 = time.time()

model = LassoLarsCV(cv=20).fit(X, y)

t_lasso_lars_cv = time.time() - t1

# 显示结果

plt.figure()

plt.semilogx(model.cv_alphas_ + EPSILON, model.mse_path_, ':')

plt.semilogx(model.cv_alphas_ + EPSILON, model.mse_path_.mean(axis=-1), 'k',

            label='Average across the folds', linewidth=2)

plt.axvline(model.alpha_, linestyle='--', color='k',

            label='alpha CV')

plt.legend()

plt.xlabel(r'$\alpha$')

plt.ylabel('Mean square error')

plt.title('每折均方误差: Lars (训练时间 : %.2fs)'

          % t_lasso_lars_cv)

plt.axis('tight')

plt.ylim(ymin, ymax)

plt.show()

如何在执行 10 倍交叉验证时在每次拆分时获得 Lasso Regression 中的系数？

【中文标题】如何在执行 10 倍交叉验证时在每次拆分时获得 Lasso Regression 中的系数？

【英文标题】：How to get the coefficients in Lasso Regression at every split while performing 10 fold cross validation?

【发布时间】：2021-08-31 03:24:08

【问题描述】：

我正在执行随机搜索 cv 以在 Lasso 回归中查找 alpha 值，并且我正在执行 10 折交叉验证。有没有办法获取每个拆分的系数值，就像我们使用 cv_results 函数获取分数一样？

【问题讨论】：

@MustafaAydın 我相信 OP 想要在每次迭代中拟合 Lasso 估计器的 coef_。这无法从cv_results_ 中检索到，因为并非所有拟合的估计器都被保存。在我的回答中提供了解决方法。

@afsharov 哦，谢谢，我现在明白了，感谢您的回答，+1。​​

【参考方案1】：

在我看来，将系数保存为附加分数比修改估计器本身（如@afsharov's answer）更巧妙。定义记分器并将其传递给搜索

def coefs_scorer(estimator, X, y):
    return estimator.coef_

rs = RandomizedSearchCV(
    ...
    scoring='r2': 'r2', 'coefs': coefs_scorer,
    refit='r2',
)

失败，因为有记录员返回单个数字的检查。所以你需要解压缩系数，我最终得到了这个：

def coefs_scorer(estimator, X, y, i):
    return estimator.coef_[i]

from functools import partial
scoring = 'r2': 'r2'
for i in range(X_train.shape[1]):
    scoring[f'coefi'] = partial(coefs_scorer, i=i)

param_distributions = 'alpha': [0.01, 0.1, 1]
rs = RandomizedSearchCV(
    Lasso(),
    param_distributions=param_distributions,
    cv=2,
    n_iter=3,
    random_state=42,
    scoring=scoring,
    refit='r2',
)

请注意，对于多个指标，您需要指定用于改装的指标。由于所有额外的工作，我不太确定这是否比自定义类更好。不过它确实有一些优点：

如果您想挑选最佳估算器，则无需打包自定义类。 分数以编程方式保存，而不仅仅是打印出来。 由于它们是分数，因此您可以获得存储在 cv_results_ 中的各个折叠系数的平均值和标准偏差（当然，自己计算它们并不困难）。

缺点：

我们必须为每个功能指定一个指标。这很丑陋，但更糟糕的是，它假设您事先知道特征的数量（如果您的估算器是具有特征选择或某些特征工程步骤的管道，它将失败）。 如果您返回火车分数，您将复制cv_results_ 中的系数。 这些实际上不是分数，所以从语义上讲这是hacky。 记分员假定coef_ 存在并且是一维的。

【参考方案2】：

没有直接的方法可以通过RandomizedSearchCV 执行此操作。但是您可以通过定义自己的类来解决这个问题，例如调用 predict 函数时将系数打印到控制台：

from sklearn.linear_model import Lasso


class MyLasso(Lasso):

    def predict(self, X):
        print(self.coef_)
        return super().predict(X)

MyLasso 的行为与Lasso 相同，可以照常使用：

from sklearn.datasets import make_regression
from sklearn.model_selection import train_test_split, RandomizedSearchCV


X, y = make_regression(n_features=5, random_state=42)
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=42)

param_distributions = 'alpha': [0.01, 0.1, 1]
rs = RandomizedSearchCV(
    MyLasso(),
    param_distributions=param_distributions,
    cv=2,
    n_iter=3,
    random_state=42
)

rs.fit(X_train, y_train)

上述示例的输出（2 折交叉验证的三次迭代得到六个结果）：

[64.57650818 98.64237403 57.07123743 60.56898095 35.59985227]
[64.57001187 98.63679695 57.06557977 60.56304163 35.59888746]
[64.43774582 98.55938568 57.01219706 60.49221968 35.51151313]
[64.37690435 98.49805298 56.95345309 60.43375789 35.5018112 ]
[63.05012223 97.72950224 56.42179336 59.72460697 34.62812171]
[62.44582912 97.11061327 55.83218634 59.14092054 34.53104869]