机器学习第八章笔记——预测数值型数据：回归_数据回归预测

目录

一、用线性回归找到最佳拟合直线

二、 局部加权线性回归 

三、预测鲍鱼的年龄 

四、缩减系数来“理解”数据 

4.1 岭回归

4.2 lasso 

4.3 向前逐步回归

一、用线性回归找到最佳拟合直线

线性回归的优缺点;

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

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

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

回归的目的是预测数值型的目标值。最直接的方法就是依据输入写出一个目标值的计算公式。

假设预测汽车的功率大小也能会这样计算：

HorsePower = 0.0015 * annualSalary - 0.99 * hoursListenIngToPublicRadio

这就是所谓回归方程，其中0.0015和-0.99称为回归系数，求这些回归系数的过程就是回归。有了这些回归系数，再给定输入，就能容易的给出预测了。具体做法是：用回归系数乘以输入值，再将结果全部加在一起，就可以得到预测值。线性回归可以将输入项分别乘以一些常量，再将结果加起来得到输出。

回归的一般方法：

        1、收集数据：采用任意方法收集数据

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

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

        4、训练算法：找到回归系数

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

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

        假设输入数据存放在矩阵x中，而回归系数存放在向量w中，对于给定的数据预测结果将会通过$Y_{1} =X_{1}^{T}w$给出。对于如何找出w，一个常用的方法是找出使误差最小的w，误差指的是预测y值和真实y值之间的差值，使用该误差的简单累加将使正差值与负差值相互抵消，因此采用平方误差。

        平方误差可以写成：

用矩阵表示还可以写作。对w求导，得到，令其等于零，解出w如下：

        w上方的小标记表示，这是当前可以估计出的w的最优解。上述公式中包含，因此需要对矩阵求逆。故需要在代码中对是否存在逆矩阵进行判断。

        最佳w求解是统计学中的常见问题，可以通过调用Numpy库中的矩阵方法。这种方法被称作OLS，即普通最小二乘法。

        对图1中的数据给出最佳拟合直线，代码为：

def loadDataSet(fileName):
 
    numFeat = len(open(fileName).readline().split('\t')) - 1  # 得到每行的列数
    xArr = []; yArr = []  # 为可能的数据和标签创建空列表
    fr = open(fileName)
    for line in fr.readlines():
        lineArr =[]
        curLine = line.strip().split('\t')
        for i in range(numFeat):
            lineArr.append(float(curLine[i]))
        xArr.append(lineArr)
        yArr.append(float(curLine[-1]))
    return xArr, yArr
 
 
def standRegres(xArr,yArr):
    xMat = mat(xArr); yMat = mat(yArr).T
    xTx = xMat.T * xMat  # 根据文中推导的公示计算回归系数
    if linalg.det(xTx) == 0.0:  # 判断行列式是否为0，若为0则不存在逆矩阵
        print("矩阵为奇异矩阵,不能求逆")
        return
    ws = xTx.I * (xMat.T*yMat)   # 求出w，.I是对矩阵求逆矩阵
    return ws

第一个函数打开文本文件并提取其中数据，第二个函数用于计算最佳拟合函数，计算方法就为上面讨论中提及的w的计算公式。 

接下来是绘制图像：

def plotRegression():
    xArr, yArr = loadDataSet('D:\\learning\\ex0.txt')                                    # 加载数据集
    ws = standRegres(xArr, yArr)                                        # 计算回归系数
    xMat = mat(xArr)                                                    # 创建xMat矩阵
    yMat = mat(yArr)                                                    # 创建yMat矩阵
    xCopy = xMat.copy()                                                    # 深拷贝xMat矩阵
    xCopy.sort(0)                                                        # 排序
    yHat = xCopy * ws                                                     # 计算对应的y值
    fig = plt.figure()
    ax = fig.add_subplot(111)                                            # 添加subplot
    ax.plot(xCopy[:, 1], yHat, c = 'red')                                # 绘制回归曲线
    ax.scatter(xMat[:,1].flatten().A[0], yMat.flatten().A[0], s = 20, c = 'blue',alpha = .5)                # 绘制样本点
    plt.title('DataSet')                                                # 绘制title
    plt.xlabel('X')
    plt.show()

调用函数，输出结果就能得到下图 

二、 局部加权线性回归 

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

        其中的一个方法是局部加权线性回归。算法中给待预测点附近的每个点赋予一定的权重，基于最小的均方差来进行普通的回归。同KNN一样，算法每次预测均需要事先选取出对应的数据子集，回归系数w的形式为：

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

        LWLR使用“核”来对附近的点赋予更高的权重，核的类型可以自由选择最常用的核就是高斯核，对应权重如下：

这里构建的是一个只含对角元素的权重矩阵w，点x与x(i)越近，w(i,i)将会越大。上述公式包含一个需要用户指定的参数k，它决定了对附近的点赋予多大的权重，也是使用LWLR时唯一需要考虑的参数。

局部加权线性回归函数:

