深度网络：反向传播与自动微分

上一篇训练 LeNet-5 时，loss.backward() 一次算出了几万个参数的梯度，训练循环再据此更新权重。当时没有展开的是：框架如何知道每个参数的梯度，又为什么能一次把它们全部算出来。

神经网络中的反向传播，可以看作反向模式自动微分在计算图上的具体应用。下面从链式法则出发推导它的计算过程，再用 Python 实现一个标量自动微分引擎，并用这个引擎训练一个小型多层感知机。理解正文只需要基础的 Python、偏导数和链式法则。

梯度计算#

训练模型，就是不断调整参数，使损失逐步减小。给定一批数据，模型算出标量损失 $L$。要更新某个参数 $\theta$，需要先知道它发生微小变化时，$L$ 会怎样变化。这个变化率是偏导数 $\partial L/\partial \theta$，梯度下降沿其反方向更新参数：

$$\theta \leftarrow \theta - \eta\,\frac{\partial L}{\partial \theta}.$$

因此，训练中的核心计算问题是：面对一个由大量简单运算复合而成、输入包含成千上万个参数、输出只有一个损失值的函数，怎样高效地求出 $L$ 对所有参数的偏导数。

符号微分与数值微分#

符号微分#

最直接的办法是写出导数的解析表达式。函数简单时这样做没有问题，但神经网络由许多层复合而成，表达式很快就会变得难以处理。

以 $L=(wx+b-y)^2$ 为例，对 $w$ 的导数是 $2(wx+b-y)\cdot x$，还可以直接写出。如果这个结果继续经过 tanh 和后续层的加权求和，展开后的表达式会越来越长，不同参数的导数中还会反复出现相同的子表达式。网络结构稍有变化，整套推导也要跟着修改。手动符号求导因而适合检查小例子，不适合承担模型训练中的梯度计算。

数值微分#

另一条路是回到导数的定义，用有限差分近似：

$$\frac{\partial L}{\partial \theta_i}\approx\frac{L(\theta_i+h)-L(\theta_i-h)}{2h}.$$

这里使用中心差分，其截断误差为 $O(h^2)$，通常比前向差分准确。数值微分不需要了解模型内部结构，只要能够计算 $L$，就能近似估计梯度，因此常用于检查反向传播的实现是否正确。

它的问题首先是计算量。每个参数都要分别扰动，并重新执行前向计算；若模型有 $n$ 个参数，计算全部梯度需要 $O(n)$ 次前向过程。对包含数百万参数的网络，这样的代价无法用于日常训练。

步长 $h$ 的选择也会影响结果。$h$ 过大时，有限差分对导数的近似不够准确；$h$ 过小时，两次接近的浮点数相减会损失有效位，舍入误差随之增大。因此，数值微分适合小规模的梯度检查，却不适合作为训练算法。

计算图与链式法则#

反向传播先把复合函数表示成计算图。每个节点保存一个输入或中间变量，边表示变量之间的计算依赖。前向传播按照依赖顺序计算各节点的值，并保留反向计算所需的中间结果，直到得到输出 $L$。

链式法则说明了梯度怎样沿图反向传递。设中间节点 $v$ 通过直接后继 $u=f(v,\dots)$ 影响最终输出 $L$，那么

$$\frac{\partial L}{\partial v}=\frac{\partial L}{\partial u}\cdot\frac{\partial u}{\partial v}.$$

其中，$\partial L/\partial u$ 是从后续计算传回的梯度，$\partial u/\partial v$ 是当前运算的局部导数。反向计算只需将二者相乘，再把结果传给 $v$。

如果 $v$ 同时通过多个后继影响 $L$，各条路径产生的梯度需要相加：

$$\frac{\partial L}{\partial v}=\sum_{u\in \text{后继}(v)}\frac{\partial L}{\partial u}\cdot\frac{\partial u}{\partial v}.$$

这就是后面的实现使用 += 累积梯度，而不直接用 = 覆盖梯度的原因。

标量计算图#

下面把 $L=(wx+b-y)^2$ 拆成计算图。取 $w=2,\ x=3,\ b=1,\ y=5$，并引入中间变量：

$$u=wx,\quad v=u+b,\quad e=v-y,\quad L=e^2.$$

前向计算依次得到 $u=6,\ v=7,\ e=2,\ L=4$。反向传播从 $\partial L/\partial L=1$ 开始。

