从零开始学Python数据分析【25】--岭回归及LASSO回归（实战部分）

作者：刘顺祥

个人微信公众号：每天进步一点点2015

前文传送门：

前言

---

      在《从零开始学Python【24】--岭回归及LASSO回归（理论部分）》一文中我们详细介绍了关于岭回归和LASSO回归的理论知识，其实质就是在线性回归的基础上添加了2范数和1范数的惩罚项。这两个模型的关键点是找到一个合理的lambda系数，来平衡模型的方差和偏差，从而得到比较符合实际的回归系数。本期是基于之前讨论的理论部分，采用Python和R语言，完成对岭回归和LASSO回归的实战。文中所利用的数据集来自R语言ISLR包中的Hitters数据集，描述的是关于棒球运动员收入的相关信息，数据集和脚本可以之文末查看链接获得。

岭回归

---

      原始数据集存在收入的缺失，我们不妨先把这样的观测删除，同时，数据集中含有离散变量，需要将这些变量转换为哑变量后方可建模，故第一步需要对原始数据集进行清洗。

# ===== Python3 =====

# 导入第三方包
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn.cross_validation import train_test_split
from sklearn.linear_model import Ridge,RidgeCV
from sklearn.linear_model import Lasso,LassoCV
from sklearn.metrics import mean_squared_error

# 读取数据
df = pd.read_csv('Hitters.csv')
# 哑变量处理
dummies = pd.get_dummies(df[['League', 'Division', 'NewLeague']])
# 将原始数据集与哑变量数据合并起来
mydf = df.join(dummies)
# 缺失值删除
mydf = mydf.dropna()
# 删除不必要的变量(字符串变量和各哑变量中的一个变量)
mydf = mydf.drop([ 'League', 'Division', 'NewLeague', 'League_N', 'Division_W', 'NewLeague_N'], axis = 1)
# 前5行展示
mydf.head()

上面的数据集清洗完毕，展现的是干净数据的前5行信息，下面要基于这个数据集进行建模。建模之前还需要将数据集拆分为两部分，一部分用于建模，另一部分用于模型的测试。

# 将数据集拆分成训练集和测试集
predictors = list(mydf.columns) predictors.remove('Salary') # Salay变量为因变量，故需要排除

X_train, X_test, y_train, y_test = train_test_split(mydf[predictors],mydf['Salary'], train_size = 0.8, random_state = 1234 )

• 基于可视化，选择lambda参数

# 通过不确定的alphas值，生成不同的岭回归模型
alphas = 10**np.linspace(-3,3,100)
ridge_cofficients = []

for alpha in alphas:
    ridge = Ridge(alpha = alpha, normalize=True)
    ridge.fit(X_train, y_train)
    ridge_cofficients.append(ridge.coef_)

# 绘制alpha的对数与回归系数的关系
# 中文乱码和坐标轴负号的处理
plt.rcParams['font.sans-serif'] = ['Microsoft YaHei']
plt.rcParams['axes.unicode_minus'] = False
# 设置绘图风格
plt.style.use('ggplot')
plt.plot(alphas, ridge_cofficients)
plt.xscale('log')
plt.axis('tight')
plt.title('alpha系数与岭回归系数的关系')
plt.xlabel('Log Alpha')
plt.ylabel('Cofficients')
plt.show()

从上面的图形结果来看，alpha在10附近时，所有的自变量系数基本趋于稳定（但也不能完全确定是这个值）。接下来，我们采用交叉验证（CV）方法确定最佳的lambda值。

基于CV选择lambda值

# 为了找到最佳的lambda值，我们采用交叉验证方法

# 岭回归模型的交叉验证
ridge_cv = RidgeCV(alphas = alphas, normalize=True, scoring='mean_squared_error', cv = 10)
ridge_cv.fit(X_train, y_train)
# 取出最佳的lambda值ridge_best_alpha = ridge_cv.alpha_
ridge_best_alpha

不出所料，得到的lambda值确实在10附近，这里最佳的lambda值为10。下面，我们要基于这个最佳的lambda值进入岭回归模型的创建和模型验证的阶段。

# 基于最佳的lambda值建模
ridge = Ridge(alpha = ridge_best_alpha, normalize=True)
ridge.fit(X_train, y_train)
# 岭回归系数
ridge.coef_

# 预测
ridge_predict = ridge.predict(X_test)
# 预测效果验证
RMSE = np.sqrt(mean_squared_error(y_test,ridge_predict))
RMSE

经过模型的验证，得到的RMSE为319.9。接下来，我们利用同样的逻辑，对比一下LASSO回归模型的效果。

LASSO回归

---

• 基于可视化，选择lambda参数

# 通过不确定的alphas值，生成不同的LASSO回归模型
alphas = 10**np.linspace(-3,3,100)
lasso_cofficients = []

for alpha in alphas:
    lasso = Lasso(alpha = alpha, normalize=True, max_iter=10000)
    lasso.fit(X_train, y_train)
    lasso_cofficients.append(lasso.coef_)

# 绘制alpha的对数与回归系数的关系# 中文乱码和坐标轴负号的处理
plt.rcParams['font.sans-serif'] = ['Microsoft YaHei']
plt.rcParams['axes.unicode_minus'] = False
# 设置绘图风格
plt.style.use('ggplot')
plt.plot(alphas, lasso_cofficients)
plt.xscale('log')
plt.axis('tight')
plt.title('alpha系数与LASSO回归系数的关系')
plt.xlabel('Log Alpha')
plt.ylabel('Cofficients')
plt.show()

