线性回归概述

一 概念

1 线性模型的基本形式

给定n个属性描述的记录$x=(x_1,x_2,...,x_n{)}^T,x$是n维列向量，线性模型试图得到一个属性的线性组合来描述预测函数，即
$$\begin{array}{r}f(x)=w_1{x}_1+w_2{x}_2+......+w_n{x}_n+b\end{array}$$
用向量表示为：
$$\begin{array}{r}f(x)=w^{T}x+b\end{array}$$
其中$w=(w_1,w_2,...,w_n{)}^T$,$w_i$是权重系数，$b$是偏移量。

2 线性回归

给定数据集D={(x1,y1),(x2,y2),...,(xn,yn)}D=\lbrace(x_1,y_1),(x_2,y_2),...,(x_n,y_n)\rbraceD={(x1​,y1​),(x2​,y2​),...,(xn​,yn​)},其中$x_i$是n维列向量，$y_i\in R$。线性回归就是试图学得一个线性模型使之能尽可能地符合实值输出标记。显然线性模型中x和y是已知量，w和b是未知量，线性回归就是求解最佳的w和b的过程。

广义地说，线性回归分为三个步骤：

• 假设线性函数
• 构造损失函数
• 求解最小损失值下的w和b

根据2、3步的不同选择，线性回归可采用多种方法，常见的有最小二乘法和梯度下降法。

二 最小二乘法实现线性回归

1 简单情况

先考虑最简单地情形：x只有一个属性。这时线性模型就是最简单地一次函数，线性回归的目标就是求解：
$$\begin{array}{r}f(x_i)=w{x}_i+b,使得f(x_i)\approx y_i\end{array}$$
最小二乘法的目标是均方误差最小化（从几何的角度理解均方误差最小化就是找到一条直线，使所有样本到直线上的欧式距离之和最小）：
$$\begin{array}{rl}定义损失函数：L(w,b)& =\sum _{i=1}^n(f(x_i)-y_i{)}^2\\ & =\sum _{i=1}^n(y_i-w{x}_i-b{)}^2\end{array}$$
$$\begin{array}{r}最优解：(w^{\ast },b^{\ast })=min\sum _{i=1}^n(y_i-w{x}_i-b{)}^2\end{array}$$
此时线性回归转化为了一个数学问题：已知二元函数$L(w,b)$和样本集$D=(X,Y)$,当$L(w,b)$取得最小值时,求w和b的估计值。

如何求？实际上，$L(w,b)$是典型的关于w和b的凸函数，具有唯一的极值点且是最值点。此时，我们对$w$和$b$分别求偏导：
$$\begin{array}{rl}\frac{\partial L(w,b)}{\partial w}& =2(w\sum _{i=1}^n{x}_i^2-\sum _{i=1}^n(y_i-b)x_i)\\ \frac{\partial L(w,b)}{\partial b}& =2(mb-\sum _{i=1}^n(y_i-w{x}_i))\end{array}$$
令$\frac{\partial L(w,b)}{\partial w}=0$和$\frac{\partial L(w,b)}{\partial b}=0$此时求得最优解：
$$\begin{array}{rl}w& =\frac{{\sum }_{i=1}^n{y}_i(x_i-\overline{x})}{{\sum }_{i=1}^n{x}_i^2-\frac{1}{n}({\sum }_{i=1}^n{x}_i{)}^2}\\ b& =\frac{1}{n}\sum _{i=1}^n(y_i-w{x}_i)\end{array}$$
上述的解看起来非常复杂，实际上这只是为了演示最小二乘法，通常我们会用矩阵的方法求。

2 一般情况

给定数据集D={(x1,y1),(x2,y2),...,(xn,yn)}D=\lbrace(x_1,y_1),(x_2,y_2),...,(x_n,y_n)\rbraceD={(x1​,y1​),(x2​,y2​),...,(xn​,yn​)},其中$x_i$是n维列向量，$y_i\in R$，$w=(w_1,w_2,...,w_n{)}^T$。线性回归试图求得：
$$\begin{array}{r}f(x_i)=w^T{x}_i+b,使得f(x_i)\approx y_i\end{array}$$
我们令$\tilde{w}=\left[\begin{array}{c}w\\ b\end{array}\right]$，$X=\left[\begin{array}{cc}{x}_1^T& 1\\ x_2^T& 1\\ :& :\\ x_n^T& 1\end{array}\right]$,$y=(y_1,y_2,...,y_n{)}^T$。此时线性模型可表示为：
$$\begin{array}{r}\left[\begin{array}{c}f(x_1)\\ f(x_2)\\ :\\ f(x_n)\end{array}\right]=\left[\begin{array}{cc}{x}_1^T& 1\\ x_2^T& 1\\ :& :\\ x_n^T& 1\end{array}\right]\ast \left[\begin{array}{c}w\\ b\end{array}\right]\end{array}$$
即$y\approx X\tilde{w}$

