The Jacobian & Multivariate Chain Rule

A vector-valued function $f: \mathbb{R}^n \to \mathbb{R}^m$ is defined by $m$ component functions $f = (f_1, f_2, \ldots, f_m)$, where each $f_i: \mathbb{R}^n \to \mathbb{R}$ is a scalar-valued function. The function maps a point $a \in \mathbb{R}^n$ to the vector

$f(a) = (f_1(a), f_2(a), \ldots, f_m(a)) \in \mathbb{R}^m.$

When $m = 1$, this reduces to the scalar-valued functions of Topic 9.

Let $f: \mathbb{R}^n \to \mathbb{R}^m$ and suppose all partial derivatives $\frac{\partial f_i}{\partial x_j}(a)$ exist. The Jacobian matrix (or simply Jacobian) of $f$ at $a$ is the $m \times n$ matrix

$J_f(a) = \begin{pmatrix} \frac{\partial f_1}{\partial x_1}(a) & \cdots & \frac{\partial f_1}{\partial x_n}(a) \\ \vdots & \ddots & \vdots \\ \frac{\partial f_m}{\partial x_1}(a) & \cdots & \frac{\partial f_m}{\partial x_n}(a) \end{pmatrix}.$

Row $i$ of $J_f(a)$ is $\nabla f_i(a)^T$ — the gradient of the $i$-th component function. Column $j$ is the vector of all partial derivatives with respect to $x_j$. Common notations include $J_f(a)$, $Df(a)$, and $\frac{\partial(f_1, \ldots, f_m)}{\partial(x_1, \ldots, x_n)}$.

Let $f(x) = Ax + b$ where $A$ is a constant $m \times n$ matrix. Then $J_f(a) = A$ for all $a$ — the Jacobian of a linear function is the matrix itself. The derivative of a linear map is the linear map. This is the multivariable analog of ”$f(x) = ax + b$ has derivative $a$.”

Let $f(r, \theta) = (r\cos\theta, r\sin\theta)$. The Jacobian is

$J_f(r, \theta) = \begin{pmatrix} \cos\theta & -r\sin\theta \\ \sin\theta & r\cos\theta \end{pmatrix}.$

At $(r, \theta) = (2, \pi/4)$, we get $J_f = \begin{pmatrix} 1/\sqrt{2} & -\sqrt{2} \\ 1/\sqrt{2} & \sqrt{2} \end{pmatrix}$. The first column tells us how $(x, y)$ changes when we increase $r$ (move radially outward); the second column tells us how $(x, y)$ changes when we increase $\theta$ (move tangentially).

A fully connected layer with activation is $f(x) = \sigma(Wx + b)$ where $\sigma$ is applied component-wise. By the chain rule (which we will prove shortly),

