数据挖掘之线性回归原理（附代码）

1.什么是线性回归

这里引用百度百科的解释：利用数理统计中的回归分析，来确定变量间相互依赖的定量关系的一种统计分析方法，其表达形式为 $y^{(i)}=w^T{x}^{(i)}+e^{(i)}$，$e$为误差（用来代替统计误差或者抽样误差），它服从均值为0的标准正态分布。例如成年人的月收入与月消费的关系，二者之间肯定存在某种关系。
在求解之前，通常会给定几个假设：

• （1）自变量$(X_i,X_j)$ 之间应相互独立；
• （2）误差项$(e_i，e_j)$ 之间应相互独立;
• （3）模型线性可加。假设因变量为Y，自变量为X1，X2，则回归分析的默认假设为$Y=b+a1X1+a2X2+e$。
  线性性 ：X1每变动一个单位，Y相应变动a1个单位，与X1的绝对数值大小无关。
  可加性：X1对Y的影响是独立于其他自变量（如X2）的。

使用python生成模拟数据，具体代码与注释如下：

import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression

font = {'family': 'Microsoft YaHei',
        'weight': 'bold',
        'size': 10}
plt.rc('font', **font) 
# 设置字体

''' 1.模拟数据x, y'''
x1 = np.arange(1, 5, 0.1)  # 0.1为步长
y = x1 * 4 + 3 + 2 * np.random.randn(len(x1))

# 2.通过API实现参数估计
model_lr = LinearRegression()
model_lr.fit(x1.reshape(-1, 1), y.reshape(-1, 1))  # 转换为列向量形
# model_lr.coef_ (1，1), model_lr.intercept_(1,)都是ndarray类型
print('API模型求解的w={}, b={}'.format(round(model_lr.coef_[0][0], 4), round(model_lr.intercept_[0], 4)))

下图是通过python接口求解到的W和b，绘制出原始数据和拟合线。

2.最大似然估计（MLE）

假设 $y=X\ast W+e$，最后的$J(w)$就是我们要优化的目标函数也称为损失函数。
$$y^{(i)}=w\cdot x^{(i)}+e^{(i)}$$

$$概率密度：p(e^{(i)})=\frac{1}{\sqrt{2\pi }\sigma }\cdot exp(-\frac{{e^{(i)}}^2}{2{\sigma }^2})$$

$$p(y^{(i)}|w,x^{(i)})=\frac{1}{\sqrt{2\pi }\sigma }\cdot exp(-\frac{[y^{(i)}-w\cdot x^{(i)}{]}^2}{2{\sigma }^2})$$

对上式求极大似然估计MLE
$$MLE=\prod _{i=1}^n\frac{1}{\sqrt{2\pi }\sigma }\cdot exp(-\frac{[y^{(i)}-w\cdot x^{(i)}{]}^2}{2{\sigma }^2})$$
对上式左右两边求对数
$$ln(MLE)=\sum _{i=1}^n-\frac{[y^{(i)}-w\cdot x^{(i)}{]}^2}{2{\sigma }^2}+constant$$
由于是想对上式求最大值，因此负号去掉等价于求它的最小值，而方差虽然未知，但是是一个定值。得到下式
$$minJ(W)=\sum _{i=1}^n-\frac{[y^{(i)}-w\cdot x^{(i)}{]}^2}{2}$$

3.利用正规方程求解（矩阵求导）

将$J(w)$写成矩阵方程后，表达式和计算过程更加简洁，其中 $X^{T}X$ 是半正定矩阵，半正定矩阵在现实数据下，多数情况是存在逆矩阵的，如果不存在可以给它加一个$\lambda I$矩阵，可以避免无逆。推导过程如下：

$$J(W)=\frac{1}{2}{(y-Xw)}^T(y-Xw)，w为列向量$$
$$J(W)=\frac{1}{2}(w^T{X}^{T}Xw-w^T{X}^{T}y-y^{T}Xw+y^{T}y)$$

