A Neural Network from Scratch

Building an automatic-differentiation engine and a neural network from nothing but Python’s arithmetic — no PyTorch, no NumPy, no libraries — and the mathematics that makes it work.

A note before you start

I wanted to understand what training a model means at the level of individual numbers, so I built the whole thing from scratch: a tiny automatic-differentiation engine, a neuron, a network of neurons, and the loop that makes them learn — using only Python’s built-in arithmetic, no numerical libraries whatsoever. This post is the account I wish I had been handed. It is written to be self-contained and mathematically explicit: by the end you should be able to derive backpropagation yourself, not merely recognise it.

Coming from econometrics, the punchline that stayed with me is that almost none of this is new mathematics. A neural network is a composite function; training it is maximum likelihood estimation; backpropagation is the chain rule applied with good bookkeeping. The novelty is not in the ideas but in their scale. Let me show you the whole machine with the panels off.

1. The neuron

Everything is built from one atom. A neuron takes a vector of inputs $x = (x_1, \dots, x_n) \in \mathbb{R}^n$, keeps a vector of weights $w = (w_1, \dots, w_n)$ and a scalar bias $b$, and produces a single number in two steps: an affine combination followed by a nonlinear activation.

Formally, the neuron computes the pre-activation

$$ z ;=; \sum_{i=1}^{n} w_i x_i + b ;=; w^{\top}x + b, $$

and then the activation

$$ a ;=; \phi(z), \qquad \phi(z) = \tanh(z) = \frac{e^{z}-e^{-z}}{e^{z}+e^{-z}}. $$

The output $a$ is a single scalar. Two remarks earn their place here.

First, the affine part $w^{\top}x + b$ is exactly a linear predictor — the same object as the right-hand side of a linear regression. If we stopped there, stacking neurons would be pointless: a composition of affine maps is itself affine, $A(Bx+c)+d = (AB)x + (Ac+d)$, so a hundred linear layers collapse into one. The nonlinearity $\phi$ is precisely what breaks that collapse and lets depth express functions a single layer cannot. It is not decoration; it is load-bearing.

Second, $\tanh$ has a derivative that will make the backward pass pleasant:

$$ \phi’(z) ;=; 1 - \tanh^2(z) ;=; 1 - a^2 . $$

The derivative is expressible through the output $a$ alone — a small convenience we will exploit.

2. A loss, and the goal of training

A model is only as good as our ability to say how wrong it is. Given a target $y$ for an input $x$, collapse the badness of the prediction $a$ into a single number with a loss function. For a scalar regression target the natural choice is squared error,

$$ L ;=; (a - y)^2 . $$

Training is the search for parameters $\theta = (w, b)$ that make $L$ small. In linear regression this search has a closed form, $\hat\beta = (X^{\top}X)^{-1}X^{\top}y$, because the model is linear and the loss is quadratic. The moment we wrap the affine term inside a nonlinearity and nest neurons, that luxury disappears — there is no formula to solve. What remains is the only general strategy available:

$$ \theta ;\leftarrow; \theta - \eta,\nabla_{\theta} L , $$

gradient descent: start from arbitrary parameters and repeatedly nudge each one against its gradient, where $\eta > 0$ is the learning rate. To take a single step we need $\nabla_\theta L$ — the sensitivity of the loss to every parameter. Producing those sensitivities efficiently is the whole game, and it is what backpropagation does.

3. The backward pass, by hand, for one neuron

Here is the part worth slowing down for. I will take one neuron with concrete numbers and compute every gradient by hand, because once you have seen it done at the level of arithmetic, the general case is only bookkeeping.

Fix a two-input neuron and one data point:

$$ x_1 = 3,\quad x_2 = 1,\qquad w_1 = 0.5,\quad w_2 = -1,\quad b = 0.5,\qquad y = 1 . $$

Forward pass. Compute left to right, storing each intermediate value:

$$ z = (0.5)(3) + (-1)(1) + 0.5 = 1.0, \qquad a = \tanh(1.0) = 0.7616, \qquad L = (0.7616 - 1)^2 = 0.0568 . $$

The prediction $0.7616$ sits below the target $1$; the loss records the gap.

Backward pass. We want $\partial L/\partial w_1$, $\partial L/\partial w_2$, $\partial L/\partial b$. The trick is to decompose the computation into a chain $w, x, b \to z \to a \to L$ and apply the chain rule one link at a time, starting from the trivial seed $\partial L/\partial L = 1$ and working backwards:

$$ \begin{aligned} \frac{\partial L}{\partial a} &= 2(a - y) = 2(-0.2384) = -0.4768, \\[4pt] \frac{\partial a}{\partial z} &= 1 - a^2 = 1 - 0.7616^2 = 0.4200, \\[4pt] \frac{\partial L}{\partial z} &= \frac{\partial L}{\partial a}\cdot\frac{\partial a}{\partial z} = (-0.4768)(0.4200) = -0.2002 . \end{aligned} $$

