Skip to article frontmatterSkip to article content

4.1. The Gradient Vector

The big theme in our course so far has been the three-step modeling recipe, of choosing a model, choosing a loss function, and then minimizing empirical risk (i.e. average loss) to find optimal model parameters.

In Chapter 1.4, we used calculus to find the slope w1∗w_1^* and intercept w0∗w_0^* that minimized mean squared error,

Rsq(w0,w1)=1n∑i=1n(yi−(w0+w1xi))2R_\text{sq}(w_0, w_1) = \frac{1}{n} \sum_{i=1}^n (y_i - (w_0 + w_1 x_i))^2

by computing ∂Rsq∂w0\frac{\partial R_\text{sq}}{\partial w_0} (the partial derivative with respect to w0w_0) ∂Rsq∂w1\frac{\partial R_\text{sq}}{\partial w_1}, setting both to zero, and solving the resulting system of equations.

Then, in Chapters 2 and 3, we focused on the multiple linear regression model and the squared loss function, which saw us minimize

Rsq(w⃗)=1n∥y⃗−Xw⃗∥2=1n∑i=1n(yi−w⃗⋅Aug(x⃗i))2R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2 = \frac{1}{n} \sum_{i=1}^n (y_i - \vec w \cdot \text{Aug}(\vec x_i))^2

where XX is the n×(d+1)n \times (d + 1) design matrix, y⃗∈Rn\vec y \in \mathbb{R}^n is the observation vector, and w⃗∈Rd+1\vec w \in \mathbb{R}^{d+1} is the parameter vector we’re trying to pick. In Chapter 2.10, we minimized Rsq(w⃗)R_\text{sq}(\vec w) by arguing that the optimal w⃗∗\vec w^* had to create an error vector, e⃗=y⃗−Xw⃗∗\vec e = \vec y - X \vec w^*, that was orthogonal to colsp(X)\text{colsp}(X), which led us to the normal equations.

It turns out that there’s a way to use our calculus-based approach from Chapter 1.4 to minimize the more general version of RsqR_\text{sq} for any dd, that doesn’t involve computing dd partial derivatives. To see how this works, we need to define a new object, the gradient vector, which we’ll do here in Chapter 4.1. After we’re familiar with how the gradient vector works, we’ll use it to build a new approach to function minimization, one that works even when there isn’t a closed-form solution for the optimal parameters: that technique is called gradient descent, which we’ll see in Chapter 4.2.

Domain and Codomain

As we saw in Chapter 2.9 when we first introduced the concept of the inverse of a matrix, the notation

f:Rd→Rnf: \mathbb{R}^d \to \mathbb{R}^n

means that ff is a function whose inputs are vectors with dd components and whose outputs are vectors with nn components. Rd\mathbb{R}^d is the domain of the function, and Rn\mathbb{R}^n is the codomain. I’ve used dd and nn to match the notation we’ve used for matrices and linear transformations. In general, if AA is an n×dn \times d matrix, then any vector x⃗\vec x multiplied by AA (on the right) must be in Rd\mathbb{R}^d and the result Ax⃗A \vec x will be in Rn\mathbb{R}^n.

Given this framing, consider the following four types of functions.

TypeDomain and CodomainExamples
Scalar-to-scalarf:R→Rf: \mathbb{R} \to \mathbb{R}f(x)=x2+sin⁡(x)f(x) = x^2 + \sin(x)
Rsq(w)=1n∑i=1n(yi−w)2R_\text{sq}(w) = \frac{1}{n} \sum_{i=1}^n (y_i - w)^2
Vector-to-scalarf:Rd→Rf: \mathbb{R}^d \to \mathbb{R}f(x⃗)=x⃗Tx⃗f(\vec x) = \vec x^T \vec x
f(x⃗)=(x1−1)2+(x1−x2)2+3f(\vec x) = (x_1 - 1)^2 + (x_1 - x_2)^2 + 3
Rsq(w⃗)=1n∥y⃗−Xw⃗∥2R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2
Scalar-to-vectorf:R→Rnf: \mathbb{R} \to \mathbb{R}^nf(x)=[121]+x[30−4]⏟parametric form of a linef(x) = \underbrace{\begin{bmatrix} 1 \\ 2 \\ 1 \end{bmatrix} + x \begin{bmatrix} 3 \\ 0 \\ -4 \end{bmatrix}}_{\text{parametric form of a line}}

f(x)=[x2+ex−3]f(x) = \begin{bmatrix} x^2 + e^{x} \\ -3 \end{bmatrix}
Vector-to-vectorf:Rd→Rnf: \mathbb{R}^d \to \mathbb{R}^nf(x⃗)=[342001]x⃗⏟linear transformation\underbrace{f(\vec x) = \begin{bmatrix} 3 & 4 \\ 2 & 0 \\ 0 & 1 \end{bmatrix} \vec x}_{\text{linear transformation}}

f(x⃗)=[x12+x1x2+cos⁡(x14)3x1x2+45x1x2/x3]f(\vec x) = \begin{bmatrix} x_1^2 + x_1x_2 + \cos(x_1^4) \\ 3x_1x_2 + 4 \\ 5x_1 \\ x_2 / x_3 \end{bmatrix}

The first two types of functions are “scalar-valued”, while the latter two are “vector-valued”. These are not the only types of functions that exist; for instance, the function f(A)=rank(A)f(A) = \text{rank}(A) is a matrix-to-scalar function.

