自研MNIST数字识别神经网络故障排查求助

两个月前我开始自研用于MNIST数字识别的神经网络，但程序无法正常运行。此后我尝试了多种修复方案，查阅了大量视频与教程，虽已理解其工作原理，但仍未定位故障所在。恳请协助我解决这一难题。

以下是我的推导计算过程，完整可运行代码已上传至GitHub仓库：

// 前向传播
[] - 一维数组
[,] - 二维数组

a0[784] - 输入层
a1[16] - 隐藏层#1
a2[16] - 隐藏层#2
a3[10] - 输出层

w1[16,784] - 隐藏层#1权重
w2[16,16] - 隐藏层#2权重
w3[10,16] - 输出层权重

b1[16] - 隐藏层#1偏置
b2[16] - 隐藏层#2偏置
b3[10] - 输出层偏置

y[10] - 预期标签

Sigmoid(x) = 1 / (1 + e^-x)
Sigmoid导数(x) = x * (1 - x)

z1[i] = 求和(a0[j] * w1[i,j]) + b1[i]
a1[i] = Sigmoid(z1[i])
z2[i] = 求和(a1[j] * w2[i,j]) + b2[i]
a2[i] = Sigmoid(z2[i])
z3[i] = 求和(a2[j] * w3[i,j]) + b3[i]
a3[i] = Sigmoid(z3[i])
cost[i] = (a3[i] - y[i])^2

// 反向传播

alpha - 学习率
(x)` - x的导数

// 输出层
dcost/da3 = (cost)` = 2 * (a3 - y)
da3/dz3 = Sigmoid导数(a3)
dz3/dw3 = (z3)` = a2
dz3/db3 = (z3)` = 1 = da3/dz3
// 隐藏层#2
dz3/da2 = (z3)` = w3
da2/dz2 = Sigmoid导数(a2)
dz2/dw2 = (z2)` = a2
dz2/db2 = (z2)` = 1 = da2/dz2
// 隐藏层#1
dz2/da1 = (z2)` = w2
da1/dz1 = Sigmoid导数(a1)
dz1/dw1 = (z1)` = a1
dz1/db1 = (z1)` = 1 = da1/dz1

// 参数更新
dcost/dz3 = [10]
dcost/dw3 = [10,16]
dcost/dz3[i] = dcost/da3 * da3/dz3 = 2(a3[i] - y[i]) * Sigmoid导数(a3[i])
dcost/dw3[i,j] = dcost/dz3 * dz3/dw3 = dcost/dz3[i] * a2[j]
dcost/db3 = dcost/dz3 * dz3/b3 = dcost/dz3
w3[i,j] = w3[i,j] - alpha * dcost/dw3[i,j]
b3[i] = b3[i] - alpha * dcost/db3[i]

dcost/dz2 = [16]
dcost/dw2 = [16,16]
dcost/dz2[i] = dcost/dz3 * dz3/da2 * da2/dz2 = 求和(dcost/dz3[j] * w3[j,i]) * Sigmoid导数(a2[i])
dcost/dw2[i,j] = dcost/dz2 * dz2/dw2 = dcost/dz2[i] * a1[j]
dcost/db2 = dcost/dz2
w2[i,j] = w2[i,j] - alpha * dcost/dw2[i,j]
b2[i] = b2[i] - alpha * dcost/db2[i]

dcost/dz1 = [16]
dcost/dw1 = [16,784]
dcost/dz1[i] = dcost/dz2 * dz2/da1 * da1/dz1 = 求和(dcost/dz2[j] * w2[j,i]) * Sigmoid导数(a1[i])
dcost/dw1[i,j] = dcost/dz1 * dz1/dw1 = dcost/dz1[i] * a0[j]
dcost/db1 = dcost/dz1
w1[i,j] = w1[i,j] - alpha * dcost/dw1[i,j]
b1[i] = b1[i] - alpha * dcost/db1[i]

内容的提问来源于stack exchange，提问作者Esper