python 最小二乘回归 高斯核_机器学习实战个人学习分享（七）：线性回归

点击蓝字 关注我们

这次的分享是线性回归的一些内容，涉及到四种不同的线性回归内容，本次采用每个理论加代码的形式编辑，另外《机器学习实战》书上有一些公式是直接给出来的，但是在这里有必要给出推导过程，方便理解，推导过程结合西华书，因为之前阅读西瓜书的过程中，每个有疑问的公式大都推导了一遍，所以这里直接附上并做简单说明，照例，这次内容结合Jack Cui的博客的内容，由于黄海安的视频就更到AdaBoost算法，所以以后的内容就主要根据博客和自己跑代码的理解，博客内容得到了Jack Cui的授权，链接如下：

Python3《机器学习实战》学习笔记(十一)：线性回归基础篇之预测鲍鱼年龄

https://blog.csdn.net/c406495762/article/details/78760239

Python3《机器学习实战》学习笔记(十二)：线性回归提高篇之乐高玩具套件二手价预测

https://blog.csdn.net/c406495762/article/details/82967529

01

普通线性回归

常规的线性回归用来做数值型的回归预测，简单的理解就是直线拟合，根据已有的数据点，以最小化平方误差来求出合适的直线方程参数，然后来预测未知的数据，下面来推导一下理论部分。

首先数据都应该是数值型的，特别是标签值需要连续，设线性回归的方程如下：

那么我们的目的是试图找到一个合适的线性方程来使得方程在每个样本点的值和实际值很接近(完全吻合最好)，这就引出了我们研究的问题，即试图寻找是均方误差(这里选择均方误差是因为怕直接使用误差存在正误差和负误差相互抵消的问题，)最小化的参数，这里选择的是均方误差是很好的几何意义，它对应了欧氏距离，以均方误差作为性能度量最小化的模型称之为“最小二乘法”，公式如下：

那么要求解的问题出来后，就是求解了，这个过程称之为最小二乘的“参数估计”，这里有两种形式，一种是代数形式，一种是向量形式，首先介绍代数形式。

1代数形式

对于代数问题，求导是常见的办法，求导如下：

然后代入：

这里西瓜书是直接给的结果，在南瓜书里有推导，我自己也推导了一下，下面把结果大致附上分析：

这里主要的转化就是上面两行长的部分，不难。

2向量形式

如果是向量形式的话，将数据集合D写成如下形式：

这是为了将参数b向量化后方便运算，向量化后，求解的问题就变成了如下：

这里依旧是求导但是涉及到矩阵求导的一些概念，用得到的两个公式如下：

那么求导可以得到：

可以看出这里可以令等式为0，求出w，但是求逆矩阵有个前提是X转置乘以X必须是满秩的，这就需要讨论。

当满秩情况下，可以直接求出w，然后代入直线方程，可以得到：

当不满秩情况下(属性多于样本数)，这里有一个伪逆的概念，之前看过知乎和博客，在这里简单说一下，参考麻省理工线性代数笔记(二十九)- 左右逆和伪逆，链接如下：

https://zhuanlan.zhihu.com/p/38346210

伪逆是在奇异值分解的基础上得到的

学过数值计算的课程，所以了解UV分解，这里也是利用这个，一个矩阵可以写成奇异值分解的形式，那么它的逆可以写成如上形式，其中如果S对角线的值小于tol，那么被视为0，然后三部分都求逆按照乘法倒序排列即可。可以看出，逆和伪逆的乘积并不是严格的单位阵，这就是因为上面所说的置零的情况。MATLAB和numpy都有函数求伪逆，可以方便使用。

那么线性回归大致就清晰了，求出直线方程的参数，这里如果引出对数变换的话，就有一个对数几率函数，它从阶跃函数而来并且弥补了阶跃函数不连续的问题，函数表达式如下：

这就可以写出如下：

这称之为对数几率回归，一般的话需要借助概率的知识，进行极大似然估计，最大化对数似然，上式可以写成：

然后问题转化成如下：