任然采用最小二乘法：
$$\begin{array}{rl}& 令损失函数:L(\tilde{w})=(y-X\tilde{w}{)}^T(y-X\tilde{w})\\ & 最优解：{\tilde{w}}^{\ast }=minL(\tilde{w})\end{array}$$
对$\tilde{w}$求导：
$$\begin{array}{r}\frac{\partial L\tilde{w}}{\partial \tilde{w}}=2X^T(X\tilde{w}-y)\end{array}$$
令$\frac{\partial L\tilde{w}}{\partial \tilde{w}}=0$,得：$X^{T}X\tilde{w}=X^{T}y$

若$X^{T}X$是满秩方阵，$\tilde{w}=(X^{T}X{)}^{-1}{X}^{T}y$

若$X^{T}X$不是，$X^{T}X\tilde{w}=X^{T}y$非齐次线性方程组有无穷多解，$\tilde{w}$有n个。此时应该根据学习算法的归纳偏好，从n个中选取一个作为输出。常见的方法是引入正则化。

三 梯度下降实现线性回归

1 梯度

梯度是微分学的重要概念，对于一元函数来说梯度就是某点斜率，对于多元函数来说梯度是某点变化最快的方向。设一个二元函数：
$$\begin{array}{r}L（x,y）=x^2+y^2\end{array}$$
其梯度函数为：$\nabla L(x,y)=(\frac{\partial L}{\partial x},\frac{\partial L}{\partial y})=(2x,2y)$ 梯度是个向量

2 梯度下降

我们令$\tilde{w}=\left[\begin{array}{c}w\\ b\end{array}\right]$，$X=\left[\begin{array}{cc}{x}_1^T& 1\\ x_2^T& 1\\ :& :\\ x_n^T& 1\end{array}\right]$,$y=(y_1,y_2,...,y_n{)}^T$，同样用均方误差来表示损失函数但略有不同:$L(\tilde{w})=\frac{1}{2m}(y-X\tilde{w}{)}^T(y-X\tilde{w})$，除以m是为了表示方差，除以1/2是为了求导约掉。此时梯度函数为$\nabla L(\tilde{w})=\frac{1}{m}{X}^T(X\tilde{w}-y)$。

上一节讲过，线性回归 <=>求解损失函数最小时的$\tilde{w}$。最小二乘法对于人来说很好计算，但是对于计算机来说不容易，而梯度下降也可以做到并且易于计算机实现。梯度下降的数学形式为：
$$\begin{array}{r}{\tilde{w}}^{k+1}=w^k-\alpha \nabla L({\tilde{w}}^k)\end{array}$$
其中$\alpha$为步长，表示每一步下降的广度。梯度下降是个迭代的过程，直到梯度趋近于0时，这时表示曲面比较平滑（意味着达到了最低点），此时${\tilde{w}}^k$是最优解。 用一个图片来描述梯度下降的过程：

2 代码实现

#梯度下降
import numpy as np
import matplotlib. pyplot as plt

#代价函数:均方误差的半
def  loss_funciton(X,Y,theta):
    diff=np.dot(X,theta)-Y
    return (1.0/(2*m))*np.dot(np.transpose(diff),diff)

#梯度函数：代价函数求导得的梯度向量
def  gradient_function(X,Y,theta):
    diff =np.dot(X,theta)-Y
    return (1.0/m)*np.dot(np.transpose(X),diff)

#梯度下降：直到梯度趋于0
def  gradient_descent(X,Y,alpha):
    m,n=X.shape
    #最后一列拼接1
    X=np.c_[X,np.ones((m,1))]
    #初始化一个theta
    theta = np.ones((n+1,1))
    # 迭代过程
    gradient = gradient_function(X, Y, theta)
    while not (np.all(np.absolute(gradient) <= 1e-5 )):

        theta = theta - alpha * gradient
        gradient = gradient_function(X, Y, theta)
    return  theta

m=20
alpha=0.01
X = np.arange(1, m+1).reshape(m, 1)
X1=np.c_[X,np.ones((m,1))]

Y= np.array([
    3, 4, 5, 5, 2, 4, 7, 8, 11, 8, 12,
    11, 13, 13, 16, 17, 18, 17, 19, 21
]).reshape(m, 1)

theta_result=gradient_descent(X,Y,alpha)
print("theta的值是：", theta_result)

#做出拟合图像
f=np.dot(X1,theta_result)
plt.plot(X,f)
plt.plot(X,Y,'ro')
plt.show()

注意，步长不能过大。过大的话可能导致一两次迭代后就使得theta越过了最低点，这样theta会在另一边无限增大，从图像来说就是一次越过了山谷到了另一边又无限往上走。

四 总结

线性回归的本质是找到一个线性模型（求解出w和b），使样本集D在该模型上具有最小的误差。然而，度量这个误差的函数是多样的，所以求解方法也是多样的，但本质都是在求最小值点。本文叙述的最小二乘法和梯度下降法，是一种很好地求解最小值点的方法。