DLSYS-HW0 记录
本作业要求实现一个基础的 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, yQ3 softmax_loss
题目已经给出了计算 loss 的公式,直接按照公式来就行:
这个函数声明为 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。算法伪代码如下(梯度的推导参照此处):
其代码实现如下:
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 函数的一个重要性质:,通过减去最大值将指数函数的值域控制在 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_W2Q6 softmax_regression_epoch_cpp
本题要求用 C++ 写 softmax 回归单轮训练代码。这个函数应该在由 X 和 y(以及大小 m、n、k)定义的数据上跑一轮,并且直接修改 theta。你的函数大概需要分配(然后再删除)一些辅助数组来存 logits 和梯度。
按照 Q3 的步骤直接实现就行。注意:
- 矩阵点乘用三重循环实现,矩阵
(i, j)位置处的表示法 - 矩阵转置对于 ,有
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;
}