Calculus & Optimization: Derivatives, Gradients, and Higher-Order Approximations

Welcome to the foundational bedrock of modern machine learning and artificial intelligence! Whether you are training a deep neural network with billions of parameters or optimizing a simple linear regression model, calculus is the engine under the hood.

In this guide, we will journey from the absolute basics of single-variable rates of change to the multidimensional matrices that power state-of-the-art AI optimization algorithms.

---

Curriculum Checklist & Milestones

Before diving in, here is our roadmap. By the end of this tutorial, you will master every single one of these core optimization building blocks:

• Functions and Limits: Understanding continuous mappings and the concept of approaching infinity/zero.
• Derivatives (Single Variable): Calculating instantaneous rates of change and tangent lines.
• Partial Derivatives: Measuring change along individual coordinate axes in multidimensional space.
• Gradient Vectors: Assembling partial derivatives into directional steepness vectors.
• Directional Derivatives: Calculating rate of change along any arbitrary vector path.
• The Chain Rule: Propagating derivatives through nested function compositions (the core of Backpropagation).
• Jacobians and Hessians: Structuring first-order systems and second-order curvature matrices.
• Taylor Series Approximation: Linearizing and quadratically approximating complex loss landscapes locally.

---

1. The Intuitive Mental Model: The Foggy Mountain Hiker

Imagine you are standing on a rugged, foggy mountain range late in the evening. You cannot see the landscape as a whole—visibility is limited to the exact spot where your boots touch the ground. Your ultimate objective is to reach the lowest possible valley (the lowest point of loss in a neural network).

Here is how each mathematical concept relates to your hike:

1. Function $f(x, y)$: The elevation map of the mountain. Input coordinates $(x, y)$ give you your altitude $z$.
2. Limit: Taking smaller and smaller micro-steps with your boots until you feel the continuous slope beneath your feet.
3. Single-Variable Derivative: The slope of the ground if you are locked into walking strictly East-West on a narrow ridge.
4. Partial Derivative: The slope under your right foot along the North-South axis while keeping your East-West position completely still.
5. Gradient Vector ($\nabla f$): A compass needle pointing in the direction of steepest ascent (uphill). If you take a step in the exact opposite direction ($-\nabla f$), you take the steepest path downhill.
6. Directional Derivative: The precise steepness you feel if you decide to walk at a custom angle, say $30^\circ$ North-East.
7. Chain Rule: How changes in your walking pace directly translate to changes in altitude over time.
8. Hessian Matrix: The curvature of the mountain under your boots. Is the ground shaped like a smooth bowl, a flat ramp, or a treacherous saddle-shaped pass?
9. Taylor Series: Drawing a local mini-map on a notepad based on your foot's current elevation, slope, and ground curvature.

---

2. Functions, Limits, and Single-Variable Derivatives

Functions and Limits

A function $f$ takes an input $x$ from a domain and maps it to a single output $y = f(x)$. In machine learning, $x$ could be model weights, and $y$ could be prediction error (loss).

A limit evaluates what value a function approaches as its input gets infinitely close to a specific point $a$:

$\lim_{x \to a} f(x) = L$

Limits allow us to deal with instantaneous states without dividing by zero.

Single-Variable Derivatives

The average rate of change of a function over an interval $[x, x + h]$ is the slope of the secant line:

$\text{Average Rate of Change} = \frac{f(x + h) - f(x)}{h}$

To make this rate instantaneous, we take the limit as the step size $h$ shrinks to zero. This limit defines the derivative $f'(x)$ or $\frac{df}{dx}$:

$f'(x) = \lim_{h \to 0} \frac{f(x + h) - f(x)}{h}$

Geometrically, $f'(x)$ is the slope of the line tangent to $f(x)$ at point $x$.

Common Single-Variable Rules

• Power Rule: $\frac{d}{dx}[x^n] = n x^{n-1}$
• Exponential Rule: $\frac{d}{dx}[e^x] = e^x$
• Sum Rule: $\frac{d}{dx}[f(x) + g(x)] = f'(x) + g'(x)$