这里可以采用梯度下降法和牛顿法进行优化求解，牛顿法迭代更新如下：

以上就是线性回归的理论部分，下面开始上代码实战：3代码实战 这个实战的数据如下： 总共有200个数据，三列，这里的后两列就是数据的横纵坐标，第一列是1，这个可以理解为x0，因为后面求出来w，会与数据相乘，所以这个x0为1，就保证了b的系数为1，即为截距。下面简单处理数据然后求出w，主函数代码如下：

if __name__ == '__main__':    xArr,yArr = loadDataSet('ex1.txt')    ws = standRegres(xArr, yArr)    xMat=mat(xArr)    yMat=mat(yArr)    yHat=xMat*ws    fig = plt.figure()    ax = fig.add_subplot(111)  # 添加subplot    ax.plot(xMat[:, 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()    print(np.corrcoef(yHat.T, yMat))

装载数据就比较简单了，还是按行读取，转化成矩阵，主要看第二个函数如下：

def standRegres(xArr, yArr):    xMat = mat(xArr)    yMat = 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)    return ws

这里函数比加简单，判断是否行列式等于0，满足矩阵可逆条件的话，就直接计算返回w，然后画出数据图和回归直线，结果如下：

看图形可能无法准确的判断模型如何，所以需要用相关系数来衡量，numpy中有函数corrcoef函数可以直接调用，相关系数代码调用如下：

print(np.corrcoef(yHat.T, yMat))

运行打印如下：

相关系数为0.9852，还比较不错(保持疑问，以后有机会介绍一下几个衡量回归型模型性能的指标，就觉得这个0.98就不是特别好了)。

02

局部加权线性回归

对于线性回归，容易出现欠拟合现象，所以引入局部加权线性回归(LWLR)的概念，给予每个点赋予一定的权重，这里使用核函数，通过调节k来确定拟合的程度，高斯核函数在这里常被使用，公式如下：

加了权重的回归系数为：

下面简单分析一下代码：

if __name__ == '__main__':    xArr, yArr = loadDataSet('ex1.txt')    # 单点预测    print(lwlr(xArr[0],xArr,yArr,1.0))    # 全部数据点测试    print(lwlrTest(xArr,xArr,yArr,0.001))    plotlwlrRegression()

主要是三个，单点预测，全部预测以及绘图。单点预测简单分析一下，全部数据点就是反复调用单点预测，代码如下：

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):  # next 2 lines create weights matrix        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("This matrix is singular, cannot do inverse")        return    ws = xTx.I * (xMat.T * (weights * yMat))    return testPoint * ws

这个就是初始化权重向量为单位阵，然后带入高斯核函数，每个点都需要和全部数据进行计算得到一个权重，最后返回ws。测试第一个点

1.000000    0.635975   4.093119

预测的结果是：[[4.06916648]]结果还不错，下面是全部数据的结果：

可以看出，选取不同的K值，拟合的程度是不同的，比较三个拟合数据与原数据的相关系数如下：

由相关系数也可以看出，K越小，拟合越好，但是也容易过拟合。

03

预测鲍鱼年龄

运用局部加权回归，预测一下鲍鱼的年龄，代码如下：

if __name__ == '__main__':    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 = standRegres(abX[0:99], abY[0:99])    yHat = np.mat(abX[100:199]) * ws    print('简单的线性回归误差大小:', rssError(abY[100:199], yHat.T.A))

这里除了数据不同，调用的函数都是上面提到的，运行代码，结果如下：

结果不必分析，一句话，所有算法离开了具体问题都是空谈。

04

岭回归

刚才的两个线性回归都是满足满秩的，但是有时候也会遇到非奇异矩阵的情况，例如当特征多于样本点时候，这个就引入了“岭回归”的概念，属于一种缩减方法。 岭回归即我们所说的 L2 正则线性回归，在一般的线性回归最小化均方误差的基础上增加了一个参数 w 的 L2 范数的罚项，从而最小化罚项残差平方和：

那么权重更新为如下：

