线性回归算法

发表于 | 分类于

线性回归算法

注：本文所有代码和数据均在个人github下click me

回归算法一般是针对预测是连续的情况下，对于预测值是离散的，采用的算法是分类算法。线性回归算法包括很多种变形，这里提到的线性回归算法是其中的几种典型算法。在实际应用中，我们采用线性算法可以预测商品的价格，预测房子的房价等等，虽然线性回归算法比较简单，但是在实际中还是有很多的使用的。

在机器学习中，我们要紧盯三件事情。第一，算法的损失函数；第二，采用什么求值算法求损失函数的最小值；第三，算法的评价指标

一般线性回归算法

一般的线性回归又叫做最小二乘算法，最小二乘是因为算法的损失函数是最小二乘的，损失函数如下：

其中$h_\theta(x^{(i)})$是算法针对数据的预测值，而$y^{(i)}$是数据的真实值，$m$表示训练数据的条数，而${\frac 12}$是为了此公式求导方便而加入的。而$\theta$是算法的参数，在这里就是线性回归的权重值。通过此公式我们可以得到，线性回归算法的损失函数就是针对每个样本计算预测值和真实值得差，然后将差求平方，之后将全体样本的差平方相加即得到损失函数。
针对损失函数，我们有两种算法可以求取损失函数的最小值。

梯度下降算法

梯度下降算法的一般形式：

这里写的是梯度下降算法的标量形式。这个公式描述了梯度下降算法是如何更新算法的参数的，其中$\alpha$是参数的更新步长。可以看到这里的关键是如何求取损失函数关于参数j的偏导数。
将损失函数带入到梯度下降算法,即公式$(1.2)$中，并且求导，可以得到下式：

我们重复的使用公式$(1.3)$直到达到收敛条件，即可以求得线性回归算法中的参数值。
从公式中，可以看到用真实值和预测值之间的差值来更新参数值，这种方式或者思想在很多的机器学习算法中可以看到，包括深度学习的后向传播算法。同时，可以看到每一次的迭代都要使用整个数据集。这种方式叫做批量梯度下降算法。还有一种方式可以求取$\theta$值，叫做随机梯度下降算法，算法如下：

从中我们可以看到，随机梯度下降算法每次采用一个训练样本取更新所有的参数值，注意那里的(for every j)。当训练样本很多时，相比于批量梯度下降算法，随机梯度下降算法能够更快的更新算法的参数值，并且能够更快的逼近损失函数的最小值。

代数解

我们将损失函数用向量表示，如下所示：

公式中的$X$表示训练数据的矩阵。因为损失函数是凸二次函数，所以只有一个最小值，所以导数为0的点就是损失函数的最小值。
具体的推导过程如下(主要利用了矩阵的迹运算)：

令导数等于0，从而得到：

回归模型的概率解释

大家想过没有，为什么在线性回归模型里面选择最小二乘作为损失函数？接下来从概率的角度来解释选择最小二乘作为损失函数的原因。
首先，假设目标变量和输入数据存入如下的关系：

这里的$\epsilon^{(i)}$是误差项，包括模型未考虑到影响目标变量的因素和随机噪声。
接下来假设，误差项相互独立同分布，并且服从高斯分布(即正态分布)

为什么要假设误差项服从高斯分布? 第一是因为采用高斯分布，在数学上处理会比较简单；第二是因为根据中心极限定理，独立的随机变量的和，其总的影响接近高斯分布

误差项的概率密度函数为：

根据公式$(1.6)$和公式$(1.7)$,我们可以得出如下结论：

在公式中，$\theta$不是随机变量，而是实际存在的值，虽然我们不知道真实值是多少。$p(y^{(i)}|x^{(i)};\theta)$的含义是给定$x^{(i)}$,参数设定为$\theta$时,的概率密度。注意公式中用的分号。
数据的概率是由$p(\vec{y} | X;\theta)$给出的，而总的概率可以看成是在固定$\theta$时，关于$\vec{y}$的函数。换个角度，我们想要将这个函数明确的看成是关于$\theta$的函数，所以我们将其称作似然函数，从而我们得到关于$\theta$的似然函数：

