Softmax Regression

作者：光火

邮箱：victor_b_zhang@163.com

$Softmax \quad Regression$ ，也称多项逻辑回归，是 $Logistic \quad Regression$ 的一般化形式，可用于解决多分类问题，是机器学习的经典模型之一

Softmax

首先，我们介绍 $Softmax$ 函数
$\Large{g(z)_i=\frac{e^{z_i}}{\sum_{j=1}^{k}e^{z_j}} = \frac{e^{z_i}}{e^{z_i} + \sum_{j\not = i} e^j} \in (0, 1)}$
多分类问题中，模型对于样本类别的打分值往往分散在整个实数域。对此，我们可以通过 $\max$ 函数获取得分最高的类别，并以此作为样本类别的预测值

$Softmax$ 函数也可以起到这样的作用，而且它还能额外给出样本属于各个类别的概率，并保证概率和为 $1$ ，因此时常被用在网络的最后一层

不过， $Softmax$ 函数也存在一定的缺陷，就比如它是 $over-parameterized$ 的（具有无穷最优解）。以及在编程实现时，由于 $e^x$ 增长过快，很容易发生溢出。我们会在下文结合实例，介绍相应的解决方案
机器学习中，我们一般将 $Softmax$ 函数写成如下形式：
$\Large P(y^{(i)}=j|x^{(i)},\theta) = \frac{exp(\theta_j^T x^{(i)})}{\sum_{k=1}^K exp(\theta_k^T x^{(i)})}$
该式就代表了将样本 $x^{(i)}$ 分为第 $j$ 类的概率，其中 $x^{(i)}$ 的上标 $(i)$ 表示的是 $x$ 在数据集抑或是当前 $mini-batch$ 中的位置。从几何角度分析， $\theta_j$ 就相当于一个超平面，负责给出样本属于第 $j$ 类的分数

接下来，我们就结合 $Softmax$ 函数，通过最大似然估计 ( $MLE$ ) ，推导出交叉熵损失函数（ $Cross-entropy$ ）

Cross-entropy

对于给定数据集 $X$ ，分布 $P(X)$ 和 $Q(X)$ 的交叉熵被定义为：
$H(P,Q) = E_p[-\log Q] = - \sum_{x \in X} P(x) \log Q(x)$
交叉熵的概念来自于信息论，它可以衡量分布 $P$ 与 $Q$ 之间的相似度，两者越相近，交叉熵就越小
我们对上式进行改写，让 $Q(x)$ 除以一个 $P(x)$ 再乘以一个 $P(x)$
$- \sum_{x \in X} P(x) \log Q(x) = - \sum_{x \in X}P(x) \log \frac{Q(x) }{P(x)} P(x) \\ - \sum_{x \in X} P(x) \log Q(x) = - \sum_{x \in X} P(x) \log P(x) - \sum_{x \in X} P(x)\log \frac{Q(x) }{P(x)}$
其中， $- \sum_{x \in X} P(x) \log P(x)$ 为概率分布 $P(x)$ 的熵，记作 $H(P)$ ，用于衡量该分布的混乱程度，并满足 $H(P) \geq 0$

$- \sum_{x \in X} P(x)\log \frac{Q(x) }{P(x)}$ 则被称为相对熵，抑或是 $KL$ 散度，记作 $D_{KL}(P||Q)$

于是， $CE = H(P) + D_{KL}(P||Q) \geq 0$
进一步，倘若 $P$ 是固定的，则最小化 $CE$ 就等价于最小化 $D_{KL}(P||Q)$ 。换句话说，最小化交叉熵就是让 $P$ 与 $Q$ 尽可能地接近。于是，我们立刻想到，如果让 $P$ 是 $label$ 的真实值， $Q$ 为我们的预测值，这不就是一个绝佳的损失函数吗？
明确了交叉熵的含义后，我们通过最大似然估计，建立它与 $Softmax$ 的关系
多分类问题中，我们认为类别的标签 $t$ 应当服从 $Multinoulli$ 分布。这是一个合理的假设，就像在回归问题中，我们一般假定输出 $y$ 服从高斯分布，并使用均方损失函数 ( $MSE$ ) 一样
具体来说，对于一个 $k$ 分类问题，在 $Multinoulli$ 分布下，有如下关系式成立（该式的含义就相当于做了 $K$ 次随机实验）
$\Large P(t|x) = \prod_{k=1}^{K} P(t_k = 1|x)^{t_k}$
其中， $t$ 代表类别的标签，是一个 $One-Hot$ 向量，即 $t = (0,..,1,..,0)^T$

$P(t_k = 1 |x)^{t_k}$ 中出现了两个 $t_k$ ，位于左侧的 $t_k$ 表示的是 $One-Hot$ 向量 $t$ 的第 $k$ 个分量，是一个随机变量。位于指数位置的 $t_k$ 则代表 $t_k$ 的具体取值（ $0$ 或者 $1$ ）

