线性回归可以预测“多少”的问题。然而除了这些，我们还对“哪种”的分类问题感兴趣。

1 softmax 回归的数学基础

1.1 独热编码

我们输入一张 2×2 像素的图像，其像素值分别为$x_1, x_2, x_3, x_4$。

假设这张图像属于类别鸡、猫、狗中的一个，我们很自然地想到将三个类别分别表示为$\{1,2,3\}$。如果类别之间有一定的自然顺序，如$\{\text{婴儿}, \text{儿童}, \text{青少年}, \text{中年人}, \text{老年人}\}$，那么这个问题可以转变为回归问题，这样编号是有意义的。

但是大多数情况下类别之间并无关联——你很难说猫和狗之间有什么自然顺序。

因此我们采用独热编码来为不同的类别编码。独热编码是一个向量，它的长度和类别数量相同。类别对应的分量设为 1，而其他的设为 0。例如在这个鸡、猫、狗的分类任务中，三种类别的独热编码分别为：

$$y \in \{(1,0,0), (0,1,0), (0,0,1)\}$$

即把四个输入$x_1, x_2, x_3, x_4$经过计算变为三个输出$o_1, o_2, o_3$。

$$\begin{pmatrix} o_1 \\ o_2 \\ o_3 \end{pmatrix} = \begin{pmatrix} w_{11} & w_{12} & w_{13} & w_{14} \\ w_{21} & w_{22} & w_{23} & w_{24} \\ w_{11} & w_{12} & w_{13} & w_{14} \\ \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \\ x_3 \\ x_4 \end{pmatrix} + \begin{pmatrix} b_1 \\ b_2 \\ b_3 \end{pmatrix}$$

像这种每个输入和每个输出都有连接的层，称为全连接层或稠密层。

1.2 softmax 函数

我们希望模型运算得到的结果$\boldsymbol{o}$为对应类别的概率。例如模型输出$\boldsymbol{o} = (0.1, 0.8, 0.1)$，我们就认为输入的是猫的图像。

但是不行，我们并不能直接把$\boldsymbol{o}$作为概率来使用。一方面，我们并未对输出分量进行归一化处理，这可能导致分量之和超过 1；另一方面，分量可能有负值，这违背了概率论公理。

于是提出 softmax 函数解决这些问题，softmax 函数被定义为：

$$\hat{\boldsymbol{y}} = {\rm softmax}({\boldsymbol{o}})$$

其中：

$$\hat{y}_i = \frac{\exp(o_i)}{\sum_j \exp(o_j)}$$

这样就解决了上面两个问题。

因此对于批量样本$\boldsymbol{X}$，softmax 回归为：

$$\begin{align*} \boldsymbol{O} &= \boldsymbol{X} \boldsymbol{W} + \boldsymbol{b} \\ \boldsymbol{\hat{Y}} &= {\rm softmax}(\boldsymbol{O}) \end{align*}$$

1.3 损失函数

上边提及 softmax 函数给出的向量$\hat{\boldsymbol{y}}$是“对给定输入$\boldsymbol{x}$的每个类的条件概率”。

例如，当给定输入$\boldsymbol{x}$条件下的判断为猫的概率为$\hat{y}_1 = P(y = \text{猫} | \boldsymbol{x})$

假设数据集$(\boldsymbol{X}, \boldsymbol{Y})$具有$n$个样本，每个样本都由特征向量$\boldsymbol{x}_i$和独热编码标签向量$\boldsymbol{y}_i$组成，每一对样本 $(\boldsymbol{x}_i, \boldsymbol{y}_i)$ 与其他样本对是相互独立的，并且它们都遵循同一个条件分布 $P(\boldsymbol{y}_i | \boldsymbol{x}_i)$。则：

$$P(\boldsymbol{Y} | \boldsymbol{X}) = P(\boldsymbol{y}_1, \cdots, \boldsymbol{y}_n | \boldsymbol{x}_1, \cdots, \boldsymbol{x}_n) = \prod_{i=1}^n P(\boldsymbol{y}_i | \boldsymbol{x}_i)$$

