Multilayer perceptrons and backpropagation

What this chapter covers

The last chapter ended at a wall. A single perceptron draws one straight line, XOR needs two, and no amount of training moves a line into a place it can't go. Minsky and Papert proved it, the field believed them, and neural networks went quiet for a decade. This chapter is the escape. Stack a second layer of units between the input and the output, bend them with a nonlinear activation, and the boundary stops being a line — it can curve, fold, wrap around a class. The exact problem that stopped the perceptron falls apart in about four lines of NumPy.

We build a one-hidden-layer network from scratch: a linear layer, a tanh nonlinearity, a second linear layer, and a sigmoid that reads out a probability. Then we train it with backpropagation, which sounds like a proper noun and is really just the chain rule run backwards through the two layers — get the error at the output, push it back through the output weights, multiply by the activation's slope to get the error at the hidden layer, and read off the gradient for every weight on the way. Gradient descent does the rest. We check the whole thing against scikit-learn's MLPClassifier on the same data.

The centerpiece is an animation of the boundary learning to bend. You watch it start as a badly placed near-straight cut and curve into shape epoch by epoch as backprop shapes the hidden units, a loss panel draining underneath it. And the headline is a single fact worth the whole chapter: one hidden layer with a nonlinear activation solves XOR exactly, at 100%, which a single perceptron provably cannot. That's the two ideas — depth and nonlinearity — that everything after this is built on.

A bit of history

The fix was imaginable in 1969 and nobody could train it. Stacking perceptrons into layers was an obvious idea; the problem was that the perceptron's learning rule is local to one unit, and once you bury a unit in the middle of a network there's no label to compare it against, so the rule has nothing to correct toward. You need a way to assign blame to a hidden unit for a mistake made at the output two layers away. That's the credit assignment problem, and it's what kept the stack untrainable.

The answer is backpropagation, and its history is a story of a good idea invented several times before it stuck. Paul Werbos wrote it down in his 1974 Harvard PhD thesis, as a way to compute derivatives through any chain of functions, and almost nobody in machine learning noticed. Others — Seppo Linnainmaa's reverse-mode automatic differentiation in 1970, work by David Parker and Yann LeCun in the early eighties — circled the same idea. What made it land was a 1986 paper in Nature by David Rumelhart, Geoffrey Hinton, and Ronald Williams, "Learning representations by back-propagating errors." They weren't the first to derive it, and they said so, but they were the ones who showed it working — a network learning useful internal representations in its hidden layer, solving problems a single layer couldn't, XOR among them. That paper restarted neural network research. Hinton got a Turing Award for the line of work it opened.

The theory caught up fast. In 1989 George Cybenko proved the universal approximation theorem: a network with a single hidden layer of sigmoidal units can approximate any continuous function on a bounded region to any accuracy you like, given enough hidden units. Kurt Hornik generalized it in 1991 — the magic isn't the sigmoid, it's the nonlinearity plus a hidden layer. This is the theorem people wave at when they say neural networks can learn anything, and the caveat is important: it says a network exists that fits the function, not that gradient descent will find it, and "enough hidden units" can be an absurd number. But it settled the representational question Minsky and Papert had raised. One hidden layer is, in principle, enough to represent anything. The rest of deep learning is about making that practical.

The intuition

Here's the whole idea, and it's worth getting before any math. A single perceptron computes one weighted sum and thresholds it — one line, one cut. A hidden layer computes several weighted sums at once, one per hidden unit, each its own line through the input space. Then the output layer takes those hidden responses as its new features and draws a line through them. A line in the space of "which side of each hidden line are you on" is a bent, piecewise shape back in the original space. Two lines become a corner; a handful become a curve; enough of them wrap a blob. The network isn't drawing one boundary — it's learning a new set of coordinates in the hidden layer where the problem is linearly separable, and then separating it there.

The nonlinearity is what makes this work, and it's not optional. If the hidden activation were linear — just pass the weighted sum through unchanged — then two stacked linear layers collapse into one linear layer, algebraically, and you're back to a single perceptron with extra steps. A product of matrices is a matrix. The tanh in the middle is what breaks that collapse; it bends each hidden unit's response so the composition can't be flattened back into a line. Depth without nonlinearity buys you nothing. Depth with nonlinearity buys you everything.