$J_f(x) = \text{diag}(\sigma'(Wx + b)) \cdot W.$

When $\sigma$ is the identity (linear layer), $J_f = W$ — recovering Example 1. When $\sigma = \text{sigmoid}$, the diagonal matrix of $\sigma'$ values modulates each row of $W$: the partial derivative $\frac{\partial f_i}{\partial x_j}$ is $\sigma'(z_i) \cdot W_{ij}$, where $z = Wx + b$ is the pre-activation vector.

When $m = 1$, the Jacobian is a $1 \times n$ row matrix: $J_f(a) = \begin{pmatrix} \frac{\partial f}{\partial x_1}(a) & \cdots & \frac{\partial f}{\partial x_n}(a) \end{pmatrix} = \nabla f(a)^T$. The gradient is the transpose of the single-row Jacobian. When $n = 1$, the Jacobian is an $m \times 1$ column matrix — just the vector of ordinary derivatives of each component. When $m = n = 1$, the Jacobian is a $1 \times 1$ matrix containing the single-variable derivative $f'(a)$ from Topic 5. The Jacobian unifies all of these cases.

A function $f: \mathbb{R}^n \to \mathbb{R}^m$ is differentiable at $a$ if there exists a linear map $Df(a): \mathbb{R}^n \to \mathbb{R}^m$ such that

$\lim_{h \to 0} \frac{\|f(a+h) - f(a) - Df(a)(h)\|}{\|h\|} = 0.$

When $Df(a)$ exists, it is unique and represented by the Jacobian matrix: $Df(a)(h) = J_f(a) \cdot h$ for all $h \in \mathbb{R}^n$. This extends Definition 4 of Topic 9 from $m = 1$ to arbitrary $m$.

$f: \mathbb{R}^n \to \mathbb{R}^m$ is differentiable at $a$ if and only if each component $f_i: \mathbb{R}^n \to \mathbb{R}$ is differentiable at $a$ (in the sense of the total derivative from Topic 9). When $f$ is differentiable, the rows of the Jacobian are the gradients of the component functions.

($\Rightarrow$) If $f$ is differentiable at $a$, then $\frac{\|f(a+h) - f(a) - J_f(a)h\|}{\|h\|} \to 0$. For each component $f_i$, the absolute value of a single component is bounded by the norm of the full vector:

$\frac{|f_i(a+h) - f_i(a) - \nabla f_i(a) \cdot h|}{\|h\|} \le \frac{\|f(a+h) - f(a) - J_f(a)h\|}{\|h\|} \to 0.$

So each $f_i$ is differentiable at $a$ with total derivative $\nabla f_i(a) \cdot h$.

($\Leftarrow$) If each $f_i$ is differentiable at $a$, then for each $i$: $|f_i(a+h) - f_i(a) - \nabla f_i(a) \cdot h| = o(\|h\|)$. We square and sum over all components:

$\|f(a+h) - f(a) - J_f(a)h\|^2 = \sum_{i=1}^m |f_i(a+h) - f_i(a) - \nabla f_i(a) \cdot h|^2.$

Each summand is $o(\|h\|^2)$ by hypothesis, so the sum is $o(\|h\|^2)$. Taking the square root: $\|f(a+h) - f(a) - J_f(a)h\| = o(\|h\|)$. $\square$

If $f: \mathbb{R}^n \to \mathbb{R}^m$ is differentiable at $a$, then

$f(a + h) \approx f(a) + J_f(a) \cdot h$

with error $o(\|h\|)$. This is the multivariable Taylor expansion to first order for vector-valued functions. Geometrically: the affine map $h \mapsto f(a) + J_f(a)h$ is the best linear approximation to $f$ near $a$.

Near $(r, \theta) = (2, \pi/4)$, the polar-to-Cartesian map is approximated by

$f\begin{pmatrix} 2 + \Delta r \\ \pi/4 + \Delta\theta \end{pmatrix} \approx \begin{pmatrix} \sqrt{2} \\ \sqrt{2} \end{pmatrix} + \begin{pmatrix} 1/\sqrt{2} & -\sqrt{2} \\ 1/\sqrt{2} & \sqrt{2} \end{pmatrix} \begin{pmatrix} \Delta r \\ \Delta\theta \end{pmatrix}.$

A small change in $(r, \theta)$ maps to a change in $(x, y)$ via matrix multiplication — the Jacobian acts on the perturbation vector. The interactive visualization below lets you see how the linear approximation (an ellipse) converges to the true image (a deformed curve) as the perturbation shrinks.

Let $g: \mathbb{R}^n \to \mathbb{R}^k$ be differentiable at $a \in \mathbb{R}^n$, and let $f: \mathbb{R}^k \to \mathbb{R}^m$ be differentiable at $g(a) \in \mathbb{R}^k$. Then the composition $f \circ g: \mathbb{R}^n \to \mathbb{R}^m$ is differentiable at $a$, and

$J_{f \circ g}(a) = J_f(g(a)) \cdot J_g(a).$

The Jacobian of the composition is the matrix product of the Jacobians, evaluated at the appropriate points. The dimensions are consistent: $J_f$ is $m \times k$, $J_g$ is $k \times n$, so $J_f \cdot J_g$ is $m \times n$ — matching the expected size for a function from $\mathbb{R}^n$ to $\mathbb{R}^m$.

This is the most important proof in this topic. We track the error terms through two layers of the total derivative.

By differentiability of $g$ at $a$, we can write

$g(a + h) = g(a) + J_g(a)h + \varepsilon_g(h)\|h\|$

where $\varepsilon_g(h) \to 0$ as $h \to 0$ (here $\varepsilon_g(h)$ is a vector in $\mathbb{R}^k$). By differentiability of $f$ at $b = g(a)$, we can write

$f(b + \ell) = f(b) + J_f(b)\ell + \varepsilon_f(\ell)\|\ell\|$

where $\varepsilon_f(\ell) \to 0$ as $\ell \to 0$. Set $\ell = g(a+h) - g(a) = J_g(a)h + \varepsilon_g(h)\|h\|$. Then

$(f \circ g)(a+h) = f(g(a) + \ell) = f(b) + J_f(b)\ell + \varepsilon_f(\ell)\|\ell\|.$

Substituting the expression for $\ell$:

$= f(b) + J_f(b)\bigl[J_g(a)h + \varepsilon_g(h)\|h\|\bigr] + \varepsilon_f(\ell)\|\ell\|$

$= f(b) + J_f(b)J_g(a)h + \underbrace{J_f(b)\varepsilon_g(h)\|h\| + \varepsilon_f(\ell)\|\ell\|}_{\text{must show this is } o(\|h\|)}.$

First error term: $\|J_f(b)\varepsilon_g(h)\|h\|\| \le \|J_f(b)\| \cdot \|\varepsilon_g(h)\| \cdot \|h\|$. The operator norm $\|J_f(b)\|$ is a fixed constant, and $\|\varepsilon_g(h)\| \to 0$ as $h \to 0$, so this entire term is $o(\|h\|)$.

Second error term: We need $\|\ell\|$. From the definition: $\|\ell\| = \|J_g(a)h + \varepsilon_g(h)\|h\|\| \le \|J_g(a)\| \cdot \|h\| + \|\varepsilon_g(h)\| \cdot \|h\|$. For $h$ sufficiently small, $\|\varepsilon_g(h)\| < 1$, so $\|\ell\| \le (\|J_g(a)\| + 1)\|h\| = C\|h\|$ for a constant $C$. As $h \to 0$, we have $\ell \to 0$, so $\varepsilon_f(\ell) \to 0$. Therefore $\|\varepsilon_f(\ell)\| \cdot \|\ell\| \le \|\varepsilon_f(\ell)\| \cdot C\|h\| = o(\|h\|)$.

Combining both error terms, we have shown

$(f \circ g)(a+h) = (f \circ g)(a) + J_f(g(a)) \cdot J_g(a) \cdot h + o(\|h\|),$

which proves that $f \circ g$ is differentiable at $a$ with $J_{f \circ g}(a) = J_f(g(a)) \cdot J_g(a)$. $\square$

Let $g: \mathbb{R}^2 \to \mathbb{R}^2$ by $g(x, y) = (x^2 + y, xy)$ and $f: \mathbb{R}^2 \to \mathbb{R}$ by $f(u, v) = u^2 + v^2$. Their Jacobians are

$J_g(x,y) = \begin{pmatrix} 2x & 1 \\ y & x \end{pmatrix}, \qquad J_f(u,v) = \begin{pmatrix} 2u & 2v \end{pmatrix}.$

At $(x,y) = (1,1)$: $g(1,1) = (2, 1)$, so $J_g(1,1) = \begin{pmatrix} 2 & 1 \\ 1 & 1 \end{pmatrix}$ and $J_f(2,1) = \begin{pmatrix} 4 & 2 \end{pmatrix}$.

Chain rule:

$J_{f \circ g}(1,1) = \begin{pmatrix} 4 & 2 \end{pmatrix} \begin{pmatrix} 2 & 1 \\ 1 & 1 \end{pmatrix} = \begin{pmatrix} 10 & 6 \end{pmatrix}.$

Verification by direct computation: $(f \circ g)(x,y) = (x^2+y)^2 + (xy)^2$. Then $\frac{\partial}{\partial x}(f \circ g) = 2(x^2+y)(2x) + 2(xy)(y) = 4x(x^2+y) + 2xy^2$. At $(1,1)$: $4(1)(2) + 2(1)(1) = 10$. Confirmed.

For a three-layer chain $h \circ g \circ f$:

$J_{h \circ g \circ f}(a) = J_h(g(f(a))) \cdot J_g(f(a)) \cdot J_f(a).$

The chain rule applies recursively — each Jacobian is evaluated at the appropriate intermediate value. This is the prototype for backpropagation through $K$ layers: the end-to-end Jacobian is $J_L \cdot J_{f_K} \cdot J_{f_{K-1}} \cdots J_{f_1}$, with each factor evaluated at the activation from the forward pass.

The chain rule is matrix multiplication because derivatives are linear maps, and composition of linear maps corresponds to matrix multiplication. The single-variable chain rule $(f \circ g)' = f' \cdot g'$ multiplies $1 \times 1$ matrices — it just happens to look like scalar multiplication. In $\mathbb{R}^n$, the matrix structure is exposed, and the product $J_f \cdot J_g$ captures how sensitivities propagate through layers of computation.

For $f: \mathbb{R}^n \to \mathbb{R}^n$ with all partial derivatives existing at $a$, the Jacobian determinant is

$\det J_f(a) = \det \begin{pmatrix} \frac{\partial f_1}{\partial x_1} & \cdots & \frac{\partial f_1}{\partial x_n} \\ \vdots & \ddots & \vdots \\ \frac{\partial f_n}{\partial x_1} & \cdots & \frac{\partial f_n}{\partial x_n} \end{pmatrix}(a).$

The Jacobian determinant is defined only when the Jacobian is square ($m = n$).

Let $f: \mathbb{R}^n \to \mathbb{R}^n$ be continuously differentiable near $a$, and let $R$ be a small rectangular region containing $a$. Then the image $f(R)$ has volume approximately

$\text{vol}(f(R)) = |\det J_f(a)| \cdot \text{vol}(R) + o(\text{vol}(R))$

as the diameter of $R$ shrinks to zero. The absolute value handles orientation: $\det J_f(a) > 0$ means $f$ preserves orientation; $\det J_f(a) < 0$ means it reverses orientation (like a reflection).

For $f(r, \theta) = (r\cos\theta, r\sin\theta)$:

$\det J_f = \det \begin{pmatrix} \cos\theta & -r\sin\theta \\ \sin\theta & r\cos\theta \end{pmatrix} = r\cos^2\theta + r\sin^2\theta = r.$

This is the $r$ in "$r\,dr\,d\theta$" for polar integration. A small $dr \times d\theta$ rectangle in polar coordinates maps to a region of area approximately $r\,dr\,d\theta$ in Cartesian coordinates. Far from the origin ($r$ large), the same angular increment covers more Cartesian area; near the origin ($r$ small), less.

For $f(x, y) = (2x, 3y)$: $J_f = \begin{pmatrix} 2 & 0 \\ 0 & 3 \end{pmatrix}$, $\det J_f = 6$. Every region’s area is multiplied by 6 — stretching by 2 in $x$ and by 3 in $y$ scales areas by $2 \times 3 = 6$. The determinant of a diagonal matrix is the product of the diagonal entries, and each entry is a stretch factor along one axis.

For $f, g: \mathbb{R}^n \to \mathbb{R}^n$,

$\det J_{f \circ g}(a) = \det J_f(g(a)) \cdot \det J_g(a).$

This follows from the chain rule (Theorem 2) and the multiplicativity of determinants:

$\det J_{f \circ g}(a) = \det\bigl(J_f(g(a)) \cdot J_g(a)\bigr) = \det J_f(g(a)) \cdot \det J_g(a).$

Geometrically: volume distortions compose multiplicatively. If $g$ scales volume by factor 3 and $f$ scales volume by factor 2, then $f \circ g$ scales volume by factor 6. $\square$

$\det J_f(a) = 0$ means the linear approximation $J_f(a)$ is not invertible — it “collapses” a neighborhood of $a$ into a lower-dimensional set. This is the boundary between the Inverse Function Theorem applying (nonzero determinant implies $f$ is locally invertible) and failing (zero determinant implies a critical value). Inverse & Implicit Function Theorems develops this connection.

The pattern we have seen — $\int_{f(U)} h(x)\,dx = \int_U h(f(u)) \cdot |\det J_f(u)|\,du$ — is formalized as the Change of Variables theorem in the Multivariable Integral Calculus track. That topic provides the rigorous justification. This topic gives you the differential machinery: the Jacobian determinant is the local volume scaling factor, and integration sums up these local contributions.