The type of function we’re most concerned with at the moment are vector-to-scalar functions, i.e. functions that take in a vector (or equivalently, multiple scalar inputs) and output a single scalar.

Rsq(w⃗)=1n∥y⃗−Xw⃗∥2R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2

is one such function, and it’s the focus of this section.

Rates of Change

Let’s think from the perspectives of rates of change, since ultimately what we’re building towards is a technique for minimizing functions. We’re most familiar with the concept of rates of change for scalar-to-scalar functions.

If

f(x)=x2sin⁡(x)f(x) = x^2 \sin(x)

then its derivative,

dfdx=2xsin⁡(x)+x2cos⁡(x)\frac{\text{d}f}{\text{d}x} = 2x \sin(x) + x^2 \cos(x)

itself is a scalar-to-scalar function, which describes how quickly ff is changing at any point xx in the domain of ff. At x=3x = 3, for instance, the instantaneous rate of change is

dfdx(3)=2⋅3sin⁡(3)+32cos⁡(3)≈−8.06\frac{\text{d}f}{\text{d}x}(3) = 2\cdot 3 \sin(3) + 3^2 \cos(3) \approx -8.06

meaning that at x=3x = 3, ff is decreasing at a rate of (approximately) 8.06 per unit change in xx. Perhaps a more intuitive way of thinking about the instantaneous rate of change is to think of it as the slope of the tangent line to ff at x=3x = 3.

Loading...

The steeper the slope, the faster ff is changing at that point; the sign of the slope tells us whether ff is increasing or decreasing at that point.

In Chapter 1.4, we saw how to compute derivatives of functions that take in multiple scalar inputs, like

f(x,y,z)=x2+2xy+3xz+4(y−z)2f(x, y, z) = x^2 + 2xy + 3xz + 4(y - z)^2

In the language of Chapter 4.1, we’d call such a function a vector-to-scalar function, and might use the notation

f(x⃗)=x12+2x1x2+3x1x3+4(x2−x3)2f(\vec x) = x_1^2 + 2x_1x_2 + 3x_1x_3 + 4(x_2 - x_3)^2

This function has three partial derivatives, each of which describes the instantaneous rate of change of ff with respect to one of its inputs, while holding the other two inputs constant. There’s a good animation of what it means to hold an input constant in Chapter 1.4 that is worth revisiting.

Here,

∂f∂x1=2x1+2x2+3x3,∂f∂x2=2x1+8x2−8x3,∂f∂x3=3x1−8x2+8x3\frac{\partial f}{\partial x_1} = 2x_1 + 2x_2 + 3x_3, \quad \frac{\partial f}{\partial x_2} = 2x_1 + 8x_2 - 8x_3, \quad \frac{\partial f}{\partial x_3} = 3x_1 - 8x_2 + 8x_3

The big idea of this section, the gradient vector, packages all of these partial derivatives into a single vector. This will allow us to think about the direction in which ff is changing, rather than just looking at its rates of change in each dimension independently.

The Gradient Vector

As usual, we’ll start with an example. Suppose x∈R2x \in \mathbb{R}^2, and let

f(x⃗)=x1e−x12−x22f(\vec x) = x_1 e^{-x_1^2 - x_2^2}
Loading...

To find ∇f(x⃗)\nabla f(\vec x), we need to compute the partial derivatives of ff with respect to each component of x⃗\vec x. The “input variables” to ff are x1x_1 and x2x_2, so we need to compute ∂f∂x1\frac{\partial f}{\partial x_1} and ∂f∂x2\frac{\partial f}{\partial x_2}, but if you’d like, replace x1x_1 and x2x_2 with xx and yy if it makes the algebra a little cleaner, and then replace xx and yy with x1x_1 and x2x_2 at the end.

f(x⃗)=x1e−x12−x22f(\vec x) = x_1 e^{-x_1^2 - x_2^2}
∂f∂x1=∂∂x1(x1e−x12−x22)=1⋅e−x12−x22+x1⋅e−x12−x22⋅(−2x1)⏟product rule=(1−2x12)e−x12−x22\frac{\partial f}{\partial x_1} = \frac{\partial}{\partial x_1} \left( x_1 e^{-x_1^2 - x_2^2} \right) = \underbrace{1 \cdot e^{-x_1^2 - x_2^2} + x_1 \cdot e^{-x_1^2 - x_2^2} \cdot (-2x_1)}_{\text{product rule}} = (1 - 2x_1^2) e^{-x_1^2 - x_2^2}
∂f∂x2=∂∂x2(x1e−x12−x22)=x1⋅e−x12−x22⋅(−2x2)⏟chain rule=−2x1x2e−x12−x22\frac{\partial f}{\partial x_2} = \frac{\partial}{\partial x_2} \left( x_1 e^{-x_1^2 - x_2^2} \right) = \underbrace{x_1 \cdot e^{-x_1^2 - x_2^2} \cdot (-2x_2)}_{\text{chain rule}} = -2x_1 x_2 e^{-x_1^2 - x_2^2}

Putting these together, we have

∇f(x⃗)=[(1−2x12)e−x12−x22−2x1x2e−x12−x22]\nabla f(\vec x) = \begin{bmatrix} (1 - 2x_1^2) e^{-x_1^2 - x_2^2} \\ -2x_1 x_2 e^{-x_1^2 - x_2^2} \end{bmatrix}

Remember, ∇f(x⃗)\nabla f(\vec x) itself is a function. If we plug in a value of x⃗\vec x, we get a new vector back.