从图形结果来看，lambda值应该在1附近，此时LASSO回归的系数也基本趋于稳定（但也不能完全确定是这个值）。同样，我们利用CV方法，来寻找最佳的lambda值。

基于CV选择lambda值

# LASSO回归模型的交叉验证
lasso_cv = LassoCV(alphas = alphas, normalize=True, cv = 10, max_iter=10000)
lasso_cv.fit(X_train, y_train)
# 取出最佳的lambda值
lasso_best_alpha = lasso_cv.alpha_
lasso_best_alpha

通过CV方法得到的lambda结果是0.23，这与与我们看图得到的1这个值还是有一点差异的。下面，我们就基于交叉验证得到的最佳lambda值重新构造LASSO回归模型。

# 基于最佳的lambda值建模
lasso = Lasso(alpha = lasso_best_alpha, normalize=True, max_iter=10000)
lasso.fit(X_train, y_train)
# 岭回归系数lasso.coef_

# 预测
lasso_predict = lasso.predict(X_test)
# 预测效果验证
RMSE = np.sqrt(mean_squared_error(y_test,lasso_predict))
RMSE

对LASSO回归模型进行验证发现得到的RMSE更小，说明LASSO回归模型的拟合效果会更贴近于Hitters数据集的原貌。

      上面的内容是基于Python工具对岭回归模型和LASSO回归模型的实战，接下来，我们再利用R语言对上面的过程再作一次复现，希望对R语言感兴趣的朋友能够有帮助。

R语言对比

---

      由于上面的逻辑我们已经通过Python进行了一一说明，这里就不再赘述，只给出R语言代码仅供参考。

• 数据清洗

# 加载第三方包
library(caret)
library(glmnet)
library(ISLR)

# 哑变量处理
dummies <- dummyVars(~League+Division+NewLeague, data = Hitters)
dummies <- predict(dummies, newdata = Hitters)
# 数据合并
Hitters_dummy <- cbind(Hitters, dummies)

# 删除缺失值
Hitters_dummy <- na.omit(Hitters_dummy)
# 删除不必要的变量
Hitters_dummy <- subset(Hitters_dummy, 
                       select = -c(League,Division,NewLeague,League.N,Division.W,NewLeague.N))
head(Hitters_dummy)

• 岭回归—可视化选lambda

# 构建训练集和测试集set.seed(1) index <- sample(1:nrow(Hitters_dummy), size = 0.8*nrow(Hitters_dummy)) train <- Hitters_dummy[index,] test <- Hitters_dummy[-index,] # 绘制lambda值与岭回归系数的关系
fit_ridge <-  glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 0) plot(fit_ridge,xvar = 'lambda',label = T)

• 岭回归—交叉验证选lambda

# 岭回归的交叉验证，确定最佳的lambda值
fit_ridge_cv <- cv.glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 0) best_lambda_ridge <- fit_ridge_cv$lambda.min best_lambda_ridge

• 岭回归—建模与验证

# 根据最佳lambda构建岭回归模型
fit_ridge <- glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 0)
coeff_ridge <- predict(fit_ridge, s = best_lambda_ridge, type = 'coefficients')
coeff_ridge

# 模型评估
pred_ridge <- predict(fit_ridge, s = best_lambda_ridge, newx = as.matrix(test[,-17]))
RMSE <- sqrt(mean((test$Salary-pred_ridge)**2))
RMSE

• LASSO—可视化选lambda

# 绘制lambda值与LASSO回归系数的关系
fit_lasso <-  glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 1) plot(fit_lasso,xvar = 'lambda',label = T)

• LASSO—交叉验证选lambda

# 岭回归的交叉验证，确定最佳的lambda值
fit_lasso_cv <- cv.glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 1) best_lambda_lasso <- fit_lasso_cv$lambda.min best_lambda_lasso

• LASSO—建模与验证

# 根据最佳lambda构建岭回归模型
fit_lasso <- glmnet(x = as.matrix(train[,-17]), y = train[,17], alpha = 1)
coeff_lasso <- predict(fit_lasso, s = best_lambda_lasso, type = 'coefficients')
coeff_lasso

# 模型评估
pred_lasso <- predict(fit_lasso, s = best_lambda_lasso, newx = as.matrix(test[,-17]))
RMSE <- sqrt(mean((test$Salary-pred_lasso)**2))
RMSE

结语

---

      OK，今天关于岭回归和LASSO回归的实战部分就介绍到这里，希望对数据挖掘或机器学习的朋友有所帮助，同时，也希望读者能够静下心来好好的复现一遍。如果你有任何问题，欢迎在公众号的留言区域表达你的疑问。同时，也欢迎各位朋友继续转发与分享文中的内容，让更多的朋友学习和进步。

• 文内代码及数据

链接: https://pan.baidu.com/s/1qYOGuc8 密码: xb3r

Python爱好者社区历史文章大合集：

福利：文末扫码立刻关注公众号，“Python爱好者社区”，开始学习Python课程：

关注后在公众号内回复“课程”即可获取：

0.小编的Python入门视频课程！！！

1.崔老师爬虫实战案例免费学习视频。

2.丘老师数据科学入门指导免费学习视频。

3.陈老师数据分析报告制作免费学习视频。

4.玩转大数据分析！Spark2.X+Python 精华实战课程免费学习视频。

5.丘老师Python网络爬虫实战免费学习视频。