Neural Networks and Backpropagation

Every modern AI system, from image classifiers to large language models, is built from the same ingredients: layers of simple units that compute weighted sums followed by non-linear functions, trained by gradient descent with gradients computed by backpropagation. Transformers add attention and convolutional networks add weight sharing, but the training engine is identical. If you understand this lesson thoroughly, including being able to write backprop by hand, you understand what every deep learning framework is doing under the hood, and you can debug it when it misbehaves.

Intuition: a committee of simple detectors

Picture a single neuron as a tiny detector with a dial for each input. It multiplies each input by its dial setting (a weight), adds them up with an offset (a bias), and then passes the total through a "squashing" or "gating" function. One detector can only draw a single straight boundary. But a layer of many detectors, each drawing a different boundary, followed by another layer that combines their verdicts, can carve out curved and patchy regions. Learning is the process of turning all the dials, a little at a time, in whichever direction reduces the mistakes. Backpropagation is the bookkeeping that tells every dial, even ones buried deep in the network, which way to turn and by how much.

The neuron

A neuron with input vector $x \in \mathbb{R}^d$, weights $w \in \mathbb{R}^d$, bias $b \in \mathbb{R}$ and activation function $\phi$ computes

$z = w^\top x + b, \qquad a = \phi(z).$

$z$ is the pre-activation and $a$ the activation. With $\phi$ the identity this is linear regression; with $\phi$ the sigmoid and a cross-entropy loss it is logistic regression. A neural network is many of these, wired together.

Activation functions

Saturation means the derivative is near zero over most of the input range: for $|z| > 5$, the sigmoid's slope is below 0.007. Gradients flowing back through many saturated units shrink toward zero (the vanishing-gradient problem, covered in the next lesson). ReLU (popularized for deep networks around 2010 to 2012) avoids saturation for positive inputs, which made deep networks much easier to train; a dead ReLU is a unit whose pre-activation is negative for every input, so it never receives gradient again. GELU (Hendrycks and Gimpel, 2016) weights the input by the probability that a standard normal variable is below it, giving a smooth, slightly non-monotone curve; GELU and related gated variants such as SwiGLU are the usual choice in Transformer feed-forward blocks.

Why non-linearity is essential. Stack two linear layers: $W_2(W_1 x + b_1) + b_2 = (W_2 W_1)x + (W_2 b_1 + b_2)$, which is just another linear map. Without a non-linear $\phi$, any depth of network collapses to a single linear model.

The multilayer perceptron

A multilayer perceptron (MLP, also called a fully connected feed-forward network) stacks layers. We use the row-vector batch convention that NumPy and PyTorch use: a batch of $n$ examples is a matrix $X \in \mathbb{R}^{n \times d}$, and a layer with $m$ units has weights $W \in \mathbb{R}^{d \times m}$ and bias $b \in \mathbb{R}^{m}$ (added to every row). A two-layer network (one hidden layer) for binary classification:

$Z_1 = X W_1 + b_1, \quad H = \phi(Z_1), \quad Z_2 = H W_2 + b_2, \quad P = \sigma(Z_2).$

Here $W_1 \in \mathbb{R}^{d\times m}$, $H \in \mathbb{R}^{n\times m}$ are the hidden activations, $W_2 \in \mathbb{R}^{m\times 1}$, and $P \in (0,1)^{n\times 1}$ are predicted probabilities. Computing these in order is the forward pass. The hidden layer learns a new representation of the input in which the classes are (closer to) linearly separable, and the output layer is a logistic regression on that representation. This is the core idea of deep learning: learned features.

Universal approximation, stated precisely

Theorem (Cybenko 1989 for sigmoidal activations; Hornik 1991; Leshno, Lin, Pinkus and Schocken 1993 in the general form below). Let $\phi: \mathbb{R}\to\mathbb{R}$ be continuous and not a polynomial. Let $K \subset \mathbb{R}^d$ be compact and $f: K \to \mathbb{R}$ continuous. Then for every $\varepsilon > 0$ there exist a width $N$ and parameters $a_j, b_j \in \mathbb{R}$, $w_j \in \mathbb{R}^d$ such that

