Calculus for Deep Learning
Training a neural network means adjusting millions of numbers so that one number — the loss, which measures how wrong the network currently is — gets smaller. A central task is working out which way to push each of those millions of numbers.
Here calculus connects changes in the parameters to changes in the loss. We will build that connection through derivatives, partial derivatives, gradients, and the chain rule. The examples use basic algebra and Python; run the code blocks in order, since later blocks reuse earlier functions. NumPy is needed for the vector calculations.
A derivative is a slope
Forget the notation for a moment and ask a practical question. You have a machine with one dial. You turn the dial a little. How much does the output move?
If turning the dial by 1 moves the output by 6, the average rate of change over that interval is 6. A derivative describes the rate at a single setting: it is the limit of these ratios as the input change approaches zero, when that limit exists.
Watch that happen. Below, the machine squares its input, and we sit at the setting 3 and try smaller and smaller nudges.
def f(a):
return a ** 2
a = 3.0
for h in (1.0, 0.1, 0.01, 0.001):
moved = f(a + h) - f(a)
print(f" nudge input by {h:<6} output moves {moved:8.5f} ratio {moved / h:7.4f}")
# nudge input by 1.0 output moves 7.00000 ratio 7.0000
# nudge input by 0.1 output moves 0.61000 ratio 6.1000
# nudge input by 0.01 output moves 0.06010 ratio 6.0100
# nudge input by 0.001 output moves 0.00600 ratio 6.0010
The ratio is settling on 6. A nudge of 1 is too coarse and reports 7, because the machine’s sensitivity has already changed by the time you finish turning. Shrink the nudge and the answer converges. The derivative is the number the ratio is heading toward — here, exactly 6.
Written out, the derivative of \(f\) at the point \(a\) is
\[f'(a)=\lim_{h\to 0}\frac{f(a+h)-f(a)}{h}\]
The fraction is the ratio we just computed: output movement over input movement. For this example, \(((3+h)^2-3^2)/h=6+h\) for nonzero \(h\), so the limit is exactly 6. The limit uses both positive and negative changes in \(h\); the table sampled only positive ones. The notations \(f'(a)\) and \(df/da\) both denote the derivative here.
A derivative of 6 has a concrete reading you can use without thinking: increase the input by \(0.001\) and the output rises by about \(0.006\). A negative derivative means the output falls when the input rises.
For numerical checks, we will use a central difference: sample at \(a-h\) and \(a+h\), then divide the output difference by the distance \(2h\) between those inputs. For \(f(a)=a^2\), the numerator expands to \(4ah\), giving exactly \(2a\) in exact arithmetic. For other smooth functions it is generally an approximation; floating-point rounding also affects the calculation.
def slope(f, a, h=1e-6):
"""Central difference: a numerical stand-in for the derivative."""
return (f(a + h) - f(a - h)) / (2 * h)
print(round(slope(lambda a: 4 * a, 2.0), 6)) # 4.0 line: slope is constant
print(round(slope(lambda a: 4 * a, 9.0), 6)) # 4.0 same slope somewhere else
print(round(slope(lambda a: a ** 2, 2.0), 6)) # 4.0 equals 2a
print(round(slope(lambda a: a ** 2, 5.0), 6)) # 10.0 same function, new point
One caution about the helper. Very small steps can amplify rounding errors, because the subtraction \(f(a+h)-f(a-h)\) loses significant digits when the two values are nearly equal. With Python floats in this run, at \(h=10^{-14}\) the slope of \(a^2\) at 3 comes out as 6.128 instead of 6, and at \(h=10^{-16}\) it comes out as exactly 0. The appropriate step size depends on the function, the input scale, and the numerical precision; the default \(h=10^{-6}\) works adequately for the examples here.
Compare the first two lines with the last two. A straight line has one slope everywhere — move along it and nothing changes. A curve does not: the same function reports 4 at one setting and 10 at another.
During gradient-based training, we recompute the gradient after each update because it generally changes with the parameters. A changing gradient does not itself rule out a direct solution: linear least squares, for example, can be solved by matrix factorization. Neural-network training usually relies on iterative optimization. Even with piecewise-linear ReLU activations, the loss gradient can change as the weights change.
The handful of rules you actually need
Derivatives of common functions are looked up. This is the table this series draws on. The third column names where each rule turns up later — softmax, cross-entropy, LSTM and GRU are architectures and losses covered in their own articles, so treat those names as signposts rather than as things to know now.
| Function | Derivative | Where it appears |
|---|---|---|
| \(c\) (constant) | \(0\) | \(\partial(wx+b)/\partial x=w\): the fixed bias contributes 0; \(\partial(wx+b)/\partial b=1\) |
| \(a^{n}\) | \(na^{n-1}\) | squared error |
| \(e^{a}\) | \(e^{a}\) | softmax, sigmoid |
| \(\log a\), \(a>0\) (natural logarithm) | \(1/a\) | cross-entropy |
| \(f+g\) | \(f’+g’\) | summed losses |
| \(fg\) | \(f’g+fg’\) | gated units (LSTM, GRU) |
ReLU, \(\max(0,a)\), uses two of these rules: below zero it is constant, so its derivative is 0; above zero it equals \(a\), so its derivative is 1. At zero the two one-sided slopes disagree and the derivative does not exist. Implementations choose a backward-pass convention there, commonly 0. The sigmoid \(\sigma(z)=1/(1+e^{-z})\) needs the chain rule as well. The next block checks its derivative at three moderate inputs; this simple exponential implementation is not intended for very large negative inputs, where it can overflow.
import math
def sigmoid(z):
return 1 / (1 + math.exp(-z))
def sigmoid_prime(z): # assembled from the table + chain rule
s = sigmoid(z)
return s * (1 - s)
for z in (-2.0, 0.0, 1.5):
print(f" z={z:5.1f} from the rules {sigmoid_prime(z):.6f}"
f" numerically {slope(sigmoid, z):.6f}")
# z= -2.0 from the rules 0.104994 numerically 0.104994
# z= 0.0 from the rules 0.250000 numerically 0.250000
# z= 1.5 from the rules 0.149146 numerically 0.149146
The numerical checks agree with the formula at these test points. The derivative \(\sigma'(z)=\sigma(z)(1-\sigma(z))\) never exceeds \(0.25\) and falls toward zero as \(z\) moves away from the origin. Along a path through sigmoid units, their local derivative factors are multiplied together and can contribute to vanishing gradients; the exercises derive it algebraically and place it alongside the other factors involved.
Many dials at once: partial derivatives
Real machines have more than one dial. A shower has a hot tap and a cold tap, and the temperature depends on both. If someone asks how sensitive the temperature is to the hot tap, the question only has an answer once you agree not to touch the cold one.
A partial derivative measures the rate of change in one input while the others are held fixed. It is written with \(\partial\), as in \(\partial f/\partial x\). To calculate it, treat the other inputs as constants. For \(f(x,y)=x^2y+3y\), differentiating with respect to \(x\) gives \(2xy\); differentiating with respect to \(y\) gives \(x^2+3\).
def f(x, y):
return x ** 2 * y + 3 * y # df/dx = 2xy, df/dy = x^2 + 3
x0, y0, h = 2.0, 5.0, 1e-6
print(round((f(x0 + h, y0) - f(x0 - h, y0)) / (2 * h), 6)) # 20.0 = 2*2*5
print(round((f(x0, y0 + h) - f(x0, y0 - h)) / (2 * h), 6)) # 7.0 = 2^2 + 3
Each line changes one argument and leaves the other fixed. At \((x,y)=(2,5)\), the partial derivatives are 20 and 7. For equally small changes in these coordinates, changing \(x\) has the larger first-order effect. Such a comparison depends on the units and scaling of the inputs.
A network with 100 million parameters has 100 million partial derivatives, and each one is defined as though the other 99,999,999 were frozen constants. Computing them all is not 100 million separate passes: backpropagation gets its efficiency by reusing intermediate results across them.
The gradient
Collect all those partial derivatives into a single list and you have the gradient, written \(\nabla f\).
\[\nabla f=\left[\frac{\partial f}{\partial x_1},\;\frac{\partial f}{\partial x_2},\;\ldots,\;\frac{\partial f}{\partial x_n}\right]^{T}\]
On a hillside, the gradient describes the local uphill slope. More precisely, for a differentiable function with a nonzero gradient, it points toward the greatest first-order increase among directions of equal Euclidean length; its negative gives the greatest first-order decrease. This describes a local rate, not a guarantee that a finite step decreases the loss. At a zero gradient there is no preferred direction from this first-order information.
Basic gradient descent uses the update:
\[\theta \leftarrow \theta-\alpha\nabla_\theta J\]
Three symbols appear there. \(J\) is the loss, the single number being made smaller. \(\theta\) is the collection of parameters being adjusted. \(\alpha\) is the learning rate, the number that scales how large each update is. Read the line as an instruction: take the current parameters \(\theta\), measure which way is uphill, and subtract \(\alpha\) times the gradient. Repeat. Gradient descent is that one line, and the loop below runs it five times on a bowl-shaped function whose minimum sits at the origin.
import numpy as np
def J(theta): # a simple bowl, minimum at (0, 0)
return theta[0] ** 2 + 3 * theta[1] ** 2
def grad_J(theta):
return np.array([2 * theta[0], 6 * theta[1]])
theta = np.array([4.0, 2.0])
for step in range(5):
theta = theta - 0.1 * grad_J(theta)
print(step, np.round(theta, 4), round(J(theta), 4))
# 0 [3.2 0.8] 12.16
# 1 [2.56 0.32] 6.8608
# 2 [2.048 0.128] 4.2435
# 3 [1.6384 0.0512] 2.6922
# 4 [1.3107 0.0205] 1.7192
We can check the first update by hand. At \(\theta=[4,2]\) the gradient is \([2\times 4,\,6\times 2]=[8,12]\), so the update is \([4,2]-0.1\times[8,12]=[3.2,0.8]\), which is what the loop prints. The loss falls to 1.72 in five steps. Later optimizers change how gradients are converted into updates, while automatic differentiation supplies the gradients.
The two coordinates move at very different rates. Each step multiplies the first coordinate by \(1-0.1\times 2 = 0.8\) and the second by \(1-0.1\times 6 = 0.4\), so the second reaches 0.02 while the first is still at 1.31. That difference comes from curvature — how fast the slope itself changes as you move. The bowl bends three times faster along the second axis than the first, and one step size is being applied to both. Momentum and adaptive optimizers change the update rule to address difficulties such as uneven curvature; their behavior still depends on their settings.
The chain rule
Backpropagation applies the chain rule through the operations of a network. Suppose a crank turns a gear, and that gear turns a wheel. If the gear makes two turns for every turn of the crank, and the wheel makes five turns for every turn of the gear, then the wheel makes ten turns per crank. Rates multiply along a chain.
Written down: if \(y=f(u)\) and \(u=g(x)\), then
\[\frac{dy}{dx}=\frac{dy}{du}\cdot\frac{du}{dx}\]
# y = (3x + 1)^2 -> outer u^2, inner 3x + 1
def y(x): return (3 * x + 1) ** 2
x = 2.0
u = 3 * x + 1 # 7
dy_du = 2 * u # 14
du_dx = 3
print(dy_du * du_dx) # 42 by the chain rule
print(round(slope(y, x), 6)) # 42.0 numerically, as a check
A network composes functions in the same way, usually with vector inputs and outputs. The next block builds a three-link chain and computes the sensitivity of the output to the input in two ways: by multiplying the local slopes, and by nudging the whole thing.
def layer1(x): return 2 * x + 1 # -> u
def layer2(u): return u ** 2 # -> v
def layer3(v): return 3 * v # -> output
def whole(x): return layer3(layer2(layer1(x)))
x = 1.5
u = layer1(x)
v = layer2(u)
d3_dv = 3.0 # slope of 3v
d2_du = 2 * u # slope of u^2 at u
d1_dx = 2.0 # slope of 2x + 1
print(f" values along the chain : x={x} u={u} v={v}")
print(f" local slopes : {d3_dv} {d2_du} {d1_dx}")
print(f" multiplied together : {d3_dv * d2_du * d1_dx}")
print(f" numerical check : {round(slope(whole, x), 6)}")
# values along the chain : x=1.5 u=4.0 v=16.0
# local slopes : 3.0 8.0 2.0
# multiplied together : 48.0
# numerical check : 48.0
The three local slopes multiply to give the derivative of the composed scalar function. In reverse-mode differentiation, we start with the final output derivative equal to 1, then propagate it backward: here the successive derivatives are 3 with respect to v, 24 with respect to u, and 48 with respect to x. Backpropagation organizes these local operations efficiently for a network, reusing forward-pass values and adding contributions wherever paths meet.
Repeated multiplication can make derivatives very small or very large. Along a scalar path, ten factors of \(0.25\) give \(0.25^{10}\approx9.54\times10^{-7}\), while ten factors of 2 give 1024. These illustrate vanishing and exploding gradients. In a network the result also depends on weight matrices and contributions from multiple paths, so activation derivatives alone do not determine the full gradient.
One extension is needed for real networks. When a variable reaches the output through several different paths, the contributions add together, which is why gradients accumulate:
\[\frac{\partial L}{\partial x}=\sum_{i}\frac{\partial L}{\partial u_i}\cdot\frac{\partial u_i}{\partial x}\]
For example, let \(u=2x\), \(v=x^2\), and \(L=u+v\). At \(x=3\), the path through \(u\) contributes \(1\times2=2\), and the path through \(v\) contributes \(1\times(2\times3)=6\). Adding the two contributions gives \(dL/dx=8\).
Reverse-mode differentiation, which is how this is organized efficiently, is derived in Derivatives and Computation Graphs for Neural Network Learning. The matrix shapes assumed throughout are introduced in Linear Algebra for Deep Learning: Only What You Actually Need.
Local minima, global minima, and saddle points
A local minimum is a point whose value is no greater than the values at nearby points; a global minimum attains the lowest value over the entire domain, and several points can share that lowest value. A bowl has one local minimum and it is also the global one; a landscape with several dips has many local minima, only some of which reach the lowest value. A function is convex if, between any two points in its convex domain, its graph lies at or below the straight segment joining their graph values. For such a function, any local minimum is also a global minimum. However, the minimum need not be unique — \(f(x,y)=x^2\) is convex with a whole line of minima — and convergence still depends on the algorithm and its step size, as the second exercise shows. Training a multilayer neural network generally gives a non-convex objective in its weights, so a local minimum need not be global.
A saddle point is a stationary point with nearby values both above and below its value. For example, \(f(x,y)=x^2-y^2\) has gradient zero at the origin but rises along the \(x\)-axis and falls along the \(y\)-axis. Plain gradient descent stays at an exact stationary point. A perturbation or minibatch noise can help it leave a saddle, but this depends on whether the perturbation reaches a descending direction. Small gradients can also make progress slow in flat regions. Useful solutions do not require a guarantee of reaching the global optimum.
Where this foundation leads
You can follow the calculations above without integration techniques, formal limit proofs, or differential equations. Later topics introduce additional tools as they are needed.
Two related objects are useful later. The Jacobian is the matrix of partial derivatives when a function maps a vector to a vector; the Hessian is the matrix of second derivatives, which describes curvature. Forming either matrix in full is usually too expensive for a real network, so frameworks compute their products with a vector directly instead. Both are introduced where they matter.
As you move to larger networks, automatic differentiation will handle most of the derivative calculations. Understanding the chain rule remains essential for diagnosing why a gradient vanishes, explodes, or stops flowing. The optimizers mentioned above are compared in Deep Learning Optimizers: Mini-Batch, Momentum, RMSprop, and Adam.
Exercises
1. Compose the rules. Differentiate \(\sigma(z)=1/(1+e^{-z})\) using only the chain rule and the derivatives in the table, then show algebraically that the result equals \(\sigma(z)(1-\sigma(z))\).
You should get: an expression with \(e^{-z}\) in it that you can rewrite in terms of \(\sigma\) alone. The peak value is \(0.25\), at \(z=0\).
Solution
Write \(\sigma=u^{-1}\) with \(u=1+e^{-z}\). The power rule gives \(d\sigma/du=-u^{-2}\), and the exponential rule with the chain rule gives \(du/dz=-e^{-z}\). Multiplying the two, the minus signs cancel:
\[\frac{d\sigma}{dz}=\frac{e^{-z}}{(1+e^{-z})^{2}}\]
Now split that fraction deliberately. One factor is \(\frac{1}{1+e^{-z}}\), which is \(\sigma\) itself. The other is \(\frac{e^{-z}}{1+e^{-z}}\), and adding and subtracting 1 in the numerator shows it equals \(1-\sigma\). So the derivative is \(\sigma(1-\sigma)\), matching the numerical check in the article.
Two consequences follow. The local sigmoid derivative can be computed from the output \(s=\sigma(z)\), without recomputing it from \(z\). To obtain a loss gradient, the backward pass also needs the incoming derivative: \(dL/dz=(dL/ds)s(1-s)\). Saving intermediate values can avoid recomputation; which values are needed depends on the operation. And since \(\sigma\) lies in \((0,1)\), the product \(\sigma(1-\sigma)\) peaks at \(0.25\) when \(\sigma=0.5\) and shrinks toward zero at both ends. Across ten sigmoid layers the activation derivatives contribute a factor of at most \(0.25^{10}\), about \(10^{-6}\). That is only one part of the product: the weight matrices contribute their own factors, and gradients arriving by several paths add. Saturated sigmoids are one contributor to vanishing gradients rather than the whole account.
2. Step size matters. Run the gradient descent loop from this article on \(J(\theta)=\theta_0^2+3\theta_1^2\) with learning rates \(0.1\), \(0.3\), and \(0.4\). Describe what happens in each case and identify the convergence interval for positive \(\alpha\), including what happens at its endpoint.
You should get: one case approaching zero without sign changes, one oscillating with decreasing magnitude, and one oscillating with increasing magnitude. The threshold involves the coefficient 3.
Solution
import numpy as np
def run(alpha, steps=20):
t = np.array([4.0, 2.0])
for _ in range(steps):
t = t - alpha * np.array([2*t[0], 6*t[1]])
return t
for a in (0.1, 0.3, 0.4):
print(a, np.round(run(a), 4))
# 0.1 [0.0461 0. ]
# 0.3 [0. 0.0231]
# 0.4 [ 0. 1673.3651]
At \(\alpha=0.1\) both coordinates approach zero. At \(0.3\) they still approach zero, but the second changes sign each step. At \(0.4\) the second coordinate alternates sign with increasing magnitude; the displayed value is positive because 20 steps is even. Entries printed as 0 are rounded: the first coordinate after 20 steps at \(\alpha=0.4\) is \(4(0.2)^{20}\approx4.19\times10^{-14}\), not exactly zero.
The threshold can be derived exactly. For a one-dimensional quadratic \(c\theta^2\) with \(c>0\) and a nonzero starting coordinate the update is \(\theta \leftarrow \theta(1-2\alpha c)\), so the coordinate converges to the minimum at zero exactly when \(|1-2\alpha c|<1\), that is for \(0<\alpha<1/c\). The endpoint is excluded: at \(\alpha=1/c\) the factor is \(-1\) and the coordinate oscillates forever without converging — here \(\alpha=1/3\) sends the second coordinate to \(-2, 2, -2, 2\) and so on. The second coordinate has \(c=3\) and so needs \(\alpha<1/3\); the first has \(c=1\) and requires \(0<\alpha<1\), also excluding its endpoint.
To make both coordinates converge from this starting point, the learning rate must satisfy the stricter limit \(0<\alpha<1/3\). This limit comes from the larger curvature, 6. The contraction factors show how quickly each coordinate then approaches zero.
3. Where the gradient points. For \(J(\theta)=\theta_0^2+3\theta_1^2\) at \(\theta=[4,2]\), compute \(\nabla J\) and compare the descent direction \(-\nabla J\) with the straight line from \(\theta\) to the minimum at the origin. Are they the same? Can any single step in this descent direction reach the origin?
Expected observation: two clearly different directions. The negative gradient does not point directly toward the minimum in this example.
Solution
\(\nabla J=[8,12]\), which normalizes to about \([0.55,0.83]\), so the descent direction \(-\nabla J\) normalizes to about \([-0.55,-0.83]\). The straight line to the origin is \(-[4,2]\), normalizing to about \([-0.89,-0.45]\). If the first coordinate is horizontal and the second vertical, the descent direction is more vertical than the direction to the origin. No single learning rate reaches the origin from this point: the first coordinate would require \(\alpha=1/2\), while the second would require \(\alpha=1/6\).
The gradient is the direction of steepest local increase, which is not the same as an arrow pointing at the optimum. On an elongated bowl it points across the valley more than along it, so progress along the shallow direction is slow. Whether the path also zigzags depends on the step size: it does so once \(\alpha\) exceeds \(1/(2c)\) for some direction, which makes that coordinate change sign each step. At \(\alpha=0.1\) here both factors, \(0.8\) and \(0.4\), stay positive and the path does not oscillate — it is simply slow.
For this bowl, the larger curvature limits the learning rate, and the contraction factors determine the remaining error in each coordinate. Momentum combines the current gradient with a decaying accumulation of past gradients. Alternating components can partly cancel while consistent components accumulate, but cancellation is not guaranteed; a poorly chosen learning rate or momentum coefficient can increase oscillations or cause divergence.
Discover more from Insightful Data Lab
Subscribe to get the latest posts sent to your email.
