此为一元线性回归第3篇。数值解由编程实现，最为深入。
1. 线性回归的概念
线性回归是利用数理统计中回归分析，来确定两种或两种以上变量间相互依赖的定量关系的一种统计分析方法 ——百度百科
线性：经过模型训练，得到自变量和因变量之间是线性关系
回归：根据已知的输入输出的到模型，根据模型进行输入，得到连续的输出
2. 线性回归模型训练的内在逻辑
假设输入、输出之间存在线性模型 y=w0+w1x，其中w0、w1为模型参数，此时我们需要对模型进行训练，获取最优的w0、w1的组合。
假设有已知 (x1, y1)，
则对应的单样本误差为 (y - y1)^2/2
总体样本误差为 Σ(y - y)^2/2
此时可得到对应的损失函数：loss = Σ(w0+w1x-y)^2/2
损失函数就是总体误差关于模型参数w0、w1的函数，属于三维数学模型，即需要找到一组w0，w1使得loss取极小值
损失函数对w0、w1的偏导函数
Σ(w0+w1x-y)^2/2 = Σ(w0^2+2(w0(w1x-y))+(w1x-y)^2)/2 = Σ(w1^2x^2+2w1x(w0-y)+(w0-y)^2)/2
loss' = Σ(w0+w1x-y)
loss' = Σ(w1x^2+x(w0-y)) = Σ(x(w1x+w0-y))
这里采用梯度下降的方式根据偏导求loss的最小值所在位置即， 求此位置对应的参数w0、w1
此方式只是线性回归的原理，即根据损失函数(偏导函数)求出最小值对应的位置，
但是线性回归底部的实际代码是根据最小二乘法处理的，内部逻辑更加复杂。
下面说一下根据损失函数求模型参数的方式。
3. 手写模型的创建、训练和测试
首先的到w0 w1 和 loss 的图像
代码示例：
import numpy as np
import matplotlib.pyplot as plt
tx = np.array(
[0.5, 0.6, 0.8, 1.1, 1.4, 1.6, 1.7, 2.0, 2.5, 2.6, 3.0, 3.1, 3.3, 3.6, 3.8, 3.9, 4.0, 4.1, 4.2])
ty = np.array(
[5.0, 5.5, 6.0, 6.8, 7.0, 11, 8.4, 9.0, 9.2, 10.5, 11, 11.1, 12.0, 13.3, 13.4, 15.7, 16, 14.8, 14.9])
n = 500
w0_grid, w1_grid = np.meshgrid(np.linspace(-3, 10, n), np.linspace(-3, 10, n))
loss_grid = 0
loss = 0
# 获取loss
for x0, y0 in zip(tx, ty):
loss += (w0_grid + w1_grid * x0 - y0) ** 2 / 2
# 画图展示w0 w1 和 loss 的关系
plt.figure("Loss-Function", facecolor='lightgray')
loss3d = plt.gca(projection='3d')
loss3d.set_xlabel('w0')
loss3d.set_ylabel('w1')
loss3d.set_zlabel('loss')
loss3d.plot_surface(w0_grid, w1_grid, loss, cstride=30, rstride=30, cmap='jet')
plt.tight_layout()
plt.show()
由三维图可以看出，最优即最小的损失函数在深蓝色区域
梯度下降求loss最小值对应的w0、w1
代码示例：
times = 1000 # 梯度下降次数
lrate = 0.01 # 没梯度下降参数的变化率
# 初始化
w0 = 1
w1 = 1
for i in range(1, times + 1):
# 输出迭代过程
loss = ((w0 + w1 * tx - ty) ** 2).sum() / 2
"""
其中 {} 都表示占位符
{:4} 表示四个占位符，目的是为了数据更加整齐；
{:.8f} 表示保留8为小数位
.format(i, w0, w1, loss) 表示将format里边的数赋值到对应的占位符
"""
print('{:4}> w0={:.8f}, w1={:.8f}, loss={:.8f}'.format(i, w0, w1, loss))
# 偏导迭代
loss_w0 = (w0 + w1 * tx - ty).sum()
loss_w1 = (tx * (w0 + w1 * tx - ty)).sum()
w0 = w0 - lrate * loss_w0
w1 = w1 - lrate * loss_w1
print("w0", w0, '\t', 'w1', w1, '\n')
# 绘制回归线，预测结果
linex = np.linspace(tx.min(), tx.max(), 10)
liney = w1 + w0 * linex
# 画图
plt.figure("Linear-Regression", facecolor='lightgray')
plt.title('Linear-Regression', fontsize=18)
plt.grid(linestyle=':')
plt.scatter(tx, ty, s=80, marker='o', color='dodgerblue', label='Samples')
plt.plot(linex, liney, color='orangered', linewidth=2, label='Linear-Regression')
plt.legend()
plt.show()
以图形的方式展示迭代过程 即 w0 w1 loss的变化曲线
代码示例：
import numpy as np
import matplotlib.pyplot as plt
tx = np.array(
[0.5, 0.6, 0.8, 1.1, 1.4, 1.6, 1.7, 2.0, 2.5, 2.6, 3.0, 3.1, 3.3, 3.6, 3.8, 3.9, 4.0, 4.1, 4.2])
ty = np.array(
[5.0, 5.5, 6.0, 6.8, 7.0, 11, 8.4, 9.0, 9.2, 10.5, 11, 11.1, 12.0, 13.3, 13.4, 15.7, 16, 14.8, 14.9])
# 以图形的方式展示迭代过程 即 w0 w1 loss的变化曲线
times = 100 # 梯度下降次数
lrate = 0.02 # 没梯度下降参数的变化率
epoches = [] # 记录每次梯度下降的索引
# 创建列表
w0 = [1]
w1 = [1]
losses = []
for i in range(1, times + 1):
loss = ((w0[-1] + w1[-1] * tx - ty) ** 2).sum() / 2
print('{:4}> w0={:.8f}, w1={:.8f}, loss={:.8f}'.format(i, w0[-1], w1[-1], loss))
losses.append(loss)
# 偏导迭代
epoches.append(i)
loss_w0 = (w0[-1] + w1[-1] * tx - ty).sum()
loss_w1 = (tx * (w1[-1]*tx + w0[-1] - ty)).sum()
# print(loss_w0, loss_w1)
w0.append(w0[-1] - lrate * loss_w0)
w1.append(w1[-1] - lrate * loss_w1)
# 画图
plt.figure("Linear-Regression", facecolor='lightgray')
plt.title('Linear-Regression', fontsize=18)
plt.subplot(311) # 311 表示3行1列第1幅图
plt.grid(linestyle=':')
plt.ylabel(r'$w_0$', fontsize=14)
plt.plot(epoches, w0[:-1], color='dodgerblue', label=r'$w_0$')
plt.subplot(312)
plt.grid(linestyle=':')
plt.ylabel(r'$w_1$', fontsize=14)
plt.plot(epoches, w1[:-1], color='dodgerblue', label=r'$w_1$')
plt.subplot(313)
plt.grid(linestyle=':')
plt.ylabel(r'$loss$', fontsize=14)
plt.plot(epoches, losses, color='orangered', label=r'$loss$')
# plt.legend()
# plt.tight_layout()
# plt.show()
# 图像显示梯度下降曲线
n = 500
w0_grid, w1_grid = np.meshgrid(np.linspace(0, 9, n), np.linspace(0, 3.5, n))
loss_grid = 0
for x0, y0 in zip(tx, ty):
loss_grid += (w0_grid + w1_grid * x0 - y0) ** 2 / 2
plt.figure("loss func", facecolor='lightgray')
ax3d = plt.gca(projection='3d')
ax3d.set_xlabel('w0')
ax3d.set_ylabel('w1')
ax3d.set_zlabel('loss')
ax3d.plot_surface(w0_grid, w1_grid, loss_grid, cstride=30, rstride=30, cmap='jet')
ax3d.plot(w0[:-1], w1[:-1], losses, 'o-', color='red') # o- 表示连点成线
# plt.legend()
plt.tight_layout()
lrate = 0.001时
lrate = 0.01时 此时变化幅度过大，会来回横跳
# 等高线的方式绘制梯度下降过程
plt.figure("Batch Gradient Descent", facecolor='lightgray')
plt.title('Batch Gradient Descent', fontsize=20)
ax3d.set_xlabel('w0', fontsize=14)
ax3d.set_ylabel('w1', fontsize=14)
plt.grid(linestyle=':')
plt.contourf(w0_grid, w1_grid, loss_grid, 10, cmap='jet')
ctf = plt.contour(w0_grid, w1_grid, loss_grid, 10, colors='black', linewidths=0.5)
plt.clabel(ctf, inline_spacing=0.1, fmt='%.2f', fontsize=8)
plt.plot(w0, w1, 'o-', c='orangered', label='BGD')
plt.legend()
plt.tight_layout()
plt.show()
lrate = 0.001时
lrate = 0.01时 此时变化幅度过大，会来回横跳
4. 线性回归api
api：sklearn.linear_model
代码示例：
import numpy as np
import matplotlib.pyplot as plt
import sklearn.linear_model as lm
import sklearn.metrics as sm
# 数据样本
tx = np.array(
[0.5, 0.6, 0.8, 1.1, 1.4, 1.6, 1.7, 2.0, 2.5, 2.6, 3.0, 3.1, 3.3, 3.6, 3.8, 3.9, 4.0, 4.1, 4.2, 4.4, 4.5, 4.7])
ty = np.array(
[5.0, 5.5, 6.0, 6.8, 7.0, 11, 8.4, 9.0, 9.2, 10.5, 11, 11.1, 12.0, 13.3, 13.4, 15.7, 16, 14.8, 14.9, 33, 30.0, 35])
# 将一维数组转换为二维数组
tx_arr = tx.reshape(-1, 1)
# 创建模型
model = lm.LinearRegression()
"""
训练模型
输入为样本矩阵
输出为每个样本的结果,一维数组
"""
model.fit(tx_arr, ty)
# 预测结果
res_y = model.predict(tx_arr)
# print(res_y)
# 画图
plt.figure('linear regression', facecolor='lightgray')
plt.title('linear regression', fontsize=18)
plt.grid(linestyle=':')
plt.scatter(tx, ty, s=80, color='dodgerblue', label='simple points')
plt.plot(tx, res_y, color='orangered', label='LR')
plt.legend()
plt.show()
使用线性回归api训练模型只需要有 输入：样本矩阵(二维数组)、输出：每个样本的结果(一维数组)，不需要尝试对梯度下降的变化率进行调试。但是从上图可以看出样本中有三个特殊样本使得整体的星星回归模型出现一定的偏差，所以需要一定的优化。
5. 岭回归
普通的线性回归使用基于梯度下降的最下二乘法，在最下损失函数的前提下，寻找最优模型参数，但是无法对异常值进行识别，使得对模型参数造成影响。岭回归模型在模型迭代过程所依据的损失函数增加了正则项，以限制模型参数对异常值的匹配程度，提高模型的拟合精度
# 岭回归
# 创建模型 sklearn.linear_model.Ridge(正则强度, fit_intercept=是否训练截距, max_iter=最大训练次数)
model2 = lm.Ridge(10, fit_intercept=True, max_iter=1000)
# 训练模型
model2.fit(tx_arr, ty)
# 获取预测值
res_y2 = model2.predict(tx_arr)
print('R2得分: ', sm.r2_score(ty, res_y2))
# 画图
plt.plot(tx, res_y2, color='blue', label='Ridge')
plt.legend()
plt.show()
由上图可以看出，在存在特殊样本的情况下，岭回归对线性回归具有一定的修正作用
6. 线性回归模型的误差评估
代码示例：
# 误差评估
# 平均误差绝对值误差 即实际值-预测值的绝对值的平均值
print('平均误差绝对值误差: ', sm.mean_absolute_error(ty, res_y))
# 平均平方误差 即 平方差
print('平均平方值误差:', sm.mean_absolute_percentage_error(ty, res_y))
# 中位绝对值误差 即即实际值-预测值的绝对值的中位数
print('中位绝对值误差: ', sm.median_absolute_error(ty, res_y))
# R2的分，[0-1),分数越高，误差越小
print('R2得分: ', sm.r2_score(ty, res_y))
一般情况下，R2的得分即可判断模型的好坏，但是特定场景也需要对其他的误差进行分析。
下边使用去除特殊样本的数组测试上边的线性回归模型和岭回归模型
代码示例：
tx = np.array(
[0.5, 0.6, 0.8, 1.1, 1.4, 1.7, 2.0, 2.5, 2.6, 3.0, 3.1, 3.3, 3.6, 3.8, 3.9, 4.0, 4.1, 4.2])
ty = np.array(
[5.0, 5.5, 6.0, 6.8, 7.0, 8.4, 9.0, 9.2, 10.5, 11, 11.1, 12.0, 13.3, 13.4, 15.7, 16, 14.8, 14.9])
tx_t = tx.reshape(-1, 1)
res_t1 = model.predict(tx_t)
res_t2 = model2.predict(tx_t)
print('LR R2得分: ', sm.r2_score(ty, res_t1), '\n', 'Ridge R2得分: ', sm.r2_score(ty, res_t2))
输出结果为：
LR R2得分: 0.09601050799783617
Ridge R2得分: 0.42383925419855595
显然岭回归的预测的更加准确一些。
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
