机器学习实战之路 —— 1 线性回归二
1. 欠拟合和过拟合问题
线性回归的一个问题是有可能出现欠拟合现象,因为它求的是最小均方误差的无偏估计。显而易见,如果模型欠拟合将不能取得最好的预测效果。所以有些方法允许在估计中引入一些偏差,从而降低预测的均方误差。
——《机器学习实战》
P
141
P_{141}
P141
处理欠拟合的方法:
*增加特征,增大模型复杂度
*增大调试参数C
处理过拟合的方法:
*减少特征,降低模型复杂度
*减小调试参数C
*增加训练数量
——《机器学习之路》
P
49
P_{49}
P49
一般而言,欠拟合就是模型能力不足,没有有效的表示出数据集内在的表达式。与欠拟合相对的就是模型的过拟合问题,在实际中,当特征维度增加到一定程度时,虽然模型在训练集的成绩越来越好,但在测试集的最终成绩反而越来越差。例举示意图如下:
总之模型越简单,模型适应力越强,越容易发生欠拟合;反之,模型越复杂,模型适应性越差,越容易发生过拟合。针对以上问题,我们可以对线性模型进行相应的优化处理。
2. 模型优化
2.1 局部加权线性回归(LWLR)
在局部加权线性回归(LocallyWeightedLinearRegression, LWLR )算法中,给待预测点附近的每个点赋予一定的权重,在这个子集上基于最小二乘法来进行普通的回归,LWLR试图求得:
arg min
ω
∑
i
=
1
m
ω
i
(
y
i
−
f
(
x
i
)
)
2
\begin{aligned} \argmin_{ω}\sum_{i=1}^{m}ω_i(y_i-f(x_i))^2 \end{aligned}
ωargmini=1∑mωi(yi−f(xi))2
此处的
ω
i
ω_i
ωi 是非负的权重(weights)值。直观看来,如果对于某个特定的 i ,它的
ω
i
ω_i
ωi 很大,那么在选择
ω
i
ω_i
ωi的时候,我们将会尽可能地使
(
y
i
−
f
(
x
i
)
)
2
(y_i-f(x_i))^2
(yi−f(xi))2 更小。如果
ω
i
ω_i
ωi 很小,那么在拟合的过程中
(
y
i
−
f
(
x
i
)
)
2
(y_i-f(x_i))^2
(yi−f(xi))2误差项就能够大大地忽略。
参考上一篇博文—线性回归(一),LWLR 回归系数推导过程和(一)类似,最终求得的回归系数ω的形式如下:
ω
^
=
(
X
T
W
X
)
−
1
X
T
W
y
\hat\omega=(X^TWX)^{-1}X^TWy
ω^=(XTWX)−1XTWy
W是一个权重矩阵,LWLR使用 “核”(与支持向量机中的核类似)来对附近的点赋予更高的权重。核的类型可以自由选择,最常用的核就是高斯核,高斯核对应的权重如下:
W
(
i
,
i
)
=
e
x
p
(
∣
x
(
i
)
−
x
∣
−
2
k
2
)
W(i,i)=exp\left({\frac{|x^{(i)}-x|}{-2k^2}}\right)
W(i,i)=exp(−2k2∣x(i)−x∣)
代码实现如下:
def lwlr(testPoint,xArr,yArr,k=1.0):
xMat = mat(xArr); yMat = mat(yArr).T
m = shape(xMat)[0]
# 创建对角矩阵
weights = mat(eye((m)))
# 创建权重矩阵,权重大小以指数级衰减
for j in range(m):
diffMat = testPoint - xMat[j,:] # xMat[j,:]为某一行的值
weights[j,j] = exp(diffMat*diffMat.T/(-2.0*k**2))
xTx = xMat.T * (weights * xMat)
if linalg.det(xTx) == 0.0:
print("矩阵为奇异矩阵,不能求逆")
return
# 计算回归系数
ws = xTx.I * (xMat.T * (weights * yMat))
return testPoint * ws
2.2 岭回归
如果特征比样本点还多(n > m ) , 也就是说输入数据的矩阵X不是满秩矩阵,非满秩矩阵在求逆时会出现问题。同时为了防止过拟合,在标准线型回归函数中增加λ扰动。
岭回归的目标函数在一般的线性回归的基础上加入了二范数正则项,在保证最佳拟合误差的同时,使得参数尽可能的“简单”,使得模型的泛化能力强。岭回归试图求得:
arg min
ω
∑
i
=
1
m
(
y
i
−
f
(
x
i
)
)
2
+
λ
∣
∣
ω
∣
∣
2
\begin{aligned} \argmin_{ω}\sum_{i=1}^{m}(y_i-f(x_i))^2 + λ||ω||^2 \end{aligned}
ωargmini=1∑m(yi−f(xi))2+λ∣∣ω∣∣2
回归系数的计算公式为:
ω
^
=
(
X
T
X
+
λ
I
)
−
1
X
T
y
\hat\omega=(X^TX+\lambda I)^{-1}X^Ty
ω^=(XTX+λI)−1XTy
代码实现如下:
def ridgeRegres(xMat,yMat,lam=0.2):
xTx = xMat.T*xMat
denom = xTx + eye(shape(xMat)[1])*lam
if linalg.det(denom) == 0.0:
print("矩阵为奇异矩阵,不能求逆")
return
# 计算回归系数
ws = denom.I * (xMat.T*yMat)
return ws
使用sklearn.linear_model中岭回归代码如下:
from sklearn.linear_model import Ridge
# sklearn的岭回归
model = Ridge()
# 0.001 ~ 100 的等比数列的20个数
alpha_can = np.logspace(-3, 2, 20)
# 超参数自动搜索模块GridSearchCV,系统地遍历多种参数组合,通过交叉验证确定最佳效果参数 5折交叉验证
lasso_model = GridSearchCV(model, param_grid={'alpha': alpha_can}, cv=5)
lasso_model.fit(x_train, y_train)
2.3 Lasso回归
区别于岭回归的二范数正则项,Lasso回归采用一范数来约束,使参数非零个数最少。Lasso回归试图求得:
arg min
ω
∑
i
=
1
m
(
y
i
−
f
(
x
i
)
)
2
+
λ
∣
ω
∣
\begin{aligned} \argmin_{ω}\sum_{i=1}^{m}(y_i-f(x_i))^2 + λ|ω| \end{aligned}
ωargmini=1∑m(yi−f(xi))2+λ∣ω∣
与岭回归相比其唯一的不同点在于,这个约束条件使用绝对值取代了平方和。虽然约束形式只是稍作变化,结果却大相径庭:在
λ
\lambda
λ足够小的时候,一些系数会因此被迫缩减到0,这个特性可以帮助我们更好地理解数据。这两个约束条件在公式上看起来相差无几,但细微的变化却极大地增加了计算复杂度(为了在这个新的约束条件下解出回归系数,需要使用二次规划算法)。
Lasso回归使用sklearn.linear_model的核心代码如下:
from sklearn.linear_model import Lasso
# sklearn的Lasso回归
model = Lasso()
# 0.001 ~ 100 的等比数列的20个数
alpha_can = np.logspace(-3, 2, 20)
# 超参数自动搜索模块GridSearchCV,系统地遍历多种参数组合,通过交叉验证确定最佳效果参数 5折交叉验证
lasso_model = GridSearchCV(model, param_grid={'alpha': alpha_can}, cv=5)
lasso_model.fit(x_train, y_train)
2.4 前向逐步回归
前向逐步回归算法可以得到与Lasso差不多的效果,但更加简单。它属于一种贪心算法,即每一步都尽可能减少误差。一开始,所有的权重都设为1,然后每一步所做的决策是对某个权重增加或减少一个很小的值。即步长每个特征。
代码实现如下:
def regularize(xMat):
'''
数据标准化
:param xMat: 数据集
:return: 标准化后的数据
'''
inMat = xMat.copy()
# 求均值
inMeans = mean(inMat,0)
# 求方差
inVar = var(inMat,0)
# 数据减均值除以方差实现标准化
inMat = (inMat - inMeans)/inVar
return inMat
def stageWise(xArr,yArr,eps=0.01,numIt=100):
'''
:param xArr: 输入数据
:param yArr: 预测变量
:param eps: 每次迭代步长
:param numIt: 迭代次数
:return: 回归系数
'''
# 数据标准化
xMat = mat(xArr); yMat=mat(yArr).T
yMean = mean(yMat,0)
yMat = yMat - yMean
xMat = regularize(xMat)
m,n=shape(xMat)
# 初始化回归系数矩阵
returnMat = zeros((numIt,n)) #testing code remove
ws = zeros((n,1)); wsTest = ws.copy(); wsMax = ws.copy()
for i in range(numIt):
print(ws.T)
lowestError = inf;
# 遍历每个特征的回归系数
for j in range(n):
for sign in [-1, 1]:
wsTest = ws.copy()
# 回归系数步进
wsTest[j] += eps * sign
yTest = xMat * wsTest
rssE = rssError(yMat.A, yTest.A)
# 若误差更小,选择最佳回归系数
if rssE < lowestError:
lowestError = rssE
wsMax = wsTest
ws = wsMax.copy()
returnMat[i, :] = ws.T
return returnMat
2.5 时间序列分析ARIMA
时间序列
时间序列(或称动态数列)是指将同一统计指标的数值按其发生的时间先后顺序排列而成的数列。时间序列分析的主要目的是根据已有的历史数据对未来进行预测。
时间序列分析
时间序列分析是根据系统观察得到的时间序列数据,通过曲线拟合和参数估计来建立数学模型的理论和方法。时间序列分析常用于国民宏观经济控制、市场潜力预测、气象预测、农作物害虫灾害预报等各个方面。
使用statsmodels.tsa.arima_model代码如下:
from statsmodels.tsa.arima_model import ARIMA
p = 2; q = 2
model = ARIMA(endog=x, order=(p, d, q)) # 自回归函数p,差分d,移动平均数q
arima = model.fit(disp=-1) # disp<0:不输出过程
prediction = arima.fittedvalues
y = prediction.cumsum() + x[0]
3. 实战:广告花费与销售额
问题描述:
已知某公司分别在电视(TV)、广播(Radio)、报纸(Newspaper)上所投入的广告花费,及对应所产生的销售额(Sales),预测在不同广告投入下所能产生的销售额。
相应的数据集如下:
NO. | TV | Radio | Newspaper | Sales |
---|---|---|---|---|
1 | 230.1 | 37.8 | 69.2 | 22.1 |
2 | 44.5 | 39.3 | 45.1 | 10.4 |
3 | 17.2 | 45.9 | 69.3 | 9.3 |
4 | 151.5 | 41.3 | 58.5 | 18.5 |
5 | 180.8 | 10.8 | 58.4 | 12.9 |
6 | 8.7 | 48.9 | 75 | 7.2 |
7 | 57.5 | 32.8 | 23.5 | 11.8 |
8 | 120.2 | 19.6 | 11.6 | 13.2 |
9 | 8.6 | 2.1 | 1 | 4.8 |
10 | 199.8 | 2.6 | 21.2 | 10.6 |
…… | …… | …… | …… | …… |
200 | 232.1 | 8.6 | 8.7 | 13.4 |
绘制该数据的示意图如下:
各数据分别显示:
后续将综合线性回归(一)和(二)的一些关键算法实现做具体展开。
3.1 线性回归
实现代码如下:
#!/usr/bin/python
# -*- coding:utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression
import sklearn.metrics as metrics
import matplotlib as mpl
if __name__ == "__main__":
# pandas读入数据(TV、Radio、Newspaper、Sales)
data = pd.read_csv('Advertising.csv')
x = data[['TV', 'Radio', 'Newspaper']]
y = data['Sales']
print("数据维度:", data.shape)
print("x数据集:\n", x)
print("y数据集:\n", y)
# 拆分数据为训练集(80%)和测试集(20%)
x_train, x_test, y_train, y_test = train_test_split(x, y, test_size=0.2, random_state=1)
print("训练数据维度:", x_train.shape, "测试数据维度:", x_test.shape)
# sklearn线性回归
linreg = LinearRegression()
# 用训练集来训练出线型回归模型
model = linreg.fit(x_train, y_train)
print(model)
print("回归系数:", linreg.coef_, " 截距:", linreg.intercept_)
# 返回测试集从小到大的索引值,沿第一轴排序(向下)
order = y_test.argsort(axis=0)
y_test = y_test.values[order]
x_test = x_test.values[order, :]
# 用测试集进行预测
y_hat = linreg.predict(x_test)
# Mean Squared Error(均方误差)
mse = np.average((y_hat - np.array(y_test)) ** 2)
# Root Mean Squared Error(均方根误差)
rmse = np.sqrt(mse)
# 用sklearn计算 Mean Squared Error(均方误差)
mse = metrics.mean_squared_error(y_test, y_hat)
# Root Mean Squared Error(均方根误差)
rmse = np.sqrt(mse)
print("======== 均方误差,均方根误差 ========")
print("MSE = ", mse)
print("RMSE = ", rmse)
# 对模型进行评价
print("======== 模型评价 R^2 ========")
print("Training R^2 = ", linreg.score(x_train, y_train))
print("Test R^2 = ", linreg.score(x_test, y_test))
# 线型回归预测图
mpl.rcParams['font.sans-serif'] = ['simHei']
mpl.rcParams['axes.unicode_minus'] = False
plt.figure(facecolor='w')
t = np.arange(len(x_test))
plt.plot(t, y_test, 'r-', linewidth=2, label='真实数据')
plt.plot(t, y_hat, 'g-', linewidth=2, label='预测数据')
plt.legend(loc='upper left')
plt.title('线性回归预测销量', fontsize=18)
plt.grid(b=True, ls=':')
plt.show()
输出结果如下:
数据维度: (200, 5)
x数据集:
TV Radio Newspaper
0 230.1 37.8 69.2
1 44.5 39.3 45.1
2 17.2 45.9 69.3
3 151.5 41.3 58.5
4 180.8 10.8 58.4
5 8.7 48.9 75.0
.. ... ... ...
199 232.1 8.6 8.7
[200 rows x 3 columns]
y数据集:
0 22.1
1 10.4
2 9.3
3 18.5
4 12.9
5 7.2
...
199 13.4
Name: Sales, Length: 200, dtype: float64
训练数据维度: (160, 3) 测试数据维度: (40, 3)
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=1, normalize=False)
回归系数: [ 0.0468431 0.17854434 0.00258619] 截距: 2.90794702082
======== 均方误差,均方根误差 ========
MSE = 1.99188555183
RMSE = 1.41134175586
======== 模型评价 R^2 ========
Training R^2 = 0.895937263233
Test R^2 = 0.892760591462
绘制的线性回归预测销量如下:
拟合的线型回归函数为:
y
=
2.9
+
0.046
∗
T
V
+
0.179
∗
R
a
d
i
o
+
0.0026
∗
N
e
w
s
p
a
p
e
r
y = 2.9 + 0.046 * TV + 0.179 * Radio + 0.0026 * Newspaper
y=2.9+0.046∗TV+0.179∗Radio+0.0026∗Newspaper
从拟合出的回归系数可知,Newspaper对销售额的影响较其他两种小一到两个数量级,不妨尝试忽略Newspaper,只选取TV、Radio两个重要特征来进行线型回归,得到结果如下:
回归系数: [ 0.04686997 0.1800065 ] 截距: 2.94751503603
======== 均方误差,均方根误差 ========
MSE = 1.95522188501
RMSE = 1.39829248908
======== 模型评价 R^2 ========
Training R^2 = 0.895852846878
Test R^2 = 0.894734495003
从结果可知,只考虑重要特征(TV、Radio),得到的结果可能会优于考虑全部特征,即特征不见得越多越好。在机器学习实践中,在深入研究应用业务的基础上,对特征的分析和恰当选取是非常重要的。
3.2 岭回归(Ridge)
主要实现代码如下:
if __name__ == "__main__":
# pandas读入数据(TV、Radio、Newspaper、Sales)
# 拆分数据为训练集(80%)和测试集(20%)
# sklearn的岭回归
model = Ridge()
# 0.001 ~ 100 的等比数列的20个数
alpha_can = np.logspace(-3, 2, 20)
np.set_printoptions(suppress=True)
print('alpha_can = ', alpha_can)
# 超参数自动搜索模块GridSearchCV,系统地遍历多种参数组合,通过交叉验证确定最佳效果参数 5折交叉验证
lasso_model = GridSearchCV(model, param_grid={'alpha': alpha_can}, cv=5)
lasso_model.fit(x_train, y_train)
print('超参数:\n', lasso_model.best_params_)
# 结果输出
输出结果如下:
超参数:
{'alpha': 16.237767391887211}
======== 均方误差,均方根误差 ========
MSE = 1.99368993021
RMSE = 1.41198085334
======== 模型评价 R^2 ========
Training R^2 = 0.895937204701
Test R^2 = 0.89266344709
3.3 Lasso
主要实现代码如下:
if __name__ == "__main__":
# pandas读入数据(TV、Radio、Newspaper、Sales)
# 拆分数据为训练集(80%)和测试集(20%)
# sklearn的Lasso回归
model = Lasso()
# 同岭回归(Ridge)
# 结果输出
输出结果如下:
超参数:
{'alpha': 0.78475997035146072}
======== 均方误差,均方根误差 ========
MSE = 2.03270682226
RMSE = 1.42573027683
======== 模型评价 R^2 ========
Training R^2 = 0.895833980993
Test R^2 = 0.890562850285
4. 参考学习的书目及论文
- 机器学习算法视频 - 邹博
- 《线性回归分析导论》机械工业出版社
- 《线性代数 - 吴传生》第3章 矩阵的运算
- 《机器学习实战》第8章 预测数值型数据:回归
- 《机器学习之路》第2章 机器学习进阶:调试模型
- 《机器学习 - 周志华》第3章 线性模型
- 《基于机器学习的纯电动汽车的行驶里程预测研究》 高航 2018硕士论文 第4章 行驶里程的多元线性回归建模
=文档信息=
本学习笔记由博主整理编辑,仅供非商用学习交流使用
由于水平有限,错误和纰漏之处在所难免,欢迎大家交流指正
如本文涉及侵权,请随时留言博主,必妥善处置
版权声明:非商用自由转载-保持署名-注明出处
署名(BY) :zhudj
文章出处:https://zhudj.blog.csdn.net/