《白话机器学习》（三）线性回归

一.什么是线性回归

在讲解线性回归之前，我们来看一组数据

姓名	房屋单价（元/平方米）	房屋面积（平方米）	房屋总价（元）
张三	6000	58	30000
李四	9000	77	55010
王五	11000	89	73542
赵六	15000	54	63201

以上例子是关于房屋价格的数据。

这个表格表包含了 房屋总价 与 房屋单价 和 房屋面积 等数据。

特征：其中 $x_1$为 房屋单价 ，$x_2$为 房屋面积 ，两者即为 特征。

标签：y为 房屋总价 ，我们称作为标签

我们可以利用一个线性函数来拟合给定的 房屋总价 与 房屋单价 和 房屋面积 数据，从而训练出模型$h_{\theta }$

训练完模型之后，如果给定我们一组新的数据，只包含 房屋单价 和 房屋面积 而不包含房屋总价，我们的目的就是预测 房屋总价 ，我们就能根据模型${ｈ}_{\theta }$来预测新的 房屋总价

$$h_{\theta }(x)={\theta }_0+{\theta }_1{x}_1+{\theta }_2{x}_2$$

参数：${\theta }_1$和${\theta }_{\mathrm{２}}$我们称作为参数,表示的是$x_1$和$x_2$两个特征对于目标函数的重要度

偏置项：${\theta }_0$表示的是偏置项

目标函数值：$h_{\theta }$是我们通过模型预测出来的房屋总价，即目标函数值。

上面的这个式子是当一个模型只有两个特征($x_1$,$x_2$)的时候的线性回归式子。

正常情况下，房屋总价和很多特征相关联，并不只是简单的这两个特征。所以我们需要把这个式子进行通用化。

假如有n个特征的话，那么式子就会变成下面的样子

$h_{\theta }(x)={\theta }_0+{\theta }_1{x}_1+{\theta }_2{x}_2+...+{\theta }_n{x}_n={\theta }_0+{\sum }_{i=1}^n{\theta }_i{x}_i$ 式①

上面的式子是一个多项求和的式子，我们将其向量化表示如下：

我们令$\theta =\left[\begin{array}{c}{\theta }_1\\ \\ {\theta }_2\\ \\ {\theta }_3\\ \\ ...\\ \\ {\theta }_n\end{array}\right]$，$x=\left[\begin{array}{c}{x}_1\\ \\ x_2\\ \\ x_3\\ \\ ...\\ \\ x_n\end{array}\right]$,

则式①表示为：

$$h_{\theta }(x)={\theta }_0+\left[\begin{array}{ccccc}{\theta }_1,& {\theta }_2,& {\theta }_3,& ...& {\theta }_n\end{array}\right]\left[\begin{array}{c}{x}_1\\ \\ x_2\\ \\ x_3\\ \\ ...\\ \\ x_n\end{array}\right]＝\sum _{i=1}^n{\theta }_i{x}_i＝{\theta }_0+{\theta }^{T}X$$

通常情况下,我们为了使得方程更为简洁，往往单独增加一个为１的特征，则我们重新令

$\theta =\left[\begin{array}{c}{\theta }_0\\ \\ {\theta }_1\\ \\ {\theta }_2\\ \\ {\theta }_3\\ \\ ...\\ \\ {\theta }_n\end{array}\right]$，$x=\left[\begin{array}{c}1\\ \\ x_1\\ \\ x_2\\ \\ x_3\\ \\ ...\\ \\ x_n\end{array}\right]$,则$h_{\theta }$则简化为$h_{\theta }(x)＝{\theta }^{T}X$　式②

因此，我们称式②为线性回归的目标函数，即我们训练出来的线性回归模型。

二.线性回归的损失函数

2.1损失函数的定义

我们通常需要一个机制去评估我们训练出来的模型$h_{\theta }$的好坏，，一般这个机制我们称为损失函数（loss　function）或者错误函数(error function)，作用是用来描述h函数不好的程度，通常用$J_{\theta }$表示，而对于线性回归模型，我们用以下式子来作为其损失函数：

$$J(\theta )=\frac{1}{2}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$$

我们的目的就是求得使得损失函数$J(\theta )$的值为最小值时候，$\theta$的所取值，从而得到我们最优的模型$h_{\theta }(x)$

2.2损失函数的推导过程（最大似然估计与最小二乘法）

