DLSYS-HW0 记录

· Tech

本作业要求实现一个基础的 softmax 回归算法,再加一个简单的两层神经网络。
本笔记用于记录学习 CMU10-414 过程,代码仅支持通过测试样例。

HW0 需要完成 6 个函数的实现,其中 Question 1 仅用于熟悉实验环境。接下来,记录其余 5 个函数的实现。

Q2 parse_mnist

这个函数声明为 parse_mnist(image_filename, label_filename) ,返回值为一个元组 (X, y),其中 X 为存储训练集图像数据的二维数组,y 为存储训练集标签的一维数组。

该函数用于读取 MNIST 训练数据集,数据集及其格式需参见 lecun。MNIST 数据集主要分为训练集和测试集,在读取图像数据与标签时,需要特别注意其数据格式以正确读取。

具体实现中,需要使用 gzip 库读取数据文件,并且使用 numpy 创建符合题意的数据结构。

其代码实现如下:

def parse_mnist(image_filename, label_filename):
    image_file_object = gzip.open(image_filename, "rb")
    label_file_object = gzip.open(label_filename, "rb")
    
    image_file_object.read(16)
    label_file_object.read(8)
    
    image_data = image_file_object.read()
    label_data = label_file_object.read()
    
    image_file_object.close()
    label_file_object.close()
    
    X = np.frombuffer(image_data, dtype=np.uint8).reshape(-1, 28 * 28).astype(np.float32)
    X = X / 255.0
    y = np.frombuffer(label_data, dtype=np.uint8)
    return X, y

Q3 softmax_loss

题目已经给出了计算 loss 的公式,直接按照公式来就行:

ℓsoftmax(z,y)=log⁡∑i=1kexp⁡zi−zy\ell_{\mathrm{softmax}}(z, y) = \log\sum_{i=1}^k \exp z_i - z_y

这个函数声明为 softmax_loss(Z, y),返回值为样本上的平均softmax损失。接收参数 Z 为一个形状为 (batch_size, num_classes) 的二维 numpy 数组,里面装着每个类别的 logit 预测值;接收参数 y 为形状为 (batch_size, ) 的一维 numpy 数组,包含每个样本的真实标签。

具体实现中,因为我们最终需要返回的是样本上的平均 softmax 损失,而公式给出的是单个样本loss的计算方式,因此我们需要逐行处理样本并取均值。

其代码实现如下:

def softmax_loss(Z, y):
    rows = np.arange(Z.shape[0])
    return np.mean(np.log(np.sum(np.exp(Z), axis=1)) - Z[rows, y])

Q4 softmax_regression_epoch

本题要求在数据上对 softmax 回归跑一轮 SGD,使用步长 lr 和指定的批量大小。这个函数应该原地修改 theta 矩阵,并且要按 X 中的顺序遍历各个批次,不要随机打乱顺序。一个 epoch 就是把整个训练集都完整跑完一遍,而整个训练集按照 batch_size 的大小被切分成各个 batch,每轮 batch 按照 SGD 的要求更新一次 theta。算法伪代码如下(梯度的推导参照此处):

for each batch:Z=XbθP=softmax(Z)P[real label]−=1G=1BXbTPθ←θ−lr⋅G\boxed{ \begin{aligned} &\text{for each batch:}\\ &\qquad Z = X_b\theta \\ &\qquad P = \text{softmax}(Z) \\ &\qquad P[\text{real label}] -= 1 \\ &\qquad G = \frac{1}{B}X_b^T P \\ &\qquad \theta \leftarrow \theta - \text{lr}\cdot G \end{aligned} }

其代码实现如下:

def softmax_regression_epoch(X, y, theta, lr = 0.1, batch=100):
    for start in range(0, X.shape[0], batch):
        # 处理最后一个batch可能发生的截断
        end = min(start + batch, X.shape[0])                        
        X_batch = X[start:end]
        y_batch = y[start:end]
        
        Z = X_batch @ theta # 计算 logits
        # 对 logits 做 softmax
        P = np.exp(Z) / np.sum(np.exp(Z), axis=1, keepdims=True)    
        rows = np.arange(X_batch.shape[0])
        P[rows, y_batch] -= 1 # i.e. P - Y                                
        
        grad = X_batch.T @ P / X_batch.shape[0]
        theta -= lr * grad