对平方运算 $L=e^2$，局部导数为 $2e=4$，所以 $\partial L/\partial e=4$。减法 $e=v-y$ 对 $v$ 和 $y$ 的局部导数分别为 $1$ 和 $-1$，因此 $\partial L/\partial v=4$，$\partial L/\partial y=-4$。

加法 $v=u+b$ 将传入的梯度分别传给两个加数，得到 $\partial L/\partial u=4$ 和 $\partial L/\partial b=4$。最后，乘法 $u=wx$ 对 $w$、$x$ 的局部导数分别是 $x=3$、$w=2$，所以 $\partial L/\partial w=12$，$\partial L/\partial x=8$。

可以再用数值微分检查这个结果。将 $w$ 分别扰动 $\pm10^{-5}$，代入 $L(w)=(3w+1-5)^2$ 做中心差分，可得 $\partial L/\partial w\approx 12.000000$，与解析导数一致。反向传播经过一次前向和一次反向便得到 $w,x,b,y$ 的全部梯度；数值微分则要逐个扰动变量。

前向模式与反向模式#

自动微分主要有前向模式和反向模式两种。前向模式沿函数的计算方向传播导数：选定一个输入后，依次计算各中间量对该输入的导数，直到输出。若要分别得到输出对 $n$ 个输入的偏导数，通常需要传播 $n$ 次。反向模式从输出出发，沿计算图反向传播，一次得到一个标量输出对所有输入的梯度。

两种模式的适用场景取决于输入和输出的数量。前向模式更适合输入少、输出多的函数；反向模式更适合输入多、输出少的函数。训练神经网络时，待求导的参数很多，而损失通常是一个标量，因此框架普遍采用反向模式。后文所说的反向传播都指这一过程。

梯度消失与梯度爆炸#

沿一条路径反向传播时，梯度会连续乘上各层的局部导数。若这些因子的绝对值长期小于 1，梯度传到较早的层时可能已经接近 0，使参数几乎无法更新，这称为梯度消失。反之，若因子的绝对值持续大于 1，梯度可能迅速增大，造成数值不稳定或训练发散，这称为梯度爆炸。

激活函数的导数是这些因子之一。sigmoid 在输入绝对值较大时进入饱和区，导数接近 0，多层相乘后容易削弱梯度。ReLU 在正区间的导数为 1，不会因为这一项继续缩小梯度，这也是它更适合深层网络的原因之一。

循环网络还会沿时间反复使用同一组权重，因此同样的问题可能表现得更明显。初始化、归一化和梯度裁剪等方法可以缓解梯度消失或爆炸，但不影响本文要讨论的反向传播机制，后续文章再分别展开。

标量自动微分引擎#

接下来实现一个只处理标量的微型引擎，整体思路参考 Andrej Karpathy 的 micrograd。核心数据结构是 Value：它既保存当前数值，也记录生成该数值的运算和前驱节点。程序执行前向运算时，计算图也随之建立。

每个 Value 包含四项主要信息：data 保存前向数值，grad 保存累积到该节点的梯度并初始化为 0，_backward 是负责向前驱传播梯度的闭包，_prev 记录前驱节点。重载后的运算符在返回结果时，也会为结果记录前驱并设置相应的 _backward。

import math


class Value:
    """标量节点：记录数值、梯度、局部反向函数和前驱节点。"""

    def __init__(self, data, _prev=(), _op=""):
        self.data = data
        self.grad = 0.0
        self._backward = lambda: None
        self._prev = set(_prev)
        self._op = _op

    def __add__(self, other):
        other = other if isinstance(other, Value) else Value(other)
        out = Value(self.data + other.data, (self, other), "+")

        def _backward():
            self.grad += out.grad
            other.grad += out.grad

        out._backward = _backward
        return out

    def __mul__(self, other):
        other = other if isinstance(other, Value) else Value(other)
        out = Value(self.data * other.data, (self, other), "*")

        def _backward():
            self.grad += other.data * out.grad
            other.grad += self.data * out.grad

        out._backward = _backward
        return out

    def __pow__(self, k):
        assert isinstance(k, (int, float)), "只支持常数幂"
        out = Value(self.data**k, (self,), f"**{k}")

        def _backward():
            self.grad += k * self.data ** (k - 1) * out.grad

        out._backward = _backward
        return out

    def tanh(self):
        t = math.tanh(self.data)
        out = Value(t, (self,), "tanh")

        def _backward():
            self.grad += (1 - t * t) * out.grad

        out._backward = _backward
        return out