def lwlr(testPoint, xArr, yArr, k = 1.0):
    xMat = mat(xArr); yMat = mat(yArr).T
    m = shape(xMat)[0]
    weights = mat(eye((m)))                   # 创建权重对角矩阵，200行
    for j in range(m):                                                  # 遍历数据集计算每个样本的权重
        diffMat = testPoint - 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
 
 
def lwlrTest(testArr, xArr, yArr, k=1.0):
    m = shape(testArr)[0]                                            # 计算测试数据集大小,200行
 
    yHat = zeros(m)
    for i in range(m):                                                    # 对每个样本点进行预测
        yHat[i] = lwlr(testArr[i], xArr, yArr, k)  # 第一个参数testArr[i]表示训练集中的每一行，每次拿一个样本进行预测
    return yHat

        给定x空间中的任意一点，计算出对应的预测值yHat。函数lwlr的开头用于读入数据并创建所需矩阵之后创建对角权重剧中weights。权重矩阵是一个方阵，阶数等于样本个数。即该矩阵为每个样本点初始化了一个权重，利用算法遍历数据集，计算每个样本点对应的权重值；随样本点与预测值点距离的递增，权重将以指数级衰减。k用于控制衰减的速度。计算完毕可以得到对回归系数w的一个估计lwlrtest则用于为每一个点调用lwlr()，方便求解k的大小。

xArr, yArr = regression.loadDataSet("D:\\learning\\ex0.txt")
result1 = regression.lwlr(xArr[0], xArr, yArr, 1.0)
result2 = regression.lwlr(xArr[0], xArr, yArr, 0.001)
print(result1)
print(result2)

得到的单点估计： 

若想得到所有点的估计，可以调用lwlrtest()函数。

绘制出估计值与原始值，以观察拟合效果：

xArr, yArr = regression.loadDataSet("D:\\learning\\ex0.txt")
result1 = regression.lwlr(xArr[0], xArr, yArr, 1.0)
result2 = regression.lwlr(xArr[0], xArr, yArr, 0.001)
yHat = regression.lwlrTest(xArr, xArr, yArr, 0.003)
 
xMat = mat(xArr)
srtInd = xMat[:, 1].argsort(0)
xSort = xMat[srtInd][:, 0, :]
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(xSort[:, 1], yHat[srtInd])
ax.scatter(xMat[:,1].flatten().A[0], mat(yArr).T.flatten().A[0], s=2, c='red')
plt.show()

最终结果为： 

当然局部加权线性回归也存在一个问题，即增加了计算量。 

三、预测鲍鱼的年龄 

        新添一个rssError()函数用于分析预测误差的大小：

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

输入：

abX,abY = regression.loadDataSet('D:\\learning\\abalone.txt')
yHat01 = regression.lwlrTest(abX[0:99], abX[0:99], abY[0:99], 0.1)
yHat1 = regression.lwlrTest(abX[0:99], abX[0:99], abY[0:99], 1)
yHat10 = regression.lwlrTest(abX[0:99], abX[0:99], abY[0:99], 10)
 
result1 = regression.rssError(abY[0:99],yHat01.T)
result2 = regression.rssError(abY[0:99], yHat1.T)
result3 = regression.rssError(abY[0:99], yHat10.T)
print(result1)
print(result2)
print(result3)

 运行后：

        从上面的数据可以看出：使用较小的核将得到较低的误差，但是在操作中并不能使用最小的核处理数据集，因为使用最小核处理将会造成过拟合，对新数据不一定能达到最好的预测效果。

新数据： 

abX,abY = regression.loadDataSet('D:\\learning\\abalone.txt')
yHat01 = regression.lwlrTest(abX[100:199], abX[0:99], abY[0:99], 0.1)
yHat1 = regression.lwlrTest(abX[100:199], abX[0:99], abY[0:99], 1)
yHat10 = regression.lwlrTest(abX[100:199], abX[0:99], abY[0:99], 10)
 
result1 = regression.rssError(abY[100:199],yHat01.T)
result2 = regression.rssError(abY[100:199], yHat1.T)
result3 = regression.rssError(abY[100:199], yHat10.T)
print(result1)
print(result2)
print(result3)

从结果上看，核大小等于10时的测试误差最小，不过其在训练集上的误差却是最大的 。

四、缩减系数来“理解”数据 

        当数据的特征比样本点还多时就不能在使用之前的方法进行预测了，否则计算时会出错。

        当特征比样本点还多时，输入数据的矩阵x就不是满秩矩阵了，非满秩矩阵在求逆矩阵时会出现错误。为解决这个问题引入岭回归的概念。

4.1 岭回归

        岭回归是在矩阵上加上一个从而使得矩阵非奇异，对求逆，x是一个m*m的单位矩阵，对角线上全是1，其他元素全是0。$\lambda$是随便定义的一个数值。这时回归系数的计算公式就为：

        岭回归通过引入$\lambda$来限制了所有w之和，通过引入该惩罚项就能够减少不重要的参数，统计学中将该技术称为缩减。

        岭回归使用了单位矩阵乘以常量$\lambda$观察其中单位矩阵I，可以看到1贯穿整个对角线。其余元素均为0.

        添加岭回归函数到.py文件中：