---

3. Multivariable Calculus: Partial Derivatives, Gradients, and Directional Derivatives

In machine learning, loss functions depend on millions of variables (weights), not just one. We must extend single-variable derivatives to multivariable space $\mathbb{R}^n \to \mathbb{R}$.

Partial Derivatives

Suppose elevation is given by $f(x, y) = x^2 + 3xy + y^2$. A partial derivative measures how $f$ changes with respect to one variable while holding all other variables constant.

• Partial derivative with respect to $x$: $\frac{\partial f}{\partial x} = \lim_{h \to 0} \frac{f(x + h, y) - f(x, y)}{h}$
• Partial derivative with respect to $y$: $\frac{\partial f}{\partial y} = \lim_{h \to 0} \frac{f(x, y + h) - f(x, y)}{h}$

Example Calculation

For $f(x, y) = x^2 + 3xy + y^2$:

• Treating $y$ as a constant: $\frac{\partial f}{\partial x} = 2x + 3y$
• Treating $x$ as a constant: $\frac{\partial f}{\partial y} = 3x + 2y$

The Gradient Vector ($\nabla f$)

The gradient of a scalar-valued function $f(\mathbf{x})$ is a vector collecting all its first-order partial derivatives:

$\nabla f(\mathbf{x}) = \begin{bmatrix} \frac{\partial f}{\partial x_1} \\ \frac{\partial f}{\partial x_2} \\ \vdots \\ \frac{\partial f}{\partial x_n} \end{bmatrix}$

Two Fundamental Properties of the Gradient Vector:

1. Direction: $\nabla f(\mathbf{x})$ points in the direction of maximum rate of increase of $f$ at point $\mathbf{x}$.
2. Magnitude: $\|\nabla f(\mathbf{x})\|$ represents the maximum rate of increase per unit distance.

Consequently, $-\nabla f(\mathbf{x})$ points in the direction of steepest descent, forming the backbone of Gradient Descent optimization!

Directional Derivatives

What if you don't move purely along the $x$-axis or $y$-axis, but along a unit direction vector $\mathbf{u} = [u_1, u_2]^T$ (where $\|\mathbf{u}\| = 1$)?

The directional derivative $D_{\mathbf{u}} f(\mathbf{x})$ measures the rate of change of $f$ along vector $\mathbf{u}$:

$D_{\mathbf{u}} f(\mathbf{x}) = \nabla f(\mathbf{x}) \cdot \mathbf{u} = \sum_{i=1}^n \frac{\partial f}{\partial x_i} u_i$

Using the dot product formula $\mathbf{a} \cdot \mathbf{b} = \|\mathbf{a}\| \|\mathbf{b}\| \cos(\theta)$:

$D_{\mathbf{u}} f(\mathbf{x}) = \|\nabla f(\mathbf{x})\| \|\mathbf{u}\| \cos(\theta) = \|\nabla f(\mathbf{x})\| \cos(\theta)$

• Maximum increase occurs when $\theta = 0^\circ$ ($\cos(\theta) = 1$), meaning $\mathbf{u}$ aligns with $\nabla f$.
• Maximum decrease occurs when $\theta = 180^\circ$ ($\cos(\theta) = -1$), meaning $\mathbf{u}$ aligns with $-\nabla f$.
• Zero change occurs when $\mathbf{u}$ is orthogonal ($\theta = 90^\circ$) to $\nabla f$ (walking along a contour line).

---

4. Function Composition and the Chain Rule

In deep learning, neural networks are stacked compositions of functions: $L(w) = f_n(f_{n-1}(\dots f_1(w)\dots))$. To calculate gradients through deep layers, we use the Chain Rule.

Single-Variable Chain Rule

If $y = g(x)$ and $z = f(y)$, then $z = f(g(x))$. The derivative of $z$ with respect to $x$ is:

$\frac{dz}{dx} = \frac{dz}{dy} \cdot \frac{dy}{dx}$

Multivariable Chain Rule