举个例子，假设有三个类别，且 $P(t = 1) = 0.2、P(t = 2) = 0.5, P(t = 3) = 0.3$

则 $P(t = (1, 0, 0)^T) =P(t_1 = 1)^1 P(t_2 = 1)^0 P(t_3 = 1)^0 = 0.2$

可见，当且仅当 $t_k = 1$ 时， $P(t_k = 1|x)^{t_k}$ 才是有效概率，否则就会因指数部分为 $0$ 而不起作用

这里之所以采用独热编码，是因为我们不希望人为地引入类别之间的数值与距离关系（类与类之间应当是相互独立的）
对于 $P(t_k = 1|x)$ ，我们指定它为 $Softmax$ 函数（注意：这是人为设置的，如果换用其他的函数，就无法导出交叉熵）这也对应了推导的一般流程，即随机实验应当采用什么分布，以及这个分布的参数应当用什么函数逼近（这里我们采用了 $Multinoulli$ 分布和 $Softmax$ 函数）

读者可以自行证明，如果是二分类问题，根据最大似然估计，利用 $Bernoulli$ 分布和 $Sigmoid$ 函数，就可以导出二分类下的 $Cross-entropy$
$Let \quad P(t_k = 1 | x) = \frac{exp(\theta^{(k)T}x)}{\sum_{j=1}^K exp(\theta^{(j)T}x)} = h_k(x) \\ \because P(t|x) = \prod_{k=1}^K P(t_k = 1|x)^{t_k} \\ \\ \therefore 对于样本间彼此独立的数据集 \quad \{(x^{(1)}, t^{(1)}),..,(x^{(N)},t^{(N)})\} \quad 有 \\ P(t^{(1)},...,t^{(N)}|X) = \prod_{n = 1}^N \prod_{k = 1}^K P(t_k^{(n)}=1|x^{(n)})^{t_k^{(n)}}$
最大化该似然函数，就等价于最小化其负对数似然函数
$E(\theta) = -\frac{1}{N} \sum_{n=1}^N \sum_{k=1}^K t_k^{(n)} \ln \frac{exp(\theta^{(k)T}x^{(n)})}{\sum_{j=1}^K exp(\theta^{(j)T}x^{(n)})} \\ E(\theta) = -\frac{1}{N} \sum_{n=1}^N \sum_{k=1}^K t_k^{(n)} \ln h_k^{(n)} \quad (*) \\ E(\theta) = -\frac{1}{N} \sum_{n = 1}^N \ln P(t_k^{(n)} = 1|x^{(n)})$
第二步中，我们用 $h_k^{(n)}$ 代表整个分式。第三步中，因为对于一个确定的样本 $x^{(n)}$ ， $t^{(n)}$ 中有且仅有一个 $t_k^{(n)}$ 为 $1$ ，所以内层的求和号就被消去了

此时，我们回看 $(*)$ 式，它就是交叉熵损失函数的形式，所有的 $t_k^{(n)}$ 就构成了真实标签的 $t^{(n)}$ 的 $One-Hot$ 向量，而 $\ln h_k$ 正是我们的预测值
总结： $H(P,Q) = - \sum_{x \in X} P(x) \log Q(x)$ ，对于 $Softmax \quad Regression$ 而言， $P(x)$ 即为真实标签（ $label$ ）的 $One-Hot$ 编码， $Q(x)$ 则为经过 $Softmax$ 函数作用后的概率值
最后，我们对损失函数求梯度，获取权重 $W$ 的更新方式：
$W = W - \alpha[Softmax(WX^T)-y^T]X \\ X = \begin{bmatrix} x_1^T \\ x_2^T \\ ... \\ x_m^T \end{bmatrix} \quad Y = \begin{bmatrix} y_1^T \\ y_2^T \\ ... \\ y_m^T \end{bmatrix} \\ 其中,x_i、y_i 均为列向量$
至此，我们已经基本熟悉 $Softmax\quad Regression$ ，可以尝试用它解决一些简单的问题了

MNIST Classification

这里，我们结合 MNIST，做一个手写体数字识别任务（下述代码只展示思路）

MNIST digits dataset is a widely used dataset for image classification in machine learning field. It contains 60,000 training examples and 10,000 testing examples. The digits have been size-normalized and centered in a fixed-size image. Each example is a 784 × 1 matrix, which is transformed from an original 28 × 28 gray-scale image. Digits in MNIST range from 0 to 9.

略去数据的读入部分

超参数设置，实际操作时应根据验证集进行调整

max_epoch = 10
batch_size = 100
learning_reate = 0.01

# L2 正则项的系数
factor = 0.01

初始化权重，并编写训练过程的代码

# weight initialization
w = np.random.randn(28 * 28, 10) * 0.001

loss_list = []
accuracy_list = []
display_frequency = 100