Here's the data the network will learn on: two interleaving half-moons, 200 points, the kind of thing no straight line comes close to splitting. One class curls into the other. Hold this picture — the whole animation later is this scatter with the network's boundary bending across it.

The math

Four ideas: the forward pass, the loss, the backward pass, the update. The backward pass is the only one with any teeth, and even that is just the chain rule applied twice. Let $X$ be the inputs, one row per example, $D$ features wide. The network has a hidden layer of $H$ units and a single output unit.

The forward pass pushes the data through two linear layers with a nonlinearity between them. First layer: a linear map, then tanh.

Second layer: another linear map, then a sigmoid to squash the score into a probability of class 1.

Here $W_1$ is $D \times H$, $W_2$ is $H \times 1$, and $A_1$ is the hidden representation — the new coordinates the network learns. Everything the model knows lives in those four parameter arrays.

The loss is binary cross-entropy, the same log loss as logistic regression, because the output unit is a logistic unit. For labels $y_i \in \{0, 1\}$ and predictions $\hat{y}_i$, averaged over $N$ examples:

Now backpropagation, which is the chain rule, output first. The reason it starts clean is a small miracle of design: sigmoid and cross-entropy were built for each other, and when you differentiate the loss through the sigmoid the messy factors cancel, leaving the bare error at the output.

That error, dotted with the hidden activations, is the gradient for the output weights. The bias gradient is just the error summed.

Then the step that names the whole method. Push the output error back through the output weights to see how much each hidden unit contributed, and multiply by the tanh derivative $1 - A_1^2$ to convert "error in the hidden activation" into "error in the hidden pre-activation." That product is the hidden-layer error.

And the hidden error, dotted with the inputs, gives the first-layer gradients — the same shape of expression as the output layer, one rung down the chain.

That is the entire algorithm. Two deltas, four gradients, and the only genuinely new idea is the middle line: error flows backward through the same weights that carried the signal forward, bent by the activation's slope at each layer. The update is plain gradient descent — nudge every parameter downhill by a learning rate $\eta$.

Do the forward pass, measure the loss, backprop the gradients, take the step, repeat. That loop, run a few thousand times, is training.

What it's good at, what it isn't

What the MLP is good at is the thing the perceptron couldn't do: it represents nonlinear boundaries, and the universal approximation theorem says one hidden layer can represent essentially any of them. Give it enough hidden units and a differentiable activation and it will fit curves, corners, disconnected regions — whatever the data needs — and it learns those shapes itself, from the gradient, with no feature engineering from you. That is a genuine break from everything earlier in the book, where you either got a linear model or hand-crafted the nonlinear features yourself. The network builds its own features in the hidden layer. That is the whole idea, and it scales: swap the two dense layers for convolutions and you have the network in the next chapter.

What it isn't is convex, cheap, or interpretable. The loss surface has many local minima and long flat regions, so the answer you get depends on the random weight initialization and the learning rate, and two honest runs can land in different places. There's no clean convergence theorem like the perceptron's — you tune, you watch the loss, you hope. It wants more data than a linear model to fill in all those parameters without overfitting, the training is slower, and when it's done you have a pile of weights that predict well and explain nothing. On a genuinely linear problem it's overkill and a logistic regression will match it with a tenth of the fuss. You reach for the MLP when the boundary is curved and you don't know its shape — which, to be fair, is most interesting problems.

The data

Two datasets, one to prove the point and one to measure it.

The first is XOR, the four-point problem that broke the perceptron: inputs at the corners of a square, $(\pm 1, \pm 1)$, labeled 1 when the coordinates have opposite signs and 0 when they match. Diagonally opposite corners share a class, so no straight line separates them. Four points, no line — the whole of Minsky and Papert's argument. We keep it because it's the cleanest possible demonstration that the extra layer changes what's representable, not just how well it fits.

The second is the moons set from the intuition section: 200 points from scikit-learn's make_moons with a noise of 0.2, two interleaving crescents, split 140 for training and 60 held out for the test. XOR is a capacity check — can the model express the answer at all — and moons is the honest accuracy measurement on data with a real curved boundary and some noise to make it earn it.