当我们在数据集中随机抽取一个样本，这个样本的标签是确定的，也就是一个独热向量。抽取到的图片只能是绝对的猫，而不是有多少多少概率的猫。而 softmax 模型并不知道真实标签，它给出的结果为一个概率向量，如$[0.1, 0.7, 0.3]$。我们的优化目标就是使得猫的概率趋近于 1，或者说使得概率向量趋于独热向量$[0, 1, 0]$。

总结就是，我们要找到一组参数$\boldsymbol{W}$，使得模型给出的概率$P(\boldsymbol{y}_i | \boldsymbol{x}_i)$尽可能大。即极大似然估计。

然而极大似然不符合直觉，因为损失函数自然是越小越好。模型输出的概率都大于 0 小于 1，因此我们可以考虑使用负对数。负对数函数单调递减，当输出概率为 0 时，其负对数为无穷大；而当输出概率为 1 时，其负对数为 0。这正符合我们对损失函数的直觉理解。

因此最大化$P(\boldsymbol{Y} | \boldsymbol{X})$就等价于最小化其负对数：

$$-\log P(\boldsymbol{Y} | \boldsymbol{X}) = -\log \prod_{i=1}^n P(\boldsymbol{y}_i | \boldsymbol{x}_i) = -\sum_{i=1}^n \log P(\boldsymbol{y}_i | \boldsymbol{x}_i)$$

接下来的问题是，$P(\boldsymbol{y}_i | \boldsymbol{x}_i)$是什么？我们知道$\boldsymbol{y}_i$是一个独热向量，怎么去求它的条件概率？这其实是一个** Multinoulli 分布**。

Multinoulli 分布

Multinoulli 分布是这样一种概率分布：

随机变量$x$只能取$k$个互斥的类别之一，即$x \in \{0, 1\}^k$，且$\sum_{i=1}^k x_i = 1$。$x$取到类别$x_i$的概率$p_i$，并满足$\sum_i p_i = 1$。

则其概率质量函数可表示为：

$$P(x_i = 1) = \prod_{i=1}^k {p_i}^{x_i}$$

在这里，独热向量$\boldsymbol{y}_i = [y_1, y_2, y_3, \cdots, y_q]$，所有的分量有且仅有一个为 1。而输入$\boldsymbol{x}$经过 softmax 运算后，输出的概率向量$\hat{\boldsymbol{y}}$中，每个分量代表取到该类的概率。将$\boldsymbol{y}$和$\hat{\boldsymbol{y}}$代入 Multinoulli 分布的概率质量公式中：

$$\log P(\boldsymbol{y}_i | \boldsymbol{x}_i) = \log \prod_{j=1}^q {\hat{y}_j}^{y_j} = \sum_{j=1}^q y_j \log \hat{y}_j$$

于是我们应该最小化的是：

$$-\log P(\boldsymbol{Y} | \boldsymbol{X}) = -\sum_{i=1}^n \log P(\boldsymbol{y}_i | \boldsymbol{x}_i) = \sum_{i=1}^n l(\boldsymbol{y}_i, \hat{\boldsymbol{y}}_i )$$

其中损失函数

$$l(\boldsymbol{y}, \hat{\boldsymbol{y}}) = -\sum_{j=1}^q y_j \log \hat{y}_j$$

通常称为交叉熵损失。

1.4 softmax 损失函数的导数

我们将 softmax 运算

$$\hat{y}_i = \frac{\exp(o_i)}{\sum_j \exp(o_j)}$$

代入交叉熵损失，即可得到：