似然(likelihood)和概率(probability)实际上是一个东西，但是似然函数是对参数$\theta$定义的，为了加以区分，使用了似然这一术语。我们可以说参数的似然，数据的概率，但不能说数据的似然，参数的概率。

极大似然估计就是选择$\theta$,使参数的似然函数最大化，也就是选择参数使得已有样本的出现概率最大。
因为$L(\theta)$是严格单调递增的，并且对数函数也是递增的，所以取对数，得到的$\theta$跟不取对象是一样的。

对数似然函数$\mathcal {l}(\theta)$:

最大化似然函数，就是最大化对数似然函数，即最小化

这个式子刚好是线性回归算法中采用的损失函数。总结一下，最小二乘回归模型刚好就是在假设了误差独立同服从正态分布后，得到的最大似然估计。同时注意到，正态分布中的方差$\sigma^2$的取值对模型并没有影响。

编程实现

代数解实现线性回归算法

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
def standard_regression(x_arr, y_arr):
    """
    这里直接采用代数解直接求取权重大小.这里有个条件是xTx的逆必须存在,
    如果不存在,则直接返回
    :param x_arr:
    :param y_arr:
    :return:
    """
    x_mat = np.mat(x_arr)
    y_mat = np.mat(y_arr).T
    xTx = x_mat.T * x_mat
    if np.linalg.det(xTx) == 0.0:
        print("矩阵是奇异的,逆不存在")
        return
    ws = xTx.I * (x_mat.T * y_mat)
    return ws

根据standard_regression函数可以直接求取算法的权重。然后根据权重就可以求得数据的预测结果，这里选择的数据在点我
画出的图像如下：

直线为拟合直线，散点是真实值。

随机梯度下降实现线性回归

由于代码篇幅比较长，所以可以直接上我的github上面看。代码

局部加权线性回归算法

每个点都对应一个不同的$\theta$，利用$\theta$计算预测值。
下面我们用刚才介绍的局部加权线性回归来拟合一下这个模型，简单回顾一下过程：
1.用高斯核函数计算出第ｉ个样本处，其它所有样本点的权重Ｗ
2.用权重ｗ对第ｉ个样本作加权线性回归，得到回归方程，即拟合的直线方程
3.用刚才得到的经验回归直线计算出$x^i$处的估计值$y^i$
4.重复一至三步，得到每个样本点的估计值

同时,从计算公式中，可以看到对于每个点的预测值，需要利用所有的数据样本进行计算。如果数据量比较大，计算量就会是一个大问题。
相比于上面提到线性回归算法(参数算法)，局部加权线性回归算法是非参数算法。

参数算法vs非参数算法

引用地址

