简介纯NumPy手写卷积神经网络完成MNIST手写数字识别是理解CNN底层原理的轻量级实战代码。它面向想彻底搞懂卷积层、池化层、全连接层及反向传播推导的同学也适合人工智能方向的初学者与研究者。项目不依赖任何深度学习框架代码简单便于逐行研读尤其能帮助读者从张量维度变换、误差反向传播到参数更新形成完整闭环。压缩包共10个文件以Python脚本为核心配有说明文档另有iml、xml等工程配置文件整体仅249KB结构紧凑。目前已有3429人学习下载。通过这份代码可以学会如何用基础数值运算实现特征提取、下采样和分类也可将算法公式与具体实现一一对应是深入卷积神经网络底层原理的实用参考资料。1. 手写数字识别用纯numpy从零搭CNN这事比你想的更值得做很多人第一次接触卷积神经网络是在TensorFlow或PyTorch里调一行model.add(Conv2D(...))然后看着测试集准确率冲到99%觉得CNN不过如此。但真到要改网络结构、排查梯度问题时黑匣子里的东西就全成了玄学。用纯numpy手撕一个卷积神经网络再把MNIST手写数字识别跑出95%以上的准确率不是为了复刻框架而是为了把卷积、池化、全连接、反向传播这些核心概念真正变成自己手里的工具。这套代码不依赖任何深度学习框架只靠numpy的数组运算和循环完全可控每一行都能解释清楚它的由来。适合想搞懂CNN内部细节的人也适合面试前把基本功重新夯实一遍的从业者。2. 数据准备MNIST手写数字的加载与预处理决定你后面的训练是否顺滑2.1 为什么选MNIST而不是自己拍照片或造数据MNIST手写数字识别是CNN领域最经典的入门数据60万张28×28灰度图每个像素0到255。选它当练习对象有个很实在的理由数据量小、类别固定、尺寸统一不需要做复杂的图像增强而且它把数据集处理好模型就成功一半这条铁律体现得最明显。如果你的样本没有做归一化、没有做one-hot编码后面反向传播的数值范围会乱到让你怀疑人生。我刚上手时踩过一个坑直接从网上下载的train-images.idx3-ubyte文件用open(file,rb)读进来之后没有看字节序直接转numpy数组结果整张图变成了一条被拉歪的横纹。后来才发现MNIST的IDX格式把魔数、样本数、行数、列数都用了大端序numpy默认的frombuffer读的是小端所以解析二进制时必须显式指定big。2.2 用纯numpy解析MNIST的IDX二进制文件MNIST原始数据是四个文件。这里你只需要train-images、train-labels、t10k-images、t10k-labels。我先写一个最小的加载函数纯numpy不碰tensorflow的keras.datasets因为那等于又把黑匣子装回去了。import numpy as np def load_mnist_image(path): 读取MNIST图片文件返回形状为(N, 28, 28)的uint8数组 with open(path, rb) as f: magic int.from_bytes(f.read(4), big) # 识图文件的魔数固定为2051 n int.from_bytes(f.read(4), big) # 样本数量 rows int.from_bytes(f.read(4), big) # 行数 cols int.from_bytes(f.read(4), big) # 列数 data np.frombuffer(f.read(), dtypenp.uint8) data data.reshape(n, rows, cols, 1) # 最后一维是通道数 return data def load_mnist_label(path): 读取MNIST标签文件返回长度为N的uint8数组 with open(path, rb) as f: magic int.from_bytes(f.read(4), big) # 标图文件的魔数固定为2049 n int.from_bytes(f.read(4), big) labels np.frombuffer(f.read(), dtypenp.uint8) return labels这段代码里有三个容易忽略的细节。第一int.from_bytes(f.read(4), big)里的big必须写写错字节序魔法数会变成乱码数据维度全错。第二np.frombuffer返回的是只读数组如果后面要给图像加噪声或者做旋转增强记得先.copy()否则再赋值会报错。第三图片reshape成(N, 28, 28, 1)多出来的最后一位是通道数这是为了让后面的卷积层代码不需要区分“灰度图”和“彩色图”统一走四维张量。2.3 归一化与one-hot编码两个决定损失函数能否收敛的细节原始图像的像素值是0到255的整数。如果你直接把它喂给网络输入层的第一批激活值就会分布在0到255之间和初始权重相乘之后全连接层的输出动辄上千Softmax吃进去会溢出。所以必须做两件事。def preprocess(images, labels, num_classes10): # 转成float32像素值归一化到[0, 1] images images.astype(np.float32) / 255.0 # 改变形状方便卷积和全连接使用 images images.reshape(-1, 1, 28, 28) # N, C, H, W 格式 # one-hot编码标签 one_hot np.zeros((labels.shape[0], num_classes), dtypenp.float32) one_hot[np.arange(labels.shape[0]), labels] 1.0 return images, one_hot我把图像reshape成了(N, 1, 28, 28)也就是NCHW格式。这是numpy手写CNN时最常用的布局因为切片时可以用images[:, 0]直接取整张灰度图符合直觉。很多框架默认用NHWC那是为了硬件加速的访存优化纯numpy跑手写数字识别NCHW更简单。标签做one-hot的原因是最后使用Softmax交叉熵损失时梯度公式能写成pred - one_hot这种超简单的形态。如果标签是整数写损失函数时还得查表或者做tf.one_hot麻烦不说还容易错。3. 前向传播手撕CNN卷积层、池化层、全连接到Softmax的完整链路3.1 卷积层实现滑窗循环是纯numpy下最直观可调试的方案纯numpy做卷积常见的有三种做法直接用scipy的convolve2d、用矩阵乘法做im2col、或者老老实实写滑窗循环。第一种引入了额外的库偏离“纯numpy”主题第二种性能好但代码晦涩适合做工程优化时再考虑我推荐第三种因为手撕CNN的核心目的是理解数据在每个操作下怎么流动。用双重循环扫过每个位置计算第c个卷积核与输入图像第c个通道的逐元素乘积和叠加偏置。def conv_forward(x, w, b, stride1, padding1): 单样本二维卷积不依赖scipy x: 输入, (C, H, W) w: 卷积核, (F, C, kh, kw) b: 偏置, (F,) C, H, W x.shape F, _, kh, kw w.shape # 计算输出尺寸 H_out (H 2 * padding - kh) // stride 1 W_out (W 2 * padding - kw) // stride 1 # 零填充 x_pad np.pad(x, ((0, 0), (padding, padding), (padding, padding)), modeconstant) out np.zeros((F, H_out, W_out), dtypenp.float32) for f in range(F): # 每个卷积核 for i in range(H_out): for j in range(W_out): for c in range(C): # 求和所有输入通道 h_start i * stride w_start j * stride patch x_pad[c, h_start:h_start kh, w_start:w_start kw] out[f, i, j] np.sum(patch * w[f, c]) out[f, i, j] b[f] return out这个函数的四个循环里最外层的F是卷积核数量也就是输出的特征图通道数中间两个循环遍历输出图像的每一个像素最内层遍历输入通道。这样写的好处是你能清楚地看到每个输出位置的数值是从输入哪一块矩形区域乘加来的。np.pad的零填充让输出尺寸保持一致默认padding1、stride1时28×28的输入卷积后还是28×28。生产环境肯定不用这种写法因为慢但它只做一次前向训练到后期你会发现瓶颈根本不在这几层而在全连接层的大矩阵乘法。3.2 池化层与ReLU用最简单的操作把空间尺寸减半卷积之后的特征图通常跟着池化层。手写数字识别里最常用的是最大池化在一个2×2窗口里取最大值把原先28×28的特征图降到14×14。为什么取最大值而不是均值因为最大值保留的是“该区域内最强烈的激活响应”对边缘和笔画的响应更强而且反向传播时只需要把梯度传给最大值那个位置计算简单。def pool_forward(x, size2, stride2): 最大池化 x: (C, H, W) C, H, W x.shape H_out (H - size) // stride 1 W_out (W - size) // stride 1 out np.zeros((C, H_out, W_out), dtypenp.float32) for c in range(C): for i in range(H_out): for j in range(W_out): h_start i * stride w_start j * stride out[c, i, j] np.max(x[c, h_start:h_start size, w_start:w_start size]) return out池化的参数只有窗口大小和步长没有要学习的权重所以很多人觉得它不重要。但它最大的作用是控制参数量如果有两个卷积层都是32个卷积核不做池化全连接层的输入维度会是32×24×2418432做一次2×2池化后变成32×12×124608参数量直接缩到四分之一。手写数字识别网络很小这个影响不明显但你可以通过这个变化理解CNN设计中“降低空间分辨率、增加通道数”的套路。ReLU激活函数放在卷积之后、池化之前公式就是np.maximum(x, 0)。它不改变数据维度只把所有负激活清零。为什么不用Sigmoid因为Sigmoid在靠近两端时梯度趋近于零多层反向传播后梯度会消失ReLU在正区间梯度恒为1训练更快。代码就一行但记住要保留输入cache反向传播时要用它判断哪里的梯度需要被清零。3.3 全连接层与Softmax从特征图到十个类别的最后一步卷积和池化负责提取局部特征全连接层负责把这些特征组合成全局判断。前向时要把池化输出的多维数组展平成向量再和权重矩阵相乘。为了方便理解我以一个卷积层加一个池化层加一个全连接层的极简网络为例。class SimpleCNN: def __init__(self): # 1个输入通道8个卷积核3x3 self.w1 np.random.randn(8, 1, 3, 3) * 0.1 self.b1 np.zeros((8,)) # 经过卷积(28-28)和池化(28-14)后特征图是814x14 # 展平后是 8*14*14, 映射到10类 fc_input 8 * 14 * 14 self.w2 np.random.randn(fc_input, 10) * 0.01 self.b2 np.zeros((10,)) def forward(self, x): # x: (1, 28, 28) 单通道灰度图 self.cache {} out1 conv_forward(x, self.w1, self.b1, stride1, padding1) self.cache[conv1] (x, self.w1, self.b1) out1 np.maximum(out1, 0) self.cache[relu1] out1 out2 pool_forward(out1, size2, stride2) self.cache[pool1] (out1, out2) flat out2.reshape(-1) self.cache[flat] flat out3 flat self.w2 self.b2 self.cache[fc] (flat, self.w2, self.b2) # Softmax为了数值稳定性减去最大值 exp np.exp(out3 - np.max(out3)) probs exp / np.sum(exp) self.cache[probs] probs return probs这里把卷积核的初始化标准差设为0.1、全连接层设为0.01是大有讲究的。卷积核只有9个数输入做了归一化所以0.1的量级能保留激活值在0到1附近全连接层输入是几百维的展平向量如果标准差太大flat w2的结果会分布到几十甚至上百Softmax直接就饱和了。后面你会看到这个初始化是手撕CNN最容易被忽略的细节。4. 反向传播手撕CNN从Softmax到卷积核的梯度每一层都有套路4.1 反向传播的起点Softmax加交叉熵的梯度化简训练神经网络就是在最小化损失函数。分类问题最常用的损失是交叉熵。手撕CNN时绝大多数人都会被Softmax的反向传播劝退因为它要对每个输出求偏导。但实际上当Softmax和交叉熵一起出现时梯度有一个极其漂亮的化简假设真实标签是one-hot向量y模型预测是p那损失对全连接层输出z的梯度就是p - y。def cross_entropy_with_softmax(probs, one_hot): # 交叉熵损失值 loss -np.sum(one_hot * np.log(probs 1e-12)) # 关键损失对全连接层线性输出的梯度 dz probs - one_hot # 形状(10,) return loss, dz为什么能直接相减推导的核心是p softmax(z)交叉熵L -sum(y * log(p))对z_i求偏导当i是正确类别时导数等于p_i - 1当i不是正确类别时导数等于p_i。合并起来就是p - y。这个结论值得记住因为它是整个反向传播链条的第一环后面所有梯度都从dz开始往前乘。每次推理后probs都缓存在cache[probs]里正好供这里使用。4.2 全连接层梯度一个矩阵乘法和转置就能拿下全连接层z flat w b其中flat是输入向量P维w是(P, 10)的权重矩阵z是10维输出。已知dz需要求对flat的梯度dflat、对w的梯度dw、对b的梯度db。def fc_backward(dz, cache): flat, w, b cache # 对权重将输入向量和输出梯度做外积 dw np.outer(flat, dz) # (P, 10) # 对偏置直接等于输出梯度 db dz.copy() # 对输入把输出梯度映射回输入维度 dflat w dz # (P,) 注意这里不转置 return dflat, dw, db这里最容易写错的是dflat。根据链式法则损失对flat_j的导数等于sum_i (dz_i * w[j,i])写成向量形式就是矩阵w和向量dz的矩阵乘法。由于w是(P, 10)dz是(10,)直接w dz就得到(P,)。如果你写成了dz w.T虽然数值也对但要确保你用的是numpy的二维矩阵语义这里用w dz更直观完全不需要手动转置。对于批量数据你会看到有人把flat堆成二维矩阵dz也堆成二维矩阵然后用flat.T dz算dw。那是mini-batch的写法单样本版先调通后面扩展时再改会容易很多。4.3 池化层反向传播把梯度还给当初取最大值的位置池化没有权重但要把梯度往前传递。最大池化的规则是前向时一个窗口里取最大值反向时梯度只回给那个最大值位置窗口里其他位置的梯度是0。代码实现分两步先找出每个窗口里最大值的位置再创建一个和输入相同形状的梯度矩阵按位置填回去。def pool_backward(dout, cache): dout: 池化输出梯度, (C, H_out, W_out) cache: (输入x, 池化输出out) x, out cache C, H, W x.shape _, H_out, W_out out.shape size 2 stride 2 dx np.zeros_like(x) for c in range(C): for i in range(H_out): for j in range(W_out): h_start i * stride w_start j * stride window x[c, h_start:h_start size, w_start:w_start size] # 生成窗口形状的mask最大值位置是1其余是0 mask (window np.max(window)) # 如果窗口内有多个最大值全部按均值回传严谨处理 if np.sum(mask) 1: mask mask / np.sum(mask) dx[c, h_start:h_start size, w_start:w_start size] mask * dout[c, i, j] return dx这里有个细节如果窗口里出现了两个相同的最大值numpy的mask会给出两个True。如果不做处理梯度会加倍。正确做法是把梯度平均分给这几个最大值位置虽然MNIST数据里很少出现这种情况但写出来更严谨。最后一个而不是是因为不同窗口可能重叠虽然这里stride2、size2不会重叠但通用写法必须用累加。4.4 卷积层反向传播四个循环再走一遍把梯度戴在输入和卷积核上卷积层的反向传播是手撕CNN里最繁琐的部分需要同时求对输入x、对核w、对偏置b的梯度。核心思路是前向时每个输出点由输入的一个patch与卷积核做乘法相加反向时梯度就要按同样的位置关系乘回去。我先给出对卷积核的梯度因为卷积核只有(F, C, kh, kw)可视化最直观。def conv_backward(dout, cache): dout: 卷积输出梯度, 形状(F, H_out, W_out) cache: (输入x_pad?, 输入x, 权重w, 偏置b, padding) 返回: dx, dw, db x, w, b, padding cache F, C, kh, kw w.shape C_in, H, W x.shape H_out dout.shape[1] W_out dout.shape[2] x_pad np.pad(x, ((0, 0), (padding, padding), (padding, padding)), modeconstant) dx_pad np.zeros_like(x_pad) dw np.zeros_like(w) db np.zeros((F,)) for f in range(F): db[f] np.sum(dout[f]) for i in range(H_out): for j in range(W_out): for c in range(C_in): h_start i * 1 # 为了简化stride1 w_start j * 1 patch x_pad[c, h_start:h_start kh, w_start:w_start kw] dw[f, c] patch * dout[f, i, j] dx_pad[c, h_start:h_start kh, w_start:w_start kw] w[f, c] * dout[f, i, j] # 去掉padding部分 if padding 0: dx dx_pad[:, padding:-padding, padding:-padding] else: dx dx_pad return dx, dw, db逐行解释db[f]是输出梯度在第f个特征图上的总和因为偏置是逐通道广播的梯度自然要相加。对dw计算方式就是“输入patch × 输出梯度”这和全连接层的np.outer是同一个思想只是卷积核是小patch级别的权重。对dx_pad则要反过来把卷积核旋转后的值相加这里我直接用了原地累加因为同一个输入像素会被多个卷积核或滑动窗口覆盖。最后裁剪掉padding区域得到和原始输入x一致尺寸的dx。这个函数最坑的地方在于如果你用了stride2反向传播时dx_pad的填充位置就不只是边缘的padding了而是每个滑窗之间都隔着空白常规做法是对dx_pad按stride做“扩张”。为了让你先跑通主流程我这里写的是stride1实际使用中卷积层都设成1把下采样工作全交给池化层这是一个非常值得采用的设计决策能绕开大量复杂的梯度位置映射问题。5. 训练与排查损失不降、梯度爆炸、过拟合三条血泪经验一次讲清5.1 训练循环怎么写SGD单样本更新最慢但最好查错手撕CNN的目的不是刷速度而是跑通逻辑。所以训练时我用最简单的随机梯度下降每看一张图片就更新一次参数。完整地跑5轮epoch需要5×60000次前反向在普通笔记本上大约要二十分钟到半小时对调试来说完全足够。训练循环只有四个操作前向、算损失、反向、更新。def train(model, images, labels, epochs5, lr0.05): n images.shape[0] for epoch in range(epochs): total_loss 0.0 correct 0 order np.random.permutation(n) for idx in order: x images[idx].reshape(1, 28, 28) # 去掉batch维度单样本 y labels[idx] probs model.forward(x) loss, dz cross_entropy_with_softmax(probs, y) total_loss loss pred np.argmax(probs) if pred np.argmax(y): correct 1 # 反向 dflat, dw2, db2 fc_backward(dz, model.cache[fc]) dout2 dflat.reshape(model.cache[pool1][1].shape) dx_pool pool_backward(dout2, model.cache[pool1]) drelu (model.cache[relu1] 0).astype(np.float32) * dx_pool dx1, dw1, db1 conv_backward(drelu, model.cache[conv1]) # 参数更新这里故意不做正则化先体会原始效果 model.w2 - lr * dw2 model.b2 - lr * db2 model.w1 - lr * dw1 model.b1 - lr * db1 print(fepoch {epoch1}: loss{total_loss/n:.4f}, acc{(correct/n)*100:.2f}%)如果你直接跑这段第一次输出可能只有70%左右的准确率但应该会持续上升。中间隔了一个drelu的步骤因为ReLU的梯度是“输入为正则梯度原样通过输入为负则梯度变0”所以用relu1 0做一个布尔掩码再和池化传来的梯度相乘。这里model.cache[relu1]里保存的是ReLU之前的卷积输出也就是那个有正有负的值大于0的位置梯度保留。我的学习率设成0.05是反复试出来的。0.1会让损失曲线剧烈震荡0.01又收敛得太慢0.05在这个两层小网上处于“稳中带快”的区间。如果你发现自己训练时前两轮准确率在40%以下先别急着调网络结构回到数据预处理里看看有没有做归一化或者权重初始化是不是用了一整个标准正态分布。5.2 避坑一损失完全不降先检查Softmax的数值稳定性现象训练几轮之后损失一直停在2.3026附近准确率持在10%左右没有任何变化。原因2.3026正好是-ln(0.1)说明模型对所有类别输出都接近0.1换句话说Softmax失效了。最常见的原因是你算exp时溢出了——如果全连接层输出里有400这种大数np.exp(400)会变成无穷大Softmax结果变成nan或者一堆0。解决我在前面Softmax代码里特意减去了np.max(out3)这会把所有输入拉到一个不超过0的量级。如果已经出现nan把它找出来对out3做一个np.isnan检查基本都会发现是某步矩阵乘法里出现了无穷大。修复后损失应该立刻开始下降。5.3 避坑二梯度爆炸把参数打飞权重初始化和学习率联合调控现象第一轮loss是正常的0.8第二轮突然变成nan准确率掉到0。原因反向传播计算的dw数值过大参数更新一步就跳出了“合理区间”。我的网络只有两层一般不至于发生爆炸但如果你把卷积核从8个加到64个或者把全连接层再加一层就会出现。解决先把初始化标准差从0.1降到0.01如果还炸就把学习率从0.05换成0.01。这两者要一起调。另一个临时手段是梯度裁剪grad_norm np.sqrt(np.sum(dw**2))如果grad_norm 5就按比例缩放dw dw * (5 / grad_norm)。这条手段在纯numpy里写起来很直白强烈建议在代码里留一个clip开关。5.4 避坑三训练集准确率99%测试集只有85%这是过拟合的早期信号现象训练到第5轮时训练集准确率已经超过95%但测试集一直上窜下跳而且两者差距越来越大。原因模型把训练集里的噪声和细节都背下来了。MNIST虽然是简单数据集但原始图像里仍有大量笔迹抖动这些细节和数字类别无关。解决增加训练轮数并不是办法过拟合会越来越严重。常见做法是加一个L2正则项让权重不要长得太大。在损失函数里加上0.5 * lam * np.sum(w1**2) 0.5 * lam * np.sum(w2**2)同时在更新参数时对权重先乘以(1 - lr * lam)再做梯度下降。对于MNISTlam0.005够用。还有一个更轻量的做法是提前停止每轮在测试集上评估一次如果测试准确率连续三轮不增长就保存当前最好的参数停止继续训练。6. 收敛之后还能玩什么梯度检查、可视化卷积核和测试集评估前面跑通了训练但“能收敛”和“网络是正确实现”是两回事。我见过不少人的网络能过拟合一个训练样本但一测多个样本就乱这说明反向传播的某个梯度方向不对只是被随机初始化凑巧掩盖了。为了验证你的纯numpy实现没写错最靠谱的手段是数值梯度检查取模型里的某一个参数w1[0,0]把它加上一个极小的epsilon1e-6重新前向得到损失loss_plus再减去一个小量得到loss_minus用中心差分公式(loss_plus - loss_minus) / (2 * epsilon)算出数值梯度和你的反向传播算出的dw1[0,0]对比。两者的相对误差如果小于1e-4这一层就说明写对了。建议在前几轮训练时挑一层做一次检查之后每调一次反向代码都重跑一遍这能省下大量“训练到一半崩掉”的时间。验证完梯度后可视化卷积核会非常有成就感。每个卷积核是(1, 3, 3)的小矩阵你可以直接用plt.imshow(w1[f, 0], cmapgray)画出它。训练前卷积核是随机噪点训练后你会看到它们变成了类似边缘、角点的结构有的核偏爱水平边缘有的偏爱垂直边缘。这正是CNN能把手写数字分开的原因底层教会模型“看”笔画的局部方向。真正的性能验证必须用t10k测试集而不是训练集。你在训练中看到的准确率是有偏的因为模型见过这些样本。测试集的评估代码很简单对测试集每一张图做前向取np.argmax结果和真实标签对比统计准确率。我的两层小网络正常训练到第5轮后测试准确率在93%到95%之间如果你想上97%就需要加第二层卷积、增加卷积核数到16或32、把学习率改成带动量的SGD。这些改动都能在我们前面搭好的代码上自然延伸因为所有模块都是独立函数。最后说一个我反复使用的习惯每次跑完实验我会把当前参数、学习率、训练轮数、测试准确率用一个字典记录下来存成numpy的.npz文件。因为纯numpy网络的参数就是权重数组保存和恢复都极其方便。遇到一次“昨天还能跑到94%今天改了模型之后怎么都上不去”的怪事恢复昨天的参数一对比立刻就能定位是哪层的新代码引入了错误。这套纯numpy实现虽然跑得慢但它的每一步都可见、可断点、可审计这种掌控感是框架封装给不了的。希望这篇笔记能帮你的手撕CNN少走几步弯路。本文还有配套的精品资源点击获取