If $z = f(u, v)$, where $u = g(x, y)$ and $v = h(x, y)$, then the partial derivative of $z$ with respect to $x$ is formed by summing all computational pathways:

$\frac{\partial z}{\partial x} = \frac{\partial z}{\partial u} \frac{\partial u}{\partial x} + \frac{\partial z}{\partial v} \frac{\partial v}{\partial x}$

This exact mechanism—accumulating partial derivatives backwards through computational graphs—is called Backpropagation in modern AI frameworks like PyTorch and TensorFlow.

---

5. Higher-Order Structures: Jacobians, Hessians, and Taylor Series

The Jacobian Matrix

When a function outputs a vector instead of a single scalar ($\mathbf{f}: \mathbb{R}^n \to \mathbb{R}^m$), its first derivatives form an $m \times n$ matrix called the Jacobian:

$\mathbf{J} = \begin{bmatrix} \frac{\partial f_1}{\partial x_1} & \frac{\partial f_1}{\partial x_2} & \dots & \frac{\partial f_1}{\partial x_n} \\ \frac{\partial f_2}{\partial x_1} & \frac{\partial f_2}{\partial x_2} & \dots & \frac{\partial f_2}{\partial x_n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial f_m}{\partial x_1} & \frac{\partial f_m}{\partial x_2} & \dots & \frac{\partial f_m}{\partial x_n} \end{bmatrix}$

The Jacobian transformed vector outputs are essential for vector-valued neural network layers, transformations, and coordinate changes.

The Hessian Matrix

The Hessian matrix $\mathbf{H}$ contains all second-order partial derivatives of a scalar-valued function $f: \mathbb{R}^n \to \mathbb{R}$:

$\mathbf{H} = \begin{bmatrix} \frac{\partial^2 f}{\partial x_1^2} & \frac{\partial^2 f}{\partial x_1 \partial x_2} & \dots & \frac{\partial^2 f}{\partial x_1 \partial x_n} \\ \frac{\partial^2 f}{\partial x_2 \partial x_1} & \frac{\partial^2 f}{\partial x_2^2} & \dots & \frac{\partial^2 f}{\partial x_2 \partial x_n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial^2 f}{\partial x_n \partial x_1} & \frac{\partial^2 f}{\partial x_n \partial x_2} & \dots & \frac{\partial^2 f}{\partial x_n^2} \end{bmatrix}$

Key Properties of the Hessian:

1. Symmetry: By Schwarz / Clairaut's theorem, if second derivatives are continuous, then $\frac{\partial^2 f}{\partial x_i \partial x_j} = \frac{\partial^2 f}{\partial x_j \partial x_i}$, so $\mathbf{H} = \mathbf{H}^T$.
2. Curvature: The eigenvalues of $\mathbf{H}$ describe local landscape curvature:
   • All positive eigenvalues ($\mathbf{H}$ is Positive Definite): Local minimum (bowl shape).
   • All negative eigenvalues ($\mathbf{H}$ is Negative Definite): Local maximum (dome shape).
   • Mixed sign eigenvalues: Saddle point (pass between peaks).

Taylor Series Approximation

The Taylor Series allows us to approximate any smooth, complex non-linear function near a local point $\mathbf{x}_0$ using polynomials.

1D Taylor Series Expansion

$f(x) \approx f(x_0) + f'(x_0)(x - x_0) + \frac{1}{2!} f''(x_0)(x - x_0)^2 + \dots + \frac{1}{k!} f^{(k)}(x_0)(x - x_0)^k$

Multivariable Second-Order Expansion

Letting $\mathbf{\Delta x} = \mathbf{x} - \mathbf{x}_0$:

$f(\mathbf{x}_0 + \mathbf{\Delta x}) \approx \underbrace{f(\mathbf{x}_0)}_{\text{0th-order (Value)}} + \underbrace{\nabla f(\mathbf{x}_0)^T \mathbf{\Delta x}}_{\text{1st-order (Linear Gradient)}} + \underbrace{\frac{1}{2} \mathbf{\Delta x}^T \mathbf{H}(\mathbf{x}_0) \mathbf{\Delta x}}_{\text{2nd-order (Quadratic Curvature)}}$

