动手实现多元线性回归

线性回归 —— 根据气温预测冰淇淋销量

深度学习力，把收集到的真实数据叫做 样本（sample），样本里的输入叫做 特征（feature），样本里希望模型预测的值叫做 标签（label）。

训练过程就是调节函数内部的 参数（Parameter），也可以叫做权重（weight），来让预测值尽可能的接近 lablel。

麻雀虽小，五脏俱全，今天来动手实现一下多元线性回归。

对于特征和 Label 之间的 非线性问题，可以通过 构造高次特征 来解决。一般的做法是先从二次项开始，逐步增加，直到达到满意的效果。

准备

训练数据

假设训练数据有以下 7 条，每一条数据包含了当天的温度、冰淇淋的价格以及冰淇淋的销量，

多元线性回归公式

已知，$x_1$ 表示温度，$x_2$ 表示价格，$y$ 表示销量。$w_0$ 表示截距，$w_1$ 表示温度权重，$w_2$ 表示价格权重，则预测销量为：
$$\hat{y}=w_0+w_1{x}_1+w_2{x}_2$$

损失函数

损失函数使用 MSE，有 7 个数据，所以 $loss$ 等于 7 个样本的预测值减去 $lable$ 值的平方累加后除以 7，
$$loss=\frac{1}{7}\sum _{i=1}^7({\hat{y}}^i-y^i{)}^2$$
将多元线性回归公式代入，
$$loss=\frac{1}{7}\sum _{i=1}^7(w_0+w_1{x}_1^i+w_2{x}_2^i-y^i{)}^2$$

梯度下降

然后用梯度下降算法逐步来更新参数 $w_0$，$w_1$ 和 $w_2$，

Step1. 计算梯度，即损失函数对每个参数求偏导；对平方项求导用链式法则：$\frac{d}{dz}{z}^2=2z$。令误差项为 $e_i={\hat{y}}_i-y_i$。

• 对 $w_0$ 的偏导
  $$\begin{array}{rl}\frac{\partial \mathrm{loss}}{\partial w_0}& =\frac{1}{7}\sum _{i=1}^{7}2e_i\cdot \frac{\partial e_i}{\partial w_0}=\frac{2}{7}\sum _{i=1}^7{e}_i\cdot 1\\ & =\frac{2}{7}\sum _{i=1}^7(w_0+w_1{x}_1^i+w_2{x}_2^i-y^i)\end{array}$$

• 对 $w_1$ 的偏导
  $$\begin{array}{rl}\frac{\partial \mathrm{loss}}{\partial w_1}& =\frac{1}{7}\sum _{i=1}^{7}2e_i\cdot \frac{\partial e_i}{\partial w_1}=\frac{2}{7}\sum _{i=1}^7{e}_i\cdot x_1^i\\ & =\frac{2}{7}\sum _{i=1}^7(w_0+w_1{x}_1^i+w_2{x}_2^i-y_i)x_1^i\end{array}$$

• 对 $w_2$ 的偏导
  $$\frac{\partial \mathrm{loss}}{\partial w_2}=\frac{2}{7}\sum _{i=1}^7(w_0+w_1{x}_1^i+w_2{x}_2^i-y_i)x_2^i$$

Step2. 参数更新（梯度下降），设学习率为 $lr$。每次迭代 按负梯度方向更新参数：

$$\begin{array}{rl}{w}_0& \leftarrow w_0-lr\cdot \frac{\partial \mathrm{loss}}{\partial w_0},\\ w_1& \leftarrow w_1-lr\cdot \frac{\partial \mathrm{loss}}{\partial w_1},\\ w_2& \leftarrow w_2-lr\cdot \frac{\partial \mathrm{loss}}{\partial w_2}.\end{array}$$

敲 Code

Python实现

$X$ 中第一个元素是 [10, 3] 表示当天温度为 10°，冰淇淋价格为 3 元，对应的 $y$ 值为 60，表示当天卖了 60 只冰淇淋。

# Feature 数据
X = [[10, 3], [20, 3], [25, 3], [28, 2.5], [30, 2], [35, 2.5], [40, 2.5]]
y = [60, 85, 100, 120, 140, 145, 163]  # Label 数据# 初始化参数
w = [0.0, 0.0, 0.0]  # w0, w1, w2
lr = 0.0001  # 学习率
num_iterations = 10000  # 迭代次数# 梯度下降
for i in range(num_iterations):y_pred = [w[0] + w[1] * x[0] + w[2] * x[1] for x in X]  # 预测值loss = sum((y_pred[j] - y[j]) ** 2 for j in range(len(y))) / len(y)  # 计算损失# 计算梯度grad_w0 = 2 * sum(y_pred[j] - y[j] for j in range(len(y))) / len(y)grad_w1 = 2 * sum((y_pred[j] - y[j]) * X[j][0] for j in range(len(y))) / len(y)grad_w2 = 2 * sum((y_pred[j] - y[j]) * X[j][1] for j in range(len(y))) / len(y)# 更新参数w[0] -= lr * grad_w0w[1] -= lr * grad_w1w[2] -= lr * grad_w2# 打印损失if i % 100 == 0:print(f"Iteration {i}: Loss = {loss}")# 输出最终参数
print(f"Final parameters: w0 = {w[0]}, w1 = {w[1]}, w2 = {w[2]}")

