Backpropagation in Convolutional Neural Networks: Understanding Gradient Flow

Introduction

Backpropagation is the algorithm that makes deep learning possible. It computes the gradient of the loss function with respect to every parameter in the network, enabling gradient descent to update weights in the direction that reduces error. While the concept is straightforward — repeated application of the chain rule from calculus — the details matter enormously when implementing it for CNNs.

Convolution layers introduce unique challenges for backpropagation. The weight-sharing structure means that each parameter influences multiple output positions, and the gradient must be summed across all of them. Max pooling layers route gradients only to the winning neuron. Padding must be handled carefully to ensure that gradient dimensions match input dimensions. And the im2col-based implementation that makes forward convolution fast requires a corresponding col2im backward pass.

This post provides a complete treatment of backpropagation for CNNs. We derive the gradients for every layer type — convolution, pooling, ReLU, batch normalisation, and fully connected — implement each one in NumPy, and trace the full backward pass through a complete network.

For an overview of CNN architecture and the forward pass, see our guide on CNN architecture from the ground up.

The Chain Rule: A Refresher

Backpropagation is the chain rule applied to neural networks. Given a function composition f(g(x)), the derivative is df/dx = df/dg * dg/dx. In a neural network, the loss L is a composition of many functions (each layer). The gradient of the loss with respect to any parameter is the product of the local gradients along the path from that parameter to the loss.

This is implemented in two phases. The forward pass computes and caches all intermediate activations and any data needed for gradient computation (inputs, weights, batch statistics, etc.). The backward pass, starting from the loss gradient, propagates gradients backward through each layer, computing both the gradient with respect to the layer parameters (for updates) and the gradient with respect to the layer input (for continued propagation).

Two conventions govern the shapes. First, gradients always have the same shape as their corresponding variable. If a weight matrix has shape (C_in, C_out), then dw also has shape (C_in, C_out). Second, gradients sum over the batch dimension for parameters (since each parameter is shared across all batch samples), but preserve the batch dimension for activations (since each sample has its own loss contribution).

Starting Point: Loss Function Gradient

The backward pass begins with the gradient of the loss with respect to the network output. For a classification network with softmax cross-entropy loss, this gradient has a particularly elegant form. The softmax function converts logits into probabilities: p_j = exp(s_j) / sum_k exp(s_k). Cross-entropy loss measures the negative log probability of the correct class: L = -log(p_y).

The gradient of this combined function with respect to the logits is dL/ds_i = p_i - 1{i = y}. For the correct class, the gradient is (probability - 1). For all other classes, it is simply the probability. Intuitively: if the model assigns probability 0.9 to the correct class, the gradient is -0.1 for that class and +0.9 for the incorrect classes, pushing probability mass toward the correct class.

import numpy as np

def softmax_cross_entropy_loss(scores, labels):
    N = scores.shape[0]
    shifted = scores - np.max(scores, axis=1, keepdims=True)
    exp_scores = np.exp(shifted)
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    correct_log_probs = -np.log(probs[np.arange(N), labels] + 1e-8)
    loss = np.sum(correct_log_probs) / N
    grad = probs.copy()
    grad[np.arange(N), labels] -= 1
    grad /= N
    return loss, grad

Backward Pass Through Fully Connected Layers

The fully connected layer is the simplest backward pass and serves as a building block. Given an input x, weights w, and bias b, the forward pass is out = x @ w + b. With shapes: x is (N, C_in), w is (C_in, C_out), b is (C_out,), and out is (N, C_out).

Applying the chain rule: dw = x^T @ dout (each element w[i, j] affects all batch samples), db = sum_n dout[n, :] (sum over batch), and dx = dout @ w^T (each input element affects all output elements through w).

def fc_backward(dout, cache):
    x, w, b = cache
    N = x.shape[0]
    dw = x.T @ dout
    db = np.sum(dout, axis=0)
    dx = dout @ w.T
    return dx, dw, db

Backward Pass Through Convolution Layers

The convolution backward pass is the most challenging part of CNN backpropagation. We need to compute three gradients: dL/dx (gradient with respect to input), dL/dw (gradient with respect to weights), and dL/db (gradient with respect to biases).

Gradient with Respect to Bias

The bias gradient is the simplest. Each bias b_f is added to every spatial position in the f-th output channel. The gradient is therefore the sum of the upstream gradient dout over all spatial positions: db_f = sum_i sum_j dout[f, i, j].

Gradient with Respect to Weights

Each weight w[f, c, m, n] connects input channel c at offset (m, n) to output channel f. During the forward pass, this weight is multiplied by every input element x[c, i*stride + m, j*stride + n] across all spatial positions (i, j) where the kernel is applied. The gradient is the sum over all positions of the upstream gradient times the corresponding input element. This is itself a convolution: it is equivalent to a cross-correlation between the input and the upstream gradient.

Gradient with Respect to Input

Each input element x_pad[c, i, j] contributes to multiple output positions through different kernel offsets. The gradient is equivalent to a full convolution (with kernel flip) of the upstream gradient by the weight tensor. In practice, we implement it by iterating over output positions and adding the contribution of each output position to its receptive field of input positions.

def conv_backward(dout, cache):
    x, w, b, stride, padding = cache
    C_in, H, W = x.shape
    F, _, kH, kW = w.shape
    F, H_out, W_out = dout.shape
    if padding > 0:
        x_pad = np.pad(x, ((0,0), (padding,padding), (padding,padding)), mode='constant')
    else:
        x_pad = x
    dx = 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 f in range(F):
        for c in range(C_in):
            for i in range(H_out):
                for j in range(W_out):
                    h_start = i * stride
                    w_start = j * stride
                    dw[f, c] += dout[f, i, j] * x_pad[c, h_start:h_start+kH, w_start:w_start+kW]
    for c in range(C_in):
        for f in range(F):
            for i in range(H_out):
                for j in range(W_out):
                    h_start = i * stride
                    w_start = j * stride
                    dx[c, h_start:h_start+kH, w_start:w_start+kW] += w[f, c] * dout[f, i, j]
    if padding > 0:
        dx = dx[:, padding:-padding, padding:-padding]
    return dx, dw, db