• 1st-order truncation: Gives us standard Gradient Descent optimization models.
• 2nd-order truncation: Gives us Newton-Raphson optimization methods, which scale much faster toward local minima by incorporating terrain curvature!

---

6. Hands-On Python Implementation from Scratch with NumPy

Let me demonstrate how to calculate every single derivative concept covered above using NumPy, combining numerical approximations with analytical formulas. This lightweight version avoids importing and initializing PyTorch, reducing startup overhead in execution environments with a short timeout. The analytical checks below are explicit formulas, not automatic differentiation.

First, we illustrate a limit and a single-variable derivative using $q(t)=t^2$:

$\lim_{t\to2}q(t)=4,\qquad q'(2)=\lim_{h\to0}\frac{(2+h)^2-4}{h}=\lim_{h\to0}(4+h)=4.$

The multivariable code uses a slightly different surface from Section 3:

$$f(x,y)=x^2+3xy+2y^2,\qquad \nabla f=\begin{bmatrix}2x+3y\\3x+4y\end{bmatrix},\qquad H=\begin{bmatrix}2&3\\3&4\end{bmatrix}.$$

For the vector function $\mathbf{F}(x,y)=[f(x,y),\,5x+\sin y]^T$, its Jacobian is:

$$J_{\mathbf{F}}(x,y)=\begin{bmatrix}2x+3y&3x+4y\\5&\cos y\end{bmatrix}.$$

Finally, an explicit chain-rule example uses $z(t)=e^{t^2}$, so $z'(t)=e^{t^2}\cdot2t$.

import numpy as np

# =====================================================================
# SECTION 1: Pure NumPy Implementation of Fundamental Calculus Concepts
# =====================================================================

def single_variable_function(t):
    return t**2


print("=== FUNCTIONS, LIMITS, AND SINGLE-VARIABLE DERIVATIVES ===")
# Sample both sides of t=2. Numerical samples illustrate, but do not
# prove, the limit. The algebra above establishes it.
for h in (0.1, 0.01, 0.001):
    left = single_variable_function(2.0 - h)
    right = single_variable_function(2.0 + h)
    quotient = (
        single_variable_function(2.0 + h)
        - single_variable_function(2.0)
    ) / h
    print(
        f"h={h:.3f}: q(2-h)={left:.6f}, q(2+h)={right:.6f}, "
        f"forward slope={quotient:.6f}"
    )
print("Exact limit q(t) as t -> 2: 4; exact derivative q'(2): 4")


def scalar_function(x):
    """
    Sample 2D Loss Surface: f(x, y) = x^2 + 3*x*y + 2*y^2
    Analytical Gradient: [2*x + 3*y, 3*x + 4*y]
    """
    return x[0]**2 + 3.0 * x[0] * x[1] + 2.0 * x[1]**2


def finite_difference_gradient(f, x, h=1e-5):
    """
    Calculates the numerical gradient vector using finite differences.
    Formula: df/dx_i ≈ (f(x + h*e_i) - f(x - h*e_i)) / (2*h)
    """
    grad = np.zeros_like(x, dtype=np.float64)
    for i in range(len(x)):
        x_forward = np.copy(x)
        x_backward = np.copy(x)

        x_forward[i] += h
        x_backward[i] -= h

        # Central difference quotient calculation
        grad[i] = (f(x_forward) - f(x_backward)) / (2.0 * h)
    return grad


def directional_derivative(f, x, direction_vector):
    """
    Calculates directional derivative: D_u f(x) = grad(f) · u
    """
    grad = finite_difference_gradient(f, x)
    # Ensure direction vector is normalized (unit length)
    unit_u = direction_vector / np.linalg.norm(direction_vector)
    return np.dot(grad, unit_u)