需要注意的点:

  • 训练集可能不会被 batch_size 整除,如何计算更新次数(如上代码或上取整)
  • 在对计算矩阵 P 时(softmax 过程),需要注意分母的归一化后的维度
    • 若无 keepdims = True ,那么分母的 shape 是 (B,)
    • 加上 keepdims = True 后,分母的 shape 是 (B, 1),可以正常应用 NumPy broadcasting,即 (B, k) / (B, 1) -> (B, k)
  • 题目要求原地更新 theta 矩阵
    • 应该写成 theta -= lr * grad
    • 而非 theta = theta - lr * grad(这只是让局部变量 theta 指向新的数组)
  • 严格情况下,指数函数对大数非常敏感(直接导致溢出污染结果),因此需要借助 softmax 函数的一个重要性质:softmax(z)=softmax(z−c)softmax(z) = softmax(z - c),通过减去最大值将指数函数的值域控制在 0 ~ 1:
Z_shift = Z - np.max(Z, axis=1, keepdims=True)
P = np.exp(Z_shift) / np.sum(np.exp(Z_shift), axis=1, keepdims=True)

Q5 nn_epoch

本题要求对由权重 W1 和 W2 定义的两层神经网络(没有偏置项)跑一轮 SGD。实现上与 Q4 基本一致,参照题目所给公式即可。需要注意:

  • ReLU 函数用 np.maximum 实现,而非直接用 max(因为我们要对矩阵逐元素做 ReLU)

其代码实现如下:

def nn_epoch(X, y, W1, W2, lr = 0.1, batch=100):

    for start in range(0, X.shape[0], batch):
        end = min(start + batch, X.shape[0])
        X_batch = X[start:end]
        y_batch = y[start:end]
        
        # compute Z1, G2, G1
        Z1 = np.maximum(X_batch @ W1, 0)
        G2 = np.exp(Z1 @ W2)
        G2 /= np.sum(G2, axis=1, keepdims=True)
        rows = np.arange(X_batch.shape[0])
        G2[rows, y_batch] -= 1
        G1 = (Z1 > 0) * (G2 @ W2.T)
        
        grad_W1 = X_batch.T @ G1 / X_batch.shape[0]
        grad_W2 = Z1.T @ G2 / X_batch.shape[0]
        
        W1 -= lr * grad_W1
        W2 -= lr * grad_W2

Q6 softmax_regression_epoch_cpp

本题要求用 C++ 写 softmax 回归单轮训练代码。这个函数应该在由 X 和 y(以及大小 m、n、k)定义的数据上跑一轮,并且直接修改 theta。你的函数大概需要分配(然后再删除)一些辅助数组来存 logits 和梯度。

按照 Q3 的步骤直接实现就行。注意:

  • 矩阵点乘用三重循环实现,矩阵 (i, j) 位置处的表示法
  • 矩阵转置对于 m×n→n×mm \times n \to n \times m ,有res[j * n + i] = mtx[i * m + j]

其代码实现如下:

void softmax_regression_epoch_cpp(const float *X, const unsigned char *y,
								  float *theta, size_t m, size_t n, size_t k,
								  float lr, size_t batch)
{

    float *Z = new float[batch * k];
    float *grad = new float[n * k];
    
    for(size_t start = 0; start < m; start += batch) {
        size_t end = std::min(start + batch, m);
        size_t batch_size = end - start;
        
        // compute Z: b × k, logits
        for(size_t i = 0; i < batch_size; i ++) {
            for(size_t j = 0; j < k; j ++) {
                // Z[i, j]
                float sum = 0;
                for(size_t p = 0; p < n; p ++) {
                    sum += X[(start + i) * n + p] * theta[p * k + j];
                }
                Z[i * k + j] = sum;
            }
        }
        
        // softmax(Z)
        for(size_t i = 0; i < batch_size; i ++) {
            float sum = 0;
            for(size_t j = 0; j < k; j ++) {
                Z[i * k + j] = std::exp(Z[i * k + j]);
                sum += Z[i * k + j];
            }
            for(size_t j = 0; j < k; j ++) {
                Z[i * k + j] /= sum;
            }
        }
        
        // P - Y i.e. softmax(Z) - Y
        for(size_t i = 0; i < batch_size; i ++) {
            Z[i * k + y[start + i]] -= 1;
        }
        
        // compute grad: n × k
        for(size_t i = 0; i < n; i ++) {
            for(size_t j = 0; j < k; j ++) {
                // grad[i, j]
                float sum = 0;
                for(size_t p = 0; p < batch_size; p ++) {
                    sum += X[(start + p) * n + i] * Z[p * k + j];
                }
                grad[i * k + j] = sum / batch_size;
            }
        }
        
        // update theta
        for(size_t i = 0; i < n; i ++) {
            for(size_t j = 0; j < k; j ++) {
                theta[i * k + j] -= lr * grad[i * k + j];
            }
        }
    }
    
    delete[] Z;
    delete[] grad;
}

参考

  1. 官方资源

    https://dlsyscourse.org/assignments/

  2. 课程主页

    https://dlsyscourse.org/