每个 _backward 都实现对应运算的局部求导规则。加法将 out.grad 传给两个加数；乘法还要乘上另一个因子的数值；幂运算使用 $\frac{d}{dx}x^k=kx^{k-1}$；tanh 的导数为 $1-\tanh^2$。梯度必须用 += 累积，因为同一个节点可能通过多条路径影响最终输出，各条路径的贡献不能互相覆盖。

减法、除法和取负可以由这些基本运算组合得到，不需要重复实现求导逻辑。

    def __neg__(self):
        return self * -1

    def __sub__(self, other):
        return self + (-other)

    def __truediv__(self, other):
        other = other if isinstance(other, Value) else Value(other)
        return self * other**-1

    def __radd__(self, other):
        return self + other

    def __rmul__(self, other):
        return self * other

拓扑排序与梯度回传#

当某个节点开始向前驱传播梯度时，它的所有后继都应处理完毕，否则该节点的 grad 可能还没有收齐。因此，节点需要按照计算图拓扑序的逆序执行。backward 先生成拓扑序，将输出节点的梯度置为 1，然后逆序调用各节点的 _backward。

    def backward(self):
        topo, visited = [], set()
        stack = [(self, False)]
        while stack:
            node, processed = stack.pop()
            if processed:
                topo.append(node)
                continue
            if node in visited:
                continue
            visited.add(node)
            stack.append((node, True))
            for child in node._prev:
                if child not in visited:
                    stack.append((child, False))

        self.grad = 1.0
        for node in reversed(topo):
            node._backward()

代码使用显式栈完成迭代式后序遍历，以免大型计算图触发 Python 的递归深度限制。每个节点分两次处理：第一次将它的前驱压栈，第二次在 processed 为真时把节点加入拓扑序。这样可以保证前驱位于当前节点之前。

用前面的标量例子检查这个实现：

w, x, b, y = Value(2.0), Value(3.0), Value(1.0), Value(5.0)
L = (w * x + b - y) ** 2
L.backward()
print(w.grad, x.grad, b.grad, y.grad)  # 12.0 8.0 4.0 -4.0

输出与前面的解析结果一致。

多层感知机#

基于 Value 可以继续定义多层感知机。每个神经元包含一组权重和一个偏置，先计算输入的加权和，再根据需要应用 tanh。多个神经元组成一层，多层依次连接；最后一层保留线性输出，作为分类分数。

import random


class Neuron:
    def __init__(self, n_in, nonlin=True):
        self.w = [Value(random.uniform(-1, 1)) for _ in range(n_in)]
        self.b = Value(0.0)
        self.nonlin = nonlin

    def __call__(self, x):
        z = sum((wi * xi for wi, xi in zip(self.w, x)), self.b)
        return z.tanh() if self.nonlin else z

    def parameters(self):
        return self.w + [self.b]


class Layer:
    def __init__(self, n_in, n_out, **kw):
        self.neurons = [Neuron(n_in, **kw) for _ in range(n_out)]

    def __call__(self, x):
        out = [n(x) for n in self.neurons]
        return out[0] if len(out) == 1 else out

    def parameters(self):
        return [p for n in self.neurons for p in n.parameters()]


class MLP:
    def __init__(self, n_in, n_outs):
        sizes = [n_in] + n_outs
        self.layers = [
            Layer(sizes[i], sizes[i + 1], nonlin=(i != len(n_outs) - 1))
            for i in range(len(n_outs))
        ]

    def __call__(self, x):
        for layer in self.layers:
            x = layer(x)
        return x

    def parameters(self):
        return [p for layer in self.layers for p in layer.parameters()]

这些模型类只描述前向计算，并不显式编写导数。Value 的运算符在执行 +、* 和 tanh 时记录计算依赖，backward 再依据这些记录计算梯度。模型结构发生变化时，只要所用运算已经定义了局部反向规则，就不需要重新手推整套导数。

双月牙分类#

下面用二维二分类任务验证这个引擎。make_moons 生成两簇相互交错的月牙形样本，线性模型无法将它们分开。实验将标签转换为 $\pm 1$，使用 2-16-16-1 的 MLP 输出分类分数，并以合页损失 $\max(0,\,1-y\cdot\text{score})$ 加 L2 正则作为训练目标。每轮训练依次完成前向计算、梯度清零、反向传播和参数更新。

from sklearn.datasets import make_moons

