多个目标岭回归如何在scikit中学习？-Java 学习之路

我正在努力理解以下内容：

Scikit-learn为Ridge Regression提供了一个多输出版本，只需移交一个2D数组[n_samples，n_targets]，但它是如何实现的？

http://scikit-learn.org/stable/modules/generated/sklearn.linear_model.Ridge.html

假设每个目标的每个回归都是独立的，这是正确的吗？在这些情况下，我如何调整它以使用每个回归的单个alpha正则化参数？如果我使用GridSeachCV，我将不得不移交一个可能的正则化参数矩阵，或者它将如何工作？

在此先感谢 - 我一直在寻找几个小时，但在这个主题上找不到任何东西 .

1 回答

我会对此进行一次拍摄，因为我已经为自己的工作调查了一下 . 我会将问题细分为部分，这样您只能查看您感兴趣的部分：

Q1：多输出岭回归中每个目标（又名输出）的回归是否独立？

A1：我认为M输出的典型多输出线性回归与M独立单输出线性回归相同 . 我认为是这种情况，因为多输出情况的普通最小二乘的表达式与M个独立的单输出情况的（总和）的表达式相同 . 为了激发这一点，让我们考虑一个没有正规化的愚蠢的双变量输出案例 .

让我们考虑两个列向量输入 x 1和 x 2，以及相应的权向量 w 1和 w 2 .

这些给出了两个单变量输出y1 = x 1 w 1T e1和y2 = x 2 w 2T e2，其中e s是独立误差 .

平方错误的总和写成：

e12 e22 =（y1 - x 1 w 1T）2（y2 - x 2 w 2T）2

我们可以看到，这只是两个独立回归的平方误差的总和 . 现在，为了优化我们在权重方面的差异并设置为零 . 由于y1不依赖于 w 2，反之亦然y2和 w 1，因此可以针对每个目标独立地执行优化 .

我在这里考虑了一个样本用于说明，但更多样本没有太多变化 . 你可以自己写出来 . 在表格中添加惩罚期限 w 1 |或者 w 2 |也不会改变这一点，因为对y2和 w 1仍然没有 w 2依赖于y1，反之亦然 .

好的，这就是那种能让你获得C-的证据（有一位慷慨的教授） . 只要这反映了sklearn，手动实现独立的回归和内置的多输出支持就会得到相同的结果 . 所以让我们用一些代码检查一下（我使用pyp.7与Jupyter）：

我们需要的东西

import numpy as np
 import matplotlib.pyplot as plt
 from sklearn import linear_model, model_selection

设置数据

## set up some test data
# N samples, K features, M outputs (aka targets)
T = 1000
K = 100
M = 500

#get the samples from independent, multivariate normal
means_X = np.zeros(K)
cov_X = np.identity(K) 
X = np.random.multivariate_normal(means_X,cov_X,T)

#Make up some random weights. 
#Here I use an exponential form since that means some would be quite small, and thus regularization is likely to help
#However for the purposes of the example it doesn't really matter

#exponential weights
W = 2.0 ** np.random.randint(-10,0,M * K) 

#shape into a weight matrix correctly
W = W.reshape((K,M))

# get the ouput - apply linear transformation
Y = np.matmul(X, W)

# add a bit of noise to the output
noise_level = 0.1
noise_means = np.zeros(M)
noise_cov = np.identity(M) 

Y_nse = Y + noise_level * np.random.multivariate_normal(noise_means,noise_cov,T)

# Start with one alpha value for all targets 
alpha = 1

使用内置多输出支持的sklearn

#%%timeit -n 1 -r 1
# you can uncomment the above to get timming but note that this runs on a seperate session so 
# the results won't be available here 
## use built in MIMO support 

built_in_MIMO = linear_model.Ridge(alpha = alpha)
built_in_MIMO.fit(X, Y_nse)

使用独立运行优化输出

# %%timeit -n 1 -r 1 -o
## manual mimo
manual_MIMO_coefs = np.zeros((K,M))

for output_index in range(M):

    manual_MIMO = linear_model.Ridge(alpha = alpha)
    manual_MIMO.fit(X, Y_nse[:,output_index]) 
    manual_MIMO_coefs[:,output_index] = manual_MIMO.coef_

比较估计（未包括的图）

manual_MIMO_coefs_T = manual_MIMO_coefs.T

## check the weights. Plot a couple
check_these_weights = [0, 10]
plt.plot(built_in_MIMO.coef_[check_these_weights[0],:],'r')
plt.plot(manual_MIMO_coefs_T[check_these_weights[0],:], 'k--')

# plt.plot(built_in_MIMO.coef_[check_these_weights[1],:],'b')
# plt.plot(manual_MIMO_coefs_T[check_these_weights[1],:], 'y--')

plt.gca().set(xlabel = 'weight index', ylabel = 'weight value' )
plt.show()

