学习笔记【机器学习重点与实战】——1 线性回归

1.线性回归概述

线性模型形式简单、易于建模，但却蕴涵着机器学习中一些重要的基本思想。许多功能更为强大的非线性模型 (nonlinear model)可在线性模型的基础上通过引入层级结构或高维映射而得。此外，由于 ω 直观表达了各属性在预测中的重要性，因此线性模型有很好的可解释’性 (comprehensibility)。

——《机器学习-周志华》P53P53

优点：结果易于理解，计算上不复杂。

缺点：对非线性的数据拟合不好。

适用数据类型：数值型和标称型数据。

——《机器学习实战》P136P136

回归的一般方法：

收集数据：采用任意方法收集数据。

准备数据：回归需要数值型数据，标称型数据将被转成二值型数据。

分析数据：绘出数据的可视化二维图将有助于对数据做出理解和分析，在采用缩减法求得新回归系数之后, 可以将新拟合线绘在图上作为对比。

训练算法：找到回归系数。

测试算法：使用R2R2或者预测值和数据的拟合度，来分析模型的效果。

使用算法：使用回归，可以在给定输入的时候预测出一个数值，这是对分类方法的提升，因为这样可以预测连续型数据而不仅仅是离散的类别标签。

——《机器学习实战》P137P137

2.线性回归常用算法实现

2.1标准线型回归

利用最小二乘法求回归系数的最优解，即求预测y值和真实y值之间的平方误差最小时的回归系数ω：

∑i=1m(yi−xTiω)2∑i=1m(yi−xiTω)2

该矩阵还可表示为(y−Xω)T(y−Xω)(y−Xω)T(y−Xω)，对ω求导得到XT(Y−Xω)XT(Y−Xω)，令该式为零，可得标准线性回归函数：

ω^=(XTX)−1XTyω^=(XTX)−1XTy

代码实现如下：

def standRegres(xArr,yArr):
xMat = mat(xArr); yMat = mat(yArr).T
xTx = xMat.T*xMat
if linalg.det(xTx) == 0.0:
print("矩阵为奇异矩阵,不能求逆")
return
# 计算回归系数
ws = xTx.I * (xMat.T*yMat)
return ws

使用sklearn.linear_model中线性回归代码如下：

from sklearn.linear_model import LinearRegression

# sklearn线性回归
linreg = LinearRegression()
# 用训练集来训练出线型回归模型
model = linreg.fit(x_train, y_train)

2.2局部加权线性回归(LWLR)

在局部加权线性回归（LocallyWeightedLinearRegression, LWLR )算法中，给待预测点附近的每个点赋予一定的权重；在这个子集上基于最小二乘法来进行普通的回归。该算法解出回归系数ω的形式如下：

ω^=(XTWX)−1XTWyω^=(XTWX)−1XTWy

其中W是一个矩阵，用来给每个数据点赋予权重。

LWLR使用 “核”（与支持向量机中的核类似）来对附近的点赋予更高的权重。核的类型可以自由选择，最常用的核就是高斯核，高斯核对应的权重如下：

W(i,i)=exp(|x(i)−x|−2k2)W(i,i)=exp(|x(i)−x|−2k2)

代码实现如下：

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.3岭回归

如果特征比样本点还多（n > m ) , 也就是说输入数据的矩阵X不是满秩矩阵，非满秩矩阵在求逆时会出现问题。同时为了防止过拟合，在标准线型回归函数中增加λ扰动。回归系数的计算公式为：

ω^=(XTX+λI)−1XTyω^=(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.4Lasso回归

与岭回归类似，Lasso回归也对回归系数做了限定，对应的约束条件如下：

∑k=1n|ωk|≤λ∑k=1n|ωk|≤λ

唯一的不同点在于，这个约束条件使用绝对值取代了平方和。虽然约束形式只是稍作变化，结果却大相径庭：在λλ足够小的时候，一些系数会因此被迫缩减到0，这个特性可以帮助我们更好地理解数据。这两个约束条件在公式上看起来相差无几，但细微的变化却极大地增加了计算复杂度（为了在这个新的约束条件下解出回归系数，需要使用二次规划算法)。

使用sklearn.linear_model中Lasso回归代码如下：

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.5前向逐步回归

前向逐步回归算法可以得到与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.6时间序列分析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），预测在不同广告投入下所能产生的销售额。

数据集如下：

TV	Radio	Newspaper	Sales
1	230.1	37.8	69.2
2	44.5	39.3	45.1
3	17.2	45.9	69.3
4	151.5	41.3	58.5
5	180.8	10.8	58.4
6	8.7	48.9	75
7	57.5	32.8	23.5
8	120.2	19.6	11.6
9	8.6	2.1	1
10	199.8	2.6	21.2
……	……	……	……
200	232.1	8.6	8.7

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, "测试数据维度：", y_train.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) 测试数据维度： (160,)
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∗TV+0.179∗Radio+0.0026∗Newspapery=2.9+0.046∗TV+0.179∗Radio+0.0026∗Newspaper

从拟合出的回归系数可知，Newspaper对销售额的影响较其他两种小一到两个数量级，不妨尝试忽略Newspaper进行线型回归，得到结果如下：

回归系数： [ 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.3Lasso

主要实现代码如下：

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.参考

机器学习升级版视频 - 邹博

《机器学习实战》第8章 预测数值型数据：回归

《机器学习 - 周志华》第3章 线性模型

===========文档信息============

学习笔记由博主整理编辑，供非商用学习交流用

如本文涉及侵权，请随时留言博主，必妥善处置

版权声明：非商用自由转载-保持署名-注明出处

署名(BY) ：dkjkls（dkj卡洛斯）

文章出处：http://blog.csdn.net/dkjkls