def ridgeRegres(xMat, yMat, lam=0.2):
    xTx = xMat.T * xMat
    denom = xTx + eye(shape(xMat)[1]) * lam
    if linalg.dot(denom) == 0.0:
        print("this matrix is singular, cannot do inverse")
        return
    ws = denom.I * (xMat.T * yMat)
    return ws
 
 
def ridgeTest(xArr, yArr):
    xMat = mat(xArr); yMat = mat(yArr).T
    yMean = mean(yMat, 0)
    yMat = yMat - yMean
    xMeans = mean(xMat, 0)
    xVar = var(xMat, 0)
    xMat = (xMat - xMeans) / xVar
    numTestPts = 30
    wMat = zeros((numTestPts, shape(xMat)[1]))
    for i in range(numTestPts):
        ws = ridgeTest(xMat, yMat, exp(i- 10))
        wMat[i, :] = ws.T
    return wMat

包含两个函数，ridgeRegres()函数以及ridgeTest()函数，前者用于计算回归系数，后者则用于在一组$\lambda$上测试结果。

        ridgeRegres()实现了给定lambda下的岭回归求解。当没指定lambda时，默认为0.2。函数首先构建了矩阵，用lam乘以单位矩阵。在普通回归方法可能产生错误时，岭回归仍然可以正常工作。为使用岭回归和缩减技术，需要对特征做标准化处理。具体做法就是将所有特征都减去各自的均值并除以方差。

        利用鲍鱼数据进行测试：

abX,abY = regression.loadDataSet('D:\\learning\\abalone.txt')
ridgeWeights = regression.ridgeTest(abX, abY)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(ridgeWeights)
plt.show()

得到数据： 

4.2 lasso 

        在增加如下约束时，普通的最小二乘法回归会得到与岭回归的一样的公式：

上式限定了所有回归系数的平方和不能大于$\lambda$。使用普通的最小二乘法回归在两个或更多的特征相关时，可能会得到一个很大的正系数和一个很大的负系数。与岭回归系数，另一个缩减方法lasso也对回归系数做了限定，对应的约束条件如下：

        两个公式的不同点在于这个约束条件使用绝对值取代了平方和。尽管约束形式只是改了一点，但结果却大相径庭。当$\lambda$足够小时，一些系数会因此被迫缩减到0。

4.3 向前逐步回归

        向前逐步回归属于一种贪心算法，每一步都尽可能的减少误差，最开始时所有的权重都设为1，之后每一步所做的决策是对某个权重增加或减少一个很小的值。

算法伪代码：

        数据标准化，使其分布满足0均值和单位方差

        在每轮迭代中：

                设置当前最小误差lowestError为正无穷

                对每个特征

                     增大或减小

                             改变一个系数得到一个新的w

                             计算新w下的误差

                             如果误差Error小于当前最小误差lowerstError；设置Wbest等于当前的W

                        将w设置为新的Wbest

实际代码为：

def stageWise(xArr, yArr, eps=0.01, numIt=100):
    # 标准化处理
    xMat = mat(xArr); yMat = mat(yArr).T
    yMean = mean(yMat,0)
    yMat = yMat - yMean
    xMat = regularize(xMat)  # 归一化xMat
    m,n = shape(xMat)
    returnMat = zeros((numIt,n))
    ws = zeros((n,1)); wsTest = ws.copy(); wsMax = ws.copy()  # 为了实现贪心算法建立两个ws的两份副本
    for i in range(numIt):
        # 每次迭代都打印w向量
        print(ws.T)
        loweastError = 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 < loweastError:
                    loweastError = rssE
                    wsMax = wsTest
        ws = wsMax.copy()
        returnMat[i,:] = ws.T
    return returnMat
 
# xMat的归一化矩阵
def regularize(xMat):
    inMat = xMat.copy()
    inMeans = mean(inMat, 0)
    inVar = var(inMat, 0)      
    inMat = (inMat - inMeans)/inVar
    return inMat

函数是一个逐步线性回归算法的实现，它与lasso做法相近但计算简单。xArr为输入数据，yArr为预测变量。还有两个参数，一个eps，表示每次迭代需要调整的步长，另一个numIt，表示迭代次数。

        数据转换并存入矩阵中，然后把特征按照均值为0方差为1，进行标准化处理。创建向量ws保存w值，并且为了实现贪心算法建立了ws的两份副本。

        贪心算法在所有特征上运行两次for循环，分别计算增加或减少该特征对误差的影响。这里使用的是平方误差，通过之前的函数rssError()得到。误差初值设为正无穷经过与所有的误差比较后取最小误差。

xArr, yArr = regression.loadDataSet('D:\\learning\\abalone.txt')
result = regression.stageWise(xArr, yArr, 0.01, 200)
print(result)

输出结果为： 

从结果中看到w1和w6都是0，标明它们不对目标造成任何影响，这些特征很可能就是不需要的。

将结果与最小二乘法进行比较，结果如下：