机器学习实战5(回归篇)

目录

1、回归

2、回归代码

3、预测鲍鱼的年龄

---

1、回归

        前面的文章介绍了很多分类算法，分类的目标变量是标称型数据，而本文将会对连续型的数据做出预测。主要讲解简单的线性回归和局部加权线性回归，并通过预测鲍鱼年龄的实例进行实战演练。

        说到回归，一般都是指线性回归（linear regression），所以本文里的回归和线性回归代表同一个意思。线性回归意味着可以将输入项分别乘以一些常量，再将结果加起来得到输出。需要说明的是，存在另一种成为非线性回归的回归模型，该模型不认同上面的做法，比如认为输出可能是输入的乘积。

        应该怎么从一大堆数据里求出回归方程呢？假定输入数据存放在矩阵X中，结果存放在向量y中：

         而回归系数存放在向量w中：

        应当怎样从一大堆数据里求出回归方程呢?假定输人数据存放在矩阵x中，而回归系数存放在向量w中。那么对于给定的数据x，预测结果将会通过给出。现在的问题是，手里有一些x和对应的y，怎样才能找到w呢？一个常用的方法就是找出使误差最小的w。这里的误差是指预测y值和真实y值之间的差值，使用该误差的简单累加将使得正差值和负差值相互抵消，所以我们采用平方误差。

         用矩阵表示还可以写做：

         为啥能这么变化，记住一个前提：若x为向量，则默认x为列向量，x^T为行向量。将上述提到的数据矩阵X和标签向量y带进去，就知道为何这么变化了。在继续推导之前，我们要先明确一个目的：找到w，使平方误差和最小。因为我们认为平方误差和越小，说明线性回归拟合效果越好。现在，我们用矩阵表示的平方误差和对w进行求导：

         得到：

         令上述公式等于0，得到：

         w上方的小标记表示，这是当前可以估计出的w的最优解。从现有数据上估计出的w可能并不是数据中的真实w值，所以这里使用了一个"帽"符号来表示它仅是w的一个最佳估计。

        值得注意的是，上述公式中包含逆矩阵，也就是说，这个方程只在逆矩阵存在的时候使用，也即是这个矩阵是一个方阵，并且其行列式不为0。

        述的最佳w求解是统计学中的常见问题，除了矩阵方法外还有很多其他方法可以解决。通过调用NumPy库里的矩阵方法，我们可以仅使用几行代码就完成所需功能。该方法也称作OLS， 意思是“普通小二乘法”（ordinary least squares）。

2、回归代码

        数据集：数据在这！

        提取码：3liq

        第一列都为1.0，即x0。第二列为x1，即x轴数据。第三列为x2，即y轴数据。

标准回归函数和数据导入函数：

        第一个函数loadDataset ()加载数据。该函数打开一个用tab键分隔的文本文件，这里仍然默认文件每行的最后一个值是目标值。第二个函数standRegres()用来计算最佳拟合直线。该函数首先读人x和y并将它们保存到矩阵中;然后计算，然后判断它的行列式是否为零，如果行列式为零，那么计算逆矩阵的时候将出现错误。NumPy提供一个线性代数的库linalg，其中包含很多有用的函数。可以直接调用linalg.det()来计算行列式。最后，如果行列式非零，计算并返回w。如果没有检查行列式是否为零就试图计算矩阵的逆，将会出现错误。

#导入数据函数
def loadDataSet(filename):
    numFeat = len(open(filename).readline().split('\t')) - 1
    dataMat = []
    labelMat = []
    fr = open(filename)
    for line in fr.readlines():
        lineArr = []
        curLine = line.strip().split('\t')
        for i in range(numFeat):
            lineArr.append(float(curLine[i]))
        dataMat.append(lineArr)
        labelMat.append(float(curLine[-1]))
    return dataMat,labelMat