∇f([−10])=[(1−2(−1)2)e−(−1)2−02−2(−1)(0)e−(−1)2−02]=[−1/e0]\nabla f\left(\begin{bmatrix} -1 \\ 0\end{bmatrix}\right) = \begin{bmatrix} (1 - 2(-1)^2) e^{-(-1)^2 - 0^2} \\ -2(-1)(0) e^{-(-1)^2 - 0^2} \end{bmatrix} = \begin{bmatrix} -1/e \\ 0 \end{bmatrix}

What does ∇f([−10])=[−1/e0]\nabla f\left(\begin{bmatrix} -1 \\ 0\end{bmatrix}\right) = \begin{bmatrix} -1/e \\ 0 \end{bmatrix} really tell us? In order to visualize it, let me introduce another way of visualizing ff, called a contour plot.

Loading...

I think of the contour plot as a birds-eye view 🦅 of ff when you look at it from above. when you look at the surface from above. Notice the correspondence between the colors in both graphs.

The circle-like traces in the contour plot are called level curves; they represent slices through the surface at a constant height. On the right, the circle labeled 0.1 represents the set of points where f(x1,x2)=0.1f(x_1, x_2) = 0.1.

Visualizing the fact that ∇f([−10])=[−1/e0]\nabla f\left(\begin{bmatrix} -1 \\ 0\end{bmatrix}\right) = \begin{bmatrix} -1/e \\ 0 \end{bmatrix} is easier to do in the contour plot, since the contour plot is 2-dimensional, like the gradient vector is. Remember that red values are high and blue values are low.

Loading...

At the point [−10]\begin{bmatrix} -1 \\ 0 \end{bmatrix}, which is at the tail of the vector drawn in gold, ff is near the global minimum, meaning there are lots of directions in which we can move to increase ff. But, the gradient vector at this point is [−1/e0]\begin{bmatrix} -1/e \\ 0 \end{bmatrix}, which points in the direction of steepest ascent starting at [−10]\begin{bmatrix} -1 \\ 0 \end{bmatrix}. The gradient describes the “quickest way up”.

As another example, consider the fact that ∇f([1.25−0.5])≈[−0.3470.204]\nabla f\left(\begin{bmatrix} 1.25 \\ -0.5 \end{bmatrix}\right) \approx \begin{bmatrix} -0.347 \\ 0.204 \end{bmatrix}.

Loading...

Again, the gradient at [1.25−0.5]\begin{bmatrix} 1.25 \\ -0.5 \end{bmatrix} gives us the direction in which ff is increasing the quickest at that very point. If we move even a little bit in any direction (in the direction of the gradient or some other direction), the gradient will change.

As you might guess, to find the critical points of a function – that is, places where it is neither increasing nor decreasing – we need to find points where the gradient is zero. Hold that thought.

Examples

More typically, the functions we’ll need to take the gradient of will themselves be defined in terms of matrix and vector operations. In all of these examples, remember that we’re working with vector-to-scalar functions.

Example: Dot Product

Let a⃗∈Rn\vec a \in \mathbb{R}^n be some fixed vector (the equivalent of a constant in this context). Let’s find the gradient of

f(x⃗)=a⃗⋅x⃗f(\vec x) = \vec a \cdot \vec x

I find it helpful to think about f(x⃗)f(\vec x) in its expanded form,

f(x⃗)=a⃗⋅x⃗=a1x1+a2x2+⋯+anxnf(\vec x) = \vec a \cdot \vec x = a_1 x_1 + a_2 x_2 + \cdots + a_n x_n

Remember, ∇f(x⃗)\nabla f(\vec x) contains all of the partial derivatives of ff, which we now need to compute.

  • What is ∂f∂x1\frac{\partial f}{\partial x_1}? To me, that looks like a1a_1, since the first term is a1x1a_1 x_1 and none of the other terms involve x1x_1.

  • Similarly, ∂f∂x2=a2\frac{\partial f}{\partial x_2} = a_2.

  • In general, ∂f∂xi=ai\frac{\partial f}{\partial x_i} = a_i.

Putting these together, we get

∇f(x⃗)=[∂f∂x1∂f∂x2⋮∂f∂xn]=[a1a2⋮an]=a⃗\nabla f(\vec x) = \begin{bmatrix} \frac{\partial f}{\partial x_1} \\ \frac{\partial f}{\partial x_2} \\ \vdots \\ \frac{\partial f}{\partial x_n} \end{bmatrix} = \begin{bmatrix} a_1 \\ a_2 \\ \vdots \\ a_n \end{bmatrix} = \vec a

Example: Norm and Chain Rule

Here’s an extremely important example that shows up everywhere in machine learning. Find the gradients of:

  1. f(x⃗)=∥x⃗∥2f(\vec x) = \lVert \vec x \rVert^2

  2. f(x⃗)=∥x⃗∥f(\vec x) = \lVert \vec x \rVert

