一元一次线性回归

抽象问题

现有一组数据，共有两列，分别为x和y，如下所示

原数据

现将数据用python作散点图，观察其变化趋势

import pandas as pd 
import numpy as np
import matplotlib.pyplot as plt

# 读取数据并转为np.array的格式
train = pd.read_csv('click.csv').values 
# 读取所有行的第0列作为训练数据的x
train_x = train[:, 0] 
# 读取所有行的第1列作为训练数据的y
train_y = train[:, 1]
# 画散点图
plt.plot(train_x, train_y, 'o')
plt.show()

散点图

我们的问题是要找到x和y的对应关系，根据此对应关系，我们就可以在给定任意x的情况下，计算得到y的值。

定义模型

观察上面的散点图可以看出来，随着x的增大，y值也在近似的不断增大。据此，我们推测，这组数据可能符合一次函数的关系，即：

此式可改写为

其中表示含有参数 ，与变量 相关的函数。该函数即为：一元一次线性回归模型。

注：统计学中，常用 来表示未知数和推测值，采用 加数字下标的形式，是为了防止当未知数的数量增加时，表达式中大量出现这样的符号。这样不但不易理解，还容易出现符号本身不够用的情况。

只要求得了上述公式中的未知参数和，即知道了x和y的关系。

现在问题转化为求未知参数和

求解方法

已知1个自变量，根据模型公式，可以求得与其对应的因变量。

那么表示实际值与模型预测值之间的误差。

我们要求解的模型函数最理想的情况就是使得，即令误差值。

实际中，很难做到让所有点的误差都等于0，因此，我们要做的就是让所有点的误差之和尽可能的小。

我们一般使用误差的平方来表征误差，一方面是去除误差正负抵消的影响，一方面是为了使之后求解微分更加方便。

所有点的误差的平方和可以用如下公式来表示：

这个表达式称为目标函数，是Error的首字母。

现在，我们的问题转化为：求出使得目标函数取最小值的参数和

梯度下降法(最速下降法)

现有函数根据导数的意义：

1. 当，即导数为正，随x的增大而增大，因此，要使得减小，需要x负移

2. 当，即导数为负，随x的增大而减小，因此，要使得减小，需要x正移

据此可知，要使得不断减小，只需要将x向与导数符号相反的方向移动即可。用数学表达式表述为：

1. ，表示通过B来定义A

2. ，表示学习率。

   η越小，更新次数越多，速度越慢

   η越大，更新次数越少，速度越快，但有可能导致不收敛。

梯度下降法的核心：根据误差函数的导数的符号来决定x移动的方向。

根据梯度下降法的原理，我们知道了怎样更新参数和，就能使目标函数的值达到最小。

参数更新方法如下

现在只要计算出和，代入上式，就可不断更新参数，求得最佳的和。

计算上式中的两个偏导数过程如下：

同理：

因此，参数更新表达式如下：

根据该式，不断迭代，直到求得最佳的参数和，即确定了模型表达式

python实现

根据上述的理论及推导，我们知道：

1. 模型表达式为
   
2. 目标函数为
   
3. 参数更新表达式为
   

在进行参数更新之前，应该将数据标准化（也称z-score规范化）。这个步骤非必须，但是做了之后，参数的收敛会更快

# 读取数据并转为np.array的格式
train = pd.read_csv('click.csv').values 
train_x = train[:, 0] 
train_y = train[:, 1]

# 对变量train_x进行标准化
mu = train_x.mean()
sigma = train_x.std()
def standardize(x):
    return (x-mu) / sigma
train_x_std = standardize(train_x)

python逐步实现

# 1,随机生成初始参数theta0 和 theta1
theta0, theta1 = np.random.rand(2)

# 2,定义模型函数f(x)
def f(x):
    y = theta0 + theta1 * x
    return y

# 3,定义目标函数E(x, y)
def E(x, y):
    e = 0.5 * np.sum((y-f(x)) ** 2)
    return e 

# 4,初始化学习率η
ETA = 0.001

# 5,初始化误差变化量diff
diff = 1

# 6,初始化更新次数count
count = 0 

# 7,迭代学习
error = E(train_x_std, train_y)    # 计算初始误差
# 开始迭代更新参数，直到误差的变化小于0.01
while diff > 0.01:
    # 更新参数theta0和theta1
    theta0 = theta0 - ETA * np.sum((f(train_x_std) - train_y))
    theta1 = theta1 - ETA * np.sum((f(train_x_std) - train_y) * train_x_std)
    # 计算更新参数后的当前误差
    error_current = E(train_x_std, train_y)
    # 计算原误差与当前误差的差值
    diff = error - error_current
    # 更新误差值
    error = error_current
    # 输出日志
    count += 1
    print(f'第{count}次：theta0 = {theta0:.3f}, theta1 = {theta1:.3f}, 差值 = {diff:.4f}')