def compute_hessian_finite_difference(f, x, h=1e-4):
    """
    Calculates the 2nd-order derivative matrix (Hessian) numerically.
    """
    n = len(x)
    hessian = np.zeros((n, n), dtype=np.float64)

    for i in range(n):
        for j in range(n):
            if i == j:
                # Second partial derivative w.r.t same variable
                x_plus = np.copy(x)
                x_minus = np.copy(x)
                x_plus[i] += h
                x_minus[i] -= h
                hessian[i, j] = (f(x_plus) - 2.0 * f(x) + f(x_minus)) / (h**2)
            else:
                # Mixed second partial derivatives
                x_pp, x_pm = np.copy(x), np.copy(x)
                x_mp, x_mm = np.copy(x), np.copy(x)

                x_pp[i] += h; x_pp[j] += h
                x_pm[i] += h; x_pm[j] -= h
                x_mp[i] -= h; x_mp[j] += h
                x_mm[i] -= h; x_mm[j] -= h

                hessian[i, j] = (f(x_pp) - f(x_pm) - f(x_mp) + f(x_mm)) / (4.0 * h**2)
    return hessian


def taylor_expansion_2nd_order(f, x0, delta_x):
    """
    Computes 2nd-order Taylor approximation at point (x0 + delta_x):
    f(x0 + Δx) ≈ f(x0) + ∇f(x0)^T Δx + 0.5 * Δx^T H(x0) Δx
    """
    val_x0 = f(x0)
    grad = finite_difference_gradient(f, x0)
    hessian = compute_hessian_finite_difference(f, x0)

    linear_term = np.dot(grad, delta_x)
    quadratic_term = 0.5 * np.dot(delta_x.T, np.dot(hessian, delta_x))

    approximation = val_x0 + linear_term + quadratic_term
    return approximation


# --- Evaluation of NumPy Implementations ---
point = np.array([2.0, 3.0])
print("\n=== NUMPY NUMERICAL CALCULUS ===")
print(f"Point of interest x0: {point}")
print(f"Exact Function Value f(x0): {scalar_function(point)}")

grad_num = finite_difference_gradient(scalar_function, point)
print(f"Numerical partial df/dx: {grad_num[0]:.6f}")
print(f"Numerical partial df/dy: {grad_num[1]:.6f}")
print(f"Numerical Gradient vector: {grad_num}")

# Arbitrary compass direction: [1.0, 1.0]
direction = np.array([1.0, 1.0])
dir_deriv = directional_derivative(scalar_function, point, direction)
print(f"Directional Derivative along [1, 1]: {dir_deriv:.4f}")

hess_num = compute_hessian_finite_difference(scalar_function, point)
print(f"Numerical Hessian Matrix:\n{hess_num}")

step = np.array([0.1, -0.05])
true_val_at_step = scalar_function(point + step)
approx_val = taylor_expansion_2nd_order(scalar_function, point, step)
print(f"\nExact f(x0 + step): {true_val_at_step:.6f}")
print(f"2nd-Order Taylor Approximation: {approx_val:.6f}")


# =====================================================================
# SECTION 2: Analytical Checks and an Explicit Chain Rule
# =====================================================================

print("\n=== ANALYTICAL CHECKS WITH NUMPY ===")

def analytical_gradient(x):
    return np.array([
        2.0 * x[0] + 3.0 * x[1],
        3.0 * x[0] + 4.0 * x[1],
    ])


# Vector function mapping R^2 -> R^2 to showcase Jacobians
def vector_function(x):
    # f1(x, y) = x^2 + 3*x*y + 2*y^2  (Scalar output 1)
    # f2(x, y) = 5*x + sin(y)         (Scalar output 2)
    f1 = scalar_function(x)
    f2 = 5.0 * x[0] + np.sin(x[1])
    return np.array([f1, f2])


def analytical_jacobian(x):
    # Each row contains the partial derivatives of one output.
    return np.array([
        [2.0 * x[0] + 3.0 * x[1], 3.0 * x[0] + 4.0 * x[1]],
        [5.0, np.cos(x[1])],
    ])


grad_exact = analytical_gradient(point)
jacobian_matrix = analytical_jacobian(point)
hessian_matrix = np.array([[2.0, 3.0], [3.0, 4.0]])

