常见线性回归算法
ordinary?least?squares	普通最小二乘法
generalized?least?squares	广义最小二乘法
other...	其他...
gradient?decent	梯度下降法

3.1 梯度下降法

有关梯度下降法的定义我暂时就不做说明，在下面的一些案例分析中，我相信大家对于梯度下降法就能够有一个清晰的了解。

3.2 模型分析

其实在线性回归中，根据我们前面所说，最重要的就是要找到一个y函数，也就是找到一个模型，使得该模型能够最大程度地符合各个散点的分布。那么这个什么模型的选取有很多，我们应该如何选取呢？或者说计算机如何判定该模型是否是一个好模型呢？

我们来看如下几个图

?在该图中我们选取了y=1这样一个函数，或者说这样一个模型。如果说一个房子不管面积多大，都是一口价一万元的话，显然是不可能的，如果有，请给我来10套。当然这是从生活角度来理解的，从图上我们也可以很轻松的看到，我们所选择的这个模型与各个散点的距离都非常的远。因此这个模型并不是一个好模型。

而在这个图中，我们选择了y=200这个模型。很明显与上面y=1模型相比，这个模型要好很多。最起码完全符合了一个散点的要求，并且两百万的价格与其他各个散点的差距就不是太大了。

?如果说起初的两个图，我们都可以用经验法或者看图法来判断模型的好坏，那么上面这两个图我们还能够准确判断出哪个模型更好吗？（是不是不行了？没事，适当承认自己的不行，也不是那么丢人）

因此我们必须引入数学的方法来计算出或者说定义出一个模型的好坏。而此时所涉及到的数学方法，或者说数学函数，被称作损失函数（lost function）

3.3 损失函数（lost function）

什么叫做损失函数？其实就是计算一个模型与各个点之间的差距到底有多大！！如果差距大，说明我们所选择的模型比较差，如果差距很小，说明这个模型不错。

那么损失函数到底该如何计算？其实就是计算最小均方差，而最小均方差的公式，如下所示

?其中的n指的是有多少个散点。猜测值指的就是我们所构建的模型猜测该点的数值。

有些人可能还是不能非常理解，我们对上面的几个模型来求一下它的最小均方差，我相信大家就能够理解了。

首先是y=1这个模型

?根据最小均方差的公式，我们不难得到，y=1这个模型的最小均方差为60132。

单看这一个数，我们没有任何意义，我们并不是特别清楚，这个模型是好还是坏。毕竟好坏是比较出来的。所以说我们再来看看其他的模型。

我们再来看一下y=200这个模型

在上图中我不仅计算出了 y=200模型下的最小均方差，我还将其与y=1模型进行了对比，很明显，y=200这个模型要更好。

?所以说到了这个时候，大家对损失函数，也就是最小均方差应该有一定的了解了。最小均方差越接近0，说明这个模型越好。但是由于数据样本的不同与实际情况的复杂性，这个最小均方差也不可能为0。而我们的目标就是尽量找到最小均方差接近0的模型。

?4.?开发算法，让程序计算出m和b(y=mx+b)

4.1 简化模型

?在明白了如何判断模型好坏之后，我们的任务自然就是将判断的过程转换为算法，并用程序来实现。

当然我们不能一口吃成胖子，我们先从简单的来。我们先来简化模型，先假设m=0

那么为了更好的理解，方便我们设计程序，我们还是先从Excel入手，来看一下这个算法应该如何设计。

?我们用Excel可以看到，最小均方差随着我们b的增加发生了变化，先变小，然后又随之增大。而从最小均方差的散点图，我相信大家也可以初步了解什么叫做梯度下降了。

4.2 代码实现没那么容易

那么有些人就会膨胀了，诶呀，这不是很简单吗，用python轻轻松松几行代码就解决了。来，看我发挥。?

mse_list = []
for m in range(0,1000):
	for b in range(0,1000):
		计算mse
		把m，b和 mse 放入mse_list[]
查找最小的mse对应的m和b

但是事实真的如此吗？这样的代码就能实现了？小伙子，你想的也太简单了。

我们来分析一下为什么不行。首先我们来看b的值，现在是用for循环，b的变换是每次变换1，此时只要我们for循环的次数足够的多，我们会发现，mse的值会出现很小的情况。但是mse到底多小才算好呢？如果你将b每次变换1变成每次变换0.5，你会发现，得到的mse变得更小了，也就是说这样子得到的模型更好了。同理，只要b的变换越小，所得到的模型会更好，那么我们如何选择我们的b得变换速度？？（当然，我们说得专业一点，b的变换速度也就是学习速率 Learning rate）同时很容易理解，b的变换速度越快，我for循环所需要的次数越少，b的变换速度越慢，所需要的次数越多。因此，我们如何选取for循环的次数？？