Build it, one function at a time

Six pieces, bottom up, in the order the forward-then-backward story wants them. Start with initialization, because for the first time in the book the starting weights matter. The perceptron began at zero and it was fine. Here, zero is fatal.

def init_params(n_in, n_hidden, n_out=1, seed=0):
    """Random small weights, zero biases, for a one-hidden-layer net.

    The weights CANNOT start at zero the way the perceptron's did. If every
    hidden unit starts identical it stays identical — they all get the same
    gradient forever, and the layer collapses to a single unit. Random init
    breaks that symmetry. The 1/sqrt(fan-in) scale (Xavier-style) keeps the
    tanh activations away from their flat saturated ends at the start, where
    gradients would be too small to move.
    """
    rng = np.random.default_rng(seed)
    return {
        "W1": rng.normal(0.0, np.sqrt(1.0 / n_in), size=(n_in, n_hidden)),
        "b1": np.zeros(n_hidden),
        "W2": rng.normal(0.0, np.sqrt(1.0 / n_hidden), size=(n_hidden, n_out)),
        "b2": np.zeros(n_out),
    }

The comment says why, and it's worth reading twice. If every hidden unit starts identical, every hidden unit receives the identical gradient, so they update in lockstep and stay identical forever — the layer collapses to a single unit and you've lost the whole point of having a layer. Random initialization breaks that symmetry so the units can specialize. The $1/\sqrt{\text{fan-in}}$ scale is the old Xavier trick: it keeps the tanh inputs in the middle of their range at the start, away from the flat saturated ends where the slope $1 - A_1^2$ is near zero and gradients die.

Now the forward pass — linear, tanh, linear, sigmoid — exactly the four equations from the math section. The sigmoid is written in its numerically stable two-branch form so a large negative score doesn't overflow the exponential.

def sigmoid(z):
    """Squash to (0, 1), in the numerically stable branch form."""
    z = np.asarray(z, dtype=float)
    out = np.empty_like(z)
    pos = z >= 0
    out[pos] = 1.0 / (1.0 + np.exp(-z[pos]))
    ez = np.exp(z[~pos])
    out[~pos] = ez / (1.0 + ez)
    return out


def forward(X, p):
    """Forward pass: linear -> tanh -> linear -> sigmoid.

        z1 = X W1 + b1        a1 = tanh(z1)
        z2 = a1 W2 + b2       yhat = sigmoid(z2)

    X is (N, D). The hidden layer a1 is (N, H); the output yhat is (N, 1),
    the probability of class 1. Returns yhat and a cache of the intermediate
    activations that backprop needs to reconstruct the gradients.
    """
    z1 = X @ p["W1"] + p["b1"]
    a1 = np.tanh(z1)
    z2 = a1 @ p["W2"] + p["b2"]
    yhat = sigmoid(z2)
    cache = {"X": X, "a1": a1, "yhat": yhat}
    return yhat, cache

Notice it returns a cache. The backward pass needs the intermediate activations to reconstruct the gradients, so we stash them on the way through rather than recompute them. That is the space-for-time trade at the heart of backprop: keep the forward pass's work around, and the backward pass is cheap.

The loss is the cross-entropy, one line of it, with a clip to keep the logarithm off zero.

def cross_entropy(yhat, y, eps=1e-12):
    """Mean binary cross-entropy. y is (N, 1) in {0, 1}.

    The same log loss as logistic regression, because the output layer IS a
    logistic unit. A confident, correct probability costs almost nothing; a
    confident, wrong one costs a lot. Clipping keeps log() off zero.
    """
    yhat = np.clip(yhat, eps, 1.0 - eps)
    return float(-np.mean(y * np.log(yhat) + (1 - y) * np.log(1 - yhat)))

Then the piece that is the whole reason this chapter exists. Backprop, output delta first, then the push back to the hidden layer, then the gradients. It is a direct transcription of the four backward equations — read it against them.