$$\begin{align*} l(\boldsymbol{y}, \hat{\boldsymbol{y}}) &= -\sum_{j=1}^q y_j \log \hat{y}_j \\ &= -\sum_{j=1}^q y_j \log \frac{\exp(o_j)}{\sum_{k=1}^q \exp(o_k)} \\ &= \sum_{j=1}^q y_j \left( \log \sum_{k=1}^q \exp{o_k} - \log \exp(o_j) \right) \\ &= \sum_{j=1}^q y_j \log \sum_{k=1}^q \exp(o_k) - \sum_{j=1}^q y_j o_j \\ &\textcolor{red}{=} \log \sum_{k=1}^q \exp(o_k) - \sum_{j=1}^q y_j o_j \end{align*}$$

提示

红色等号能成立的原因是，在第一项中，除了$y_j$之外不含任何与$j$有关的变量，因此单独计算$\sum_{j=1}^q y_j$，显然其和为 1。

对上式求关于$o_j$的梯度：

$$\begin{align*} \frac{\partial}{\partial o_j} l(\boldsymbol{y}, \hat{\boldsymbol{y}}) &= \frac{\partial}{\partial o_j} \left(\log \sum_{k=1}^q \exp(o_k) - \sum_{j=1}^q y_j o_j \right) \\ &= \frac{\exp{o_j}}{\sum_{k=1}^q \exp(o_k)} - y_j \\ &= {\rm softmax}(\boldsymbol{o}_j) - y_j \end{align*}$$

提示

上述推导过程省略了对第一项$f(\boldsymbol{o}) = \log \sum_{k=1}^q \exp(o_k)$关于$o_j$求导的步骤，在这里补充之。

这是一个复合函数求导，令$z = \sum_{k=1}^q \exp(o_k)$，则$f(\boldsymbol{o}) = \log z$。

根据链式法则：

$$\begin{align*} \frac{\partial}{\partial o_j} f &= \frac{{\rm d} \log z}{{\rm d} z} \cdot \frac{\partial z}{\partial o_j} \\ &= \frac{1}{z} \cdot \frac{\partial}{\partial o_j} \sum_{k=1}^q \exp(o_k) \end{align*}$$

而只有$k = j$的那一项与$o_j$有关，其余是常数：

$$\frac{\partial}{\partial o_j} \sum_{k=1}^q \exp(o_k) = \frac{\partial}{\partial o_j} \exp(o_j) = \exp(o_j)$$

因此

$$\frac{\partial}{\partial o_j} f = \frac{\exp(o_j)}{z} = \frac{\exp{o_j}}{\sum_{k=1}^q \exp(o_k)} = {\rm softmax}(\boldsymbol{o}_j)$$

2 使用框架实现 softmax 回归

2.1 读取数据集

我们使用 Fashion-MNIST 数据集。该数据集由 10 个类别的图像组成，每个类别由训练集的 6000 张图像和测试集的 1000 张图像组成。因此 Fashion-MINST 数据集共包含训练集 60000 张图片和测试集 10000 张图片。Fashion-MINST 数据集中的每个图像均为 28×28 的单通道灰度图像。

我们使用torchvison模块来下载数据集，使用 PyTorch 的DataLoader模块来加载数据集。

from torch.utils import data
import torchvision
from torchvision import transforms

def load_data_fashion_mnist(batch_size):
    # 将 Fashion-MNIST 数据集转换为 Tensor
    trans = [transforms.ToTensor()]
    trans = transforms.Compose(trans)

    minst_train = torchvision.datasets.FashionMNIST(
        root='data/FashionMNIST',
        train=True,
        transform=trans,
        download=True
    )

    minst_test = torchvision.datasets.FashionMNIST(
        root='data/FashionMNIST',
        train=False,
        transform=trans,
        download=True
    )

    train_iter = data.DataLoader(
        minst_train,
        batch_size,
        shuffle=True
    )

    test_iter = data.DataLoader(
        minst_test,
        batch_size,
        shuffle=False
    )

    return train_iter, test_iter

2.2 训练

与线性神经网络类似，基本流程为定义神经网络、初始化参数、定义损失函数、定义优化器，然后开始训练。

需要注意的是，我们输入的是一个图片，就是二维张量。然而线性神经网络接收一维张量，因此要引入nn.Flatten()将二维张量展平。