The pre-activation $z = w_1x_1 + w_2x_2 + b$ is linear in the parameters, so its local derivatives are immediate — $\partial z/\partial w_i = x_i$ and $\partial z/\partial b = 1$ — and one more application of the chain rule finishes each parameter:

$$ \begin{aligned} \frac{\partial L}{\partial w_1} &= \frac{\partial L}{\partial z}\cdot x_1 = (-0.2002)(3) = -0.6007, \\[4pt] \frac{\partial L}{\partial w_2} &= \frac{\partial L}{\partial z}\cdot x_2 = (-0.2002)(1) = -0.2002, \\[4pt] \frac{\partial L}{\partial b} &= \frac{\partial L}{\partial z}\cdot 1 = -0.2002 . \end{aligned} $$

The whole computation, forward values and backward gradients together, lives on a single graph:

Two sanity checks, both of which I now run reflexively. Every parameter gradient is negative, so gradient descent will increase each parameter — which is right, since the prediction is below target and raising the parameters raises it. And $\partial L/\partial w_1$ is exactly three times $\partial L/\partial w_2$, because $x_1 = 3x_2$: a parameter’s gradient scales with the input it multiplies. The mathematics agrees with common sense, which is the first thing to demand of it.

One step of descent. With $\eta = 0.1$,

$$ \begin{aligned} w_1 &\leftarrow 0.5 - 0.1(-0.6007) = 0.5601, \\ w_2 &\leftarrow -1 - 0.1(-0.2002) = -0.9800, \\ b &\leftarrow 0.5 - 0.1(-0.2002) = 0.5200 . \end{aligned} $$

Recomputing the forward pass with the updated parameters gives $z = 1.2203$, $a = 0.8398$, and

$$ L = (0.8398 - 1)^2 = 0.0257 , $$

down from $0.0568$ — a $55%$ reduction in the loss from a single step. Repeat this a few thousand times and you have training.

4. From one neuron to a network

Nothing above changes when we scale up; it only gets wider and longer. A layer is a collection of neurons reading the same input in parallel. Stacking their weight vectors as rows of a matrix $W \in \mathbb{R}^{d_{\text{out}} \times d_{\text{in}}}$ and their biases into $b \in \mathbb{R}^{d_{\text{out}}}$, the entire layer is one matrix–vector product followed by an elementwise activation:

$$ z = Wx + b \in \mathbb{R}^{d_{\text{out}}}, \qquad a = \phi(z) . $$

A network (a multilayer perceptron) is a composition of such layers, each feeding the next. Writing $a^{(0)} = x$ for the input and $L$ for the number of layers,

$$ z^{(l)} = W^{(l)} a^{(l-1)} + b^{(l)}, \qquad a^{(l)} = \phi!\left(z^{(l)}\right), \qquad l = 1, \dots, L, $$

with the prediction read off the final layer, $\hat y = a^{(L)}$.

Now the point I find genuinely beautiful, and the reason building this by hand pays off. The backward pass does not change. A network is a longer composite function, and backpropagation is nothing but the derivative of a composition, evaluated by the same seed-and-reverse-walk we used for one neuron. Concretely, defining the sensitivity of the loss to a layer’s pre-activation as $\delta^{(l)} = \partial L / \partial z^{(l)}$, the chain rule gives the backpropagation recursion

$$ \delta^{(L)} = \nabla_{a^{(L)}} L ;\odot; \phi’!\left(z^{(L)}\right), \qquad \delta^{(l)} = \left(W^{(l+1)\top}\delta^{(l+1)}\right) \odot \phi’!\left(z^{(l)}\right), $$

where $\odot$ is the elementwise product, and the parameter gradients fall straight out:

$$ \frac{\partial L}{\partial W^{(l)}} = \delta^{(l)}, a^{(l-1)\top}, \qquad \frac{\partial L}{\partial b^{(l)}} = \delta^{(l)} . $$

Look at what these formulas are. The seed $\nabla_{a^{(L)}}L$ is the same $2(a-y)$ from the single neuron. The factor $\phi’(z^{(l)}) = 1 - (a^{(l)})^2$ is the same local $\tanh$ derivative. The term $W^{(l+1)\top}\delta^{(l+1)}$ is the same “multiply by the weight on the way back” rule, now bundled into a matrix transpose because a layer has many neurons at once. It is one neuron’s arithmetic, vectorised and iterated. When I extended my own code from a single neuron to a full network, I did not change a single line of the backward routine — it operates on the graph of elementary operations and is blind to whether those operations spell out one neuron or ten thousand. That invariance is what the word deep in deep learning is trading on: the chain rule scales for free; the ideas do not have to be reinvented at each layer.

5. The training loop

With gradients in hand for arbitrary networks, learning is a loop. For a dataset ${(x_i, y_i)}_{i=1}^{N}$ define the objective as the mean loss,