def backward(cache, y, p):
    """Backprop: the chain rule, layer by layer, output first.

    The output delta is the bare error, because sigmoid + cross-entropy were
    built for each other and their derivatives cancel:

        d2 = yhat - y                          (N, 1)

    Output-layer gradients — the hidden activations dotted with that error:

        dW2 = a1^T d2 / N        db2 = mean(d2)

    Now the key step. Push the output error back through the output weights,
    then multiply by the tanh derivative (1 - a1^2) to get the error AT the
    hidden layer. That product is the whole of backpropagation:

        d1 = (d2 W2^T) * (1 - a1^2)             (N, H)

    Hidden-layer gradients — the inputs dotted with the hidden error:

        dW1 = X^T d1 / N         db1 = mean(d1)
    """
    X, a1, yhat = cache["X"], cache["a1"], cache["yhat"]
    N = X.shape[0]
    d2 = (yhat - y) / N                 # fold the 1/N in once, at the source
    dW2 = a1.T @ d2
    db2 = d2.sum(axis=0)
    d1 = (d2 @ p["W2"].T) * (1.0 - a1**2)
    dW1 = X.T @ d1
    db1 = d1.sum(axis=0)
    return {"W1": dW1, "b1": db1, "W2": dW2, "b2": db2}

The line that does the work is d1 = (d2 @ p["W2"].T) * (1.0 - a1**2). The first factor sends the output error back through the output weights; the second bends it by the tanh slope. Everything else is a dot product turning an error into a gradient. I fold the $1/N$ into d2 once, at the source, so it flows into every downstream gradient without repeating.

The update is the plainest gradient-descent step there is — walk every parameter downhill.

def sgd_update(p, grads, lr):
    """One gradient-descent step: walk every parameter downhill."""
    for k in p:
        p[k] = p[k] - lr * grads[k]
    return p

And the loop that runs it: initialize, then each epoch do a forward pass, measure the loss, backprop, step. Full-batch, because these datasets are tiny and every epoch can afford to see all the data at once. The record_every hook snapshots the weights on a schedule so we can replay the boundary bending into shape.

def fit(X, y, n_hidden=8, lr=0.5, n_epochs=4000, seed=0, record_every=0):
    """Full-batch gradient descent with backprop.

    Init the params, then each epoch: forward pass, measure the loss,
    backprop the gradients, take one step against them. Full-batch here
    because the datasets are tiny — every epoch sees all the data at once.
    With record_every > 0 we snapshot (epoch, params, loss) on a schedule
    so the chapter can replay the boundary bending into shape.
    """
    y = np.asarray(y, dtype=float).reshape(-1, 1)
    p = init_params(X.shape[1], n_hidden, 1, seed=seed)
    history = []
    for epoch in range(n_epochs):
        yhat, cache = forward(X, p)
        loss = cross_entropy(yhat, y)
        if record_every and epoch % record_every == 0:
            history.append((epoch, {k: v.copy() for k, v in p.items()}, loss))
        grads = backward(cache, y, p)
        p = sgd_update(p, grads, lr)
    if record_every:
        final_loss = cross_entropy(forward(X, p)[0], y)
        history.append((n_epochs, {k: v.copy() for k, v in p.items()}, final_loss))
        return p, history
    return p

Watch it work

This is the network learning the moons, and it's the chapter in one picture. The top panel is the training data; the shaded grid behind it is the network's current prediction everywhere in the plane — cyan where it predicts class 0, purple where it predicts class 1, and the seam between them is the decision boundary. The bottom panel is the cross-entropy loss falling as training runs. Every frame is a real snapshot of the actual weights at that epoch — no smoothing, no staging. Reset and replay it as many times as you like; the seed is fixed, so it's the same run every time.

Watch the boundary at the start. Epoch 0 is random weights, and the boundary is a single near-straight cut sitting in the wrong place — a perceptron's-eye view of the world, one line, badly aimed. Then backprop goes to work. The line develops a kink, then a bend, then it curls to follow the gap between the crescents. That curl is the hidden units specializing: each one is claiming a region of the input, and the output layer is stitching their responses into a shape no single line could make. By the last frame the boundary has wrapped the interleaving moons and the training accuracy is about 0.96, the loss down near 0.07. A single perceptron, given forever, could not draw that curve. This one drew it in a few thousand cheap steps.