print('Average diff between manual and built in weights is %f ' % ((built_in_MIMO.coef_.flatten()-manual_MIMO_coefs_T.flatten()) ** 2).mean())


## FYI, our estimate are pretty close to the orignal too, 
plt.clf()
plt.plot(built_in_MIMO.coef_[check_these_weights[1],:],'b')
plt.plot(W[:,check_these_weights[1]], 'y--')
plt.gca().set(xlabel = 'weight index', ylabel = 'weight value' )
plt.legend(['Estimated', 'True'])

plt.show()

print('Average diff between manual and built in weights is %f ' % ((built_in_MIMO.coef_.T.flatten()-W.flatten()) ** 2).mean())

输出是（不包括这里的情节）

Average diff between manual and built in weights is 0.000000 

Average diff between manual and built in weights is 0.000011

所以我们看到内置的sklearn估计与我们的手册相同 . 然而，内置的一个更快，因为它使用矩阵代数一次解决整个事物，而不是我在这里使用的循环（对于脊正则化的矩阵公式，请参阅Tikhonov正则化的维基） . 您可以通过取消注释上面的%% timeit魔法来自行检查）

Q2：我们如何为每个回归使用单独的alpha正则化参数？

A2：sklearn Ridge接受每个输出（目标）的不同正则化 . 例如，继续上面的代码，为每个输出使用不同的alpha：

# now try different alphas for each target.
# Simply randomly assign them between min and max range 
min_alpha = 0
max_alpha = 50
alphas = 2.0 ** np.random.randint(min_alpha,max_alpha,M)
built_in_MIMO = linear_model.Ridge(alpha = alphas)
built_in_MIMO.fit(X, Y_nse)

如果将其与M个独立回归的手动实现进行比较，每个回归都有自己的alpha：

manual_MIMO_coefs = np.zeros((K,M))

for output_index in range(M):

    manual_MIMO = linear_model.Ridge(alpha = alphas[output_index])
    manual_MIMO.fit(X, Y_nse[:,output_index]) 
    manual_MIMO_coefs[:,output_index] = manual_MIMO.coef_

你得到相同的结果：

manual_MIMO_coefs_T = manual_MIMO_coefs.T

## check the weights. 
print('Average diff between manual and built in weights is %f ' % ((built_in_MIMO.coef_.flatten()-manual_MIMO_coefs_T.flatten()) ** 2).mean())

Average diff between manual and built in weights is 0.000000

所以这些都是一样的 .

However ，在这种情况下，性能很大程度上取决于求解器（由@Vivek Kumar直觉） .

默认情况下，Ridge.fit（）使用Cholesky分解（至少对于非稀疏数据），查看github上的代码（https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/linear_model/ridge.py中的_solve_cholesky），我看到当为每个目标单独选择alpha时，sklearn实际上确实适合他们 . 我不是后者 .

但是，对于不同的求解器，例如基于SVD（_solve_svd（）），代码似乎已经将不同的alpha合并到问题的矩阵公式中，所以整个东西只运行一次 . 这意味着当为每个输出单独选择alpha时，并且当有许多输出时，svd求解器可以快得多 .

问题3：我如何使用GridSeachCV？我是否交出了可能的正则化参数矩阵？

A3：我没有使用内置网格搜索，因为它不太适合我的问题 . 但是，通过上述说明，可以直接实现这一点;只需使用sklearn.model_selection.KFold（）或类似内容获得一些CV折叠，然后使用不同的alpha训练每个折叠：

from sklearn import metrics, model_selection
# just two folds for now
n_splits = 2
#logarithmic grid
alphas = 2.0 ** np.arange(0,10) 
kf = model_selection.KFold(n_splits=n_splits)

# generates some folds
kf.get_n_splits(X)

# we will keep track of the performance of each alpha here
scores = np.zeros((n_splits,alphas.shape[0],M))

#loop over alphas and folds
for j,(train_index, test_index) in enumerate(kf.split(X)):

    for i,alpha in enumerate(alphas):

        cv_MIMO = linear_model.Ridge(alpha = alpha)
        cv_MIMO.fit(X[train_index,:], Y_nse[train_index,:]) 
        cv_preds = cv_MIMO.predict(X[test_index,:])
        scores[j,i,:] = metrics.r2_score(Y_nse[test_index,:], cv_preds, multioutput='raw_values')

## mean CV score  
mean_CV_score = scores.mean(axis = 0)
# best alpha for each target
best_alpha_for_target = alphas[np.argmax(mean_CV_score,axis = 0)]

我匆匆写了这个，所以仔细检查一下 . 请注意，我们需要使用度量模块，因为内置的Ridge.score（）平均所有目标，这是我们不想要的 . 通过使用multioutput ='raw_values'选项，我们保留每个目标的原始值 .

希望有所帮助！

回复于 2024-05-05T01:25:22+08:00

多个目标岭回归如何在scikit中学习？

1 回答

相关问题