$\sup_{x \in K} \left| f(x) - \sum_{j=1}^{N} a_j \, \phi(w_j^\top x + b_j) \right| < \varepsilon .$

In words: a network with one hidden layer and enough units can approximate any continuous function on a bounded, closed region as closely as you like, in the uniform (worst-case) sense. ReLU, tanh, sigmoid and GELU all qualify. (Leshno et al. also cover a broader class of activations that are not necessarily continuous; the statement above is a simplified version.)

What the theorem does not say:

• Nothing about how large $N$ is. For some functions the required width grows exponentially with $d$. Depth can be exponentially more efficient than width for certain function classes, which is one theoretical motivation for deep networks.
• Nothing about learning. Existence of good weights does not mean gradient descent will find them.
• Nothing about generalization. Fitting a function on the training set says nothing about unseen data.
• Only compact domains. Outside the region seen in training, behaviour is unconstrained, which is one reason networks extrapolate poorly.

The loss

For binary classification with labels $y_i \in \{0,1\}$ and predicted probabilities $p_i$, the binary cross-entropy (negative log-likelihood of a Bernoulli model) is

$\mathcal{L} = -\frac{1}{n}\sum_{i=1}^{n}\big[y_i \log p_i + (1-y_i)\log(1-p_i)\big].$

For regression the usual choice is mean squared error, and for multiclass classification softmax with categorical cross-entropy. Training means finding parameters $\theta = (W_1, b_1, W_2, b_2)$ that minimize $\mathcal{L}$ by gradient descent: $\theta \leftarrow \theta - \eta \nabla_\theta \mathcal{L}$ with learning rate $\eta$.

Backpropagation: the derivation

Backpropagation is the chain rule applied systematically from the loss backwards, reusing intermediate results so that computing all gradients costs about the same as two or three forward passes. (The method has several independent origins; it was popularized for neural networks by Rumelhart, Hinton and Williams in 1986.) Write $\delta$ for the gradient of the loss with respect to a pre-activation.

Step 1: output layer. For a single example, with $p = \sigma(z_2)$ and BCE loss $\ell = -[y\log p + (1-y)\log(1-p)]$:

$\frac{\partial \ell}{\partial p} = -\frac{y}{p} + \frac{1-y}{1-p}, \qquad \frac{\partial p}{\partial z_2} = p(1-p).$

Multiplying, $\frac{\partial \ell}{\partial z_2} = -y(1-p) + (1-y)p = p - y$. This clean result is why sigmoid is paired with cross-entropy (and softmax with categorical cross-entropy, which gives the same form). For the batch mean loss:

$\Delta_2 = \frac{\partial \mathcal{L}}{\partial Z_2} = \frac{1}{n}(P - Y) \in \mathbb{R}^{n\times 1}.$

Step 2: output weights. Since $Z_2 = H W_2 + b_2$, each entry of $W_2$ multiplies hidden activations, so

$\frac{\partial \mathcal{L}}{\partial W_2} = H^\top \Delta_2, \qquad \frac{\partial \mathcal{L}}{\partial b_2} = \sum_{i} (\Delta_2)_{i} .$

Step 3: back into the hidden layer. The gradient with respect to the hidden activations is $\frac{\partial \mathcal{L}}{\partial H} = \Delta_2 W_2^\top$. Passing through the element-wise activation multiplies by its derivative:

$\Delta_1 = \frac{\partial \mathcal{L}}{\partial Z_1} = (\Delta_2 W_2^\top) \odot \phi'(Z_1),$

where $\odot$ is element-wise multiplication.

Step 4: first-layer weights.

$\frac{\partial \mathcal{L}}{\partial W_1} = X^\top \Delta_1, \qquad \frac{\partial \mathcal{L}}{\partial b_1} = \sum_i (\Delta_1)_{i,:} .$

The pattern generalizes to any depth: each layer receives an upstream gradient, multiplies by its local derivative, produces weight gradients as (layer input)$^\top \times \delta$, and passes $\delta W^\top$ further back. A quick sanity check: every gradient has the same shape as the parameter it belongs to.

