---
author:
- "[Eugeniy E. Mikhailov](http://physics.wm.edu/~evmik/) and [Greg
  Bentsen](https://physics.wm.edu/~gbentsen/)"
title: "Numerical Optimization III: Gradient Descent"
---

## Logistics and Agenda

Homework 5 due Monday Oct 26 at 11:59pm

\

Midterm 2

-   Take-home posted Monday Oct 26
-   In-class exam Monday Nov 2

\

Agenda

-   This week:  Optimization in multiple dimensions
-   Next week:  Advanced optimization, combinatoric optimization,
    machine learning

## Today

Optimization in multiple dimensions

-   Gradient descent optimization
-   Multi-start optimization + AI demo
-   Advanced techniques
    -   Nelder-Mead (derivative-free)
    -   Simulated annealing

\

**Goal:** Efficiently navigate the global energy landscape.

Hunt for global minimum (but often get stuck in local minima instead).

------------------------------------------------------------------------

## Where we left off

Given adjustable parameters $\vec{\theta}=(\theta_1,\ldots,\theta_d)$,
minimize the **objective function** $E(\vec{\theta})$

::: columns
::: {.column width="35%"}
Two examples from last time:

-   **Linear regression:**
    -   parameters $(m,b)$
    -   objective $\chi^2(m,b)$
-   **Four-bead protein:**
    -   parameters $(\theta_1,\theta_2)$
    -   objective $E(\theta_1,\theta_2)$
:::

::: {.column width="55%"}
::: {style="text-align: left;"}
![](./pics/linear_regression_gui.png){width="75%"}
![](./pics/protein_gui.png){width="75%"}
:::
:::
:::

**Question:** How do we find the bottom of a landscape without scanning
every point?

# Gradient-Based Optimization

## The gradient: which way is uphill?

For $E(\theta_1,\theta_2)$, define the **gradient:**

$$
\boxed{\vec\nabla E=\begin{pmatrix}
\partial E/\partial\theta_1\\[3pt]
\partial E/\partial\theta_2
\end{pmatrix}}
$$

-   $\vec\nabla E$ points **uphill** (steepest local increase)
-   $-\vec\nabla E$ points **downhill** (steepest local decrease)
-   At a local minimum $\vec\nabla E=\vec 0$ (but so do maxima and
    saddle points!)

::: {style="text-align: center;"}
![Source: generated with ChatGPT](./pics/grad_descent.png){width="55%"}
:::

## Gradient descent is an iterated map

Start with $\vec{\theta}_0$, then repeatedly step downhill:

$$
\boxed{\vec\theta_{n+1}=\vec\theta_n-\eta \ \vec\nabla E(\vec\theta_n)}
$$

-   $\eta>0$: **learning rate** (step size)
-   Evaluate gradient $\rightarrow$ step downhill $\rightarrow$ repeat
-   Stop when energy / parameters stop changing appreciably

This is an *iterated map*, just like our earlier dynamical-systems
examples.

## Gradient descent in Python

Start with $\vec{\theta}_0$, then repeatedly step downhill:

$$
\boxed{\vec\theta_{n+1}=\vec\theta_n-\eta \ \vec\nabla E(\vec\theta_n)}
$$

``` {.python .numberLines}
import numpy as np

# grad_energy(theta) takes a length-2 array and returns the gradient (also length-2)

theta = np.array([1.0, -2.0])   # initial guess (radians)
eta = 0.05                      # learning rate
tol = 1e-6                      # convergence tolerance
max_steps = 1e4                 # maximum number of steps

grad = grad_energy(theta)
step = 0

while (np.linalg.norm(grad) > tol) & (step < max_steps):
    theta -= eta * grad
    grad = grad_energy(theta)
    step = step + 1
```

**In practice:** $\vec\nabla E$ may be derived analytically or
approximated using finite differences.

## Demo: gradient descent for linear regression

**First:** linear-regression loss

$$
\chi^2(m,b)=\sum_i[y_i-(mx_i+b)]^2
$$

**Demo:** `linear_fit_gradient_descent.py`

**Observe:** Which minimum do we reach? How do the lines / protein
change during descent?

## Adjusting the learning rate

::: {style="text-align: center;"}
![](./pics/learning_rate.png){width="75%"}
:::

-   Small $\eta$: many tiny steps, slow convergence
-   Moderate $\eta$: efficient progress
-   Large $\eta$: overshooting, oscillation, possible instability

**Takeaway:** a larger step is not always faster!

## Too slow, too fast, or just right?

Return to 1D and use a simple objective $f(x)=x^2$

Newton-Raphson:

$$
\frac{df}{dx}=2x
\qquad\Rightarrow\qquad
\boxed{x_{n+1}=(1-2\eta)x_n}
$$

Fixed point analysis: $|F'| = |1-2\eta|$

-   $0<\eta<\tfrac12$: approaches minimum **without changing sign**
-   $\tfrac12<\eta<1$: approaches minimum **while oscillating**
-   $\eta>1$: **diverges** (unless started at $x=0$)

At $\eta=1$, the iterates alternate with constant magnitude.

## Demo: gradient descent for protein model

**What changes in a more complicated landscape?**

**Demo:** `folding_gradient_descent.py`

**Different initial guesses:** may converge to *different* minima

**Question:** How can we be confident we found the *global* minimum?

# Other Ways to Navigate the Landscape

What if the gradient $\vec \nabla E$ isn't available?

## Nelder--Mead: gradient-free optimization **\[optional\]**

Some objective functions come from expensive simulations, making
gradients inconvenient.

**Idea:** move a triangle (a **simplex** in two dimensions) across the
landscape.

::: columns
::: {.column width="45%"}
-   Evaluate $f$ at its three vertices
-   Reflect / expand / contract to seek lower values
-   Repeat, shrinking around a candidate minimum
:::

::: {.column width="35%"}
::: {style="text-align: left;"}
![Source:
https://www.slideserve.com/nay/nelder-mead](./pics/nelder_mead.png){width="85%"}
:::
:::
:::

**Like gradient descent, Nelder--Mead is generally a local optimizer.**

# Pitfalls and Takeaways

## Common optimization traps

-   **Local $\ne$ global minima:** following downhill is not enough: a
    low valley may not be the lowest valley
-   **Initial guesses:** different $\vec{\theta_0}$ lead to different
    minima
-   **Step size (learning rate):** a poor choice can cause either slow
    or unstable convergence
-   **Dimensionality:** $m$ grid samples per coordinate require $m^d$
    evaluations!

**Practical strategy:** run the optimizer from *many* starting
configurations.

# Avoiding Local Minima

## A downhill algorithm can get trapped

-   Gradient descent is a **local** method
-   The final answer can depend on the **initial guess**
-   A vanishing gradient does **not** guarantee the lowest energy
    everywhere

::: {style="text-align: center;"}
![Source:
https://commons.wikimedia.org/wiki/File%3AGradient_Descent_in_2D.webm](./pics/local_minima.png){width="30%"}
:::

Movie:  [Wikimedia Commons: Gradient Descent in
2D](https://commons.wikimedia.org/wiki/File%3AGradient_Descent_in_2D.webm)

**Question:** How could we make the search more reliable?

## Multistart optimization: try many initial guesses

1.  Draw many starting configurations $\boldsymbol\theta_0^{(k)}$
2.  Run gradient descent from each
3.  Compare final energies and **keep the lowest**

$$
\boxed{k_* = \mathop{\mathrm{argmin}}_k E\!\left(\boldsymbol\theta_{\mathrm{final}}^{(k)}\right)}
$$

::: {style="text-align: center;"}
![](./pics/multistart_optimization.png){width="70%"}
:::

**Still no mathematical guarantee of a global minimum.**

## AI-coding demo: automate multistart searches

**Prompt a local AI model:**

> Starting from our four-bead protein `energy(theta1, theta2, ...)`
> function located in `folding_model.py`, write Python code that samples
> 100 random initial angle configurations
> $theta1, theta2 \in [-\pi,\pi)^2$, runs gradient descent from each,
> and plots the trajectories on an energy contour map. Highlight the
> lowest final energy. Do not change any existing code, and write this
> Python code in a new file named `multistart_protein.py`. Do **not**
> examine any files other than `folding_model.py` and do **not** leave
> the directory `Lecture 14`.

**Verify its output:**

-   Does it actually launch **100 independent searches**?
-   Are the angles treated as periodic?
-   Are the final energies recalculated and compared?
-   Does changing the random seed change which minima are discovered?

## Simulated annealing: sometimes go uphill!

**Problem:** strictly downhill steps cannot cross an energy barrier.

**New strategy:** propose a random change and sometimes accept an uphill
move:

$$
\boxed{P(\mathrm{accept})=\min\left(1,e^{-\Delta E/T}\right)}
$$

with

$$
\Delta E=E_{\mathrm{new}}-E_{\mathrm{old}}
$$

\

-   Downhill moves: always accepted
-   Uphill moves: accept with probability $e^{-\Delta E/T}$

Uphill moves are especially probable at **high temperature** $T$, but
occur less frequently at low temperature.

\

**Gradually lower $T$:** explore first, settle later.

## Demo: hot and cold searches

::: {style="text-align: center;"}
![](./pics/annealing_protein.png){width="50%"}
:::

**Same energy landscape, different temperatures:**

-   **High $T$:** many uphill moves; explore broadly
-   **Low $T$:** uphill moves rare; become trapped more easily
-   **Cooling schedule:** seek a compromise between exploration and
    convergence

**Connection to physics:** thermal fluctuations and Boltzmann weights.

## Which method should we use?

  -----------------------------------------------------------------------
  Method            Needs gradients?  Main strength     Main limitation
  ----------------- ----------------- ----------------- -----------------
  Gradient descent  Yes (or finite    Simple, efficient Initial-guess and
                    differences)      local descent     step-size
                                                        sensitivity

  Multistart        Yes (repeatedly)  Explores multiple More
  descent                             basins            computational
                                                        work

  Nelder--Mead      No                Uses function     Can stagnate;
                                      values only       local search

  Simulated         No                Can cross energy  Requires many
  annealing                           barriers          evaluations and
                                                        temperature
                                                        tuning
  -----------------------------------------------------------------------

## Takeaways

-   **Gradient descent** moves against the local gradient
-   **Learning rate** controls convergence and stability
-   A **local minimum** is not necessarily the global minimum
-   **Multistart** and **annealing** improve exploration in different
    ways
-   The best method depends on the landscape and what we can compute

**Next:** combinatorial optimization --- when our choices are discrete
rather than continuous.