有些人会很天真的说，没事啊，我只要把我的b变换取得足够小，for循环的次数取得足够多，我不就可以得到非常不错的模型了吗？？对的，你很棒，说得很对。但就是不知道你需要等多长时间，你的算法才能够跑完。

懂我的意思了吗？learning rate 我们无法知道选择多少合适，循环的次数我们也不知道选择多少合适。因此这样粗糙的代码是无法满足我们机器学习的需要的。

4.3 引入新的方法，mse对b进行求导

所以我们必须想办法引入新的方法来解决这个问题。这个新的方法，就叫做求导（微分、斜率）

废话不多说我们来看一下图

你觉得图中三个点，哪个地方的mse好？答案显然易见，最好的是橙色的那个点。那么我们是如何判断的呢？其实就是我们高中或者说高数里面学过的求微分，也就是求导。很显然，红色点与绿色点的斜率，也就是该点的导数要比橙色点大很多的。当斜率越大时，说明该点离我们所需的最小mse越远，当斜率越小时，越接近我们所需要的最小mse。

所以说，我们应该计算出最小均方差的导数

求导的方法相信大家都知道，因为此时是在我们简化的模型下进行研究的，所以此时我们只需要求mse对b的导数即可。如下图所示

?为了避免有些同学忘记了如何求导，我把计算的过程放到下图。

?4.3.1 用Excel验证mse对b的求导

我们可以通过表格轻易的看到，与我们之前预估的一样，mse对b的导数的绝对值越小（越接近0），那么最小均方差mse就越小。?

除此以外，我们还可以看出，当mse对b的导数为负时，说明此时的b还小，我们所需的最优最小均方差的点还在它的右边，也就是说b还需要继续增大。反之，mse对b的导数为正时，说明此时的b大了，需要调小一点。.

所以说，至此，我们可以发现，b的选取，与mse对b导数的的正负有关，也与该导数绝对值得大小有关。因此，我们对b的操作，不应该再是简单的每次加一或者加十，而应该采用算式，使其与mse对b的导数有关。而根据我们前面所分析的内容，我们不难得到b的公式应该为。

b=上一次的b-mse对b的导数*learning rate

learning rate就是学习速率，这个概念咱们之前提过，我们再从数学的角度来考虑。比如说上图中的第一个导数，-488，是一个负数很大的值，它值越大，说明离最优mse点越远，那么它乘以learning rate 也就越大，b就快点变化，使其尽快到达最优mse点附近。而导数值越小，说明离最优点近了，此时导数乘以learning rate也就自然小了。正好符合要求。

在此再附上一张Excel验证的图形

再总结一下这个简化模型下梯度下降算法的流程图，如下图所示

4.4 真实模型?

4.4.1 偏导数分别求解m和b的导数

在4.1中我提到了，为了方便理解，我们将模型进行简化，使得m=0，只观察b。当经过了前面的讲解之后，我们自然要回到真实的模型，对m和b一起研究。那么原先只有b的时候只需要对b进行求导即可。当有m和b两个变量以后，我们就需要分别求解mse对m和b的偏导。我相信这个地方大家应该很好理解，就是我们高数里面学的多元函数求偏导的过程。

求偏导的结果我放到了下图

?如果有同学对如何求偏导有些遗忘，我把求解的过程放到了下图，大家可以进行参考

4.4.3 用Excel验证m和b下的真实模型?

废话不多说，我只放一张图，不做详细说明?

4.5 python代码实现梯度下降算法

为了提高代码的可读性与逻辑性，我们需要分模块去写不同的函数。本次代码所需要的函数如下图所示。

? ? ? ?4?

当然，为了方便理解，我们还是先选取假数据（也就是刚刚Excel里面的房屋面积与房价的数据）来进行演示与验证。后面再使用真实的数据。

废话不多说，上代码

import numpy as np

data = np.array([
    [80, 200],
    [95, 230],
    [104, 245],
    [112, 274],
    [125, 259],
    [135, 262]
])

m=1
b=1
# 将data数据的两列进行拆分
xarray = data[:, 0]
yreal = data[:, -1]
learningrate = 0.00001