Worked example with numbers

A tiny regression network: input $x = (1, 2)$, two ReLU hidden units, one linear output, loss $\ell = \frac{1}{2}(\hat y - y)^2$ with target $y = 1$.

Parameters: $W_1 = \begin{pmatrix} 0.5 & -1 \\ 0.25 & 0.5 \end{pmatrix}$, $b_1 = (0, 0.5)$, $w_2 = (2, -1)$, $b_2 = 0.5$.

Forward.

• $z_1 = x W_1 + b_1 = (1\cdot0.5 + 2\cdot0.25,\; 1\cdot(-1) + 2\cdot0.5) + (0, 0.5) = (1.0, 0.0) + (0, 0.5) = (1.0, 0.5)$.
• $h = \text{ReLU}(z_1) = (1.0, 0.5)$.
• $\hat y = h \cdot w_2 + b_2 = 2.0 - 0.5 + 0.5 = 2.0$.
• $\ell = \frac{1}{2}(2.0 - 1)^2 = 0.5$.

Backward.

• $\partial \ell/\partial \hat y = \hat y - y = 1.0$.
• $\partial \ell/\partial w_2 = 1.0 \cdot h = (1.0, 0.5)$; $\partial \ell/\partial b_2 = 1.0$.
• $\partial \ell/\partial h = 1.0 \cdot w_2 = (2, -1)$.
• Both pre-activations are positive, so ReLU passes the gradient through: $\delta_1 = (2, -1)$.
• $\partial \ell/\partial W_1 = x^\top \delta_1 = \begin{pmatrix} 1\cdot 2 & 1\cdot(-1) \\ 2\cdot 2 & 2\cdot(-1)\end{pmatrix} = \begin{pmatrix} 2 & -1 \\ 4 & -2 \end{pmatrix}$; $\partial \ell/\partial b_1 = (2, -1)$.

Gradient check. Nudge $W_1[1,0]$ (row 2, column 1) from 0.25 to $0.25 + 10^{-6}$: the first hidden pre-activation grows by $x_2 \times 10^{-6} = 2\times10^{-6}$, $\hat y$ by $4\times10^{-6}$, and the loss by about $4\times 10^{-6}$, so the numerical derivative is 4, matching the analytic value.

Update with $\eta = 0.05$: $W_1 \to \begin{pmatrix} 0.4 & -0.95 \\ 0.05 & 0.6\end{pmatrix}$, $b_1 \to (-0.1, 0.55)$, $w_2 \to (1.95, -1.025)$, $b_2 \to 0.45$. New forward pass: $z_1 = (0.4, 0.8)$, $\hat y = 1.95(0.4) - 1.025(0.8) + 0.45 = 0.41$, loss $= \frac12(0.59)^2 \approx 0.174$, down from 0.5. (With $\eta = 0.1$ this example overshoots to $\hat y \approx -0.75$ and loss about 1.54: a first taste of why the learning rate matters.)

A from-scratch NumPy MLP

The code below trains a 2-16-1 network with tanh hidden units and a sigmoid output on an XOR-style problem (label 1 when both coordinates have the same sign), which no linear model can solve.

import numpy as np
rng = np.random.default_rng(0)
X = rng.uniform(-1, 1, size=(400, 2))
y = (X[:, 0] * X[:, 1] > 0).astype(float).reshape(-1, 1)

def sigmoid(z): return 1 / (1 + np.exp(-z))

d, m = 2, 16
W1 = rng.normal(0, np.sqrt(1 / d), (d, m)); b1 = np.zeros(m)
W2 = rng.normal(0, np.sqrt(1 / m), (m, 1)); b2 = np.zeros(1)
lr, n = 0.5, len(X)

