2025, Nov 14 00:03

Почему ваш MNIST на NumPy застрял на 10% и как починить оси градиентов

Разбор типичной ошибки в двухслойной сети MNIST на NumPy: неверные оси градиентов и суммирование смещений. Пошаговое исправление, устойчивый softmax и код.

Обучение небольшого классификатора MNIST «с нуля» в NumPy часто срывается по банальной причине: градиенты считаются по неверным осям. Модель выглядит исправно, лосс меняется, но точность упрямо держится на уровне случайного угадывания. В двухслойной сети достаточно небольшой ошибки с осями в одном‑двух местах, чтобы испортить каждый шаг обновления. Ниже — минимальный, воспроизводимый путь от проблемы к рабочему исправлению.

Постановка задачи и симптомы

Сеть проста: скрытый слой с ReLU и выходной слой с softmax. Данные развёрнуты до формы (образцы, признаки). Прямой проход этой договорённости следует, а вот обратный смешивает ориентации матриц и «сжимает» градиенты по смещениям вдоль неверной оси. В результате обновления параметров используют тензоры неправильной формы; благодаря неявному broadcasting NumPy код работает, но веса перезаписываются так, что точность топчется на уровне 10–14%, а затем откатывается к случайности.

Код, в котором проявляется ошибка

Приведённый ниже фрагмент воспроизводит неверный поток градиентов. Названия переменных краткие, но логика полностью повторяет сломанную версию, включая проблемные оси и транспонирования.

import numpy as np
w_in = np.random.randn(784, 10)
b_in = np.random.randn(10)
w_out = np.random.randn(10, 10)
b_out = np.random.randn(10)
def act_relu(z, deriv=False):
    if deriv:
        return z > 0
    return np.maximum(z, 0)
def act_softmax(y):
    ex = np.exp(y - np.max(y))
    return ex / ex.sum(axis=1, keepdims=True)
def to_onehot(y):
    z = np.zeros((len(y), 10))
    for i in range(len(y)):
        z[i][y[i]] = 1
    return z
def forward_pass(x):
    z1 = x.dot(w_in) + b_in
    a1 = act_relu(z1)
    z2 = a1.dot(w_out) + b_out
    a2 = act_softmax(z2)
    return z1, a1, z2, a2
def accuracy(probs, y_true):
    preds = probs.argmax(axis=1)
    total = 0
    i = 0
    while i < len(preds):
        if preds[i] == y_true[i]:
            total += 1
        i += 1
    return (total / len(preds)) * 100
def apply_step(lr, d_w1, d_b1, d_w2, d_b2):
    global w_in, b_in, w_out, b_out
    w_in = w_in - lr * d_w1.T
    b_in = b_in - lr * d_b1
    w_out = w_out - lr * d_w2
    b_out = b_out - lr * d_b2
def fit_loop(x_data, y_data, epochs=40, lr=0.1):
    for ep in range(epochs):
        y_hot = to_onehot(y_data)
        m = len(y_hot)
        z1, a1, z2, a2 = forward_pass(x_data)
        d_z2 = a2 - y_hot
        d_w2 = (1 / m) * d_z2.T.dot(a1)
        d_b2 = (1 / m) * np.sum(d_z2, axis=1)
        d_z1 = w_out.dot(d_z2.T).T * act_relu(z1, deriv=True)
        d_w1 = (1 / m) * d_z1.T.dot(x_data)
        d_b1 = (1 / m) * np.sum(d_z1)
        apply_step(lr, d_w1, d_b1, d_w2, d_b2)
        probs = forward_pass(x_data)[-1]
        print("epoch", ep + 1, "acc", accuracy(probs, y_data))

Что именно идёт не так

Тензоры градиентов ориентированы не в соответствии с раскладкой данных. Для градиентов весов допущены лишние транспонирования, а оба градиента по смещениям суммируются по неправильной оси, из‑за чего их формы не совпадают с векторами параметров. Поскольку NumPy скрывает эти промахи за счёт broadcasting, код исполняется, но на каждом шаге обновления параметры калечатся. Поэтому модель застревает около 10%.

Исправление: правильные оси и формы в обратном проходе

Решение несложное: считать градиенты в согласованных формах матриц, суммировать смещения по оси классов и убрать компенсирующие транспонирования при обновлении. Помогает и численно устойчивый softmax. Ниже — исправленный проход с тем же интерфейсом, что и выше.