参数算法：
如果我们对所要学习的问题有足够的认识，具备一定的先验知识，此时我们一般会假定要学习的目标函数f(x)或分布P(y|x)的具体形式。然后，通过训练数据集，基于ERM、SRM、MLE、MAP等学习策略，可以估计出f(x)或P(y|x)中含有的未知参数。一旦未知参数估计完毕，训练数据一般来说，就失去其作用了，因为这些估计出来的参数就是训练数据的浓缩。通过这种方式建立起来的模型就是参数模型。参数模型的一个很重要的特点是，如果对于模型的假设正确，那么只需要很少的训练数据就可以从假设空间中学出一个很好的模型。但是，如果模型的假设错误，那么无论训练的数据量有多大，甚至趋于无穷大，学出的模型都会与实际模型出现不可磨灭的偏差。感知机、逻辑斯特回归、高斯判别分析、朴素贝叶斯、线性支持向量机都属于参数模型。对于神经网络来说，当固定了隐层的数目以及每一层神经元的个数，它也属于参数模型。但由于隐层数目与每一层神经元个数的不确定性，很多时候，神经网络都被归类为半参数模型。
非参数算法：
当我们对所要学习的问题知之甚少，此时我们一般不会对潜在的模型做过多的假设。在面对预测任务的时候，我们通常会用上所有的训练数据。例如简单的核密度估计(KDE)的表达式中，就带有所有训练数据的信息。通过这种方式建立的模型就是非参数模型。非参数模型的一个很重要的特点就是：let the data speak for itself. 正因为如此，非参数模型的存储开销、计算开销都会比参数模型大的多。但是，由于不存在模型的错误假定问题，可以证明，当训练数据量趋于无穷大的时候，非参数模型可以逼近任意复杂的真实模型。这正是非参数模型诱人的一点。另外需要说明的一点是，非参数模型之所以叫做非参数，并不是因为模型中没有参数。实际上，非参数模型中一般会含有一个或多个超参数，外加无穷多个普通的参数。k近邻就是典型的非参数模型。

总结就是所谓参数学习算法它有固定的明确的参数，参数 一旦确定，就不会改变了，我们不需要在保留训练集中的训练样本。而非参数学习算法，每进行一次预测，就需要重新学习一组 ， 是变化的，所以需要一直保留训练样本。

实际使用

这里采用的数据跟线性回归算法采用的相同。选取不同的k值，得到的拟合曲线如下图：

从图中可以看到，当k值较小时，对于训练数据有很好的拟合效果；当k值较大时(比如取1)，对于训练数据拟合的效果不是很好。当然这都是针对训练数据，在实际使用中，我们更关注的是对于测试数据的效果如何。在实际中，我们可以将数据进行十倍交叉验证，选择最合适的k值。
从高斯公式中，我们从两个方面看。首先，我们固定k值，那么距离$x^{(i)}$较远的样本，其对应的权重相对较小；距离$x^{(i)}$较近的样本，其对应的权重相对较大，当x取$x^{(i)}$时，对应的权重为1。而如果我们固定x(即只看某一特定样本点)，k取值较大时，此样本点对应的权重相对较大；k取值较小时，此样本点对应的权重相对较大。所以k可以控制算法选择$x^{(i)}$点附近有多少样本参与计算。

最上面的图是样本点，剩下三幅图都是针对样本点x=0.5,根据不同的k值，画出的x附近样本的权重变化。可以看到，当k取0.5时，当计算x=0.5时的预测值时，几乎所有的样本点都包括了；而当k=0.01时，仅仅取了x=0.5附近的点参与计算。如果k值取无穷大，那么w对于所有的点的权重都是1，那么局部加权线性回归就变成了普通的线性回归算法。

岭回归算法(Ridge Regression)

上面介绍的线性回归算法里面可以看到需要计算$X^TX$的逆，那么如果逆不存在呢？首先思考第一个问题，什么情况下，$X^TX$的逆不存在呢？

  1. X本身存在相关关系，即非满秩矩阵。比如其中两列是具有线性关系
  2. 如果特征列多余样本数量，那么$X^TX$也是非满秩的。

对于逆不存在的情况下，我们需要将原来的线性回归算法进行处理，使原先无法求逆的$X^TX$变成非奇异的。可以通过缩减的方式来处理这类问题，比如岭回归算法和LASSO算法。同时由于能够调整算法中权重的大小，能够防止线性回归算法的过拟合问题。

岭回归算法损失函数

通过改变$\lambda$的值，可以使得算法的方差和偏差之间达到平衡，增加$\lambda$,模型的方差减小而偏差增加。
对损失函数求取一阶导数，并另导数等于0，求得权重如下式：