for step in range(5001):
    Z1 = X @ W1 + b1; H = np.tanh(Z1)               # forward
    Z2 = H @ W2 + b2; P = sigmoid(Z2)
    loss = -np.mean(y * np.log(P + 1e-12) + (1 - y) * np.log(1 - P + 1e-12))
    dZ2 = (P - y) / n                                # backward: Delta_2
    dW2 = H.T @ dZ2; db2 = dZ2.sum(0)
    dZ1 = (dZ2 @ W2.T) * (1 - H**2)                  # Delta_1, tanh' = 1 - tanh^2
    dW1 = X.T @ dZ1; db1 = dZ1.sum(0)
    W1 -= lr * dW1; b1 -= lr * db1; W2 -= lr * dW2; b2 -= lr * db2
    if step % 1000 == 0:
        print(f"step {step:4d} loss {loss:.4f} acc {((P > 0.5) == y).mean():.3f}")

Each line of the backward pass is one equation from the derivation. Running it, the loss falls from about 0.70 (chance) to below 0.05 and training accuracy reaches about 99%. Weights are initialized randomly with variance $1/\text{fan-in}$; if all weights started equal, every hidden unit would compute the same function and receive the same gradient forever (the symmetry problem). Initialization is covered in depth in the next lesson.

The same model in PyTorch

PyTorch builds a computation graph as the forward pass runs and computes all gradients with loss.backward() (reverse-mode automatic differentiation, which is backpropagation generalized to any differentiable program).

import torch
from torch import nn

torch.manual_seed(0)
X = torch.rand(400, 2) * 2 - 1
y = (X[:, 0] * X[:, 1] > 0).float().unsqueeze(1)

model = nn.Sequential(nn.Linear(2, 16), nn.Tanh(), nn.Linear(16, 1))
loss_fn = nn.BCEWithLogitsLoss()          # sigmoid + BCE, numerically stable
opt = torch.optim.SGD(model.parameters(), lr=0.5)

for step in range(5001):
    logits = model(X)                     # forward pass
    loss = loss_fn(logits, y)
    opt.zero_grad()                       # clear old gradients
    loss.backward()                       # backprop: fills p.grad for every parameter
    opt.step()                            # p -= lr * p.grad
    if step % 1000 == 0:
        acc = ((logits > 0).float() == y).float().mean()
        print(f"step {step:4d} loss {loss.item():.4f} acc {acc.item():.3f}")

Notes: BCEWithLogitsLoss takes raw scores (logits) and applies the sigmoid internally using a numerically stable formula, so the model has no final sigmoid; a logit above 0 corresponds to probability above 0.5. zero_grad() is required because PyTorch accumulates gradients across backward() calls. nn.Linear stores its weight as (out, in) and computes $xW^\top + b$, the transpose of our notation, with the same result.

Common pitfalls

• Forgetting the non-linearity, or putting a ReLU after the final layer of a regression model so it can never predict negative values.
• Applying sigmoid twice: a sigmoid layer in the model plus BCEWithLogitsLoss. Use logits with the "with logits" losses.
• Not zeroing gradients in PyTorch, so each step uses the sum of all past gradients.
• Shape bugs that broadcast silently, such as a target of shape (n,) against predictions of shape (n, 1), producing an $n\times n$ loss matrix. Check shapes explicitly.
• Zero or identical initialization, which leaves hidden units symmetric.
• Learning rate too high or too low: divergence or loss that barely moves. Try a few values spaced by factors of 3 to 10.
• Hand-written backprop without a gradient check: compare against finite differences on a tiny network before trusting it.

In practice

• Nobody writes backprop by hand in production; frameworks (PyTorch, JAX, TensorFlow) do it. But engineers who understand it can read error messages, spot exploding or vanishing gradients, write custom layers, and reason about memory (the forward activations must be stored for the backward pass, which is why training uses far more memory than inference).
• For tabular data, MLPs are often beaten by gradient-boosted trees; their strength is learned representations of images, audio, text and other raw signals, and they appear inside every larger architecture (the feed-forward block of a Transformer is a two-layer MLP).
• Use ReLU or GELU in hidden layers, a linear output with an appropriate "with logits" loss, and the Adam or AdamW optimizer as a default; the next lesson covers making deep versions train reliably.
• When a model will not learn, shrink the problem: a tiny network, a handful of examples, a gradient check. If it cannot overfit ten examples, there is a bug.