Linear Algebra for ML · 1.17

Matrix Calculus: Gradients, Jacobians & Backprop

The chain rule in matrix form — which is precisely what backpropagation is.

On this page

In plain English — beginner to advanced

New to derivatives? A derivative just measures "if I nudge this input a tiny bit, how much does the output change, and in which direction?" A partial derivative is the same question asked about one input at a time, while every other input is frozen in place — "if I only nudge x (not y), how does the output move?" That's the entire prerequisite for this lesson; nothing more advanced is assumed.

Beginner: the gradient of a function with respect to a vector is just a vector holding every one of those partial derivatives — one number per input, all collected together. The gradient always points in the direction of steepest increase, which is exactly why gradient descent moves in the opposite direction.

Intermediate: the Jacobian generalizes this to functions that output a whole vector, not just a single number: it's a matrix where row i is the gradient of output i with respect to every input. The Hessian (section 1.16) is the Jacobian of the gradient itself — the matrix of all second derivatives.

Advanced: the chain rule, written in matrix form, says the Jacobian of a composition of functions is the product of their individual Jacobians. This one sentence is, precisely and without exaggeration, what backpropagation is: applying the chain rule layer by layer through a neural network, where "layer by layer" means "Jacobian by Jacobian."

Formula

A common, load-bearing identity in ML:

\nabla_x (x^TAx) = (A + A^T)x \ \xrightarrow{A \text{ symmetric}} \ 2Ax

This is exactly the gradient used in the descent animation above and in every derivation of linear/ridge regression's closed-form solution.

Derivation: ∇(xᵀAx) = (A+Aᵀ)x

Write the quadratic form out entirely in components:

f(\vec{x}) = \vec{x}^TA\vec{x} = \sum_i \sum_j A_{ij}x_ix_j

Differentiate with respect to one single coordinate xₖ, using the ordinary product rule on every term that contains xₖ. A term Aᵢⱼxᵢxⱼ contains xₖ either when i=k (contributing Akjxj) or when j=k (contributing Aikxi) — both cases can happen at once only when i=j=k, but the sum below already handles that correctly:

\frac{\partial f}{\partial x_k} = \sum_j A_{kj}x_j + \sum_i A_{ik}x_i = (A\vec{x})_k + (A^T\vec{x})_k

Collecting this across every k gives the full gradient vector, ∇f = Ax + Aᵀx = (A+Aᵀ)x. When A is symmetric (A = Aᵀ, section 1.7's expert note), this simplifies immediately to 2Ax.

Where this is used: setting this gradient to zero and solving is exactly how the closed-form solution to linear and ridge regression (sections 1.6, 1.9) is derived — the "normal equation" is nothing more than this identity, solved for x.

Gradient descent, literally rolling downhill

Drag the starting point, then press Descend — each step moves opposite the gradient 2Ax. Notice the zig-zag: that's the cost of a poorly conditioned Hessian.

Practical example — verifying a gradient with finite differences

"Gradient checking" — comparing an analytic gradient formula against a numeric finite-difference approximation — is a standard debugging technique whenever you implement backpropagation by hand.

python

import numpy as np

A = np.array([[2., 0.], [0., 5.]])

def f(x):
    return x @ A @ x

def analytic_grad(x):
    return 2 * A @ x   # the identity above, since A is symmetric

def numeric_grad(x, eps=1e-6):
    g = np.zeros_like(x)
    for i in range(len(x)):
        dx = np.zeros_like(x); dx[i] = eps
        g[i] = (f(x + dx) - f(x - dx)) / (2 * eps)
    return g

x = np.array([1.5, -0.7])
print(analytic_grad(x))
print(numeric_grad(x))   # should match closely -> "gradient checking"

Real-world examples

Every deep learning framework's autograd engine (PyTorch's autograd, TensorFlow's GradientTape) is a system for automatically applying the chain rule of Jacobians through an arbitrary computation graph.
Normalizing flows in generative modeling need the Jacobian determinant (section 1.12) of each transformation to correctly track how probability density changes — they're specifically engineered so that Jacobian is triangular and cheap to compute (section 1.35).
Second-order optimizers (Newton's method, L-BFGS) use the Hessian or an approximation of it to converge faster than plain gradient descent near a solution — in practice via Hessian-vector products (section 1.32) rather than the full matrix.

Common mistakes

Using ∇(xᵀAx) = 2Ax when A isn't symmetric — the correct general form is (A + Aᵀ)x, which only simplifies to 2Ax when A = Aᵀ.
Forming a full Jacobian matrix explicitly when only a Jacobian-vector product is needed — for a network with millions of parameters, the full Jacobian would be far too large to fit in memory; real autograd never materializes it.

Going deeper

Reverse-mode automatic differentiation (what backpropagation actually is) computes gradients by propagating Jacobian-vector products backward through the computation graph, without ever forming any full Jacobian matrix — this is precisely why it's efficient enough to train networks with billions of parameters.

At the master level: forward-mode automatic differentiation computes the same chain rule in the opposite order (Jacobian-vector products propagated forward) — it's more efficient when a function has few inputs and many outputs, the mirror image of the typical neural-network case (many inputs/parameters, one scalar loss output), which is exactly why reverse-mode dominates deep learning specifically.

PreviousQuadratic Forms, Convexity & the Hessian NextEinsum, Tensor Contractions & Convolution

Matrix Calculus: Gradients, Jacobians & Backprop

Stay in the loop

Newsletter