对化简后的$J(W)$进行求导

$$\nabla J(W)=\frac{1}{2}(2X^{T}Xw-X^{T}y-(y^{T}x{)}^T)$$
$$\nabla J(W)=X^{T}Xw-X^{T}y$$
$$求驻点：X^{T}Xw-X^{T}y=0$$

$$得出：w={(X^{T}X)}^{-1}{X}^{T}y$$

if you 不太了解矩阵求导，这边给出一种情况的简单推导过程，其他的可以自己尝试。

更多求导知识大家可以参考—>矩阵求导简介
附上正规方程求解方式的代码

''' 3.通过线性代数矩阵求导，也叫正规方程 X*W=y'''
X = np.array([x1])
X = np.insert(X, 1, values=[1] * len(X), axis=0)  
# 给array插入一行为1的元素 2xn矩阵， 为了最后把b也看成一个w
X = X.T  ''' X = nx2矩阵, 我们运算所需要的'''
inverse = np.linalg.inv(np.dot(X.T, X))
# 计算$X^TX的逆矩阵$
W = np.dot(np.dot(inverse, X.T), y)
print('正规方程求解的w={}，b={}'.format(round(W[0], 4), round(W[1], 4)))

上面的代码亲测有效管用，可以直接运行，求解过程或最后推导的一样，事实上API也是这样求解的。

4.梯度下降

梯度下降，顾名思义是根据梯度来不断逼近全局最小值（针对于凸函数），只不过是按照负梯度来的，因为正梯度是增大的方向。首先需要求出它的梯度。用原来的参数，减去它的梯度，为了防止学习步长过长，会跳过极小值点，因此通常会加一个小于1的学习率。

$$grad=\nabla J(W)=X^T(XW-y)$$

$$w=w-\lambda \ast grad$$

建议大家在做线性回归时，可以把损失函数写成均方误差$loss=(y_{p}redict-y{)}^2/m$，原因是降低它的梯度，可以减少迭代，可以从学习率和迭代次数上调整。最后付代码

# 4.通过梯度下降求解参数w,b  y=w*x+b
def grad(y):
    w = np.random.randn()
    # 随机生成一个w, b
    b = np.random.randn()
    iteration = 10000
    learn_rating = 0.0001
    for i in range(iteration):
        y_hat = w * x1 + b
        gradient_b = 2 / len(x1) * np.sum(y_hat - y)
        gradient_w = 2 / len(x1) * np.sum((y_hat - y) * x1)  # 每一个残差与x的积再作和
        w = w - learn_rating * gradient_w
        b = b - learn_rating * gradient_b
    # return round(w,4), round(b,4)

    # 下面是通过矩阵梯度下降,更简洁一点
    weight = np.random.randn(2, 1)  # 生成一个2x1的列向量w和b
    for i in range(iteration):
        gradient = X.T.dot(X.dot(weight) - y.reshape(-1, 1))
        weight -= learn_rating * gradient
    return weight[0][0], weight[1][0]


print('梯度下降求解的w={}, b={}'.format(grad(y)[0], grad(y)[1]))
# 和API求解的参数是一样的
MSE = np.mean((W[0] * x1 + W[1] - y) ** 2)
RSS = np.sum((W[0] * x1 + W[1] - y) ** 2)  # Residual Sum of Square! y_hat-y
ESS = np.sum((W[0] * x1 + W[1] - np.mean(y)) ** 2)  # Explain Sum of Square!
TSS = np.sum((y - np.mean(y)) ** 2)  # yi-y_bar
print('R-Square:', ESS / TSS)
plt.scatter(x1.reshape(-1, 1), y.reshape(-1, 1))
plt.plot(x1, W[0] * x1 + W[1], color='r', linestyle='--', label='一元一次线性拟合')
plt.legend()
plt.show()

好了，到此简易的线性回归原理就告一段落了，觉的有帮助的读者可以点个赞，有错误的地方也欢迎批评指正，never too late to learn。
资料参考：机器学习之线性回归