def train(lr, num_epochs, train_iter, test_iter):
    # 定义神经网络
    net = nn.Sequential(
        nn.Flatten(),
        nn.Linear(784, 10)
    )

    # 仅为线性层初始化模型参数
    def init_weights(m):
        if type(m) == nn.Linear:
            nn.init.normal_(m.weight, std=0.01)

    net.apply(init_weights)

    # 交叉熵损失
    loss = nn.CrossEntropyLoss()

    # 小批量随机梯度下降
    optimizer = torch.optim.SGD(net.parameters(), lr)

    # 训练循环
    for epoch in range(num_epochs):
        # 训练损失之和、训练准确率之和，计算平均损失和准确率
        train_loss_sum = 0.0
        train_acc_sum = 0.0

        # 样本数量
        n = 6

        for X, y in train_iter:
            y_hat = net(X)
            l = loss(y_hat, y)

            # 梯度清零
            optimizer.zero_grad()
            l.backward()
            optimizer.step()

            train_loss_sum += l.item()
            train_acc_sum += (y_hat.argmax(dim=1) == y).sum().item()
            n += y.shape[0]

        # 测试集评估
        test_acc_sum, test_n = 0.0, 0
        with torch.no_grad():
            for X, y in test_iter:
                test_acc_sum += (net(X).argmax(dim=1) == y).sum().item()
                test_n += y.shape[0]

        # 打印训练进度
        print(f'epoch {epoch + 1}, '
              f'train loss {train_loss_sum / len(train_iter):.3f}, '
              f'train acc {train_acc_sum / n:.3f}, '
              f'test acc {test_acc_sum / test_n:.3f}')

    return net

2.3 预测

为了方便对比预测结果和真实结果，我们使用matplotlib包实现 GUI 展示。

import matplotlib.pyplot as plt

def predict(net, test_iter, n=6):
    # Fashion-MNIST标签映射
    labels = ['t-shirt', 'trouser', 'pullover', 'dress', 'coat',
              'sandal', 'shirt', 'sneaker', 'bag', 'ankle boot']

    for X, y in test_iter:
        images = X[0:n]
        labels_true = y[0:n]
        with torch.no_grad():
            predictions = net(images).argmax(axis=1)

        titles = [f'true: {labels[y]}\npred: {labels[pred]}'
                 for y, pred in zip(labels_true, predictions)]
        show_predict_result(images, titles)
        break

# 使用 GUI 展示预测结果
def show_predict_result(images, labels):
    _, figs = plt.subplots(2, 3, figsize=(12, 8))
    for f, img, lbl in zip(figs.flat, images, labels):
        f.imshow(img.reshape((28, 28)).detach().numpy())
        f.set_title(lbl)
        f.axes.get_xaxis().set_visible(False)
        f.axes.get_yaxis().set_visible(False)
    plt.show()

2.4 启动主函数

if __name__ == '__main__':
    # 设置超参数
    batch_size = 256
    lr = 0.1

    # 读取数据集
    train_iter, test_iter = load_data_fashion_mnist(batch_size)

    # 训练模型
    net = train(lr=0.1, num_epochs=10, train_iter=train_iter, test_iter=test_iter)

    # 预测
    predict(net, test_iter)

训练结果如下：

epoch 1, train loss 0.784, train acc 0.749, test acc 0.789
epoch 2, train loss 0.570, train acc 0.814, test acc 0.801
epoch 3, train loss 0.526, train acc 0.824, test acc 0.815
epoch 4, train loss 0.502, train acc 0.832, test acc 0.814
epoch 5, train loss 0.485, train acc 0.836, test acc 0.825
epoch 6, train loss 0.474, train acc 0.840, test acc 0.819
epoch 7, train loss 0.465, train acc 0.842, test acc 0.826
epoch 8, train loss 0.458, train acc 0.845, test acc 0.830
epoch 9, train loss 0.452, train acc 0.847, test acc 0.832
epoch 10, train loss 0.447, train acc 0.848, test acc 0.830

预测结果如下：