在岭回归算法中$\lambda$的选择对于算法有很大的影响。下图展示了不同的$\lambda$的取值对于权重的影响，因为数据有八个特征，所以这里有八个权重。当$\lambda$较小时，权重的值跟采用常规的线性回归差不多；而当$\lambda$较大时，权重的值都会被调解的较小。

为了找到合适的$\lambda$值，我们在实践中往往会采用交叉测试来找到合适的$\lambda$值。

lasso回归算法(Least Absolute Shrinkage and Selection Operator)

从岭回归算法中，我们可以看到，算法防止过拟合主要是在损失函数中添加惩罚项。在岭回归中，惩罚项如下所示:

而在lasso回归算法中，惩罚项变成下式：

即将权重的平方和小于$\lambda$，替换为权重的绝对值和小于$\lambda$.进行了这个变化后，能够将权重缩小到0，而岭回归中无法将权重值缩小到0，只能接近0.

lasso回归算法损失函数

我们也可以结合岭回归算法和lasso的损失函数，构建新的损失函数。这就是弹性网络(ElasticNet)

逐步向前回归

LASSO计算复杂度相对较高，本部分稍微介绍一下逐步向前回归，他属于一种贪心算法，给定初始系数向量，然后不断迭代遍历每个系数，增加或减小一个很小的数，看看代价函数是否变小，如果变小就保留，如果变大就舍弃，然后不断迭代直到回归系数达到稳定。代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
def run_stage():
    x_arr, y_arr = load_data_set('./data/abalone.txt')
    all_w_001 = stage_wise(x_arr, y_arr, 0.001, 5000)
    print(all_w_001)
    all_w_01 = stage_wise(x_arr, y_arr, 0.01, 200)
    print(all_w_01)
  
def stage_wise(x_arr, y_arr, eps=0.01, num_iter=100):
    """
    forward stagewise regression算法（前向梯度算法）
    是一种近似的 lasso算法
    :param x_arr:
    :param y_arr:
    :param eps:每次特征权重的变化步长
    :param num_iter: 迭代次数
    :return:
    """
    x_mat = np.mat(x_arr)
    y_mat = np.mat(y_arr).T
    y_mean = np.mean(y_mat, 0)
    y_mat = y_mat - y_mean
    x_mean = np.mean(x_mat, 0)
    x_var = np.var(x_mat, 0)
    x_mat = (x_mat - x_mean) / x_var
    m, n = np.shape(x_mat)
    ws = np.zeros((n, 1))
    ws_best = ws.copy()
    return_mat = np.zeros((num_iter, n))  # 保存每次迭代最好的权重值
    for i in range(num_iter):
        # print(ws.T)
        lowest_error = np.inf
        for j in range(n):
            for sign in [-1, 1]:
                ws_test = ws.copy()
                ws_test[j] += eps * sign
                y_test = x_mat * ws_test
                rss_err = rss_error(y_mat.A, y_test.A)  # 将矩阵转为数组
                if rss_err < lowest_error:
                    lowest_error = rss_err
                    ws_best = ws_test
        ws = ws_best.copy()
        return_mat[i, :] = ws.T
    return return_mat

逐步回归算法的主要有点在于他可以帮助人们理解现有的模型并作出改进。当构建了一个模型后，可以运行逐步回归算法找出重要的特征，即使停止那些不重要特征的收集。

总结

上面提到的这些线性回归算法，我们应该怎么选择呢？一般来说有一点正则项的表现更好，因此通常你应该避免使用简单的线性回归。岭回归是一个很好的首选项，但是如果你的特征仅有少数是真正有用的，你应该选择Lasso和弹性网络。就像我们讨论的那样，它两能够将无用特征的权重降为零。一般来说，弹性网络的表现要比Lasso好，因为当特征数量比样例的数量大的时候，或者特征之间有很强的相关性时，Lasso可能会表现的不规律。

参考文档

  1. 《吴恩达cs229 课程讲义》
  2. 《机器学习实战》
  3. 机器学习实战_线性回归&逻辑回归（二）