import numpy as np
w_in = np.random.randn(784, 10)
b_in = np.random.randn(10)
w_out = np.random.randn(10, 10)
b_out = np.random.randn(10)
def act_relu(z, deriv=False):
    if deriv:
        return z > 0
    return np.maximum(z, 0)
def softmax_stable(z):
    z = z - z.max(axis=1, keepdims=True)
    e = np.exp(z)
    return e / e.sum(axis=1, keepdims=True)
def to_onehot(y):
    z = np.zeros((len(y), 10))
    for i in range(len(y)):
        z[i][y[i]] = 1
    return z
def forward_pass(x):
    z1 = x @ w_in + b_in              # (m,784) @ (784,10) -> (m,10)
    a1 = act_relu(z1)
    z2 = a1 @ w_out + b_out           # (m,10)
    a2 = softmax_stable(z2)           # (m,10)
    return z1, a1, z2, a2
def accuracy(probs, y_true):
    preds = probs.argmax(axis=1)
    total = 0
    i = 0
    while i < len(preds):
        if preds[i] == y_true[i]:
            total += 1
        i += 1
    return (total / len(preds)) * 100
def fit_loop(x_data, y_data, epochs=40, lr=0.1):
    for ep in range(epochs):
        y_hot = to_onehot(y_data)
        m = x_data.shape[0]
        z1, a1, z2, a2 = forward_pass(x_data)
        d_z2 = a2 - y_hot              # (m,10)
        d_w_out = a1.T @ d_z2 / m      # (10,10)
        d_b_out = d_z2.sum(0) / m      # (10,)
        d_z1 = (d_z2 @ w_out.T) * (z1 > 0)   # (m,10)
        d_w_in = x_data.T @ d_z1 / m         # (784,10)
        d_b_in = d_z1.sum(0) / m             # (10,)
        # Шаг SGD с согласованными формами, без лишних транспонирований
        global w_in, b_in, w_out, b_out
        w_out -= lr * d_w_out
        b_out -= lr * d_b_out
        w_in  -= lr * d_w_in
        b_in  -= lr * d_b_in
        probs = a2
        print("epoch", ep + 1, "acc", accuracy(probs, y_data))

Обратите внимание, что в каждом произведении матриц транспонирование стоит у левого операнда — там, где ему естественно находиться для накопления градиента. Смещения суммируются по оси 0, то есть по размерности классов, образуя векторы, совпадающие по форме с параметрами. С этими поправками модель перестаёт портить свои веса и начинает учиться.

Устойчивый softmax

Стабильный softmax вычитает максимум по строке перед возведением в экспоненту, что предотвращает переполнение и улучшает численное поведение на протяжении обучения. Функция softmax_stable выше реализует ровно этот приём и выдаёт результат формы (m, 10).

Почему это важно усвоить

Backprop в NumPy суров: стоит промахнуться с одной осью, и NumPy покорно выполнит вещание и посчитает что‑то другое. Модель будет запускаться, выглядеть стабильной и зависнет на уровне случайности. Исправление осей согласует математику с раскладкой данных и убирает необходимость компенсирующих транспонирований в обновлении параметров.

Когда формы исправлены и используется численно безопасный softmax, та же двухслойная архитектура обучается как ожидается. С этими исправлениями плюс инициализацией Хе и меньшим шагом обучения, например 0,01, такая схема достигает около 93% на MNIST за 15–20 эпох. На практике результаты стабильны в этом диапазоне. В одном отчётном запуске с инициализацией Хе получено 91% при скорости обучения 0,25 за 200 эпох. В другом эксперименте — около 60% при скорости 0,5 за 300 эпох.

Практические выводы

Если ваш классификатор на NumPy держится на уровне 10–14% и не хочет улучшаться, в первую очередь проверьте формы градиентов. Убедитесь, что dW формируются как input.T @ delta и hidden.T @ delta для соответствующих слоёв, а оба db суммируются по оси 0, совпадая с числом единиц. Уберите транспонирования в обновлении, которые компенсировали ошибки «выше по потоку», и используйте устойчивый softmax. Когда оси согласованы, та же минимальная сеть становится крепкой базовой линией для MNIST.