#标准回归函数
def standRegress(xArr,yArr):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    xTx = xMat.T*xMat
    if linalg.det(xTx) ==0.0:#矩阵为奇异矩阵,不能求逆
        print("This matrix is singular,cannot do inverse")
        return
    ws = xTx.I*(xMat.T*yMat)
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'blue',alpha = .5)
    xCopy = xMat.copy()
    xCopy.sort(0)
    yHat = xCopy*ws #预测值
    ax.plot(xCopy[:,1],yHat)
    plt.show()
    return ws

        绘画出来的数据如图：

         输出的ws的值：

 局部加权线性回归：

        线性回归的一个问题是有可能出现欠拟合现象，因为它求的是具有最小均方误差的无偏估计。显而易见，如果模型欠拟合将不能取得最好的预测效果。所以有些方法允许在估计中引入一些偏差，从而降低预测的均方误差。

        其中的一个方法是局部加权线性回归(Locally Weighted Linear Regression，LWLR)。在该算法中，我们给待预测点附近的每个点赋予一定的权重;在这个子集上基于最小均方差来进行普通的回归。与kNN一样,这种算法每次预测均需要事先选取出对应的数据子集。该算法解出回归系数w的形式如下:

         其中W是一个矩阵，这个公式跟我们上面推导的公式的区别就在于W，它用来给每个店赋予权重。

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

         这样我们就可以根据上述公式，编写局部加权线性回归，我们通过改变k的值，可以调节回归效果，编写代码如下：

#局部加权线性回归
def lwlr(testpoint,xArr,yArr,k):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    m = shape(xMat)[0]
    weights = np.mat(eye((m))) #创建权重对角矩阵
    for j in range(m): #遍历数据集计算每个样本的权重
        diffMat = testpoint - xMat[j,:]
        weights[j,j] = exp(diffMat*diffMat.T/(-2*k**2))
    xTx = xMat.T*(weights*xMat)
    if linalg.det(xTx) ==0.0:
        print("this matrix is singular,cannot do inverser")
        return
    ws = xTx.I*(xMat.T*(weights*yMat))
    return testpoint*ws

def lwlrTest(testArr,xArr,yArr,k):
    m = shape(testArr)[0]
    yHat = zeros(m)
    for i in range(m):
        yHat[i] = lwlr(testArr[i],xArr,yArr,k)  #为每个数据集调用lwlr。这有助于求解k 的大小

    #画图
    xMat = np.mat(xArr)
    strInd = xMat[:,1].argsort(0)
    xSort = xMat[strInd][:,0,:]
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(xSort[:,1],yHat[strInd])
    ax.scatter(xMat[:,1].flatten().A[0], np.mat(yArr).flatten().A[0], s = 20, c = 'red',alpha = .5)
    plt.show()
    return yHat

xArr, yArr = loadDataSet('ex0.txt')
yHat = lwlrTest(xArr,xArr,yArr,1.0) #记得修改k值

        我们可以不断修改k的值：

        当k=1.0时：

         当k=0.01时：

 当k=0.003时：

         可以看到，当k越小，拟合效果越好。但是当k过小，会出现过拟合的情况，例如k等于0.003的时候。

3、预测鲍鱼的年龄

        接下来，我们将回归用于真实数据。在data目录下有一份来自UCI数据集合的数据，记录了鲍鱼（一种介壳类水生动物）的年龄。鲍鱼年龄可以从鲍鱼壳的层数推算得到。

        可以看到，数据集是多维的，所以我们很难画出它的分布情况。每个维度数据的代表的含义没有给出，不过没有关系，我们只要知道最后一列的数据是y值就可以了，最后一列代表的是鲍鱼的真实年龄，前面几列的数据是一些鲍鱼的特征，例如鲍鱼壳的层数等。我们不做数据清理，直接用上所有特征，测试下我们的局部加权回归。

        新建abalone.py文件，添加rssError函数，用于评价最后回归结果。

def rssError(yArr,yHatArr):
    return ((yArr-yHatArr)**2).sum()

from numpy import *
import numpy as np
import matplotlib.pyplot as plt
def loadDataSet(filename):
    numFeat = len(open(filename).readline().split('\t')) - 1
    dataMat = []
    labelMat = []
    fr = open(filename)
    for line in fr.readlines():
        lineArr = []
        curLine = line.strip().split('\t')
        for i in range(numFeat):
            lineArr.append(float(curLine[i]))
        dataMat.append(lineArr)
        labelMat.append(float(curLine[-1]))
    return dataMat,labelMat

def rssError(yArr,yHatArr):
    return ((yArr-yHatArr)**2).sum()