可以尝试调整学习率，迭代次数。学习率太大的话，训练过程不会收敛，$loss$ 值可能会越来越大，直到程序出错。

向量化表示（更紧凑、便于实现）

向量化形式：更便于用 NumPy/框架实现，且一次性计算所有样本的梯度（更高效）。
向量化关键在于把带偏置项的设计矩阵 $X_{design}\in R_{n\times 3}$ 构造出来：

令 $n=7$，构造 输入矩阵 与 参数向量：

$$X_{design}=\left[\begin{array}{ccc}1& x_{11}& x_{21}\\ 1& x_{12}& x_{22}\\ ⋮& ⋮& ⋮\\ 1& x_{17}& x_{27}\end{array}\right]\in {\mathbb{R}}^{7\times 3},\mathbf{w}=\left[\begin{array}{c}{w}_0\\ w_1\\ w_2\end{array}\right],\mathbf{y}=\left[\begin{array}{c}{y}_1\\ ⋮\\ y_7\end{array}\right].$$

模型预测为 $\hat{\mathbf{y}}=X_{design}\mathbf{w}$。损失 $\mathrm{loss}=\frac{1}{n}\Vert X_{design}\mathbf{w}-\mathbf{y}{\Vert }^2$。梯度则为

$${\nabla }_{\mathbf{w}}\mathrm{loss}=\frac{2}{n}{X}_{design}^{\top }(X_{design}\mathbf{w}-\mathbf{y}).$$

参数更新：$\mathbf{w}\leftarrow \mathbf{w}-lr\cdot {\nabla }_{\mathbf{w}}\mathrm{loss}$。

在深度学习框架（PyTorch、TF）里通常直接用框架的自动微分，并使用向量化表达式或内置损失（例如 MSELoss）。

向量化表示的 Python 代码

import numpy as np# 原始数据
X = np.array([[10, 3], [20, 3], [25, 3], [28, 2.5], [30, 2], [35, 2.5], [40, 2.5]])
y = np.array([60, 85, 100, 120, 140, 145, 163], dtype=float)# 构造带常数项的一般设计矩阵 (n x 3)
n = X.shape[0]
X_design = np.hstack([np.ones((n, 1)), X])  # 第一列为 1，对应 w0# 初始化参数、超参
w = np.zeros(3, dtype=float)   # [w0, w1, w2]
lr = 0.0001
num_iterations = 10000# 训练（向量化梯度下降）
for i in range(num_iterations):y_pred = X_design.dot(w)            # (n,)error = y_pred - y                 # (n,)loss = np.mean(error ** 2)         # MSE# 向量化梯度： (2/n) * X^T (Xw - y)grad = (2.0 / n) * X_design.T.dot(error)  # (3,)# 参数更新w -= lr * gradif i % 1000 == 0:  # 每 1000 次打印一次print(f"Iteration {i}: Loss = {loss:.6f}  w = {w}")# 训练结束，输出最终结果
print("Final Loss =", loss)
print("Final parameters: w0 =", w[0], "w1 =", w[1], "w2 =", w[2])

使用 PyTorch 实现

import torch
import torch.nn as nn
import torch.optim as optim# X 每行是 [temperature, price]
X = torch.tensor([[10.0, 3.0],[20.0, 3.0],[25.0, 3.0],[28.0, 2.5],[30.0, 2.0],[35.0, 2.5],[40.0, 2.5],
], dtype=torch.float32)  # shape (7,2)y = torch.tensor([[60.0],[85.0],[100.0],[120.0],[140.0],[145.0],[163.0],
], dtype=torch.float32)      # shape (7,1)# ----------------------------
# 模型：线性回归 y = w0 + w1*x1 + w2*x2
model = nn.Linear(in_features=2, out_features=1, bias=True)# 可选：把参数初始化成 0
nn.init.zeros_(model.weight)
nn.init.zeros_(model.bias)# 损失和优化器
criterion = nn.MSELoss()            # 默认是平均 MSE（mean）
lr = 0.0001
optimizer = optim.SGD(model.parameters(), lr=lr)# 训练超参
num_iterations = 10000# ----------------------------
# 训练循环
for i in range(num_iterations):optimizer.zero_grad()          # 清零梯度y_pred = model(X)              # 前向 (7,1)loss = criterion(y_pred, y)    # MSEloss.backward()                # 反向optimizer.step()               # 参数更新if i % 1000 == 0:print(f"Iteration {i}: Loss = {loss.item():.6f}")# ----------------------------
# 输出最终参数（对应之前的 w0,w1,w2）
with torch.no_grad():weight = model.weight[0]       # shape (2,)bias = model.bias.item()       # scalarprint("\nFinal parameters (PyTorch):")
print(f"w0 (bias) = {bias:.6f}")
print(f"w1 = {weight[0].item():.6f}")
print(f"w2 = {weight[1].item():.6f}")# 打印每个样本的预测 vs 真实值
with torch.no_grad():preds = model(X).squeeze().numpy()  # shape (7,)X_np = X.numpy()y_np = y.squeeze().numpy()print("\nPredictions vs Ground Truth:")
for xi, p, t in zip(X_np, preds, y_np):print(f"X = {xi.tolist()}, pred = {p:.3f}, true = {t:.3f}")