The loss panel tells the same story as a number going down. It drops fast at first — the early gradients are large because the errors are large — then flattens into the long tail where the network is refining a boundary that's already roughly right. That shape, steep then flat, is what almost every neural network's loss curve looks like, and reading it is half of knowing whether training is going well.

Now the payoff, the thing this whole architecture was built to do. Here's the decision surface the same network learns on XOR — the four points, and the region the network assigns to each class shaded behind them.

Look at what the network did. It carved the plane into a diagonal band of one class cutting across two corners of the other — not one line but a pair of them, bent into an X, exactly the shape XOR needs and exactly the shape a single perceptron cannot make. Each of the four points sits in the right-colored region. This is the picture the field waited fifteen years for. The extra layer didn't just fit better; it changed the set of shapes the model could draw at all, and that change is the entire difference between a line and a brain-in-training.

One more view of training, the loss on its own, for both datasets on one axis. This is the number the animation's bottom panel was tracking, plotted end to end.

Both curves have the same character — a steep early drop while the errors are large, then a long flat tail of refinement. XOR falls further and faster because four clean points are an easier fit than 140 noisy ones; moons levels off higher because the noise sets a floor no boundary can beat. Neither curve is monotone- smooth if you zoom in, and that's fine — full-batch gradient descent on a nonconvex loss wobbles. What matters is the trend, and the trend is down.

The full implementation

The whole file, top to bottom — initialization, forward, loss, backward, update, the training loop, and the prediction helpers the traces use. This is the exact code the animation ran:

"""A multilayer perceptron, built from scratch.

One hidden layer is all it takes. Push the inputs through a linear layer,
bend them with a nonlinear activation (tanh), push that through a second
linear layer, and squash the result to a probability with the sigmoid.
That stack can draw a curved decision boundary — which a single perceptron,
the previous chapter, provably cannot.

Training is backpropagation, which is nothing more than the chain rule run
backwards through the two layers: get the error at the output, push it back
through the output weights, multiply by the activation's derivative to get
the error at the hidden layer, and read off the gradient for every weight
on the way. Then take a gradient-descent step and do it again.

Pure NumPy. Every function below appears in the chapter one step at a time;
the `# region:` markers are what the book's include directives pull in.
Labels live in {0, 1}; the output is P(class = 1).
"""

import numpy as np


# region: init_params
def init_params(n_in, n_hidden, n_out=1, seed=0):
    """Random small weights, zero biases, for a one-hidden-layer net.

    The weights CANNOT start at zero the way the perceptron's did. If every
    hidden unit starts identical it stays identical — they all get the same
    gradient forever, and the layer collapses to a single unit. Random init
    breaks that symmetry. The 1/sqrt(fan-in) scale (Xavier-style) keeps the
    tanh activations away from their flat saturated ends at the start, where
    gradients would be too small to move.
    """
    rng = np.random.default_rng(seed)
    return {
        "W1": rng.normal(0.0, np.sqrt(1.0 / n_in), size=(n_in, n_hidden)),
        "b1": np.zeros(n_hidden),
        "W2": rng.normal(0.0, np.sqrt(1.0 / n_hidden), size=(n_hidden, n_out)),
        "b2": np.zeros(n_out),
    }
# endregion


# region: forward
def sigmoid(z):
    """Squash to (0, 1), in the numerically stable branch form."""
    z = np.asarray(z, dtype=float)
    out = np.empty_like(z)
    pos = z >= 0
    out[pos] = 1.0 / (1.0 + np.exp(-z[pos]))
    ez = np.exp(z[~pos])
    out[~pos] = ez / (1.0 + ez)
    return out


def forward(X, p):
    """Forward pass: linear -> tanh -> linear -> sigmoid.

        z1 = X W1 + b1        a1 = tanh(z1)
        z2 = a1 W2 + b2       yhat = sigmoid(z2)

    X is (N, D). The hidden layer a1 is (N, H); the output yhat is (N, 1),
    the probability of class 1. Returns yhat and a cache of the intermediate
    activations that backprop needs to reconstruct the gradients.
    """
    z1 = X @ p["W1"] + p["b1"]
    a1 = np.tanh(z1)
    z2 = a1 @ p["W2"] + p["b2"]
    yhat = sigmoid(z2)
    cache = {"X": X, "a1": a1, "yhat": yhat}
    return yhat, cache