Backward Pass Through Pooling Layers

Pooling layers have no learnable parameters, but they must still propagate gradients. During the forward pass, max pooling selects the maximum value within each pooling window and discards the rest. During the backward pass, the gradient must be routed exclusively to the winning element. The implementation uses an argmax mask — a binary matrix with 1 at the position of the maximum and 0 elsewhere.

Average pooling distributes the gradient evenly across all elements in the pooling window. If the window size is k x k, each element in the window receives dout / (k^2).

def max_pool_backward(dout, cache):
    x, pool_size, stride = cache
    C, H, W = x.shape
    C, H_out, W_out = dout.shape
    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+pool_size, w_start:w_start+pool_size]
                mask = (window == np.max(window))
                dx[c, h_start:h_start+pool_size, w_start:w_start+pool_size] += mask * dout[c, i, j]
    return dx

Backward Pass Through Batch Normalisation

Batch normalisation is the most complex standard layer to differentiate. The forward pass computes mu = mean(x, axis=0), var = variance(x, axis=0), x_hat = (x - mu) / sqrt(var + eps), and y = gamma * x_hat + beta. The backward pass must differentiate through the mean, variance, normalisation, and affine transformation.

The key insight is that x_hat affects the output through mu and var, so the gradient must account for this dependency. Let D be the total number of elements being normalised per channel (N * H * W for a convolutional layer). The gradient computation involves several intermediate terms.

def batchnorm_backward(dout, cache):
    x, gamma, beta, mu, var, eps, x_hat = cache
    N, C, H, W = dout.shape
    D = N * H * W
    dgamma = np.sum(dout * x_hat, axis=(0, 2, 3))
    dbeta = np.sum(dout, axis=(0, 2, 3))
    dx_hat = dout * gamma.reshape(1, C, 1, 1)
    dvar = np.sum(dx_hat * (x - mu.reshape(1, C, 1, 1)) * -0.5 * (var.reshape(1, C, 1, 1) + eps) ** -1.5, axis=(0, 2, 3))
    dmu = np.sum(dx_hat * -1.0 / np.sqrt(var.reshape(1, C, 1, 1) + eps), axis=(0, 2, 3)) + dvar * np.sum(-2 * (x - mu.reshape(1, C, 1, 1)), axis=(0, 2, 3)) / D
    dx = dx_hat / np.sqrt(var.reshape(1, C, 1, 1) + eps) + dvar.reshape(1, C, 1, 1) * 2 * (x - mu.reshape(1, C, 1, 1)) / D + dmu.reshape(1, C, 1, 1) / D
    return dx, dgamma, dbeta

Tracing the Full Backward Pass

Let us trace the backward pass through a complete three-layer CNN: Conv1 -> ReLU -> MaxPool -> Conv2 -> ReLU -> MaxPool -> FC -> Softmax. After the forward pass, we have cached all intermediate activations. The backward pass proceeds in reverse:

1. dl/dscores = softmax_cross_entropy_backward(scores, labels)
2. dl/dw3, dl/db3, dl/dpool2 = fc_backward(dl/dscores, (pool2_out, w3, b3))
3. dl/drelu2 = max_pool_backward(dl/dpool2, (conv2_out_relu, 2, 2))
4. dl/dconv2 = relu_backward(dl/drelu2, (conv2_out,))
5. dl/dw2, dl/db2, dl/dpool1 = conv_backward(dl/dconv2, (pool1_out, w2, b2, 1, 1))
6. dl/drelu1 = max_pool_backward(dl/dpool1, (conv1_out_relu, 2, 2))
7. dl/dconv1 = relu_backward(dl/drelu1, (conv1_out,))
8. dl/dw1, dl/db1, dl/dx = conv_backward(dl/dconv1, (x, w1, b1, 1, 1))

At each step, the shapes must match the corresponding forward shapes. A common debugging technique is to verify this shape consistency: every backward call should produce dx with the same shape as the x that was passed to the forward call.

Vanishing and Exploding Gradients

Even with correct backpropagation, deep CNNs can suffer from gradient instability. In deep networks, the gradient magnitude can decay exponentially (vanishing) or grow exponentially (exploding) as it propagates backward through many layers. The root cause is the multiplication of many Jacobian matrices during backpropagation.

Batch normalisation and residual connections are the primary solutions: batch norm keeps activations in a well-behaved range, and skip connections provide a gradient shortcut that bypasses the problematic multiplications. Proper weight initialisation is also critical — He initialisation (for ReLU networks) sets the variance of weights to 2 / (fan_in), ensuring that the variance of activations is preserved through the forward pass.

Gradient Checking: Verifying Your Implementation

When implementing custom CNN layers, gradient checking is an essential debugging tool. The idea is to approximate the gradient using the definition of the derivative: df/dx approx (f(x + h) - f(x - h)) / (2h). By comparing this numerical gradient against your analytic gradient, you can verify correctness.

• Use a small h (1e-5 typically works well)
• Check the relative error: |analytic - numeric| / max(|analytic|, |numeric|, 1e-8)
• The relative error should be less than 1e-7 for well-implemented layers
• Errors of 1e-5 may indicate minor issues (like missing the 1/N factor in the loss)

Gradient checking is computationally expensive (each parameter requires two forward passes), so it is only practical for small networks or individual layers. Most development workflows check each layer type separately before assembling the full network.