# 定义梯度下降程序
def grandentdecent():
    b_slop = 0
    for index, x in enumerate(xarray):
        b_slop = b_slop + m * x + b - yreal[index]
    b_slop = b_slop * 2 / len(xarray)
#    print("mse对b求导={}".format(b_slop))

    m_slop = 0
    for index, x in enumerate(xarray):
        m_slop = m_slop + (m * x + b - yreal[index])*x
    m_slop = m_slop * 2 / len(xarray)
#    print("mse对m求导={}".format(m_slop))

    return (b_slop, m_slop)

# 定义训练函数
def train():
    for i in range(1, 10000000):
        b_slop, m_slop = grandentdecent()
        global m
        m= m - m_slop*learningrate
        global b
        b= b - b_slop*learningrate
        if(abs(m_slop)<0.5 and abs(b_slop)<0.5):
            break

    print("m={},b={}".format(m, b))

if __name__ == '__main__':
    train()

代码跑完的结果放到了下图，可以看到与Excel就算的结果基本一致。但是由于Python的精度比Excel精度高一些，所以二者的结果略有不同。

?好了，假数据的代码，我们已经跑完了，现在已经基本理解了线性回归代码的流程，接下来我们便要真正的跑一个真实数据的线性回归代码，也就是梯度下降算法。那么首先，我们再来看一下整体的算法流程加深记忆。

这次真实代码的实战是：梯度下降，线性回归，探索汽车马力和油耗之间的关系。

油耗=m*（汽车马力）+b

求解m和b

废话不多说，继续上代码

import numpy as np

data = np.loadtxt("cars.csv", delimiter=",", skiprows=1, usecols=(4, 1))

m = 1
b = 1
# 将data数据的两列进行拆分
xarray = data[:, 0]
yreal = data[:, -1]
learningrate = 0.00001

# 定义梯度下降程序
def grandentdecent():
    b_slop = 0
    for index, x in enumerate(xarray):
        b_slop = b_slop + m * x + b - yreal[index]
    b_slop = b_slop * 2 / len(xarray)
    print("mse对b求导={}".format(b_slop))

    m_slop = 0
    for index, x in enumerate(xarray):
        m_slop = m_slop + (m * x + b - yreal[index])*x
    m_slop = m_slop * 2 / len(xarray)
    print("mse对m求导={}".format(m_slop))

    return (b_slop, m_slop)

# 定义训练函数
def train():
    for i in range(1, 10000000):
        b_slop, m_slop = grandentdecent()
        global m
        m= m - m_slop*learningrate
        global b
        b= b - b_slop*learningrate
        if(abs(m_slop)<0.5 and abs(b_slop)<0.5):
            break

    print("m={},b={}".format(m, b))

if __name__ == '__main__':
    train()

?代码的运行结果为

?当然，这样子的结果，可能大家看不出来这个模型到底好不好。所以，要想验证这个模型的好坏，我们还需要将数据data拆分成训练集和测试集，然后验证准确率，这样才能清楚地知道模型的好坏。

?拆分训练集和测试集进行验证的过程我就不再演示，大家可以自己动手实验。

补充：KNN之多维度数据标准化

在之前学习的机器学习（一）中，KNN我们是用的二维KNN，在标准化的时候就比较容易，直接编写就可以。但是如果KNN的feature纬度比较多的话，用普通的方法就会特别复杂。所以接下来补充一种基于矩阵转置的多维度数据标准化处理方法。代码如下。

# 预测房价，经度，纬度，房屋使用面积，房价

import numpy as np

data = np.array([
    [47.5112, -122.257, 1180, 221900],
    [48.5112, -123.257, 1280, 231900],
    [46.5112, -124.257, 1380, 241900],
    [45.5112, -125.257, 1480, 251900],
    [43.5112, -126.257, 1580, 261900],
])

for i in range(0, len(data.T)-1):
    feature = (data.T[i])
    mean = (np.mean(feature))
    std = np.std(feature)
    data.T[i] = ((feature-mean) / std)

print(data)

最近几天的学习记录就如上所示，接下来可能会不定期更新，下面涉及的内容会和很多的高等数学以及线性代数中的矩阵相关。用高等数学的方法使得我们机器学习的代码更加简便以及更有拓展性。

本文中所涉及到的代码以及数据资料如有需要可以私信我获取。

本人水平有限，所写内容仅为个人学习记录，欢迎同行共同交流讨论。