# endregion


# region: loss
def cross_entropy(yhat, y, eps=1e-12):
    """Mean binary cross-entropy. y is (N, 1) in {0, 1}.

    The same log loss as logistic regression, because the output layer IS a
    logistic unit. A confident, correct probability costs almost nothing; a
    confident, wrong one costs a lot. Clipping keeps log() off zero.
    """
    yhat = np.clip(yhat, eps, 1.0 - eps)
    return float(-np.mean(y * np.log(yhat) + (1 - y) * np.log(1 - yhat)))
# endregion


# region: backward
def backward(cache, y, p):
    """Backprop: the chain rule, layer by layer, output first.

    The output delta is the bare error, because sigmoid + cross-entropy were
    built for each other and their derivatives cancel:

        d2 = yhat - y                          (N, 1)

    Output-layer gradients — the hidden activations dotted with that error:

        dW2 = a1^T d2 / N        db2 = mean(d2)

    Now the key step. Push the output error back through the output weights,
    then multiply by the tanh derivative (1 - a1^2) to get the error AT the
    hidden layer. That product is the whole of backpropagation:

        d1 = (d2 W2^T) * (1 - a1^2)             (N, H)

    Hidden-layer gradients — the inputs dotted with the hidden error:

        dW1 = X^T d1 / N         db1 = mean(d1)
    """
    X, a1, yhat = cache["X"], cache["a1"], cache["yhat"]
    N = X.shape[0]
    d2 = (yhat - y) / N                 # fold the 1/N in once, at the source
    dW2 = a1.T @ d2
    db2 = d2.sum(axis=0)
    d1 = (d2 @ p["W2"].T) * (1.0 - a1**2)
    dW1 = X.T @ d1
    db1 = d1.sum(axis=0)
    return {"W1": dW1, "b1": db1, "W2": dW2, "b2": db2}
# endregion


# region: update
def sgd_update(p, grads, lr):
    """One gradient-descent step: walk every parameter downhill."""
    for k in p:
        p[k] = p[k] - lr * grads[k]
    return p
# endregion


# region: fit
def fit(X, y, n_hidden=8, lr=0.5, n_epochs=4000, seed=0, record_every=0):
    """Full-batch gradient descent with backprop.

    Init the params, then each epoch: forward pass, measure the loss,
    backprop the gradients, take one step against them. Full-batch here
    because the datasets are tiny — every epoch sees all the data at once.
    With record_every > 0 we snapshot (epoch, params, loss) on a schedule
    so the chapter can replay the boundary bending into shape.
    """
    y = np.asarray(y, dtype=float).reshape(-1, 1)
    p = init_params(X.shape[1], n_hidden, 1, seed=seed)
    history = []
    for epoch in range(n_epochs):
        yhat, cache = forward(X, p)
        loss = cross_entropy(yhat, y)
        if record_every and epoch % record_every == 0:
            history.append((epoch, {k: v.copy() for k, v in p.items()}, loss))
        grads = backward(cache, y, p)
        p = sgd_update(p, grads, lr)
    if record_every:
        final_loss = cross_entropy(forward(X, p)[0], y)
        history.append((n_epochs, {k: v.copy() for k, v in p.items()}, final_loss))
        return p, history
    return p
# endregion


def predict_proba(X, p):
    """P(class = 1) for every row."""
    return forward(X, p)[0].ravel()


def predict(X, p, threshold=0.5):
    """0/1 labels at a decision threshold."""
    return (predict_proba(X, p) >= threshold).astype(int)


def accuracy(X, y, p, threshold=0.5):
    """Fraction of rows the net labels correctly."""
    return float((predict(X, p, threshold) == np.asarray(y)).mean())

The library version

You would not hand-roll this in production; scikit-learn ships it as MLPClassifier. We give it the same architecture — one hidden layer of eight tanh units, a logistic output, cross-entropy loss — so the comparison is honest:

