NumPy Gradient Calculation for Multivariable Functions

1. The Gradient: A Vector of Partial Derivatives

For a function with a single variable, like $f(x)=x^2$, the derivative $f^{\prime }(x)$ is a scalar that tells us the slope at any point $x$. But our loss functions will almost always depend on multiple parameters (like m and b in linear regression, or thousands of weights in a neural network).

When a function has multiple input variables, its "derivative" is called the gradient. The gradient is not a single number; it's a vector, where each component is the partial derivative with respect to one of the variables.

So, what is a partial derivative? It's simply the derivative with respect to one variable, while treating all other variables as if they were constants.

The complete gradient vector, $\nabla f$, points in the direction of the steepest increase of the function. To find the minimum, we'll want to move in the opposite direction, $-\nabla f$.

2. Analytical Gradient Calculation

The most direct way to find the gradient is to use the rules of calculus to find the formula for each partial derivative. This is known as the analytical method. It's fast and precise, but it requires you to figure out the derivatives beforehand.

Let's watch a video that demonstrates this for a multivariable function. The author defines a function $z=f(x,y)$ and then provides the functions for its partial derivatives, which together form the gradient.

As you saw, for the function:
$z(x,y)=\frac{\sin (5x)\cos (5y)}{5}$
The gradient vector is:
$\nabla z=\left[\begin{array}{c}\frac{\partial z}{\partial x}\\ \frac{\partial z}{\partial y}\end{array}\right]=\left[\begin{array}{c}\cos (5x)\cos (5y)\\ -\sin (5x)\sin (5y)\end{array}\right]$

Here is how you would implement that calculate_gradient function in Python using NumPy. It takes a vector of inputs and returns a vector representing the gradient.

import numpy as np

def calculate_gradient_analytical(v):
    """
    Calculates the analytical gradient for the function z = sin(5x)cos(5y)/5.
    
    Args:
        v (np.array): A numpy array [x, y].
        
    Returns:
        np.array: The gradient vector [dz/dx, dz/dy].
    """
    x, y = v[0], v[1]
    
    # Partial derivative with respect to x: cos(5x)cos(5y)
    dz_dx = np.cos(5 * x) * np.cos(5 * y)
    
    # Partial derivative with respect to y: -sin(5x)sin(5y)
    dz_dy = -np.sin(5 * x) * np.sin(5 * y)
    
    return np.array([dz_dx, dz_dy])

# Example usage:
point = np.array([0.1, 0.2])
gradient_at_point = calculate_gradient_analytical(point)
print(f"The gradient at {point} is {gradient_at_point}")
# Expected output: The gradient at [0.1 0.2] is [0.8357326  -0.09531189]

3. Numerical Gradient Calculation: The "From Scratch" Approach

What if a function is too complex to differentiate analytically, or what if we want a general method that works for any function? We can estimate the gradient numerically. This approach aligns perfectly with your preference for "from scratch" implementations because it builds the concept from its mathematical definition.

The definition of a derivative is:
$f^{\prime }(x)={\lim }_{h\to 0}\frac{f(x+h)-f(x)}{h}$
We can approximate this by using a very small, non-zero value for $h$:
$f^{\prime }(x)\approx \frac{f(x+h)-f(x)}{h}$
We can do the exact same thing for partial derivatives. To find the partial derivative with respect to a variable $x_i$, we add a small $h$ only to $x_i$ and see how much the function's output changes.

Let's look at a resource that shows a clean Python implementation of this idea.

Here is that code, ready for you to use. This is a powerful tool you will keep in your toolkit.

import numpy as np

def partial_difference_quotient(f, v, i, h):
    """
    Computes the ith partial difference quotient of f at v.
    
    Args:
        f (callable): A multivariable function.
        v (np.array): The point (vector) at which to compute the quotient.
        i (int): The index of the variable to differentiate with respect to.
        h (float): A very small number.
        
    Returns:
        float: The partial difference quotient.
    """
    w = v.copy()
    w[i] = w[i] + h
    return (f(w) - f(v)) / h

def estimate_gradient(f, v, h=0.00001):
    """
    Estimates the gradient of f at v.
    
    Args:
        f (callable): A multivariable function.
        v (np.array): The point (vector) at which to estimate the gradient.
        h (float): A very small number.
        
    Returns:
        np.array: The estimated gradient vector.
    """
    return np.array([partial_difference_quotient(f, v, i, h) for i in range(len(v))])

# --- Example Usage ---
# Let's test it on a simple function: f(x, y) = x^2 + y^2
# The analytical gradient is [2x, 2y]

def my_function(v):
    x, y = v[0], v[1]
    return x**2 + y**2

point_to_test = np.array([3.0, 4.0])

# Analytical gradient should be [2*3, 2*4] = [6, 8]
analytical_grad = np.array([2 * point_to_test[0], 2 * point_to_test[1]])

# Numerical gradient
numerical_grad = estimate_gradient(my_function, point_to_test)

print(f"Point: {point_to_test}")
print(f"Analytical gradient: {analytical_grad}")
print(f"Estimated gradient:  {numerical_grad}")

As you can see from running the code, the estimated gradient is extremely close to the true analytical gradient. This numerical method is a general-purpose way to find the "downhill" direction for any function, which is exactly what we need for our optimizer.

import numpy as np

# (Helper functions from the lesson)
def partial_difference_quotient(f, v, i, h):
    w = v.copy()
    w[i] = w[i] + h
    return (f(w) - f(v)) / h

def estimate_gradient(f, v, h=0.00001):
    return np.array([partial_difference_quotient(f, v, i, h) for i in range(len(v))])

# Point of interest
point = np.array([1.0, 2.0])

# 1. Analytical Method
def gradient_analytical(v):
    x, y = v[0], v[1]
    df_dx = 2*x + 3*y
    df_dy = 3*x + 3*(y**2)
    return np.array([df_dx, df_dy])

analytical_result = gradient_analytical(point)
print(f"Analytical Gradient at {point}: {analytical_result}")

# 2. Numerical Method
def f_xy(v):
    x, y = v[0], v[1]
    return x**2 + 3*x*y + y**3

numerical_result = estimate_gradient(f_xy, point)
print(f"Numerical Gradient at {point}: {numerical_result}")

# 3. Comparison
# Analytical: df/dx = 2(1) + 3(2) = 8. df/dy = 3(1) + 3(2^2) = 3 + 12 = 15.
# The analytical result is [8, 15].
# The numerical result should be very close, e.g., [8.00003, 15.00006]. They are indeed nearly identical.

Conclusion

You've now successfully implemented the core mechanism for finding the direction of steepest descent. This is a crucial building block for the optimization algorithms we'll construct next.

Key Takeaways:

• The gradient is a vector of partial derivatives that points in the direction of a function's fastest increase. For optimization, we use its negative, $-\nabla f$, to move "downhill".
• We can calculate gradients analytically by finding the derivative formulas. This is fast and exact.
• We can also numerically estimate gradients by applying the limit definition of a derivative. This estimate_gradient function is a powerful, general tool that works for any function you can write in code.

Preview of the Next Lesson:
We have the "what" (minimize a loss function) and the "how to find the direction" (calculate the gradient). In the next lesson, you will put these pieces together to implement the full batch gradient descent algorithm. You will use your new estimate_gradient function to iteratively take steps downhill and find the minimum of a simple function.