$$ \mathcal{L}(\theta) ;=; \frac{1}{N}\sum_{i=1}^{N} \ell!\left(f(x_i;\theta),, y_i\right), $$

and iterate three steps until the loss stops falling:

$$ \underbrace{f(x_i;\theta)}{\text{forward}} ;\longrightarrow; \underbrace{\nabla\theta \mathcal{L}}{\text{backward}} ;\longrightarrow; \underbrace{\theta \leftarrow \theta - \eta,\nabla\theta \mathcal{L}}_{\text{update}} . $$

One forward pass to predict, one backward pass to obtain every gradient at once, one descent step, repeat. Averaging over examples matters: a single point gives a noisy gradient — in the worked neuron above, a point with $x_2 = 0$ would have produced $\partial L/\partial w_2 = 0$ and taught $w_2$ nothing — so we average the signal over a batch before stepping.

The learning rate deserves the respect it demands in practice. Consider the residual $e = a - y$ for the linear model $a = wx$: one descent step transforms it into

$$ e^{(t+1)} = \bigl(1 - 2\eta x^2\bigr), e^{(t)} , $$

a first-order linear recursion — an AR(1) process in the residual. It converges geometrically precisely when $\lvert 1 - 2\eta x^2 \rvert < 1$, i.e. $\eta < 1/x^2$; at the boundary it oscillates without decaying, and beyond it the error doubles every step and the run diverges. The stability of gradient descent is a unit-root condition. For an econometrician this is oddly comforting: the failure mode of a diverging training run is the same object as an explosive time series.

There is a deeper identity hiding in the choice of loss. Squared error is not arbitrary — it is the negative log-likelihood of a Gaussian observation model, and for classification the natural loss is the cross-entropy, which is exactly the negative log-likelihood of a categorical one. Minimising it by gradient descent is therefore maximum likelihood estimation, carried out iteratively because no closed form exists. The estimator I had used for years reappears here, stripped of its formula and handed to an optimiser. Pretraining a large language model is this same MLE, run over a very long document.

6. Three mistakes that taught me the most

Building without libraries means the bugs are yours alone, and mine were more instructive than any success. Three in particular reshaped how I think.

The method that was never called. I implemented multiplication correctly and named it __Multiply__. Every test that used the * operator ignored it completely, because Python translates a * b into a call to the exact method name __mul__ — a fixed contract, not a descriptive label of my choosing. A perfectly correct function under the wrong name is not a function at all; it is dead code the interpreter never visits. The lesson generalises past Python: an interface is a promise about names, and being right about the mathematics does not exempt you from being right about the contract.

The comma that makes a tuple. To record a node’s parents I wrote (self) where I meant (self,). In Python the parentheses do nothing — it is the comma that builds a tuple, so (self) is just self, while (self,) is a one-element tuple. The mistake surfaced far from its cause, as a TypeError: 'Value' object is not iterable thrown deep inside the constructor. Reading the traceback from the bottom up — last line first — took me straight to it. Small syntax can carry large meaning, and the discipline of reading errors backwards is worth more than any amount of staring at the line you think is wrong.

The bug the tests could not see. This one matters most. My backward rule for the exponential read, in effect,

$$ \texttt{grad} \mathrel{+}= e^{z} - \texttt{out.data}, $$

which, since out.data already equals $e^{z}$, deposits $e^{z} - e^{z} = 0$ — always — and silently omits the incoming upstream gradient entirely. It was wrong in two independent ways, and yet every test passed. The reason is simple and chastening: not one of my tests exercised the exponential’s backward path, so the referee never looked there. A green test suite does not certify that code is correct; it certifies that the code behaves on the cases you thought to check, and is silent about everything else. I found it only by deliberately computing what the gradient should be and comparing — after which I wrote the test that would have caught it, so it can never return unseen. The habit this instilled — trust a result only against an independent check, and treat the absence of a failing test as the absence of evidence rather than evidence of absence — is, I suspect, the single most valuable thing the whole exercise gave me.

Closing

What surprised me, in the end, was how little machinery there is. An automatic-differentiation engine fits on a page: a number that remembers how it was computed, a local-derivative rule for each elementary operation, and a reverse walk over the computation graph that seeds $\partial L/\partial L = 1$ and lets every node accumulate its share of the gradient from those it feeds. A neuron is an affine map through a nonlinearity. A network is that, composed. Training is gradient descent on a loss that turns out to be a negative log-likelihood. Everything grander — the deep networks, the models with billions of parameters — is this same page of ideas, run wider and longer.

That is the reassurance I took from doing it by hand, and the reason I would recommend it to anyone with a mathematical background who uses these tools without quite trusting them: the frontier is not built from unfamiliar mathematics. It is built from the chain rule and maximum likelihood, executed at a scale that makes them feel like something else. Take the panels off, and it is machinery you already understand.