为了探索LASSO回归的工作原理，我编写了一小段代码，应该LASSO通过选择最佳 alpha 参数来优化回归。
我无法弄清楚为什么LASSO在交叉验证后回归给我的 alpha 参数结果如此不稳定。
这是我的 Python 代码：
from sklearn.linear_model import Lasso
from sklearn.cross_validation import KFold
from matplotlib import pyplot as plt
# generate some sparse data to play with
import numpy as np
import pandas as pd
from scipy.stats import norm
from scipy.stats import uniform
### generate your own data here
n = 1000
x1x2corr = 1.1
x1x3corr = 1.0
x1 = range(n) + norm.rvs(0, 1, n) + 50
x2 = map(lambda aval: aval*x1x2corr, x1) + norm.rvs(0, 2, n) + 500
y = x1 + x2 #+ norm.rvs(0,10, n)
Xdf = pd.DataFrame()
Xdf['x1'] = x1
Xdf['x2'] = x2
X = Xdf.as_matrix()
# Split data in train set and test set
n_samples = X.shape[0]
X_train, y_train = X[:n_samples / 2], y[:n_samples / 2]
X_test, y_test = X[n_samples / 2:], y[n_samples / 2:]
kf = KFold(X_train.shape[0], n_folds = 10, )
alphas = np.logspace(-16, 8, num = 1000, base = 2)
e_alphas = list()
e_alphas_r = list() # holds average r2 error
for alpha in alphas:
lasso = Lasso(alpha=alpha, tol=0.004)
err = list()
err_2 = list()
for tr_idx, tt_idx in kf:
X_tr, X_tt = X_train[tr_idx], X_test[tt_idx]
y_tr, y_tt = y_train[tr_idx], y_test[tt_idx]
lasso.fit(X_tr, y_tr)
y_hat = lasso.predict(X_tt)
# returns the coefficient of determination (R^2 value)
err_2.append(lasso.score(X_tt, y_tt))
# returns MSE
err.append(np.average((y_hat - y_tt)**2))
e_alphas.append(np.average(err))
e_alphas_r.append(np.average(err_2))
## print out the alpha that gives the minimum error
print 'the minimum value of error is ', e_alphas[e_alphas.index(min(e_alphas))]
print ' the minimizer is ', alphas[e_alphas.index(min(e_alphas))]
## <<< plotting alphas against error >>>
plt.figsize = (15, 15)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(alphas, e_alphas, 'b-')
ax.plot(alphas, e_alphas_r, 'g--')
ax.set_ylim(min(e_alphas),max(e_alphas))
ax.set_xlim(min(alphas),max(alphas))
ax.set_xlabel("alpha")
plt.show()
如果您重复运行此代码，它会为 alpha 提供截然不同的结果：
>>>
the minimum value of error is 3.99254192539
the minimizer is 1.52587890625e-05
>>> ================================ RESTART ================================
>>>
the minimum value of error is 4.07412455842
the minimizer is 6.45622425334
>>> ================================ RESTART ================================
>>>
the minimum value of error is 4.25898253597
the minimizer is 1.52587890625e-05
>>> ================================ RESTART ================================
>>>
the minimum value of error is 3.79392968781
the minimizer is 28.8971008254
>>>
为什么 alpha 值不能正确收敛？我知道我的数据是合成的，但分布是一样的。x1此外，和的变化非常小x2。
是什么导致这种情况如此不稳定？
用 R 编写的相同内容给出了不同的结果 - 它始终返回 alpha 的最高可能值作为“optimal_alpha”。
我也在R中写了这个，这给了我一个稍微不同的答案，我不知道为什么？
library(glmnet)
library(lars)
library(pracma)
set.seed(1)
k = 2 # number of features selected
n = 1000
x1x2corr = 1.1
x1 = seq(n) + rnorm(n, 0, 1) + 50
x2 = x1*x1x2corr + rnorm(n, 0, 2) + 500
y = x1 + x2
filter_out_label <- function(col) {col!="y"}
alphas = logspace(-5, 6, 100)
for (alpha in alphas){
k = 10
optimal_alpha = NULL
folds <- cut(seq(1, nrow(df)), breaks=k, labels=FALSE)
total_mse = 0
min_mse = 10000000
for(i in 1:k){
# Segement your data by fold using the which() function
testIndexes <- which(folds==i, arr.ind=TRUE)
testData <- df[testIndexes, ]
trainData <- df[-testIndexes, ]
fit <- lars(as.matrix(trainData[Filter(filter_out_label, names(df))]),
trainData$y,
type="lasso")
# predict
y_preds <- predict(fit, as.matrix(testData[Filter(filter_out_label, names(df))]),
s=alpha, type="fit", mode="lambda")$fit # default mode="step"
y_true = testData$y
residuals = (y_true - y_preds)
mse=sum(residuals^2)
total_mse = total_mse + mse
}
if (total_mse < min_mse){
min_mse = total_mse
optimal_alpha = alpha
}
}
print(paste("the optimal alpha is ", optimal_alpha))
上面 R 代码的输出是：
> source('~.....')
[1] "the optimal alpha is 1e+06"
事实上，无论我为“ alphas = logspace(-5, 6, 100)”行设置什么，我总是取回 alpha 的最高值。
我想这里实际上有两个不同的问题：
为什么用 Python 编写的版本的 alpha 值如此不稳定？
为什么用 R 编写的版本会给我不同的结果？（我意识到该
logspace函数与Rto不同python，但编写的版本R总是为我alpha提供最佳 alpha 值的最大值，而 python 版本没有）。
知道这些就太好了……