#局部加权线性回归
def lwlr(testpoint,xArr,yArr,k):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    m = shape(xMat)[0]
    weights = np.mat(eye((m))) #创建权重对角矩阵
    for j in range(m): #遍历数据集计算每个样本的权重
        diffMat = testpoint - xMat[j,:]
        weights[j,j] = exp(diffMat*diffMat.T/(-2*k**2))
    xTx = xMat.T*(weights*xMat)
    if linalg.det(xTx) ==0.0:
        print("this matrix is singular,cannot do inverser")
        return
    ws = xTx.I*(xMat.T*(weights*yMat))
    return testpoint*ws

def lwlrTest(testArr,xArr,yArr,k):
    m = shape(testArr)[0]
    yHat = zeros(m)
    for i in range(m):
        yHat[i] = lwlr(testArr[i],xArr,yArr,k)  #为每个数据集调用lwlr。这有助于求解k 的大小

    #画图
    xMat = np.mat(xArr)
    strInd = xMat[:,1].argsort(0)
    xSort = xMat[strInd][:,0,:]
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(xSort[:,1],yHat[strInd])
    ax.scatter(xMat[:,1].flatten().A[0], np.mat(yArr).flatten().A[0], s = 20, c = 'red',alpha = .5)
    plt.show()
    return yHat

def standRegress(xArr,yArr):
    xMat = np.mat(xArr)
    yMat = np.mat(yArr).T
    xTx = xMat.T*xMat
    if linalg.det(xTx) ==0.0:#矩阵为奇异矩阵,不能求逆
        print("This matrix is singular,cannot do inverse")
        return
    ws = xTx.I*(xMat.T*yMat)
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'blue',alpha = .5)
    xCopy = xMat.copy()
    xCopy.sort(0)
    yHat = xCopy*ws #预测值
    ax.plot(xCopy[:,1],yHat)
    plt.show()
    return ws


abX, abY = loadDataSet('abalone.txt')
print('训练集与测试集相同:局部加权线性回归,核k的大小对预测的影响:')
yHat01 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 0.1)
yHat1 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 1)
yHat10 = lwlrTest(abX[0:99], abX[0:99], abY[0:99], 10)
print('k=0.1时,误差大小为:',rssError(abY[0:99], yHat01.T))
print('k=1  时,误差大小为:',rssError(abY[0:99], yHat1.T))
print('k=10 时,误差大小为:',rssError(abY[0:99], yHat10.T))
print('')
print('训练集与测试集不同:局部加权线性回归,核k的大小是越小越好吗？更换数据集,测试结果如下:')
yHat01 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 0.1)
yHat1 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 1)
yHat10 = lwlrTest(abX[100:199], abX[0:99], abY[0:99], 10)
print('k=0.1时,误差大小为:',rssError(abY[100:199], yHat01.T))
print('k=1  时,误差大小为:',rssError(abY[100:199], yHat1.T))
print('k=10 时,误差大小为:',rssError(abY[100:199], yHat10.T))
print('')
print('训练集与测试集不同:简单的线性归回与k=1时的局部加权线性回归对比:')
print('k=1时,误差大小为:', rssError(abY[100:199], yHat1.T))
ws = standRegress(abX[0:99], abY[0:99])
yHat = np.mat(abX[100:199]) * ws
print('简单的线性回归误差大小:', rssError(abY[100:199], yHat.T.A))

        输出为：

 

     可以看到，当k=0.1时，训练集误差小，但是应用于新的数据集之后，误差反而变大了。这就是经常说道的过拟合现象。我们训练的模型，我们要保证测试集准确率高，这样训练出的模型才可以应用于新的数据，也就是要加强模型的普适性。可以看到，当k=1时，局部加权线性回归和简单的线性回归得到的效果差不多。这也表明一点，必须在未知数据上比较效果才能选取到最佳模型。那么最佳的核大小是10吗？或许是，但如果想得到更好的效果，应该用10个不同的样本集做10次测试来比较结果。

        本示例展示了如何使用局部加权线性回归来构建模型，可以得到比普通线性回归更好的效果。局部加权线性回归的问题在于，每次必须在整个数据集上运行。也就是说为了做出预测，必须保存所有的训练数据。