print(f"Analytical Gradient: {grad_exact}")
print(f"Vector Function Value: {vector_function(point)}")
print(f"Analytical Jacobian Matrix:\n{jacobian_matrix}")
print(f"Analytical Hessian Matrix:\n{hessian_matrix}")

np.testing.assert_allclose(grad_num, grad_exact, atol=1e-6)
np.testing.assert_allclose(hess_num, hessian_matrix, atol=1e-4)
np.testing.assert_allclose(approx_val, true_val_at_step, atol=1e-6)

# Explicit chain rule: u(t) = t^2, z(u) = exp(u)
# dz/dt = (dz/du) * (du/dt) = exp(t^2) * 2*t
t = 0.5
u = t**2
dz_du = np.exp(u)
du_dt = 2.0 * t
chain_derivative = dz_du * du_dt

h = 1e-5
numerical_chain_derivative = (
    np.exp((t + h)**2) - np.exp((t - h)**2)
) / (2.0 * h)

print(f"\nChain-rule derivative of exp(t^2) at t={t}: {chain_derivative:.6f}")
print(f"Finite-difference check: {numerical_chain_derivative:.6f}")
np.testing.assert_allclose(
    chain_derivative, numerical_chain_derivative, atol=1e-6
)

Python 3.12 (Pyodide Wasm)

At $(2,3)$, the gradient is $[13,18]^T$, and the directional derivative along the normalized direction $[1,1]^T$ is $31/\sqrt{2}\approx21.9203$. Both the true value and the second-order Taylor calculation at the chosen step are approximately 40.400000. Because this surface is quadratic, its second-order Taylor polynomial is exact; any small discrepancy here comes from numerical differentiation and floating-point arithmetic.

---

7. Mathematical Reference Summary Table

To solidify your mental model, here is a quick reference table comparing all the fundamental calculus concepts:

Concept	Mathematical Symbol	Input Dimension	Output Dimension	Core Intuition / Role
Derivative	$\frac{df}{dx}$	$\mathbb{R}^1$	$\mathbb{R}^1$	Instantaneous rate of change / Slope of single variable
Partial Derivative	$\frac{\partial f}{\partial x_i}$	$\mathbb{R}^n$	$\mathbb{R}^1$	Slope along one specific coordinate axis
Gradient Vector	$\nabla f$	$\mathbb{R}^n$	$\mathbb{R}^n$	Direction and magnitude of steepest ascent
Directional Derivative	$D_{\mathbf{u}}f$	$\mathbb{R}^n$	$\mathbb{R}^1$	Instantaneous rate of change in custom direction $\mathbf{u}$
Chain Rule	$\frac{dz}{dx} = \frac{dz}{dy}\frac{dy}{dx}$	Variable	Variable	Backpropagating rates of change through nested functions
Jacobian Matrix	$\mathbf{J}$	$\mathbb{R}^n$	$\mathbb{R}^{m \times n}$	Matrix of all 1st-order partial derivatives for vector functions
Hessian Matrix	$\mathbf{H}$	$\mathbb{R}^n$	$\mathbb{R}^{n \times n}$	Matrix of all 2nd-order partial derivatives (terrain curvature)
Taylor Series	$P_2(\mathbf{x})$	$\mathbb{R}^n$	$\mathbb{R}^1$	Polynomial approximation of local surface terrain

---

Conclusion & Next Steps

Congratulations! You have covered the full mathematical progression of calculus for optimization:

1. Started with single-variable limits and derivatives.
2. Expanded to multivariable landscapes via partial derivatives, gradient vectors, and directional derivatives.
3. Learned how the Chain Rule connects compound functions across deep computational graphs.
4. Explored matrix structures like the Jacobian and Hessian to measure multidimensional systems and local curvature.
5. Approximated non-linear surfaces using Taylor Series expansions.

You now possess the foundational calculus knowledge required to master Gradient Descent, Newton-Raphson Optimization, and Neural Network Backpropagation!