# training
for epoch in range(0, max_epoch):
    counts = train_size // batch_size
    for count in range(0, counts):
        # 获取下一个 batch 
        x, label = train_set.next_batch(batch_size)
        # 调用 softmax 分类器
        loss, gradient, prediction = softmax_classifier(w, x, label, factor)
        # 计算准确率
        label = np.argmax(label, axis=1)
        accuracy = sum(prediction == label) / float(len(label))
        # 收集数据，方便可视化展示
        loss_list.append(loss)
        accuracy_list.append(accuracy)
        # 更新权重
        w -= learning_rate * gradient
        if count % display_frequency == 0:
            print("Epoch [{}][{}]\t Batch [{}][{}]\t Training Loss {:.4f}\t Accuracy {:.4f}".format(
                epoch, max_epoch, count, counts, loss, accuracy))

实现 $Softmax$ 分类器
```
def softmax_classifier(w, x, label, factor):
    """
      Softmax Classifier

      Inputs have dimension D, there are C classes, a mini-batch have N examples.

      Inputs:
      - w: A numpy array of shape (D, C) containing weights.
      - x: A numpy array of shape (N, D) containing a mini-batch of data.
      - label: A numpy array of shape (N, C) containing labels, label[i] is a
        one-hot vector, label[i][j]=1 means i-th example belong to j-th class.
      - factor: regularization strength, which is a hyper-parameter.

      Returns:
      - loss: a single float number represents the average loss over the mini-batch.
      - gradient: shape (D, C), represents the gradient with respect to weights W.
      - prediction: shape (N,), prediction[i]=c means i-th example belong to c-th class.
    """
    grade = softmax(np.matmul(x, w))
    prediction = np.argmax(grade, axis=1)
    gradient = np.matmul(x.T, (grade - label)) + 2 * factor * w
    loss = - np.sum(label * np.log(grade + 1e-5)) + factor * np.square(w).sum()
    
    return loss, gradient, prediction
```
在书写代码时，我们应当尽可能地使用矩阵乘法来代替 for 循环，因为前者在效率上要优越得多

具体到本例， $x$ 是由 $N$ 个样本组成的矩阵，它的每一行就是 $mini-batch$ 中的一个样本。设样本的维度为 $D$ ，则 $x$ 的形状为 $N \times D$

$w$ 是权重矩阵，它的维度为 $D \times C$ ，其中 $C$ 为类别数。因为我们是对手写数字 $0 - 9$ 进行分类，所以 $C = 10$ 。 $w$ 的第 $i$ 列就对应了第 $i$ 类的各个参数

因此，我们通过 np.matmul 将 $x$ 与 $w$ 相乘，得到一个 $N\times C$ 的矩阵。该矩阵的每一行就对应了该样本在各个类别上的打分值。随后，我们将 $x * w$ 送入 $softmax$ 函数，将得分值转化为概率

np.argmax(grade, axis=1) 用于取出矩阵 $grade$ 每行的最大值所对应的索引，它就是我们对样本类别的预测值

$gradient$ 、 $loss$ 分别对应在这个 $mini-batch$ 上计算得到的梯度和损失函数值。我们引入 $L2$ 正则项，避免 $Softmax$ 的多解性。这里还有一个技巧，就是 np.log(grade + 1e-5)。为了避免 np.exp(x) 发生溢出，我们会在必要时，让 $x * w$ 整体减去其最大值（见 $Softmax$ 函数代码），这就会导致矩阵中有 $0$ 值出现，进而无法取对数，于是我们就引入一个微小的正数
```
def softmax(x):
    m = x.max()
    if m > 20:
        x -= x.max()
    return np.exp(x) / np.sum(np.exp(x), axis=1, keepdims=True)
```
关于 np.sum 的参数：

axis = 0 代表压缩行，即将每一列的元素相加，将矩阵压缩为一行

axis = 1 代表压缩列，即将每一行的元素相加，将矩阵压缩为一列

keepidims = True 结果保持原来的维数，默认为 False

进行测试

correct = 0
counts = test_size // batch_size

# Test process
for count in range(0, counts):
    data, label = test_set.next_batch(batch_size)
    
    # We only need prediction results in testing
    _, _, prediction = softmax_classifier(w, data , label, factor)
    label = np.argmax(label, axis=1)
    
    correct += sum(prediction == label)
    
accuracy = correct * 1.0 / test_size
print('Test Accuracy: ', accuracy)

绘制曲线

# training loss curvep
lt.figure()
plt.plot(loss_list, 'b--')
plt.xlabel('iteration')
plt.ylabel('loss')

# training accuracy curve
plt.figure()
plt.plot(accuracy_list, 'r--')
plt.xlabel('iteration')
plt.ylabel('accuracy')

Softmax Regression：理解与应用

Softmax Regression

Softmax

Cross-entropy

MNIST Classification