这个地方的变化就是多了一个系数乘以单位阵，使得矩阵非奇异，从残差平方和的式子来看，引入了lambda使得限制l了w，在统计学中称之为缩减。对角元素为1，类似于零，即是“岭回归”的来源。下面上代码实践：

if __name__ == '__main__':    xArr, yArr = loadDataSet('abalone.txt')    print(ridgeTest(xArr,yArr))

这是主函数，下面是核心函数：

def ridgeTest(xArr, yArr):    xMat = mat(xArr)    yMat = mat(yArr).T    yMean = mean(yMat, 0)    yMat = yMat - yMean  # to eliminate X0 take mean off of Y    # regularize X's    xMeans = mean(xMat, 0)  # calc mean then subtract it off    xVar = var(xMat, 0)  # calc variance of Xi then divide by it    xMat = (xMat - xMeans) / xVar    numTestPts = 30    wMat = zeros((numTestPts, shape(xMat)[1]))    for i in range(numTestPts):        ws = ridgeRegres(xMat, yMat, exp(i - 10))        wMat[i, :] = ws.T    return wMat

这里在数据处理的时候需要注意，为了使用岭回归，这里对数据进行了标准化处理，使得每个特征的重要性均等。这里调用了一个函数去计算回归系数，代码如下：

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("This matrix is singular, cannot do inverse")        return    ws = denom.I * (xMat.T * yMat)    return ws

代码比较简单，简单的计算和判断，为了防止lambda设置为0，所以还需要再次判别矩阵是否非奇异。设置lamnda为不同值，运行结果如下：

05

前向逐步回归

前向逐步回归属于一种贪心算法，每一步都尽可能的减小误差，回归系数的计算不再是通过公式，而是每次不断调整，然后计算预测误差，结合代码分析：

if __name__ == '__main__':    xArr, yArr = loadDataSet('abalone.txt')    print(stageWise(xArr,yArr,0.01,200))

这是主函数，下面看核心函数：

def stageWise(xArr, yArr, eps=0.01, numIt=100):    xMat = mat(xArr)    yMat = mat(yArr).T    yMean = mean(yMat, 0)    yMat = yMat - yMean  # can also regularize ys but will get smaller coef    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 # 一次加个sign，一次减个sign，达到L1范数的效果                yTest = xMat * wsTest                rssE = rssError(yMat.A, yTest.A)                if rssE < lowestError:                    lowestError = rssE                    wsMax = wsTest        ws = wsMax.copy()

这里简单分析，数据标准化后，对每个特征遍历，每次加上一个sign或者减去一个sign，达到L1范数的效果，计算新的误差，如果误差小于当前误差就更新误差，直至找到最好的。运行代码，绘图如下：

def plotstageWiseMat():    """    函数说明:绘制岭回归系数矩阵    Website:        http://www.cuijiahua.com/    Modify:        2017-11-20    """    font = FontProperties(fname=r"c:\windows\fonts\simsun.ttc", size=14)    xArr, yArr = loadDataSet('abalone.txt')    returnMat = stageWise(xArr, yArr, 0.005, 1000)    fig = plt.figure()    ax = fig.add_subplot(111)    ax.plot(returnMat)    ax_title_text = ax.set_title(u'前向逐步回归:迭代次数与回归系数的关系', FontProperties=font)    ax_xlabel_text = ax.set_xlabel(u'迭代次数', FontProperties=font)    ax_ylabel_text = ax.set_ylabel(u'回归系数', FontProperties=font)    plt.setp(ax_title_text, size=15, weight='bold', color='red')    plt.setp(ax_xlabel_text, size=10, weight='bold', color='black')    plt.setp(ax_ylabel_text, size=10, weight='bold', color='black')    plt.show()

结果如下：

这就是线性回归的主要内容，由于篇幅，还有个预测乐高玩具的例子就不写了，这里有涉及选取最优特征的，有空可以跑跑代码，sklearn的代码也有空统一更新在sklearn章节吧。

不

断

前

进

好文！必须在看