部分输出日志如下

注：若执行多次，会发现迭代次数和误差的差值在每次执行时都不一样，这是因为随机初始化的参数和的不同而导致的。

此时，最终参数为：

print(f'theta0 = {theta0:.3f}')
print(f'theta1 = {theta1:.3f}')

作图查看拟合结果如下

plt.plot(train_x_std, train_y, 'o')     # 标准化后散点图
x = np.linspace(-3, 3, 100)
plt.plot(x, f(x))                       # 拟合直线
plt.show()

线性拟合结果

sklearn方法

在实际使用中，机器学习库sklearn为我们提供了更模块化的方式来进行线性回归，如下所示

from sklearn.linear_model import LinearRegression
# 1, 定义线性回归模型
lr = LinearRegression()
# 2, 训练模型
lr.fit(train_x_std.reshape(-1,1), train_y)

# 查看tehta0 和 theta1
print(f'theta0 = {lr.intercept_}')
print(f'theta1 = {lr.coef_}')

结果如下所示

和使用python逐步计算得到的参数一致。

np.polyfit(X, y, n)方法

Numpy库也提供了线性拟合方法：Numpy.polyfit(X, y, n)，使用方法如下所示：

# 一次线性拟合
z1 = np.polyfit(train_x_std, train_y, 1)
# 线性表示
p1 = np.poly1d(z1)
print(z1)
print(p1)

结果如下所示

和使用前两种方法得到的结果一致。

一元多次线性回归-多项式回归

观察前面做出的图像，发现一元一次线性模型的拟合效果并不好。

实际上，对于给定的数据，曲线比直线拟合地更好。因此，我们重新定义模型函数为
该式为一元多次线性回归，即多项式拟合模型表达式。

如果使用更高次数的表达式，如下所示，就能表示更复杂的曲线了次数越高，拟合的越好，但可能会出现过拟合问题

对于要解决的问题，在找出合适的表达式之前，需要不断地去尝试。

这里以二次拟合函数为例，我们增加了参数，此时：

1. 模型表达式为
   
2. 目标函数为
   

使用同样的方式，对进行微分，得到参数更新表达式如下：

即使增加参数（增加次数），比如有，等，依然可以用同样的方法来求出参数更新表达式。

像这样增加函数中多项式的次数，然后再使用函数的分析方法被称为多项式回归

多元线性回归

前述的线性回归都是一元回归，即只有一个变量。

但是，实际中要解决的很多问题是变量超过2个的复杂问题。

例如：有三个变量，分别为，，。此时，模型表达式为：

将此式推广到n个变量的情况，此时

模型表达式为：

我们可以把参数和变量看作向量，用黑体表示

其中 ，则：

因此：

对目标函数 求的偏导如下：

综上：

1. 模型表达式为：
   

2. 目标函数为：
   

3. 参数更新表达式为：
   

像这样包含了多个变量的回归称为多重回归。

采用矩阵求解参数

由于训练集数据有很多，所以我们把1行数据当作一个训练数据，以矩阵的形式来处理更方便。

n个变量，n个参数，n个训练数据的模型表达式如下：

对于参数更新表达式：

当时，求和项可展开为：

令：

则：

Python逐步实现

对于多项式回归：

# 1，创建训练数据的矩阵X
def to_matrix(x):
    return np.vstack([np.ones(x.shape[0]), x, x**2]).T
X = to_matrix(train_x_std)

# 2，随机初始化参数Theta
theta = np.random.rand(3)

# 3，定义预测函数
def f(x):
    return np.dot(x, theta)

# 4, 初始化参数
diff = 1
ETA = 0.001

# 5, 迭代更新参数
error = E(X, train_y)
while diff > 0.01:
    theta = theta - ETA * np.dot(f(X)-train_y, X)
    current_error = E(X, train_y)
    diff = error - current_error
    error = current_error

运行结束后，得到参数如下

将结果绘图展示如下：

x = np.linspace(-3, 3, 100)
plt.plot(train_x_std, train_y, 'o')   # 原数据散点图
plt.plot(x, f(to_matrix(x)))             # 多项式拟合曲线
plt.show()

np.polyfit(X, y, n)方法

z2 = np.polyfit(train_x_std, train_y, 2)
p2 = np.poly1d(z2)
print(z2)
print(p2)

结果如下所示

与python逐步实现结果一致