X, y = make_moons(n_samples=120, noise=0.1, random_state=42)
y = y * 2 - 1  # 标签变为 -1 / +1

model = MLP(2, [16, 16, 1])  # 337 个参数

for epoch in range(100):
    # 前向：对每个样本算分，合页损失加 L2 正则
    scores = [model([Value(a), Value(b)]) for a, b in X]
    margins = [1 + Value(-yi) * si for yi, si in zip(y, scores)]
    losses = [m if m.data > 0 else m * 0 for m in margins]
    data_loss = sum(losses) * (1.0 / len(losses))
    reg_loss = sum((p * p for p in model.parameters()), Value(0.0)) * 1e-4
    loss = data_loss + reg_loss

    # 反向：先清零再回传
    for p in model.parameters():
        p.grad = 0.0
    loss.backward()

    # 更新：学习率随轮次线性衰减
    lr = 1.0 - 0.9 * epoch / 100
    for p in model.parameters():
        p.data -= lr * p.grad

每轮反向传播前都要将参数的 grad 清零，因为 backward 会累积梯度。若省略这一步，上一轮的梯度会继续保留。实际运行中，损失从 0.47 降至约 0.012，训练准确率在第 40 轮左右达到 100%。将模型在二维平面上的输出画出来，可以看到决策边界沿两簇月牙之间弯曲。

至此，这个标量引擎已经完成了一次完整训练。它的运行速度很慢：每个数值都是独立的 Python 对象，每次运算还会创建新的图节点，因此无法直接扩展到 MNIST 等更大规模的数据。它保留了反向模式自动微分的核心过程，但省略了张量运算、设备执行和计算图内存管理等工程机制。

PyTorch autograd#

PyTorch autograd 采用相同的反向模式自动微分原理，不过基本数据结构是张量。带有 requires_grad=True 的张量参与运算时，PyTorch 会记录生成结果所需的操作。调用 loss.backward() 后，系统沿动态图反向计算，并将结果累积到叶子张量的 .grad 中。

因此，PyTorch 的训练循环也要在每次更新前清除旧梯度，通常通过 optimizer.zero_grad() 完成。与前面的标量实现相比，PyTorch 还提供批量张量运算、优化后的矩阵内核、设备调度和更完整的计算图管理，使同一套求导方法能够用于大规模模型。

前面的 MLP 可以用 PyTorch 的张量和模块重新表达，并扩展到上一篇使用的 MNIST。每张 $28\times28$ 的图像先展平为 784 维向量，再送入一个 784-128-64-10 的 MLP。其训练步骤仍然是清除梯度、执行前向计算、反向传播和更新参数。

import torch
from torch import nn
from torchvision import datasets, transforms

model = nn.Sequential(
    nn.Flatten(),
    nn.Linear(784, 128), nn.Tanh(),
    nn.Linear(128, 64), nn.Tanh(),
    nn.Linear(64, 10),
)
optimizer = torch.optim.SGD(model.parameters(), lr=0.1)
criterion = nn.CrossEntropyLoss()

for images, labels in train_loader:
    optimizer.zero_grad()      # 对应手写循环里的 p.grad = 0.0
    loss = criterion(model(images), labels)
    loss.backward()            # 对应 loss.backward()
    optimizer.step()           # 对应 p.data -= lr * p.grad

在标准归一化、批大小为 128 的设置下，这个 MLP 训练 5 轮后，在 MNIST 测试集上的准确率约为 97%。由于输入图像被直接展平，它没有利用像素的空间结构，效果低于上一篇的 LeNet-5。不过，两段代码中的反向计算遵循相同的链式法则：框架负责记录计算依赖，并在调用 loss.backward() 时按反向拓扑顺序累积梯度。

参考文献

1. A. Karpathy. micrograd: A tiny scalar-valued autograd engine. GitHub, 2020
2. D. E. Rumelhart, G. E. Hinton, and R. J. Williams. Learning representations by back-propagating errors. Nature, 1986. DOI: 10.1038/323533a0
3. A. G. Baydin, B. A. Pearlmutter, A. A. Radul, and J. M. Siskind. Automatic Differentiation in Machine Learning: a Survey. Journal of Machine Learning Research, 2018
4. I. Goodfellow, Y. Bengio, and A. Courville. Deep Learning, Chapter 6.5: Back-Propagation and Other Differentiation Algorithms. MIT Press, 2016
5. PyTorch Contributors. Autograd mechanics. PyTorch Documentation