def sklearn_mlp(Xtr, ytr, Xte, yte, n_hidden=8, seed=0, max_iter=4000):
    """Fit sklearn's MLPClassifier and report train + test accuracy.

    Same architecture as the scratch net: one hidden layer of `n_hidden`
    tanh units, a logistic output, cross-entropy loss. Returns the fitted
    model, its train accuracy, test accuracy, and the epochs it ran.
    """
    clf = MLPClassifier(
        hidden_layer_sizes=(n_hidden,),
        activation="tanh",
        solver="adam",
        max_iter=max_iter,
        random_state=seed,
    )
    clf.fit(Xtr, ytr)
    return (
        clf,
        float(clf.score(Xtr, ytr)),
        float(clf.score(Xte, yte)),
        int(clf.n_iter_),
    )

The differences are all under the hood, and they're the differences you'd want. Where our version runs plain full-batch gradient descent at a fixed learning rate, sklearn defaults to the adam optimizer — an adaptive per-parameter step size that usually trains faster and needs less tuning — and it adds a small L2 penalty (the alpha parameter) that regularizes the weights, which ours omits. It also stops itself when the loss stops improving rather than running a fixed number of epochs. We keep those defaults on purpose. The point of the face-off is to compare against the library used the way you'd actually use it, not a hobbled clone of our own code.

Scratch versus library

Two questions, two datasets. On XOR the question is binary: can the model express the answer? On moons the question is quantitative: how accurately does it label held-out points?

On XOR both models score 1.0 — a perfect fit, all four points on the right side of a curved boundary. That's the headline cashed out as a number: the single perceptron of the last chapter scored 0.5 here, a coin flip, and one hidden layer takes it to 1.0. The extra layer didn't improve XOR, it unlocked it. If you want one comparison that justifies the rest of a deep learning course, it's the jump from 0.5 to 1.0 on a problem a child solves by drawing two lines.

On moons, measured on the 60 held-out test points, our from-scratch net scores 0.983 and sklearn's scores 0.867. Ours edges it out here, which deserves an honest footnote rather than a victory lap. The test set is small — 0.983 versus 0.867 is about fifty-nine correct against fifty-two out of sixty — and the two models made different bets. Our full-batch gradient descent ran the full four thousand epochs at a fairly aggressive learning rate and fit the training set tightly (training accuracy 0.964); sklearn's adam stopped itself after 513 iterations and carries that L2 penalty, both of which pull it toward a smoother, more conservative boundary. On this particular small, low-noise 2-D problem, fitting harder happened to pay. That is not a general claim that hand-rolled beats the library — turn sklearn's regularization down and let it run and it closes the gap — it's a reminder that on neural networks the hyperparameters are the model, and two reasonable choices give two different boundaries. Both draw a curve through the moons. That they draw slightly different curves is the whole nonconvex story in miniature.

Takeaways

The multilayer perceptron is where machine learning becomes deep learning, and the reason is two ideas you can hold in one hand. Depth: put a layer of units between the input and the output so the network can build its own features instead of using the ones you gave it. Nonlinearity: bend those units with something like tanh so the layers don't algebraically collapse back into one. With both, the boundary can curve, and the universal approximation theorem says it can curve into essentially any shape. Without either, you have a perceptron. That's the whole conceptual jump, and everything fancier — more layers, ReLU, dropout, convolutions, attention — is engineering on top of these two moves.

Three things to carry forward. First, backpropagation is not magic and it is not a separate algorithm from calculus — it is the chain rule, applied layer by layer, reusing the forward pass's cached activations so the whole gradient costs about as much as one more forward pass. Once you've written the two deltas by hand you'll never be mystified by a framework's .backward() again; it's doing exactly this, just for more layers. Second, the guarantees you had are gone. The perceptron converged provably on separable data; the MLP's loss surface is a nonconvex landscape, and initialization, learning rate, and when you stop all move the answer. You trade a guarantee for a capability, and you manage the capability by watching the loss curve — which is why we drew it. Third, this is the foundation the next chapter stands on. A convolutional network is this same forward-loss- backward-update loop with the dense first layer replaced by weight-shared convolutions; the training story does not change, only the layer does. Get the MLP cold — the forward pass, the cross-entropy, the two deltas, the step — and you have read the skeleton of every neural network there is. XOR is where the perceptron ended and this began. This is where the deep networks begin.