Solution

  1. As we did in the previous example, we can expand f(x⃗)=∥x⃗∥2f(\vec x) = \lVert \vec x \rVert^2 to get

    f(x⃗)=x⃗⋅x⃗=x12+x22+⋯+xn2f(\vec x) = \vec x \cdot \vec x = x_1^2 + x_2^2 + \cdots + x_n^2

    For each ii, ∂f∂xi=2xi\frac{\partial f}{\partial x_i} = 2x_i. So,

    ∇f(x⃗)=[2x12x2⋮2xn]=2x⃗\nabla f(\vec x) = \begin{bmatrix} 2x_1 \\ 2x_2 \\ \vdots \\ 2x_n \end{bmatrix} = \boxed{2\vec x}

    Think of this as the equivalent of the “power rule” for vectors.

  2. There are two ways to find the gradient of f(x⃗)=∥x⃗∥f(\vec x) = \lVert \vec x \rVert: directly, or by using the chain rule. It’s not immediately obvious how the chain rule should work here, so we’ll start with the direct method and reason about how the chain rule may arise.

    Direct method: Let’s start by expanding f(x⃗)=∥x⃗∥f(\vec x) = \lVert \vec x \rVert like we did above.

    f(x⃗)=x⃗⋅x⃗=(x⃗⋅x⃗)1/2=(x12+x22+⋯+xn2)1/2f(\vec x) = \sqrt{\vec x \cdot \vec x} = (\vec x \cdot \vec x)^{1/2} = (x_1^2 + x_2^2 + \cdots + x_n^2)^{1/2}

    For each ii, the (regular, scalar-to-scalar) chain rule tells us that

    ∂f∂xi=12(x12+x22+⋯+xn2)−1/2⋅2xi=xix12+x22+⋯+xn2=xi∥x⃗∥\frac{\partial f}{\partial x_i} = \frac{1}{2}(x_1^2 + x_2^2 + \cdots + x_n^2)^{-1/2} \cdot 2x_i = \frac{x_i}{\sqrt{x_1^2 + x_2^2 + \cdots + x_n^2}} = \frac{x_i}{\lVert \vec x \rVert}

    So,

    ∇f(x⃗)=[x1∥x⃗∥x2∥x⃗∥⋮xn∥x⃗∥]=x⃗∥x⃗∥\nabla f(\vec x) = \begin{bmatrix} \frac{x_1}{\lVert \vec x \rVert} \\ \frac{x_2}{\lVert \vec x \rVert} \\ \vdots \\ \frac{x_n}{\lVert \vec x \rVert} \end{bmatrix} = \boxed{\frac{\vec x}{\lVert \vec x \rVert}}

    Chain rule method: Let me start by writing f(x⃗)f(\vec x) in terms of a composition of two functions.

    f(x⃗)=∥x⃗∥=∥x⃗∥2=h(g(x⃗))f(\vec x) = \lVert \vec x \rVert = \sqrt{\lVert \vec x \rVert^2} = h(g(\vec x))

    where g(x⃗)=∥x⃗∥2g(\vec x) = \lVert \vec x \rVert^2 and h(x)=xh(x) = \sqrt{x}. Note that g:Rn→Rg: \mathbb{R}^n \to \mathbb{R} is the vector-to-scalar function we found the gradient of above, and h:R→Rh: \mathbb{R} \to \mathbb{R} is a scalar-to-scalar function.

    Then, generalizing the calculation we did with the first method, we have a “chain rule” for a function h(g(x⃗))h(g(\vec x)) (where hh is scalar-to-scalar and gg is vector-to-scalar):

    ∇f(x⃗)=(dhdx(g(x⃗)))⏟h′(g(x⃗))∇g(x⃗)\nabla f(\vec x) = \underbrace{\left(\frac{\text{d}h}{\text{d}x}(g (\vec x))\right)}_{h' (g(\vec x))} \nabla g(\vec x)

    Remember that h(x)=xh(x) = \sqrt{x}, so dhdx(x)=12x\frac{\text{d}h}{\text{d}x}(x) = \frac{1}{2\sqrt{x}} and dhdx(g(x⃗))=12g(x⃗)=12∥x⃗∥\frac{\text{d}h}{\text{d}x}(g(\vec x)) = \frac{1}{2\sqrt{g(\vec x)}} = \frac{1}{2\lVert \vec x \rVert}. This means

    f(x⃗)=(dhdx(g(x⃗)))∇g(x⃗)=(12∥x⃗∥)2x⃗=x⃗∥x⃗∥f(\vec x) = \left(\frac{\text{d}h}{\text{d}x}(g (\vec x))\right) \nabla g(\vec x) = \left( \frac{1}{2\lVert \vec x \rVert} \right) 2 \vec x = \frac{\vec x}{\lVert \vec x \rVert}

    which is what we found earlier.

Example: Norm to an Exponent

Find the gradient of f(x⃗)=∥x⃗∥pf(\vec x) = \lVert \vec x \rVert^p, where pp is some real number.

Solution

We can treat this as a composition of two functions, g(x⃗)=∥x⃗∥g(\vec x) = \lVert \vec x \rVert and h(x)=xph(x) = x^p, and use the chain rule introduced in the solution to the previous example.

dhdx(x)=pxp−1\frac{\text{d}h}{\text{d}x}(x) = p x^{p-1} and ∇g(x⃗)=x⃗∥x⃗∥\nabla g(\vec x) = \frac{\vec x}{\lVert \vec x \rVert}. Putting these together yields

∇f(x⃗)=(dhdx(g(x⃗)))∇g(x⃗)=p g(x⃗)p−1x⃗∥x⃗∥=p ∥x⃗∥p−1x⃗∥x⃗∥=p ∥x⃗∥p−2x⃗\begin{align*} \nabla f(\vec x) &= \left( \frac{\text{d}h}{\text{d}x}(g (\vec x)) \right) \nabla g(\vec x) \\ &= p\, g(\vec x)^{p-1} \frac{\vec x}{\lVert \vec x \rVert} \\ &= p\, \lVert \vec x \rVert^{p-1} \frac{\vec x}{\lVert \vec x \rVert} \\ &= p\, \lVert \vec x \rVert^{p-2} \vec x \end{align*}

Example: Log Sum Exp

If x⃗∈Rn\vec x \in \mathbb{R}^n, we can define the log sum exp function as

f(x⃗)=log⁡(∑i=1nexi)f(\vec x) = \log \left( \sum_{i=1}^n e^{x_i} \right)

What is ∇f(x⃗)\nabla f(\vec x)? (The answer is called the softmax function, and comes up all the time in machine learning, when we want our models to output predicted probabilities in a classification problem.)

Solution

Let’s look at the partial derivatives with respect to each xix_i.

∂f∂xi=∂∂xi(log⁡(∑j=1nexj))=(1∑j=1nexj)∂∂xi(∑j=1nexj)=exi∑j=1nexj\frac{\partial f}{\partial x_i} = \frac{\partial}{\partial x_i} \left( \log \left( \sum_{j=1}^n e^{x_j} \right) \right) = \left(\frac{1}{\sum_{j=1}^n e^{x_j}} \right) \frac{\partial}{\partial x_i} \left( \sum_{j=1}^n e^{x_j} \right) = \frac{e^{x_i}}{\sum_{j=1}^n e^{x_j}}

Then,

∇f(x⃗)=[ex1∑j=1nexjex2∑j=1nexj⋮exn∑j=1nexj]=1∑j=1nexj[ex1ex2⋮exn]\nabla f(\vec x) = \begin{bmatrix} \frac{e^{x_1}}{\sum_{j=1}^n e^{x_j}} \\ \frac{e^{x_2}}{\sum_{j=1}^n e^{x_j}} \\ \vdots \\ \frac{e^{x_n}}{\sum_{j=1}^n e^{x_j}} \end{bmatrix} = \frac{1}{\sum_{j=1}^n e^{x_j}} \begin{bmatrix} e^{x_1} \\ e^{x_2} \\ \vdots \\ e^{x_n} \end{bmatrix}

There isn’t really a way to simplify the expression using matrix-vector operations, so I’ll leave it as-is. As mentioned above, the gradient we’re looking at is called the softmax function. The softmax function maps Rn→Rn\mathbb{R}^n \to \mathbb{R}^n, meaning its a vector-to-vector function.

Let’s suppose we have the matrix [35−1]\begin{bmatrix} 3 \\ 5 \\ -1 \end{bmatrix}. What does passing it through the softmax function yield?

softmax([35−1])=1e3+e5+e−1[e3e5e−1]≈[0.1190.8790.0002]\text{softmax}\left(\begin{bmatrix} 3 \\ {\color{orange}5} \\ -1 \end{bmatrix}\right) = \frac{1}{e^3 + e^{5} + e^{-1}}\begin{bmatrix} e^3 \\ e^5 \\ e^{-1} \end{bmatrix} \approx \begin{bmatrix} 0.119 \\ {\color{orange}0.879} \\ 0.0002 \end{bmatrix}

The output vector has the same number of elements as the input vector, but each element is between 0 and 1, and the sum of elements is 1, meaning that we can interpret the outputted vector as a probability distribution. Larger values in the output correspond to larger values in the input, and almost all of the “mass” is concentrated at the maximum element of the input vector (position 2), hence the name “soft” max. (The “hard” max might be [010]\begin{bmatrix} 0 \\ 1 \\ 0 \end{bmatrix} in this case.)

Example: Quadratic Forms

Suppose x∈Rnx \in \mathbb{R}^n and AA is an n×nn \times n matrix. The function

f(x⃗)=x⃗TAx⃗f(\vec x) = \vec x^T A \vec x

is called a quadratic form, and its gradient is given by

∇f(x⃗)=(A+AT)x⃗\nabla f(\vec x) = (A + A^T) \vec x

We won’t directly cover the proof of this formula here; one place to find it is here. Instead, we’ll focus our energy on understanding how it works, since it’s extremely important.

  1. Let A=[abcd]A = \begin{bmatrix} a & b \\ c & d \end{bmatrix}. Expand out f(x⃗)=x⃗TAx⃗f(\vec x) = \vec x^T A \vec x and compute ∇f(x⃗)\nabla f(\vec x) directly by computing partial derivatives, and verify that the result you get matches the formula above.

  2. In quadratic forms, we typically assume that AA is symmetric, meaning that A=ATA = A^T. Why do you think this assumption is made (what does it help with)?

    • Hint: Let A=[3261]A = \begin{bmatrix} 3 & 2 \\ 6 & 1 \end{bmatrix} and B=[3441]B = \begin{bmatrix} 3 & 4 \\ 4 & 1 \end{bmatrix}. Compute ∇(x⃗TAx⃗)\nabla (\vec x^T A \vec x) and ∇(x⃗TBx⃗)\nabla (\vec x^T B \vec x).

  3. If AA is any symmetric n×nn \times n matrix, what is ∇f(x⃗)\nabla f(\vec x)?

  4. Suppose AA is symmetric and n×nn \times n, b⃗∈Rn\vec b \in \mathbb{R}^n, and c∈Rc \in \mathbb{R}. Find the gradient of

    f(x⃗)=x⃗TAx⃗+b⃗⋅x⃗+cf(\vec x) = \vec x^T A \vec x + \vec b \cdot \vec x + c
Solution

  1. If A=[abcd]A = \begin{bmatrix} a & b \\ c & d \end{bmatrix}, then

    f(x⃗)=x⃗TAx⃗=[x1x2][abcd][x1x2]=[x1x2][ax1+bx2cx1+dx2]=ax12+(b+c)x1x2+dx22\begin{align*} f(\vec x) &= \vec x^T A \vec x \\ &= \begin{bmatrix} x_1 & x_2 \end{bmatrix} \begin{bmatrix} a & b \\ c & d \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} \\ &= \begin{bmatrix} x_1 & x_2 \end{bmatrix} \begin{bmatrix} a x_1 + b x_2 \\ c x_1 + d x_2 \end{bmatrix} \\ &= a x_1^2 + (b + c) x_1 x_2 + d x_2^2 \end{align*}

    Then,

    ∂f∂x1=2ax1+(b+c)x2,∂f∂x2=(b+c)x1+2dx2\frac{\partial f}{\partial x_1} = 2a {\color{#3d81f6}x_1} + (b + c) {\color{#3d81f6}x_2}, \quad \frac{\partial f}{\partial x_2} = (b + c) {\color{#3d81f6}x_1} + 2d {\color{#3d81f6}x_2}
    ∇f(x⃗)=[2ax1+(b+c)x2(b+c)x1+2dx2]=[2ab+cb+c2d][x1x2]=(A+AT)x⃗\nabla f(\vec x) = \begin{bmatrix} 2a {\color{#3d81f6}x_1} + (b + c) {\color{#3d81f6}x_2} \\ (b + c) {\color{#3d81f6}x_1} + 2d {\color{#3d81f6}x_2} \end{bmatrix} = \begin{bmatrix} 2a & b + c \\ b + c & 2d \end{bmatrix} \begin{bmatrix} {\color{#3d81f6}x_1} \\ {\color{#3d81f6}x_2} \end{bmatrix} = (A + A^T) \vec x

    since AT=[acbd]A^T = \begin{bmatrix} a & c \\ b & d \end{bmatrix}, meaning A+AT=[2ab+cb+c2d]A + A^T = \begin{bmatrix} 2a & b + c \\ b + c & 2d \end{bmatrix}.

  2. For a particular quadratic form, there are infinitely many choices of matrices AA that represent it. To illustrate, let’s look at A=[3261]A = \begin{bmatrix} 3 & 2 \\ 6 & 1 \end{bmatrix} and B=[3441]B = \begin{bmatrix} 3 & 4 \\ 4 & 1 \end{bmatrix} as provided in the hint.

    x⃗TAx⃗=[x1x2][3261][x1x2]=3x12+(2+6)x1x2+x22\vec x^T A \vec x = \begin{bmatrix} x_1 & x_2 \end{bmatrix} \begin{bmatrix} 3 & 2 \\ 6 & 1 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} = 3 x_1^2 + (2 + 6) x_1 x_2 + x_2^2
    x⃗TBx⃗=[x1x2][3441][x1x2]=3x12+(4+4)x1x2+x22\vec x^T B \vec x = \begin{bmatrix} x_1 & x_2 \end{bmatrix} \begin{bmatrix} 3 & 4 \\ 4 & 1 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} = 3 x_1^2 + (4 + 4) x_1 x_2 + x_2^2

    Note that both x⃗TAx⃗\vec x^T A \vec x and x⃗TBx⃗\vec x^T B \vec x are equal to the expression 3x12+8x1x2+x223x_1^2 + 8x_1x_2 + x_2^2. In fact, any matrix of the form [3bc1]\begin{bmatrix} 3 & b \\ c & 1 \end{bmatrix} where b+c=8b + c = 8 would produce the same quadratic form.

    So, to avoid this issue of having infinitely many choices of the matrix AA, we pick the symmetric matrix AA, where A=ATA = A^T. As we’re about to see, this choice of AA simplifies the calculation of the gradient.

  3. If AA is any symmetric n×nn \times n matrix, then A=ATA = A^T, and A+AT=2AA + A^T = 2A. So,

    ∇(x⃗TAx⃗)=(A+AT)x⃗=2Ax⃗\nabla (\vec x^T A \vec x) = (A + A^T) \vec x = 2A \vec x

    This is also an important rule; don’t forget it.

  4. Think of f(x⃗)=x⃗TAx⃗+b⃗⋅x⃗+cf(\vec x) = \vec x^T A \vec x + \vec b \cdot \vec x + c as the matrix-vector equivalent of a quadratic function, ax2+bx+cax^2 + bx + c. The derivative of ax2+bx+cax^2 + bx + c is 2ax+b2ax + b. Check out what the gradient of f(x⃗)f(\vec x) ends up being!

f(x⃗)=x⃗TAx⃗+b⃗⋅x⃗+c∇f(x⃗)=(A+AT)x⃗+b⃗=2Ax⃗+b⃗(since A is symmetric)\begin{align*} f(\vec x) &= \vec x^T A \vec x + \vec b \cdot \vec x + c \\ \nabla f(\vec x) &= (A + A^T) \vec x + \vec b \\ &= 2A \vec x + \vec b \qquad \text{(since $A$ is symmetric)} \end{align*}

Summary of Important Gradient Rules

These are the core rules you need to know moving forward, not just because we’re about to use them in an important proof, but because they’ll come up repeatedly in your future machine learning work.

FunctionNameGradient
f(x⃗)=a⃗⋅x⃗f(\vec x) = \vec a \cdot \vec xdot product∇f(x⃗)=a⃗\nabla f(\vec x) = \vec a
f(x⃗)=∥x⃗∥2f(\vec x) = \lVert \vec x \rVert^2squared norm∇f(x⃗)=2x⃗\nabla f(\vec x) = 2\vec x
f(x⃗)=x⃗TAx⃗f(\vec x) = \vec x^T A \vec xquadratic form∇f(x⃗)=(A+AT)x⃗\nabla f(\vec x) = (A + A^T) \vec x
if AA is symmetric, ∇f(x⃗)=2Ax⃗\nabla f(\vec x) = 2A \vec x

Optimization

In the calculus of scalar-to-scalar functions, we have a well-understood procedure for finding the extrema of a function. The general strategy is to take the derivative, set it to zero, and solve for the inputs (called critical points) that satisfy that condition. To be thorough, we’d perform a second derivative test to check whether each critical point is a maximum, minimum, or neither.

In the land of vector-to-scalar functions, the equivalent is to solve for where the gradient is zero, which corresponds to finding where all partial derivatives are zero. Assessing whether we’ve arrived at a maximum or minimum is more difficult to do in the vector-to-scalar case, and we will save a discussion of this for Chapter 4.2.

As an example, consider

f(x⃗)=x⃗T[3441]x⃗+[12]⋅x⃗+3f(\vec x) = \vec x^T \begin{bmatrix} 3 & 4 \\ 4 & 1 \end{bmatrix} \vec x + \begin{bmatrix} 1 \\ 2 \end{bmatrix} \cdot \vec x + 3

As we computed earlier, the gradient of f(x⃗)=x⃗TAx⃗+b⃗⋅x⃗+cf(\vec x) = \vec x^T A \vec x + \vec b \cdot \vec x + c is ∇f(x⃗)=2Ax⃗+b⃗\nabla f(\vec x) = 2A \vec x + \vec b for symmetric AA. So,

∇f(x⃗)=2[3441]x⃗+[12]=[6x1+8x2+18x1+2x2+2]\nabla f(\vec x) = 2 \begin{bmatrix} 3 & 4 \\ 4 & 1 \end{bmatrix} \vec x + \begin{bmatrix} 1 \\ 2 \end{bmatrix} = \begin{bmatrix} 6x_1 + 8x_2 + 1 \\ 8x_1 + 2x_2 + 2 \end{bmatrix}

To find the critical points, we set the gradient to zero and solve the resulting system. We can also accomplish this by using the inverse of AA, if we happen to have it:

∇f(x⃗)=0  ⟹  2Ax⃗+b⃗=0  ⟹  x⃗∗=−12A−1b⃗\nabla f(\vec x) = 0 \implies 2 A \vec x + \vec b = 0 \implies \vec x^* = -\frac{1}{2}A^{-1} \vec b

Either way, we find that x⃗∗=[−7/261/13]\vec x^* = \begin{bmatrix} -7/26 \\ 1/13 \end{bmatrix} satisfies ∇f(x⃗∗)=0\nabla f(\vec x^*) = 0, which corresponds to a local minimum.

Loading...

Minimizing Mean Squared Error

Remember, the goal of this section is to minimize mean squared error,

Rsq(w⃗)=1n∥y⃗−Xw⃗∥2R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2

In the general case, XX is an n×(d+1)n \times (d + 1) matrix, y∈Rny \in \mathbb{R}^n, and w⃗∈Rd+1\vec w \in \mathbb{R}^{d+1}.

We’re now equipped with the tools to minimize Rsq(w⃗)R_\text{sq}(\vec w) by taking its gradient and setting it to zero. Hopefully, we end up with the same conditions on w⃗∗\vec w^* that we derived in Chapter 2.10.

In the most recent example we saw, the optimal vector x⃗∗\vec x^* corresponded to a local minimum. We know that we won’t run into such an issue here since Rsq(w⃗)R_\text{sq}(\vec w) cannot output a negative number (it is the average of squared losses), so its minimum possible output is 0, meaning that there will be some global minimizer w⃗∗\vec w^*.

Let’s start by rewriting the squared norm as a dot product and eventually matrix multiplication.

Rsq(w⃗)=1n∥y⃗−Xw⃗∥2=1n(y⃗−Xw⃗)⋅(y⃗−Xw⃗)=1n(y⃗−Xx⃗)T(y⃗−Xw⃗)⏟since u⃗⋅v⃗=u⃗Tv⃗=1n(y⃗T−(Xw⃗)T)(y⃗−Xw⃗)=1n(y⃗Ty⃗−y⃗TXw⃗−(Xw⃗)Ty⃗+(Xw⃗)TXw⃗)\begin{align*}R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2 &= \frac{1}{n} (\vec y - X \vec w) \cdot (\vec y - X \vec w) \\ &= \underbrace{\frac{1}{n} (\vec y - X \vec x)^T (\vec y - X \vec w)}_{\text{since } \vec u \cdot \vec v = \vec u^T \vec v} \\ &= \frac{1}{n} \left( \vec y^T - (X \vec w)^T \right) (\vec y - X \vec w) \\ &= \frac{1}{n} \left( \vec y^T \vec y - {\color{orange}\vec y^T X \vec w} - {\color{orange}(X \vec w)^T \vec y} + (X \vec w)^T X \vec w \right)\end{align*}

Let’s focus on the two terms in orange. They are both equal: they are both the dot product of y⃗\vec y and Xw⃗X \vec w. Ideally, I want to express each term as a dot product of w⃗\vec w with something, since I’m taking the gradient with respect to w⃗\vec w. Remember, the dot product is a scalar, and the transpose of a scalar is just that same scalar. So,

y⃗TXw⃗=(y⃗TXw⃗)T=w⃗TXTy⃗=w⃗T(XTy⃗)\vec y^T X \vec w = (\vec y^T X \vec w)^T = \vec w^T X^T \vec y = \vec w^T (X^T \vec y)

so, performing this substitution in for both orange terms gives us

Rsq(w⃗)=1n(y⃗Ty⃗−w⃗T(XTy⃗)−w⃗TXTy⃗+w⃗T(XTX)w⃗)=1n(y⃗Ty⃗−2w⃗T(XTy⃗)+w⃗T(XTX)w⃗)\begin{align*}R_\text{sq}(\vec w) &= \frac{1}{n} \left( \vec y^T \vec y - {\color{orange}\vec w^T (X^T \vec y)} - {\color{orange}\vec w^T X^T \vec y} + \vec w^T (X^T X) \vec w \right) \\ &= \frac{1}{n} \left( \vec y^T \vec y - 2 \vec w^T (X^T \vec y) + \vec w^T (X^T X) \vec w \right)\end{align*}

Now, we’re ready to take the gradient, which we’ll do term by term.

  • ∇(y⃗Ty⃗)=0⃗\nabla \left( \vec y^T \vec y \right) = \vec 0, since y⃗Ty⃗\vec y^T \vec y is a constant with respect to w⃗\vec w

  • ∇(2w⃗T(XTy⃗))=2XTy⃗\nabla \left( 2 \vec w^T (X^T \vec y) \right) = 2 X^T \vec y using the dot product rule, since this is the dot product between 2XTy⃗2X^T \vec y (a vector) and w⃗\vec w (a vector)

  • ∇(w⃗T(XTX)w⃗)=2XTXw⃗\nabla \left( \vec w^T (X^T X) \vec w \right) = 2X^T X \vec w, using the quadratic form rule, since XTXX^T X is a symmetric matrix

Plugging these terms in gives us

Rsq(w⃗)=1n(y⃗Ty⃗−2w⃗T(XTy⃗)+w⃗T(XTX)w⃗)∇Rsq(w⃗)=1n(∇(y⃗Ty⃗)−∇(2w⃗T(XTy⃗))+∇(w⃗T(XTX)w⃗))=1n(0−2XTy⃗+2XTXw⃗)=2n(XTXw⃗−XTy⃗)\begin{align*}R_\text{sq}(\vec w) &= \frac{1}{n} \left( \vec y^T \vec y - 2 \vec w^T (X^T \vec y) + \vec w^T (X^T X) \vec w \right) \\ \nabla R_\text{sq}(\vec w) &= \frac{1}{n} \left( \nabla \left(\vec y^T \vec y \right) - \nabla \left( 2 \vec w^T (X^T \vec y) \right) + \nabla \left( \vec w^T (X^T X) \vec w \right) \right) \\ &= \frac{1}{n} \left( 0 - 2 X^T \vec y + 2X^T X \vec w \right) \\ &= \boxed{\frac{2}{n} (X^T X \vec w - X^T \vec y)} \end{align*}

Finally, to find the minimizer w⃗∗\vec w^*, we set the gradient to zero and solve.

2n(XTXw⃗∗−XTy⃗)=0  ⟹  XTXw⃗∗=XTy⃗\frac{2}{n} (X^T X \vec w^* - X^T \vec y) = 0 \implies X^TX \vec w^* = X^T \vec y

Stop me if this feels familiar... these are the normal equations once again! It shouldn’t be a surprise that we ended up with the same conditions on w⃗∗\vec w^* that we derived in Chapter 2.10, since we were solving the same problem.

We’ve now shown that the minimizer of

Rsq(w⃗)=1n∥y⃗−Xw⃗∥2R_\text{sq}(\vec w) = \frac{1}{n} \lVert \vec y - X \vec w \rVert^2

is given by solving XTXw⃗∗=XTy⃗X^TX \vec w^* = X^T \vec y. These equations have a unique solution if XTXX^TX is invertible, and infinitely many solutions otherwise. If w⃗∗\vec w^* satisfies the normal equations, then Xw⃗∗X \vec w^* is the vector in colsp(X)\text{colsp}(X) that is closest to y⃗\vec y. All of that interpretation from Chapter 2.10 and Chapter 3 carry over; we’ve just introduced a new way of finding the solution.

Heads up: In Homework 9, you’ll follow similar steps to minimize a new objective function, that resembles Rsq(w⃗)R_\text{sq}(\vec w) but involves another term. There, you’ll minimize

Rridge(w⃗)=∥y⃗−Xw⃗∥2+λ∥w⃗∥2R_\text{ridge}(\vec w) = \lVert \vec y - X \vec w \rVert^2 + \lambda \lVert \vec w \rVert^2

where λ>0\lambda > 0 is a constant, called the regularization hyperparameter. (Notice the missing 1n\frac{1}{n}.) A good way to practice what you’ve learned (and to get a head start on the homework) is to compute the gradient of Rridge(w⃗)R_\text{ridge}(\vec w) and set it to zero. We’ll walk through what the significance of Rridge(w⃗)R_\text{ridge}(\vec w) is in the homework.

EECS 245 Course Notes
Chapter 4: Gradients
EECS 245 Course Notes
Chapter 5: Eigenvalues and Eigenvectors