我们刚才已经提到$J(\theta )=\frac{1}{2}{\sum }_{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$为我们的损失函数，那么，它是如何得到的呢?

接下来，就让我们利用极大似然估计推导来对损失函数进行推导

在推导之前让我们来看下面这一张图：

图中的横坐标X1 和 X2 分别代表着 两个特征(房屋单价、房屋面积) 。纵坐标Y代表标签值(房屋总价)。

其中红点代表的就是实际的标签值(房屋总价).而平面上和红点竖向相交的点代表着我们根据线性回归模型得到的值。

也就是说实际的标签值和根据模型预估的值之间是有一定误差的，这个就是误差项，我们记为${\epsilon }^{(i)}$。因此我们有

$y^{(i)}=h_{\theta }(x^{(i)})+{\epsilon }^{(i)}＝{\theta }^T{x}^{(i)}+{\epsilon }^{(i)}$ 式③

因为误差项是真实值和误差值之间的一个差距。那么肯定我们希望误差项越小越好。

在这里我们假设误差项${\epsilon }^{(i)}$是独立的并且具有相同的分布的，并且服从均值为０，标准偏差为$\sigma$的高斯分布

解释：

1.独立：张三和李四购买的房子，房子总价之间互不影响

２.具有相同分布：张三和李四是买的的是同一座城市的房子

３.高斯分布：绝大多数的情况下，在一个的空间内浮动不大

下面是高斯分布的图像，忘记的同学们可以回忆下。

由于我们假设${\epsilon }^{(i)}$是符合高斯分布的，因此，我们关于${\epsilon }^{(i)}$的概率密度函数为：

$p\left({\epsilon }^{(i)}\right)=\frac{1}{\sqrt{2\pi }\sigma }\exp \left(-\frac{{\left({\epsilon }^{(i)}\right)}^2}{2{\sigma }^2}\right)$　式④

由　式③　和　式④
$p(y^{(i)}|x^{(i)};\theta )=\frac{1}{\sqrt{2\pi }\sigma }exp(-\frac{(y^{(i)}-{\theta }^T{x}^{(i)}{)}^2}{2{\sigma }^2})$

$p\left(y^{(i)}|x^{(i)};\theta \right)$表示的意思是，在$\theta$已知的情况下，在$x^{(i)}$发生的条件下，$y^{(i)}$发生的条件概率

误差项肯定是越小越好了，那么接下来要讨论的就是什么样的参数和特征的组合能够让误差项最小呢？

这里就引入了似然函数。

似然函数的作用就是要根据样本来求什么样的参数和特征的组成能够最接近真实值。越接近真实值则误差越小。

因此，似然函数为

$$\begin{array}{rl}L(\theta )& =\prod _{i=1}^{m}p\left(y^{(i)}|x^{(i)};\theta \right)\\ & =\prod _{i=1}^m\frac{1}{\sqrt{2\pi }\sigma }\exp \left(-\frac{{\left(y^{(i)}-{\theta }^T{x}^{(i)}\right)}^2}{2{\sigma }^2}\right)\end{array}$$

由极大似然估计的定义，我们需要$L(\theta )$最大，那么我们怎么才能使的这个值最大呢？
两边取对数把多个数相乘，转化成多个数相加的形式，对这个表达式进行化简如下：
$$\begin{array}{rl}\log L(\theta )& =\log \prod _{i=1}^{m}p\left(y^{(i)}|x^{(i)};\theta \right)\\ & =\sum _{i=1}^m\log \frac{1}{\sqrt{2\pi }\sigma }\exp \left(-\frac{{\left(y^{(i)}-{\theta }^T{x}^{(i)}\right)}^2}{2{\sigma }^2}\right)\\ & =\mathrm{mlog}\frac{1}{\sqrt{2\pi }\sigma }-\frac{1}{2{\sigma }^2}\sum _{i=1}^m{\left(y^{(i)}-{\theta }^T{x}^{(i)}\right)}^2\end{array}$$

记$$l(\theta )=logL(\theta )$$

则$$l(\theta )=mlog\frac{1}{\sqrt{2\pi }\sigma }-\frac{1}{2{\sigma }^2}\sum _{i=1}^m(y^{(i)}-{\theta }^T{x}^{(i)}{)}^2$$

要使得$l(\theta )$最大，由于

$mlog\frac{1}{\sqrt{2\pi }\sigma }$是常数，我们可以忽略不计，

只要使得$-\frac{1}{2{\sigma }^2}{\sum }_{i=1}^m(y^{(i)}-{\theta }^T{x}^{(i)}{)}^2$最大即可

也就是使得$-\frac{1}{2}{\sum }_{i=1}^m(y^{(i)}-{\theta }^T{x}^{(i)}{)}^2$最大即可，即

$\frac{1}{2}{\sum }_{i=1}^m(y^{(i)}-{\theta }^T{x}^{(i)}{)}^2$最小即可

我们记
$J(\theta )=\frac{1}{2}{\sum }_{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$　(最小二乘法)

所以，我们最后由极大似然估计推导得出，我们希望 $J(\theta )$最小，而这刚好就是最小二乘法做的工作。

而回过头来我们发现，之所以最小二乘法有道理，是因为我们之前假设误差项服从高斯分布，假如我们假设它服从别的分布，那么最后的目标函数的形式也会相应变化。

好了，上边我们得到了有极大似然估计或者最小二乘法。

我们的模型求解需要最小化目标函数$J(\theta )$。

那么我们的$\theta$到底怎么求解呢？有没有一个解析式可以表示$\theta$？

接下来，我们就对$\theta$进行求解

三.正规方程求解$\theta$

在上一节，我们推导得到线性回归模型的损失函数为:
$$\begin{array}{rl}J(\theta )& =\frac{1}{2}\sum _{i=1}^m{\left(h_{\theta }\left(x^{(i)}\right)-y^{(i)}\right)}^2\\ & =\frac{1}{2}(X\theta -Y{)}^T(X\theta -Y)\\ & =\frac{1}{2}\left((X\theta {)}^T-Y^T\right)(X\theta -Y)\\ & =\frac{1}{2}\left({\theta }^T{X}^T-Y^T\right)(X\theta -Y)\\ & =\frac{1}{2}\left({\theta }^T{X}^{T}X\theta -{\theta }^T{X}^{T}Y-Y^{T}X\theta +Y^{T}Y\right)\end{array}$$

要求$J(\theta )$的最小值所对应的$\theta$值，即求驻点，我们让$\frac{\partial J(\theta )}{\partial \theta }=0$即可

因为$$J(\theta )=\frac{1}{2}\left({\theta }^T{X}^{T}X\theta -{\theta }^T{X}^{T}Y-Y^{T}X\theta +Y^{T}Y\right)$$
所以
$$\begin{array}{rl}\frac{\partial J(\theta )}{\partial \theta }& =\frac{1}{2}\left(2X^{T}X\theta -X^{T}Y-{\left(Y^{T}X\right)}^T\right)\\ & =\frac{1}{2}\left(2X^{T}X\theta -X^{T}Y-X^{T}Y\right)\\ & =\frac{1}{2}\left(2X^{T}X\theta -2X^{T}Y\right)\\ & =X^{T}X\theta -X^{T}Y\end{array}$$

因为$\frac{\partial J(\theta )}{\partial \theta }=0$即$X^{T}X\theta -X^{T}Y=0$

则$$X^{T}X\theta ＝X^{T}Y$$

因此
$$\theta ＝(X^{T}X{)}^{-1}{X}^{T}Y$$

至此，我们就可以求得到$\theta$的值了

关于向量矩阵的求导可参看【白话机器学习】机器学习算法推导过程中的数学基础知识

上述最后一步有一些问题，假如 $X^{T}X$不可逆呢？

我们知道 $X^{T}X$ 是一个半正定矩阵，所以若$X^{T}X$ 不可逆或为了防止过拟合，我们增加$\lambda$扰动，得到

$$\theta ＝(X^{T}X+\lambda I{)}^{-1}{X}^{T}Y$$

三.梯度下降法求解$\theta$

在上边我们用正规方程求解出了参数$\theta$,实际python 的 numpy库就是使用的这种方法。

但是，这仅仅适用于参数维度不大的情况下，通常我们在项目工程中，参数的维度是相当大的，使用正规方程求解就不适用了，这里我们用梯度下降算法求解$\theta$。

梯度下降算法是一个最优化算法，它是沿梯度下降的方向求解损失函数$J(\theta )$极小值。

$$损失函数：J(\theta )=\frac{1}{2m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$$

这个损失函数是我们之前推导出来的最小二乘的公式。只不过1/2后面多了个m。m是总的数据量。意味着是求多个数据之后的平均值。

我们的目标是
$$minimiz{e}_{\theta }J(\theta )$$

3.1梯度下降的通俗理解：

想象一下你正站立在这座红色山上。

在梯度下降算法中，我们要做的就是旋转360度，看看我们的周围，并问自己要在某个方向上，用小碎步尽快下山。这些小碎步需要朝什么方向？如果我们站在山坡上的一点，你看一下周围，你会发现最佳的下山方向。

你再看看周围，然后再一次想想，我应该从什么方向迈着小碎步下山？然后你按照自己的判断又迈出一步。

重复上面的步骤，直到你接近局部最低点的位置。

3.2梯度下降法的具体过程：

梯度下降法的具体过程如下：

1）对$\theta$进行初始化。

我们可以随机初始化，或者初始化为全零的向量。

图片中红圈上的点，从上面的图可以看出：初始点不同，获得的最小值也不同，因此梯度下降求得的只是局部最小值；

2）确定学习率(Learning rate)$\alpha$

学习率可以理解为下山时每一步迈的大小。

步子迈得太大有可能不收敛，无法到达山的最低处，即损失函数的最小值处。

步子迈的太小下山速度太慢。

工程经验学习率我们可以先使用0.1试下，如果不收敛或者收敛慢，再试0.01、0.001。

小窍门学习率并不一定全程不变，可以刚开始的时候大一些，后期的时候不断调小。

３）沿着负梯度方向更新$\theta$。

采用负梯度方向的原因是负梯度方向是函数值减小的最快的方向

$$\theta :=\theta -\alpha \cdot \frac{\partial J(\theta )}{\partial \theta }$$

在线性回归中,对于单个样本

$$\begin{array}{rl}\frac{\partial J(\theta )}{{\theta }_j}& =\frac{\partial }{\partial {\theta }_j}\frac{1}{2}(h_{\theta }(x^{(i)})-y^{(i)}{)}^2\\ & \\ & =\frac{1}{2}\cdot 2\cdot (h_{\theta }(x^{(i)})-y^{(i)})\frac{\partial }{\partial {\theta }_j}(h_{\theta }(x^{(i)})-y^{(i)})\\ & \\ & =(h_{\theta }(x^{(i)})-y^{(i)})\frac{\partial }{\partial {\theta }_j}(\sum _{i=0}^n{\theta }_i{x}_i-y^{(i)})\\ & \\ & =(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}\\ & \end{array}$$

因此对于ｍ个样本来说，

$${\theta }_j:={\theta }_j-\alpha \sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

４）训练终止。

不断重复步骤３）

损失函数达到预设值的一个值，或者收敛不明显时，可以终止训练。

最后得到一个局部最优解，不一定是全局最优，但是是堪用的。

此时所对应的$\theta$的值就是所求的$\theta$

3.3梯度下降法的分类：

梯度下降法作为机器学习中较常使用的优化算法，其有着三种不同的形式：批量梯度下降（Batch Gradient Descent）、随机梯度下降（Stochastic Gradient Descent）以及小批量梯度下降（Mini-Batch Gradient Descent）。其中小批量梯度下降法也常用在深度学习中进行模型的训练。接下来，我们将对这三种不同的梯度下降法进行理解。

为了便于理解，这里我们将使用只含有一个特征的线性回归来展开。此时线性回归的假设函数为：

$$h_{\theta }(x^{(i)})={\theta }_1{x}^{(i)}+{\theta }_0$$

其中，$i=1,2,...,m$表示样本数。

对应的**目标函数（代价函数）**即为：

$$J({\theta }_0,{\theta }_1)=\frac{1}{2m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$$

下图为$J({\theta }_0,{\theta }_1)$ 与参数 ${\theta }_0,{\theta }_1$的关系的图：

3.3.1、批量梯度下降（Batch Gradient Descent，BGD）

批量梯度下降法是最原始的形式，它是指在每一次迭代时使用所有样本来进行梯度的更新。从数学上理解如下：

（1）对目标函数求偏导：

$$\frac{\Delta J({\theta }_0,{\theta }_1)}{\Delta {\theta }_j}=\frac{1}{m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

其中$i=1,2,...,m$ 表示样本数， $j=0,1$ 表示特征数，这里我们使用了 偏置项 $x_0^{(i)}=1$

2）每次迭代对参数进行更新：

$${\theta }_j:={\theta }_j-\alpha \frac{1}{m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

注意这里更新时存在一个求和函数，即为对所有样本进行计算处理，可与下文SGD法进行比较。

伪代码形式为：

repeat{

$${\theta }_j:={\theta }_j-\alpha \frac{1}{m}\sum _{i=1}^m(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}(forj=0,1)$$

}

优点：

（1）一次迭代是对所有样本进行计算，此时利用矩阵进行操作，实现了并行。

（2）由全数据集确定的方向能够更好地代表样本总体，从而更准确地朝向极值所在的方向。当目标函数为凸函数时，BGD一定能够得到全局最优。

缺点：

（1）当样本数目 m 很大时，每迭代一步都需要对所有样本计算，训练过程会很慢。

从迭代的次数上来看，BGD迭代的次数相对较少。其迭代的收敛曲线示意图可以表示如下：

3.3.2随机梯度下降（Stochastic Gradient Descent，SGD）

随机梯度下降法不同于批量梯度下降，随机梯度下降是每次迭代使用一个样本来对参数进行更新。使得训练速度加快。

对于一个样本的目标函数为：

$$J^{(i)}({\theta }_0,{\theta }_1)=\frac{1}{2}(h_{\theta }(x^{(i)})-y^{(i)}{)}^2$$

（1）对目标函数求偏导：

$$\frac{\Delta J^{(i)}({\theta }_0,{\theta }_1)}{{\theta }_j}=(h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

（2）参数更新：

$${\theta }_j:={\theta }_j-\alpha (h_{\theta }(x^{(i)})-y^{(i)})x_j^{(i)}$$

注意，这里不再有求和符号

伪代码形式：

repeat{
    for i=1,...,m{
      $$\theta_j := \theta_j -\alpha (h_{\theta}(x^{(i)})-y^{(i)})x_j^{(i)}$$
      (for j =0,1)
    }
  }

优点:

（1）由于不是在全部训练数据上的损失函数，而是在每轮迭代中，随机优化某一条训练数据上的损失函数，这样每一轮参数的更新速度大大加快。

缺点：

（1）准确度下降。由于即使在目标函数为强凸函数的情况下，SGD仍旧无法做到线性收敛。

（2）可能会收敛到局部最优，由于单个样本并不能代表全体样本的趋势。

（3）不易于并行实现。

解释一下为什么SGD收敛速度比BGD要快：

答：这里我们假设有30W个样本，对于BGD而言，每次迭代需要计算30W个样本才能对参数进行一次更新，需要求得最小值可能需要多次迭代（假设这里是10）；而对于SGD，每次更新参数只需要一个样本，因此若使用这30W个样本进行参数更新，则参数会被更新（迭代）30W次，而这期间，SGD就能保证能够收敛到一个合适的最小值上了。也就是说，在收敛时，BGD计算了 10×30W 次，而SGD只计算了 1×30W 次。

从迭代的次数上来看，SGD迭代的次数较多，在解空间的搜索过程看起来很盲目。其迭代的收敛曲线示意图可以表示如下：

3.3.3 mini-batch梯度下降法

小批量梯度下降，是对批量梯度下降以及随机梯度下降的一个折中办法。其思想是：每次迭代 使用 batch_size 个样本来对参数进行更新。

这里我们假设 $batchsize=10$ ，样本数 $m=1000$ 。

伪代码形式为：

  repeat{
    for i=1,11,21,31,...,991{
     $$\theta_j := \theta_j - \alpha \frac{1}{10} \sum_{k=i}^{(i+9)}(h_{\theta}(x^{(k)})-y^{(k)})x_j^{(k)}$$
      (for j =0,1)
    }
  }

优点：

（1）通过矩阵运算，每次在一个batch上优化神经网络参数并不会比单个数据慢太多。

（2）每次使用一个batch可以大大减小收敛所需要的迭代次数，同时可以使收敛到的结果更加接近梯度下降的效果。(比如上例中的30W，设置batch_size=100时，需要迭代3000次，远小于SGD的30W次)

（3）可实现并行化。

缺点：

（1）batch_size的不当选择可能会带来一些问题。

batcha_size的选择带来的影响：

（1）在合理地范围内，增大batch_size的好处：

a. 内存利用率提高了，大矩阵乘法的并行化效率提高。

b. 跑完一次 epoch（全数据集）所需的迭代次数减少，对于相同数据量的处理速度进一步加快。

c. 在一定范围内，一般来说 Batch_Size 越大，其确定的下降方向越准，引起训练震荡越小。

（2）盲目增大batch_size的坏处：

a. 内存利用率提高了，但是内存容量可能撑不住了。

b. 跑完一次 epoch（全数据集）所需的迭代次数减少，要想达到相同的精度，其所花费的时间大大增加了，从而对参数的修正也就显得更加缓慢。

c. Batch_Size 增大到一定程度，其确定的下降方向已经基本不再变化。

下图显示了三种梯度下降算法的收敛过程：