pset 4 solutions
This commit is contained in:
Steven G. Johnson committed
1 parent
da64ac497b
commit
b0649a5e44
2 files changed
+710
-1
No files matched your search
@@ -290,7 +290,7 @@ Described the KKT conditions for a (local) optimum/extremum (Boyd, section 5.5.3
|
||||
|
||||
* Quick review of augmented Lagrangians and ADMM from last lecture. Indicator function example from section III.4.
|
||||
* CCSA interior-point algorithm
|
||||
* pset 4 solutions: coming soon
|
||||
* [pset 4 solutions](psets/pset4sol.ipynb)
|
||||
* pset 5: coming soon, due 4/21
|
||||
|
||||
Went over very different example of a nonlinear optimization scheme, solving a fairly general inequality-constrained nonlinear-programming problem: the CCSA algorithm(s), as described by Svanberg (2002). This is a surprisingly simple algorithm (the [NLopt](http://ab-initio.mit.edu/nlopt) implementation is only 300 lines of C code), but is robust and provably convergent, and illustrates a number of important ideas in optimization: optimizing an approximation to update the parameters **x**, guarding the approximation with trust regions and penalty terms, and optimizing via the dual function (Lagrange multipliers). Like many optimization algorithms, the general ideas are very straightforward, but getting the details right can be delicate!
|
||||
|
||||
@@ -0,0 +1,709 @@
|
||||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "73fdabd6",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"# 18.065 Problem Set 4 Solutions"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "0fb69bf7",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 1 (5+5 points)\n",
|
||||
"\n",
|
||||
"A \"convex set\" $S$ is one such that for any two points $x,y$ in $S$, the straight line connecting them is also in $S$. That is, $\\alpha x + (1-\\alpha)y \\in S$ for all $\\alpha \\in [0,1]$.\n",
|
||||
"\n",
|
||||
"**(a)** If $f_i(x)$ is a [convex function](https://en.wikipedia.org/wiki/Convex_function) (as defined in e.g. lecture 16), explain why the constraint $f_i(x) \\le 0$ defines a convex set of feasible points.\n",
|
||||
"\n",
|
||||
"**(b)** Explain why the intersection of two convex sets is a convex set. (Hence, the feasible set for many convex constraints is a convex set.)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "c5d653a2",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Solution:\n",
|
||||
"\n",
|
||||
"**(a)** Suppose $x$ and $y$ are feasible, i.e. $f_i(x) \\le 0$ and $f_i(y) \\le 0$. Then, by the definition from class of convex functions, we have:\n",
|
||||
"$$\n",
|
||||
"f_i(\\alpha x + (1-\\alpha) y) \\le \\alpha f_i(x) + (1-\\alpha) f_i(y) \\le 0\n",
|
||||
"$$\n",
|
||||
"for any $\\alpha \\in [0,1]$, which means that $\\alpha x + (1-\\alpha) y$ is feasible. So, the feasible set is convex.\n",
|
||||
"\n",
|
||||
"**(b)** If $S$ is convex and $T$ is convex, then consider their intersection $S \\cap T$. $x, y \\in S \\cap T$ means that $x, y \\in S$ and $x,y \\in T$, by the definition of intersection. But since $S,T$ are convex, this means that $\\alpha x + (1-\\alpha) y$ is $\\in S$ and $\\in T$, so it's also in $S \\cap T$. Hence $S \\cap T$ is convex."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "188e8b83",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 2 (10 points)\n",
|
||||
"\n",
|
||||
"Suppose we are solving the $\\ell^1$-regularized least-square problem:\n",
|
||||
"$$\n",
|
||||
"\\min_{x\\in \\mathbb{R}^n} \\left( \\Vert b - Ax \\Vert_2^2 + \\lambda \\Vert x \\Vert_1 \\right)\n",
|
||||
"$$\n",
|
||||
"where $\\lambda > 0$ is some regularization parameter and $ \\Vert x \\Vert_1 = \\sum_i |x_i|$ is the $\\ell^1$ norm.\n",
|
||||
"\n",
|
||||
"Similar to the examples in 16, show that this can be converted into a an equivalent [quadratic programming (QP)](https://en.wikipedia.org/wiki/Quadratic_programming) problem — a convex quadratic objective with affine constraints — by introducing one or more dummy variables (an \"epigraph\" formulation). This replaces the non-differentiable $\\ell^1$ problem with an equivalent *differentiable* problem with *differentiable constraints*."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "d0aa3f44",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Solution:\n",
|
||||
"\n",
|
||||
"We just need to apply the epigraph trick to *every term*. That is, suppose we introduce $n$ dummy variables $t \\in \\mathbb{R}^n$, and solve:\n",
|
||||
"$$\n",
|
||||
"\\min_{x,t\\in \\mathbb{R}^n} \\left( \\Vert b - Ax \\Vert_2^2 + \\sum_{i=1}^n t_i \\right)\n",
|
||||
"$$\n",
|
||||
"subject to the $2n$ affine constraints:\n",
|
||||
"$$\n",
|
||||
"t_i \\ge \\lambda x_i \\, \\qquad t_i \\ge -\\lambda x_i\n",
|
||||
"$$\n",
|
||||
"which will force $t_i = \\lambda |x_i|$ at the minimum. This is a [QP: a convex quadratic objective function + affine constraints](https://en.wikipedia.org/wiki/Quadratic_programming)."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "73427817",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 3 (5 + 5 points)\n",
|
||||
"\n",
|
||||
"Let\n",
|
||||
"$$\n",
|
||||
"A(p) = A_0 + \\underbrace{\\begin{pmatrix} p_1 & & & \\\\ & p_2 & & \\\\ & & \\ddots & \\\\ & & & p_m \\end{pmatrix}}_{\\text{diagm(p)}}\n",
|
||||
"$$\n",
|
||||
"where $A_0$ is some $m \\times m$ matrix and $p \\in \\mathbb{R}^m$ are $m$ parameters.\n",
|
||||
"\n",
|
||||
"Now, suppose we compute\n",
|
||||
"$$\n",
|
||||
"f(p) = g(\\underbrace{A(p)^{-1} b}_{x(p)})\n",
|
||||
"$$\n",
|
||||
"where $g(x) = x^T G x$, $b \\in \\mathbb{R}^m$ is some vector, and $G = G^T$ is some symmetric $m \\times m$ matrix.\n",
|
||||
"\n",
|
||||
"**(a)** What is $\\nabla f$? Explain how you can compute *all m* components of $\\nabla f$ by solving only *two* $m \\times m$ systems followed by $\\Theta(m)$ additional work in *total*. *Clarification*: What I really want is for you to show that, in addition to the cost of computing $f(p)$ (which involves solving one system $Ax=b$), the additional work to compute $\\nabla f$ only involves solving one more system plus $\\Theta(m)$ arithmetic.\n",
|
||||
"\n",
|
||||
"**(b)** Implement your solution from (a) by filling in the function `∇f(p)` below, and check that it correctly predicts $df = f(p + dp) - f(p)$ (*approximately*) for a random small $dp$.\n",
|
||||
"\n",
|
||||
"### Solution:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 31,
|
||||
"id": "d16dcac8",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"-8.034026689349616"
|
||||
]
|
||||
},
|
||||
"execution_count": 31,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"using LinearAlgebra\n",
|
||||
"\n",
|
||||
"# some A₀, b, G for m=5\n",
|
||||
"m = 5\n",
|
||||
"A₀ = randn(m,m)\n",
|
||||
"b = randn(m)\n",
|
||||
"G = randn(m,m); G = G' + G\n",
|
||||
"\n",
|
||||
"# our functions\n",
|
||||
"g(x) = x' * G * x\n",
|
||||
"f(p) = g((A₀ + Diagonal(p)) \\ b)\n",
|
||||
"\n",
|
||||
"p = randn(m) # some random parameters\n",
|
||||
"f(p) # make sure it gives a number out"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "6a093e50",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"**(a)**\n",
|
||||
"\n",
|
||||
"As explained in class (lecture 18) for a very similar problem, by the chain rule and the formula $d(A^{-1}) = -A^{-1} dA \\, A^{-1}$, we have:\n",
|
||||
"$$\n",
|
||||
"df = \\nabla g^T dx = \\underbrace{-\\nabla g^T A^{-1}}_{v^T} dA \\, \\underbrace{A^{-1} b}_x\n",
|
||||
"$$\n",
|
||||
"and so\n",
|
||||
"$$\n",
|
||||
"\\frac{\\partial f}{\\partial p_i} = v^T \\frac{\\partial A}{\\partial p_i} x\n",
|
||||
"$$\n",
|
||||
"where $v$ solves $A^T v = -\\nabla g$ (the \"adjoint\" equation), and $\\frac{\\partial A}{\\partial p_i}$ is just the matrix of all 0's except for a single 1 in the $i,i$ diagonal entry.\n",
|
||||
"\n",
|
||||
"So, to compute $\\nabla f$, we just need to do:\n",
|
||||
"\n",
|
||||
"1. Compute $x$ by solving $Ax=b$ and $\\nabla g = 2Gx$. (We need both $x$ and $Gx$ to compute $f$ in the first place, so this is \"free\"). \n",
|
||||
"2. Solve $A^T v = -\\nabla g = -2Gx$: this is *one* additional system of the same size.\n",
|
||||
"3. Compute the dot products $v^T \\frac{\\partial A}{\\partial p_i} x = v_i x_i$ for $i = 1,\\ldots,m$. This is $\\Theta(m)$ work: in Julia notation, we are just computing the elementwise product $\\nabla f = v \\mbox{ .* } x$.\n",
|
||||
"\n",
|
||||
"**(b)** The above algorithm is implemented as follows.\n",
|
||||
"\n",
|
||||
"To check that our $\\nabla f$ is correct, we simply want to check that our finite-difference `df` (for a small random `dp`) is $\\approx \\nabla f^T dp$. (The match won't be exact because `dp` is not infinitesimal.)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 33,
|
||||
"id": "12609b9b",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"df = 1.404071513277927e-6\n",
|
||||
"(∇f(p))' * dp = 1.4040715903119016e-6\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"1.4040715903119016e-6"
|
||||
]
|
||||
},
|
||||
"execution_count": 33,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"# your solution to (b):\n",
|
||||
"\n",
|
||||
"function ∇f(p)\n",
|
||||
" # in a practical setting, we would re-use this\n",
|
||||
" # from the f(p) computation:\n",
|
||||
" A = A₀ + Diagonal(p)\n",
|
||||
" x = A \\ b\n",
|
||||
" ∇g½ = G*x\n",
|
||||
" \n",
|
||||
" # adjoint equation:\n",
|
||||
" # (in practice you might re-use the LU factorization of A!)\n",
|
||||
" v = A' \\ (-2∇g½)\n",
|
||||
" \n",
|
||||
" ∇f = v .* x # elementwise products, from vᵀ(∂A/∂pᵢ)x\n",
|
||||
"end\n",
|
||||
"\n",
|
||||
"dp = randn(m) * 1e-8 # a random small dp\n",
|
||||
"df = f(p + dp) - f(p)\n",
|
||||
"\n",
|
||||
"# check: ∇f approximately predicts df\n",
|
||||
"@show df # finite-difference approximation\n",
|
||||
"@show ∇f(p)'dp # gradient prediction"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "43a63240",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Yup, they match to about 8 significant digits, which is about as much as we can hope for given a finite $dp \\sim 10^{-8}$."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "84c4326e",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 4 (5+5+5 points)\n",
|
||||
"\n",
|
||||
"In class, we considered steepest descent for $f(x) = \\kappa x_1^2 + x_2^2$ in $\\mathbb{R}^2$, and argued that for an arbitrary $x = [x_1,x_2]$ starting point and $\\kappa \\gg 1$, the steepest-descent step $x \\leftarrow x - sz$ is *approximately* $sz \\approx [x_1, x_2/\\kappa]$.\n",
|
||||
"\n",
|
||||
"### Part (a):\n",
|
||||
"\n",
|
||||
"**(a)** If $sz = [x_1, x_2/\\kappa]$ exactly, then on the next step we would have the new $x \\leftarrow x - sz = [0, (1-\\frac{1}{\\kappa}) x_2)]$. However, explain why a more careful calculation shows that the new $x - sz \\approx [O(1/\\kappa^2), x_2-\\frac{1}{\\kappa} x_2 + O(1/\\kappa^3)]$, i.e. the first component is proportional to $1/\\kappa^2$ to leading order in $1/\\kappa$.\n",
|
||||
"\n",
|
||||
"### Solution:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "946cd3d5",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"From class, each step is $x \\to x - sz$ where $z = Ax = [\\kappa x_1, x_2]$ and\n",
|
||||
"$$\n",
|
||||
"s = \\frac{z^T z}{z^T A z} = \\frac{\\kappa^2 x_1^2 + x_2^2}{\\kappa^3 x_1^2 + x_2^2} = \\kappa^{-1} \\frac{ 1 + (x_2/x_1)^2 \\kappa^{-2}}{1 + (x_2/x_1)^2 \\kappa^{-3}}\n",
|
||||
"$$\n",
|
||||
"Taylor expanding in $\\kappa^{-1}$, we obtain:\n",
|
||||
"$$\n",
|
||||
"s \\approx \\kappa^{-1} + (x_2/x_1)^2 \\kappa^{-3} + O(\\kappa^{-4})\n",
|
||||
"$$\n",
|
||||
"Hence\n",
|
||||
"$$\n",
|
||||
"x - sz \\approx \\begin{pmatrix} x_1 - (1 + (x_2/x_1)^2 \\kappa^{-2} + O(\\kappa^{-3}))x_1 \\\\ x_2 - \\kappa^{-1} x_2 + O(\\kappa^{-3}) \\end{pmatrix}\n",
|
||||
"= \\begin{pmatrix} - (x_2^2/x_1) \\kappa^{-2} + O(\\kappa^{-3}) \\\\ (1 - \\kappa^{-1}) x_2 + O(\\kappa^{-3}) \\end{pmatrix}\n",
|
||||
"$$\n",
|
||||
"which is of the suggested form."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "07b0921e",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Part (b):\n",
|
||||
"\n",
|
||||
"**(b)** If you start with an $x = [\\#/\\kappa^2, x_2]$, i.e. where the first component is proportional to $1/\\kappa^2$ and $\\#$ is some number of the same order of magnitude as $x_2$, show that after one steepest-descent step (for $\\kappa \\gg 1$) the $x_1$ component is *still* roughly order $1/\\kappa^2$ but of the opposite sign, and $x_2$ again subtracts a term roughly proportional to $x_2/\\kappa$.\n",
|
||||
"\n",
|
||||
"### Solution:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "faf3c6e5",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"Conversely, if $x_1 = \\# / \\kappa$, we have $x_2/x_1 = \\kappa (x_2 / \\#)$ and $s$ becomes:\n",
|
||||
"$$\n",
|
||||
"s = \\kappa^{-1} \\frac{ 1 + (x_2/x_1)^2 \\kappa^{-2}}{1 + (x_2/x_1)^2 \\kappa^{-3}} = \\kappa^{-1} \\frac{ 1 + (x_2 / \\#)^2}{1 + (x_2 / \\#)^2 \\kappa^{-1}} \\approx \\kappa^{-1} (1 + (x_2 / \\#)^2) + O(\\kappa^{-2})\n",
|
||||
"$$\n",
|
||||
"So our steepest-descent step becomes, to leading order in $1/\\kappa$:\n",
|
||||
"$$\n",
|
||||
"x - sz \\approx \\begin{pmatrix} -(x_2 / \\#)^2 \\kappa^{-1} \\\\ x_2 \\left[1 - \\kappa^{-1} (1 + (x_2 / \\#)^2)\\right] \\end{pmatrix}\n",
|
||||
"$$\n",
|
||||
"i.e. the first term is still of order $1/\\kappa$ but with the opposite sign, and the second term has shrunk by a term of $O(x_2 / \\kappa)$. \n",
|
||||
"\n",
|
||||
"It is **zig-zagging slowly towards the solution!**\n",
|
||||
"\n",
|
||||
"**Optional**:\n",
|
||||
"\n",
|
||||
"We can also analyze the case that we get after one step, in part (a):\n",
|
||||
"\n",
|
||||
"If $x_1 = \\# / \\kappa^2$, we have to re-do our calculation of $s$, since $x_2/x_1 = \\kappa^2 (x_2 / \\#)$:\n",
|
||||
"$$\n",
|
||||
"s = \\kappa^{-1} \\frac{ 1 + (x_2/x_1)^2 \\kappa^{-2}}{1 + (x_2/x_1)^2 \\kappa^{-3}} = \\kappa^{-1} \\frac{ 1 + (x_2 / \\#)^2 \\kappa^{2}}{1 + (x_2 / \\#)^2 \\kappa} \\approx 1 - (\\# / x_2)^2 \\kappa^{-1} + O(\\kappa^{-2})\n",
|
||||
"$$\n",
|
||||
"So our steepest-descent step becomes, to leading order in $1/\\kappa$:\n",
|
||||
"$$\n",
|
||||
"x - sz \\approx \\kappa^{-1} \\begin{pmatrix} -\\# \\\\ +x_2(\\# / x_2)^2 \\end{pmatrix}\n",
|
||||
"$$\n",
|
||||
"i.e. *both* terms are of order $\\kappa^{-1}$, with the first term having an opposite sign from what it started with.\n",
|
||||
"\n",
|
||||
"(So, if this is the second iteration, then we are roughly back where we started but we gained a factor of $\\kappa^{-1}$. That means, for this starting point it actually converges *much* faster, gaining about a factor of $\\kappa^{-1}$. every 2 iterations. Not every starting point zig-zags! But we are mainly interested in the *worst-case* convergence since we can't count on getting lucky with our initial guess.)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "770529a8",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Part (c):\n",
|
||||
"\n",
|
||||
"**(c)** Implement this steepest-descent process numerically for $\\kappa = 100$ and a starting point $x_1 = 0.01234567$ and $x_2 = 0.8910$. Plot $100x_1$ and $x_2$ for 100 iterations of steepest descent. A more careful analysis would show a convergence proportional to $\\left( \\frac{\\kappa - 1}{\\kappa + 1} \\right)^k$, where $k$ is the iteration number, following equation (4) in the Strang book — include this function for comparison on your plot.\n",
|
||||
"\n",
|
||||
"### Solution:\n",
|
||||
"\n",
|
||||
"Since $x_1$ is of order $1/\\kappa$, we expect something similar to the second case in part (b): from one iteration to the next, it should stay of roughly the same order but flip sign, and both components should asymptotically shrink slowly following equation (4):"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 28,
|
||||
"id": "bef5c3ee",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAHFCAYAAAD/kYOsAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAADJWElEQVR4nOxdeXgUVfY9lZ0QEgiQBRI2FxYRkAQUFFkUGHBBIgIqm4rKgIIwjoqMCs7PQR1UXFBcUBSRoBDXQRGVTQU0CbggKFtYA2ENJJC9fn9UXter6qrqV91VnQD3fF9/3fXq9quq7k73ybnn3ifJsiyDQCAQCAQC4RxCSE2fAIFAIBAIBILTIIJDIBAIBALhnAMRHAKBQCAQCOcciOAQCAQCgUA450AEh0AgEAgEwjkHIjgEAoFAIBDOORDBIRAIBAKBcM6BCA6BQCAQCIRzDkRwCAQCgUAgnHMggkMgEGoFVq1aBUmSsGTJEp+x06dPhyRJQTgr/3Do0CE0bNhQ+HoIBILzIIJDIBAIDmPChAmIioqq6dMgEM5rEMEhEAiu4cyZMzV9CkHH0qVLsXz5csycObOmT4VAOK9BBIdAIJiCpYI2btyIjIwMxMbGIi4uDiNGjMDhw4c1sS1atMD111+PrKwsXHbZZYiKisKMGTMAAL///jsGDRqEBg0aICoqCp06dcK7775reMySkhJMmTIFSUlJqFOnDnr27ImNGzcKne/ixYvRrVs31K1bFzExMejfv7/Xc8eMGYOYmBhs3boV/fv3R926dZGcnIynn34aALB+/XpcddVVqFu3Li6++GLT8zTCsWPHMGHCBDz11FNo1qyZadzWrVtx6623IjExEZGRkWjWrBlGjRqF0tJS4WMRCARrEMEhEAg+MXjwYFx44YVYsmQJpk+fjk8++QT9+/dHeXm5Ji43Nxf//Oc/MXHiRHz11Ve4+eab8eeff6J79+7YvHkzXnrpJWRlZaFdu3YYM2YMnn32Wa9jPfroo9i5cyfeeustvPXWWzhw4AB69eqFnTt3Wp7jf/7zH9x6661o164dPvzwQyxYsACnTp1Cjx498Mcff2hiy8vLkZGRgeuuuw6ffvopBgwYgKlTp+LRRx/F6NGjceedd+Ljjz9G69atMWbMGOTk5Ai9ThMnTkTLli1x3333mcb88ssv6NKlC9avX48nn3wSX375JWbOnInS0lKUlZUJHYdAIAhAJhAIBBM88cQTMgB58uTJmvGFCxfKAOT333/fM9a8eXM5NDRU/vPPPzWxw4cPlyMjI+U9e/ZoxgcMGCBHR0fLJ06ckGVZlleuXCkDkDt37ixXVVV54vLy8uTw8HB57NixXufFsGfPHjksLEy+//77Ncc4deqUnJSUJA8dOtQzNnr0aBmAvHTpUs9YeXm53LhxYxmAnJub6xk/evSoHBoaKk+ZMsXna/XFF1/I4eHh8m+//aa5no8++kgT16dPH7l+/fpyQUGBzzkJBIL/IAWHQCD4xO23367ZHjp0KMLCwrBy5UrNeIcOHXDxxRdrxr777jtcc801SE1N1YyPGTMGp0+fxrp16zTjt912m6ZCqnnz5ujevbvXsXgsX74cFRUVGDVqFCoqKjy3qKgo9OzZE6tWrdLES5KEgQMHerbDwsJw4YUXIjk5GZdddplnPD4+HgkJCdi9e7fpsQGgsLAQ9957Lx5++GG0b9/eNO706dNYvXo1hg4disaNG1vOSSAQAkNYTZ8AgUCo/UhKStJsh4WFoWHDhjh69KhmPDk52eu5R48eNRxv0qSJZ7/VsdjYL7/8Ynp+hw4dAgB06dLFcH9IiPZ/uejoaK8qp4iICMTHx3s9NyIiAiUlJabHBoBp06YhPDwc9913H06cOAEAKCoqAqCQmhMnTiAuLg7Hjx9HZWUlUlJSLOcjEAiBgwgOgUDwiYMHD6Jp06ae7YqKChw9ehQNGzbUxBn1pmnYsCHy8/O9xg8cOAAAaNSokdexjI6vPxYPNseSJUvQvHlziytxB7///jvy8vIMydno0aMBAMePH0d8fDxCQ0Oxb9++YJ8igXDegQgOgUDwiYULFyItLc2z/eGHH6KiogK9evXy+dxrrrkGH3/8MQ4cOOBRbQDgvffeQ3R0NK644gpN/KJFizBlyhQPWdq9ezd+/PFHjBo1yvQY/fv3R1hYGHbs2IGbb77Z5tUFjtmzZ3uUG4ZNmzZh8uTJmD59Onr27ImYmBiEhYWhZ8+e+Oijj/DUU095kTsCgeAciOAQCASfyMrKQlhYGPr27YvNmzfjscceQ8eOHTF06FCfz33iiSfwxRdfoHfv3nj88ccRHx+PhQsX4n//+x+effZZxMXFaeILCgowePBg3H333SgsLMQTTzyBqKgoTJ061fQYLVq0wJNPPolp06Zh586d+Nvf/oYGDRrg0KFD+Omnn1C3bl1Pybob6NSpk+m+Sy65REMEn3/+eVx11VW4/PLL8cgjj+DCCy/EoUOH8Nlnn+H1119HvXr1XDtPAuF8AhEcAoHgE1lZWZg+fTpee+01SJKEG264AbNnz0ZERITP57Zu3Ro//vgjHn30UUyYMAFnzpxB27Zt8c4772DMmDFe8f/5z3/w888/44477sDJkyfRtWtXZGZm4oILLrA8ztSpU9GuXTu8+OKLWLRoEUpLS5GUlIQuXbpg3Lhx/l664+jYsSN++uknPPHEE5g6dSpOnTqFpKQk9OnTR+j1JBAIYpBkWZZr+iQIBELtxPTp0zFjxgwcPnyY0ikEAuGsApWJEwgEAoFAOOdABIdAIBAIBMI5B0pREQgEAoFAOOdACg6BQCAQCIRzDkRwCAQCgUAgnHMggkMgEAgEAuGcw3nZB6eqqgoHDhxAvXr1DFvLEwgEAoFAqH2QZRmnTp1CkyZNvNaY0+O8JDgHDhzwWtmYQCAQCATC2YG9e/f6XLT2vCQ4rBX63r17ERsbW8NnQyAQCAQCQQQnT55Eamqq0JIm5yXBYWmp2NhYIjgEAoFAIJxlELGXkMmYQCAQCATCOQciOAQCgUAgEM45EMEhEAgEAoFwzuG89OAQCAQCwV1UVlaivLy8pk+DcBYiIiLCZwm4CIjgEAgEAsExyLKMgwcP4sSJEzV9KoSzFCEhIWjZsiUiIiICmocIDoFAIBAcAyM3CQkJiI6OpmaqBFtgjXjz8/PRrFmzgD4/RHAIBAKB4AgqKys95KZhw4Y1fTqEsxSNGzfGgQMHUFFRgfDwcL/ncdVkvGbNGtxwww1o0qQJJEnCJ598YhmflZWFvn37onHjxoiNjUW3bt2wfPlyTcz8+fMhSZLXraSkxMUrIRAIBIIvMM9NdHR0DZ8J4WwGS01VVlYGNI+rBKe4uBgdO3bEK6+8IhS/Zs0a9O3bF8uWLUNOTg569+6NG264ARs3btTExcbGIj8/X3OLiopy4xIIBAKBYBOUliIEAqc+P66mqAYMGIABAwYIx8+ePVuz/Z///AeffvopPv/8c1x22WWecUmSkJSU5NRpOofKSmDtWiA/H0hOBnr0AEJDa/qsCAQCgUA471Cr++BUVVXh1KlTiI+P14wXFRWhefPmSElJwfXXX++l8OhRWlqKkydPam6OIysLaNEC6N0buO025b5FC2WcQCAQCARCUFGrCc5zzz2H4uJiDB061DPWpk0bzJ8/H5999hkWLVqEqKgoXHnlldi2bZvpPDNnzkRcXJzn5vhK4llZwJAhwL592vH9+5VxIjkEAoEgjspKYNUqYNEi5T5ALwbh/EStJTiLFi3C9OnTsXjxYiQkJHjGr7jiCowYMQIdO3ZEjx498OGHH+Liiy/Gyy+/bDrX1KlTUVhY6Lnt3bvXuROtrAQmTQJk2XsfG3vgAfoDJRAIBBHUkBrO/mHev3+/Z2zs2LHo0KEDCgsLXT02wR3USoKzePFi3HXXXfjwww9x7bXXWsaGhISgS5culgpOZGSkZ+Vwx1cQX7vWW7nhIcvA3r1KHIFAIBDMUYNq+PDhw9G6dWvMnDkTADBjxgwsX74cX375JeLi4lw7LsE91Lo+OIsWLcKdd96JRYsW4brrrvMZL8syNm3ahEsvvTQIZ2eA/Hxn4wgEAuFcgSwDp0+LxVZWAhMnmqvhkqSo5ddeK1a8ER2tPEcQkiThqaeewpAhQ9CkSRO8+OKLWLt2LZo2bQoAGDx4MFatWoVrrrkGS5YsEZ6XUHNwleAUFRVh+/btnu1du3Zh06ZNiI+PR7NmzTB16lTs378f7733HgCF3IwaNQovvvgirrjiChw8eBAAUKdOHQ+DnjFjBq644gpcdNFFOHnyJF566SVs2rQJc+bMcfNSzJGc7GwcgUAgnCs4fRqIiXFmLllWlB1RNaWoCKhb19Yhrr/+erRr1w4zZszA119/jUsuucSzb+LEibjzzjvx7rvv2pqTUHNwNUWVnZ2Nyy67zFPiPWXKFFx22WV4/PHHAQD5+fnYs2ePJ/71119HRUUFJkyYgOTkZM9t0qRJnpgTJ07gnnvuQdu2bdGvXz/s378fa9asQdeuXd28FHP06AGkpJj/pyBJQGqqEkcgEAiEWovly5dj69atqKysRGJiomZf7969Ua9evRo6M4I/kGTZSA88t3Hy5EnExcWhsLDQGT8OyxsDWnmVkZ4lS4CMjMCPQyAQCLUYJSUl2LVrF1q2bKk0X7WTolqzBhg40HfcsmXA1Vf7jrOZosrNzUWvXr0wZ84cZGZmIjo6Gh999JEmZtWqVXjllVcoReUyvD5HHOz8ftc6D85ZiYwMhcTcd5/Wa5OSAsyeTeSGQCCcn5Ak8TRRv37Kd+b+/cY+HElS9vfr53gD1by8PFx33XV45JFHMHLkSLRr1w5dunRBTk4O0tLSHD0WIXiolVVUZyUyMgC+4eA33wC7dnmTG+rvQCAQCN4IDQVefFF5rFde2Pbs2Y6Tm2PHjmHAgAG48cYb8eijjwIA0tLScMMNN2DatGmOHosQXJCC4yQLine truncated
|
||||
"text/plain": [
|
||||
"Figure(PyObject <Figure size 640x480 with 1 Axes>)"
|
||||
]
|
||||
},
|
||||
"metadata": {},
|
||||
"output_type": "display_data"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"PyObject <matplotlib.legend.Legend object at 0x296e625b0>"
|
||||
]
|
||||
},
|
||||
"execution_count": 28,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"x = [0.01234567, 0.8910] # starting guess\n",
|
||||
"κ = 100\n",
|
||||
"A = [ κ 0\n",
|
||||
" 0 1 ]\n",
|
||||
"xvals = [ x ] # list of x values\n",
|
||||
"for k = 2:100\n",
|
||||
" z = A*x\n",
|
||||
" s = z'z / (z'*A*z) # = (κ^2*x[1]^2 + x[2]^2) / (κ^3*x[1]^2 + x[2]^2)\n",
|
||||
" x = x - s*z\n",
|
||||
" push!(xvals, x)\n",
|
||||
"end\n",
|
||||
"\n",
|
||||
"using PyPlot\n",
|
||||
"plot(1:100, getindex.(xvals,1) * κ, \"ro-\")\n",
|
||||
"plot(1:100, getindex.(xvals,2), \"b*-\")\n",
|
||||
"plot(1:100, ((κ-1)/(κ+1)).^(1:100), \"k--\")\n",
|
||||
"xlabel(\"iteration\")\n",
|
||||
"title(\"problem 4c\")\n",
|
||||
"legend([L\"x_1\", L\"x_2\", \"Strang eq. (4)\"])"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "6e9798d4",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 5 (10 points)\n",
|
||||
"\n",
|
||||
"In the Adam algorithm for stochastic optimization, we keep a \"running average\" (\"1st moment\") $m^{(k)}$ of the estimated gradients. If $\\tilde{\\nabla}f^{(k)}$ is our estimated (mini-batch) gradient at $x^{(k)}$ (the k-th step, for $k=1,2,3,\\ldots$), then the running average is computed as:\n",
|
||||
"$$\n",
|
||||
"m^{(k)} = \\beta_1 m^{(k-1)} + (1-\\beta_1) \\tilde{\\nabla}f^{(k)} \\, ,\n",
|
||||
"$$\n",
|
||||
"where $0 \\le \\beta_1 < 1$ (typically $\\beta_1 = 0.9$) is a given weight, and we define $m^{(0)} = \\vec{0}$.\n",
|
||||
"\n",
|
||||
"However, before we use the gradient to take a step, we \"de-bias\" it to:\n",
|
||||
"$$\n",
|
||||
"\\hat{m}^{(k)} = \\frac{m^{(k)}}{1 - \\beta_1^k}\n",
|
||||
"$$\n",
|
||||
"**Explain** why $\\hat{m}^{(k)}$ is exactly a weighted average of the estimates $\\tilde{\\nabla}f^{(j)}$ for $j = 1,\\ldots,k$ — that is, explain why dividing by $1 - \\beta_1^k$ exactly corresponds to dividing by the sum of the weights of each $\\tilde{\\nabla}f^{(j)}$ term. (The formulas for the [sum of a geometric series](https://en.wikipedia.org/wiki/Geometric_series) will be helpful.)\n",
|
||||
"\n",
|
||||
"### Solution:\n",
|
||||
"\n",
|
||||
"Just expanding out the formula (with $m^{(0)} = 0$ — there is no momentum to start with), we get:\n",
|
||||
"$$\n",
|
||||
"m^{(k)} = \\beta_1 m^{(k-1)} + (1-\\beta_1) \\tilde{\\nabla}f^{(k)} \\\\\n",
|
||||
"= \\beta_1^2 m^{(k-2)} + \\beta_1 (1-\\beta_1) \\tilde{\\nabla}f^{(k-1)} + (1-\\beta_1) \\tilde{\\nabla}f^{(k)} \\\\\n",
|
||||
"= \\beta_1^3 m^{(k-3)} + \\beta_1^2 (1-\\beta_1) \\tilde{\\nabla}f^{(k-2)} + \\beta_1 (1-\\beta_1) \\tilde{\\nabla}f^{(k-1)} + (1-\\beta_1) \\tilde{\\nabla}f^{(k)}\n",
|
||||
"= \\cdots \\\\\n",
|
||||
"= (1-\\beta_1) \\sum_{i=0}^{k-1} \\beta_1^i \\tilde{\\nabla}f^{(k-i)} \\, .\n",
|
||||
"$$\n",
|
||||
"The sum of all the coefficients is, using the usual formula for a finite geometric series:\n",
|
||||
"$$\n",
|
||||
"(1-\\beta_1) \\sum_{i=0}^{k-1} \\beta_1^i = (1-\\beta_1) \\frac{1 - \\beta_1^k}{1-\\beta_1} = 1 - \\beta_1^k \\, .\n",
|
||||
"$$\n",
|
||||
"So, $1 - \\beta_1^k$ is *exactly* the what we would divide by to make $\\hat{m}^{(k)}$ correspond to simply a weighted average of all the $\\tilde{\\nabla}f^{(k-i)}$ values. "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "c7f174a6",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Problem 6 (10 points)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "a3e6a871",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The Julia notebook from class (lecture 21) implemented the Adam algorithm for stochastic optimization. This code is reproduced below.\n",
|
||||
"\n",
|
||||
"A more recent algorithm, that has better convergence guarantees than Adam, is called [YOGI (Zaheer et al., 2018)](https://papers.nips.cc/paper_files/paper/2018/hash/90365351ccc7437a1309dc64e4db32a3-Abstract.html). Make a copy of the code below and **modify it to implement the YOGI algorithm** as described in algorithm 2 of the linked paper. (This algorithm omits the \"de-biasing\" step of problem 5, but you should include that in your implementation; the de-biasing is identical to that of Adam.) **Plot the convergence** of YOGI for this problem, with the same hyper-parameters, along with the original Adam algorithms (using different colors superimposed on a single plot)."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 64,
|
||||
"id": "2181180f",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA04AAAGICAYAAAB2jrABAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd1hT5xcH8G8IEPbeG9youHCLgnujiKjVVuuq1ZZaV2tta5e1vzqqto7W1bonjroXKO69xa3svXdyc39/hFCBQG7CvbkJvJ/n4Xk03HEIIbnnvuc9r4CmaRoEQRAEQRAEQRBEtfT4DoAgCIIgCIIgCELbkcSJIAiCIAiCIAhCCZI4EQRBEARBEARBKEESJ4IgCIIgCIIgCCVI4kQQBEEQBEEQBKEESZwIgiAIgiAIgiCUIIkTQRAEQRAEQRCEEiRxIgiCIAiCIAiCUIIkTgRBEARBEARBEEqQxImoNS8vL0yYMKH8/1FRURAIBIiKiuL0vGvWrMHff/9d6+O8efMGAoGAlWNpkx07dmDFihWMtw8MDERgYGD5/wsLC/Hdd99x/ntURlviIAiC0EVeXl747rvv+A5DZUzj/u677yAQCLgPiCAA6PMdAFH3tG3bFleuXIGvry+n51mzZg3s7OwqJG3Ef3bs2IGHDx9i5syZjLZfs2ZNhf8XFhbi+++/B4AKCZWmaUscBEEQhPaZPHky+vfvz3cYRD1BEqd6iqZpFBcXw9jYmPVjW1hYoFOnTqwfl+AW14munFgshkAggL4+efshCILQNYWFhTAxManyOBvXFUVFRTAyMlJpBMnNzQ1ubm5qn5MgVEFK9XTcoUOH4OfnB5FIBB8fH6xcuVLhsLVAIMAnn3yCdevWoVmzZhCJRPjnn38AAN9//z06duwIGxsbWFhYoG3btti4cSNomq5wDLFYjHnz5sHJyQkmJibo1q0brl+/XiWm6kr1bt68iaFDh8LGxgZGRkZo06YN9uzZU2Gbv//+GwKBAJGRkfj4449hZ2cHW1tbhISEIDExsXw7Ly8vPHr0COfPn4dAIIBAIICXl5fS5ysxMRFhYWEwNzeHpaUlRo0aheTkZIXbMom3sLAQc+bMgbe3N4yMjGBjYwN/f3/s3LmzwnbXrl3DkCFDYGtrCyMjIzRo0KDKSNDz58/x3nvvwcHBASKRCM2aNcPq1asVPrc7d+7EggUL4OLiAgsLC/Tu3RtPnz4t3y4wMBBHjx7F27dvy58fZR9E75bqvXnzBvb29gBkrw/5/u+O7qkS79atWzF79my4urpCJBLhxYsXSEtLw/Tp0+Hr6wszMzM4ODigZ8+eiI6OLt+frTikUil++uknNGnSBMbGxrCysoKfnx9WrlxZ43Mybdo0GBkZ4datWxWO1atXLzg6OiIpKanG/QmCILTV7t270blzZ5iamsLMzAz9+vXDnTt3KmwzYcIEmJmZ4cGDB+jbty/Mzc3Rq1cvADVfV1y8eBG9evWCubk5TExM0KVLFxw9erTCseWf96dOncLEiRNhb28PExMTlJSUqPRzKLrm8fLywuDBg3HixAm0bdsWxsbGaNq0KTZt2qT0eL/88gv09PTw77//VnkuTExM8ODBA5XiI+oWcstXh504cQIhISHo3r07du/eDYlEgqVLlyIlJUXh9gcPHkR0dDS+/fZbODk5wcHBAYDs4vSjjz6Ch4cHAODq1av49NNPkZCQgG+//bZ8/ylTpmDLli2YM2cO+vTpg4cPHyIkJAR5eXlKY42MjET//v3RsWNHrFu3DpaWlti1axdGjRqFwsLCKuV2kydPxqBBg7Bjxw7ExcVh7ty5GDduHM6dOwcAOHDgAEJDQ2FpaVleYiYSiWqMoaioCL1790ZiYiIWL16Mxo0b4+jRoxg1apTa8c6aNQtbt27FTz/9hDZt2qCgoAAPHz5ERkZG+bFOnjyJIUOGoFmzZli+fDk8PDzw5s0bnDp1qnybx48fo0uXLvDw8MCyZcvg5OSEkydPIjw8HOnp6Vi4cGGF+L766it07doVGzZsQG5uLr744gsMGTIET548gVAoxJo1azB16lS8fPkSBw4cUPr7qczZ2RknTpxA//79MWnSJEyePBkAypMYVeOdP38+OnfujHXr1kFPTw8ODg5IS0sDACxcuBBOTk7Iz8/HgQMHEBgYiLNnzyIwMJC1OH799Vd89913+Prrr9G9e3eIxWLExMQgOzu7xudhxYoVuHbtGsLCwnDr1i1YWVnh+++/R1RUFE6cOAFnZ2eVn1uCIAi+/fzzz/j666/x4Ycf4uuvv0ZpaSmWLFmCgIAAXL9+vUIFQmlpKYYOHYqPPvoIX375JSQSSfn3FF1XnD9/Hn369IGfnx82btwIkUiENWvWYMiQIdi5c2eVz9yJEydi0KBB2Lp1KwoKCmBgYMDKz3jv3j3Mnj0bX375JRwdHbFhwwZMmjQJDRs2RPfu3avd74svvkB0dDTGjx+PO3fuwNPTE5s3b8Y///yDDRs2oGXLlqzER+gomtBZ7du3p93d3emSkpLyx/Ly8mhbW1u68q8WAG1paUlnZmbWeEyKomixWEz/8MMPtK2tLS2VSmmapuknT57QAOjPP/+8wvbbt2+nAdDjx48vfywyMpIGQEdGRpY/1rRpU7pNmza0WCyusP/gwYNpZ2dnmqIomqZpevPmzTQAevr06RW2+/XXX2kAdFJSUvljzZs3p3v06FHjz/OutWvX0gDoQ4cOVXh8ypQpNAB68+bNKsfbokULetiwYTWet0GDBnSDBg3ooqKiarfp168f7ebmRufk5FR4/JNPPqGNjIzKf2/y53bgwIEVttuzZw8NgL5y5Ur5Y4MGDaI9PT1rjO1dPXr0qPB8pqWl0QDohQsX1jre7t27Kz2/RCKhxWIx3atXL3r48OGsxjF48GC6devWSmNQ5Pnz57SFhQU9bNgw+syZM7Senh799ddfq3UsgiAITfP09Kzw/hkbG0vr6+vTn376aYXt8vLyaCcnJzosLKz8sfHjx9MA6E2bNlU5bnXXFZ06daIdHBzovLy88sckEgndokUL2s3Nrfy6Qv55/8EHHzCKuzoLFy6scs3j6elJGxkZ0W/fvi1/rKioiLaxsaE/+ugjpcdMT0+n3dzc6A4dOtC3b9+mTUxM6HHjxindj6j7SKmejiooKMDNmzcxbNgwGBoalj9uZmaGIUOGKNynZ8+esLa2rvL4uXPn0Lt3b1haWkIoFMLAwADffvstMjIykJqaCkA2AgMAY8eOrbBvWFiY0rkqL168QExMTPm+Eomk/GvgwIFISkqqUGYGAEOHDq3wfz8/PwDA27dvazwXAFAUVeEcUqm0/GcwNzevcuz33ntP7Xg7dOiA48eP48svv0RUVBSKiooqHOvZs2d4+fIlJk2aBCMjI4XxFhcX4+zZsxg+fDhMTEyqnK+4uBhXr15l7fmpLXXiHTFihMJjrVu3Dm3btoWRkRH09fVhYGCAs2fP4smTJ6zG0aFDB9y7dw/Tp0/HyZMnkZuby/jnbdiwIdavX4+DBw9i8ODBCAgI0MkOVQRBEICsCkIikeCDDz6o8L5pZGSEHj16KOxgWt17eOXrioKCAly7dg2hoaEwMzMrf1woFOL9999HfHx8lc/76o5dW61bty6vpAEAIyMjNG7cmNHnpK2tLXbv3o3bt2+XVzWsW7eOkzgJ3UISJx2VlZUFmqbh6OhY5XuKHgOgsKzo+vXr6Nu3LwBg/fr1uHTpEm7cuIEFCxYAQHkiIC89c3JyqrC/vr4+bG1ta4xVXjo4Z84cGBgYVPiaPn06ACA9Pb3CPpWPKS/Dq5yYKNKgQYMK5/jhhx/KfwZFz03ln0mVeFetWoUvvvgCBw8eRFBQEGxsbDBs2DA8f/4cAMrL0WqauJqRkQGJRILff/+9yvkGDhzI+vNTW+rEq+i1t3z5cnz88cfo2LEj9u/fj6tXr+LGjRvo378/o59DlTjmz5+PpUuX4urVqxgwYABsbW3Rq1cv3Lx5k9HPPGjQIDg6OqK4uBizZs2CUChktB9BEIS2kX/GtW/fvsp75+7du6u8f5uYmMDCwkLhsSq/t8uvTRS957u4uABAhVJ2Rcdgi6JrE5FIxPhzsmPHjmjevDmKi4vx8ccfw9TUlO0QCR1E5jjpKGtrawgEAoXzmaprdqCoOcCuXbtgYGCAI0eOVBgROXjwYIXt5G9AycnJcHV1LX9cIpFUeROszM7ODoDs4jUkJEThNk2aNKnxGKr4999/K0wulb9Z29raKmxmUfn5UiVeU1NTfP/99/j++++RkpJSPvo0ZMgQxMTElM/FiY+PrzZea2vr8rtxM2bMULiNt7d3tftrmjrxKnrtbdu2DYGBgVi7dm2Fx5nMmVM1Dn19fcyaNQuzZs1CdnY2zpw5g6+++gr9+vVDXFycwg5R75o2bRry8vLQvHlzhIeHIyAgQOHoLUEQhLaTf8bt27cPnp6eSrevqbFQ5e9ZW1tDT09PYeMceYMn+fmZHJ9PCxcuxIMHD9CuXTt8++23GDx4MHx8fPgOi+AZSZx0lKmpKfz9/XHw4EEsXbq0vFwvPz8fR44cYXwceVvod++gFxUVYevWrRW2k3db2759O9q1a1f++J49eypMFFWkSZMmaNSoEe7du4eff/6ZcWzKVHfnqLqJm0FBQdizZw8OHz5codRtx44drMTr6OiICRMm4N69e1ixYgUKCwvRuHFjNGjQAJs2bcKsWbMUNrAwMTFBUFAQ7ty5Az8/vwqll7Whyp216vYHqo5isRWvQCCo8nzcv38fV65cgbu7O2dxWFlZITQ0FAkJCZg5cybevHlTYyv2DRs2YNu2bdi0aRN69OiBtm3b4sMPP6xyc4EgCEIX9OvXD/r6+nj58iXrZXKmpqbo2LEjIiIisHTp0vLW5FKpFNu2bYObmxsaN27M6jm5cPr0aSxevBhff/01Zs6cidatW2PUqFG4dOkSa5/RhG4Line truncated
|
||||
"text/plain": [
|
||||
"Figure(PyObject <Figure size 1000x400 with 2 Axes>)"
|
||||
]
|
||||
},
|
||||
"metadata": {},
|
||||
"output_type": "display_data"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"PyObject Text(0.5, 1.0, '|error| in x')"
|
||||
]
|
||||
},
|
||||
"execution_count": 64,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"using LinearAlgebra, PyPlot, Random\n",
|
||||
"\n",
|
||||
"# typical Adam hyper-parameters\n",
|
||||
"ϵ = 1e-8\n",
|
||||
"β₁ = 0.9\n",
|
||||
"β₂ = 0.999\n",
|
||||
"\n",
|
||||
"# more parameters\n",
|
||||
"log₁₀α = -1\n",
|
||||
"α = exp10(log₁₀α) # approx. maximum stepsize |Δx|\n",
|
||||
"M = 1 # mini-batch size\n",
|
||||
"iterations = 300\n",
|
||||
"\n",
|
||||
"Random.seed!(314159)\n",
|
||||
"\n",
|
||||
"x̂₀ = [1,2] # ground-truth model\n",
|
||||
"x = [0.0,0.0]\n",
|
||||
"X = [x]\n",
|
||||
"ΔX = [norm(x - x̂₀)]\n",
|
||||
"m = zero(x) # initial 1st-moment (mean gradient)\n",
|
||||
"v = zero(x) # initial 2nd-moment (mean gradient.^2)\n",
|
||||
"for t = 1:iterations\n",
|
||||
" # generate a mini-batch of size M:\n",
|
||||
" ã = randn(M)\n",
|
||||
" b̃ = @. 1 + 2ã + randn()\n",
|
||||
" Ã = [one.(ã) ã]\n",
|
||||
"\n",
|
||||
" # approximate gradient from this mini-batch\n",
|
||||
" ∇f̃ = (2/M)*(Ã'*(Ã*x - b̃))\n",
|
||||
"\n",
|
||||
" # Adam update\n",
|
||||
" m = β₁*m + (1-β₁)*∇f̃\n",
|
||||
" v = β₂*v + (1-β₂)*(∇f̃.^2)\n",
|
||||
" m̂ = m / (1 - β₁^t) # de-bias: normalize by total weights\n",
|
||||
" v̂ = v / (1 - β₂^t) # de-bias: normalize by total weights\n",
|
||||
" x = @. x - α * m̂ / (√v̂ + ϵ)\n",
|
||||
"\n",
|
||||
" push!(X, x)\n",
|
||||
" push!(ΔX, norm(x - x̂₀))\n",
|
||||
"end\n",
|
||||
"\n",
|
||||
"# objective function for contour plot:\n",
|
||||
"a = randn(1000)\n",
|
||||
"b = @. 1 + 2a + randn()\n",
|
||||
"A = [one.(a) a]\n",
|
||||
"f(x) = norm(b-A*x)^2 / length(b)\n",
|
||||
"f(x1,x2) = f([x1,x2])\n",
|
||||
"x1 = range(-2,3, length=100)\n",
|
||||
"x2 = range(0,4,length=100)\n",
|
||||
"fvals = f.(x1',x2)\n",
|
||||
"\n",
|
||||
"# plot progress\n",
|
||||
"\n",
|
||||
"figure(figsize=(10,4))\n",
|
||||
"\n",
|
||||
"subplot(1,2,1)\n",
|
||||
"plt.subplots_adjust(wspace=0.4)\n",
|
||||
"contour(x1', x2, fvals)\n",
|
||||
"plot(x̂₀[1], x̂₀[2], \"r*\")\n",
|
||||
"plot(getindex.(X,1), getindex.(X,2), \"ko-\")\n",
|
||||
"xlabel(L\"x_1\")\n",
|
||||
"ylabel(L\"x_2\")\n",
|
||||
"title(\"gradient-descent iterates x\")\n",
|
||||
"\n",
|
||||
"subplot(1,2,2)\n",
|
||||
"semilogy(ΔX, \"ko-\")\n",
|
||||
"ylim(1e-2, 2)\n",
|
||||
"xlabel(\"iteration\")\n",
|
||||
"ylabel(L\"\\Vert x - \\hat{x} \\Vert\")\n",
|
||||
"title(\"|error| in x\")"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "f238f2e7",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Solution:\n",
|
||||
"\n",
|
||||
"The only difference from Adam is the update of the \"variance\" estimate,\n",
|
||||
"which is now:\n",
|
||||
"$$\n",
|
||||
"v^{(k)} = \\beta_2 v^{(k-1)} - (1 - \\beta_2) \\underbrace{\\text{sign}.(v^{(k-1)} - g^{(k)}\\mbox{.^2})}_\\mbox{new} \\mbox{ .* }g^{(k)}\\mbox{.^2}\n",
|
||||
"$$\n",
|
||||
"where $g^{(k)} = \\tilde{\\nabla} f^{(k)}$ is the estimated (mini-batch) gradient, the [sign function](https://docs.julialang.org/en/v1/base/math/#Base.sign) is ±1 depending on the sign of each element, and all the multiplications are taken elementwise (denoted with a `.` in Julia notation). Note also the minus sign in front of the second term.\n",
|
||||
"\n",
|
||||
"(In Julia, we can simply make *everything* elementwise in a line of code by putting `@.` before the calculation.)\n",
|
||||
"\n",
|
||||
"This is implemented below. In order to compare it quantitatively to Adam, we want to do things:\n",
|
||||
"\n",
|
||||
"1. Plot the Adam and Yogi results on the same plot.\n",
|
||||
"2. To make the comparison more exact, use the *same \"random\" numbers* for both (i.e. the same mini-batches). (Otherwise we will get iterations that are qualitatively similar but quantitatively quite different.)\n",
|
||||
"\n",
|
||||
"One way to do the latter is to merge the loops — perform both the Adam and Yogi updates in the same loop, and save the iterates in different variables. An even more slick way to do this is to use the [`Random.seed!`](https://docs.julialang.org/en/v1/stdlib/Random/#Random.seed!) function to set a deterministic [seed](https://en.wikipedia.org/wiki/Random_seed) for the random-number generator. This allows us to generate the same \"random\" numbers both times. We can pick any seed we want, but I used `Random.seed!(314159)`:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 75,
|
||||
"id": "e024300e",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA04AAAGbCAYAAADzzMA0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd3hT1RvA8W+a7r0HHZRN2RuK7LL3EhQUEAQZioo4EARFEf2BCCpLGcpesqcIFMouU7aArJaWLrpXxv39EVqabiBp0nI+z9NHe3Nz7puQ5t73nnPeI5MkSUIQBEEQBEEQBEEokImhAxAEQRAEQRAEQTB2InESBEEQBEEQBEEogkicBEEQBEEQBEEQiiASJ0EQBEEQBEEQhCKIxEkQBEEQBEEQBKEIInESBEEQBEEQBEEogkicBEEQBEEQBEEQiiASJ0EQBEEQBEEQhCKIxEkQBEEQBEEQBKEIInESBKHE/f7778hksuwfU1NTfHx8eOuttwgPD9fpsWQyGe+++26R+wUHByOTyQgODtbp8V+UQqFgzpw51K5dGysrKxwdHWnevDnHjx9/4bYvXbqETCbDzMyMiIiIYj/v7t27yGQyfv/99xeOobQpzmsfN24cZmZmnDt3Ls9jmZmZ1K5dm8qVK5OSkqLT2Nq0aUObNm2eeb/U1FS+/PJLg3/2C4sj6zvj7t27JR6XIAhCFlNDByAIwstr+fLlVK9enbS0NI4cOcLMmTM5fPgwly5dwsbGxtDhGZxKpaJPnz4cPXqUTz75hObNm5OSksLZs2d1ctG9ZMkSAJRKJStWrODTTz994TYFmDVrFvv372fo0KGcPXsWc3Pz7Me+/PJLrl69ypEjR3T+GV+wYMFz7ZeamspXX30FUKzES18Ki6Nbt26cOHECLy8vA0QmCIKgIRInQRAMplatWjRq1AiAtm3bolKp+Prrr9m6dSuDBw/O9zmpqalYW1uXZJgG8/PPP7Nnzx6OHTtGs2bNsrd369bthdvOyMhg9erV1K1bl5iYGJYtWyYSJx2xtrbmjz/+oGXLlkybNo2ZM2cCEBoayv/+9z8mTpzIK6+8ovPj1qhRQ6f7vSiFQpHdo/yi3NzccHNz00FUgiAIz08M1RMEwWhkJQf37t0DYNiwYdja2nLp0iU6duyInZ0dQUFBAMTFxTF27Fi8vb0xNzenYsWKTJ48mYyMjHzbXrx4MVWrVsXCwoIaNWqwbt26YsV05swZevbsibOzM5aWltSvX58NGzZo7ZM1jOjgwYOMHDkSFxcX7O3tGTJkCCkpKURGRjJgwAAcHR3x8vJi4sSJKBSKIo89b948WrVqpZU06crWrVuJjY3l7bffZujQofz7778cPXo0z34PHz5kwIAB2NnZ4eDgwMCBA4mMjMyz35kzZ3jttdfw9/fHysoKf39/Xn/99ex/yyz6eq/Wr19Px44d8fLywsrKioCAAD777LM8PXNZn6lbt27RtWtXbG1t8fX15aOPPsrz2Snua89PYGAgH3/8MbNmzeLUqVNkZGQwbNgwAgICmD59OgBHjx4lKCgIOzs7rK2tad68Obt27crT1tGjRwkMDMTS0hJvb2+++OILlixZkmfo2vMM1bt79252QvLVV19lD58dNmxY9v43b95k0KBBuLu7Y2FhQUBAAPPnz9dqM2uo68qVK/noo4/w9vbGwsKCW7duER0dzdixY6lRowa2tra4u7vTrl07QkJCsp9fVBwFDdVbtmwZdevWxdLSEmdnZ/r06cO1a9e09nmWf/OFCxdSt25dbG1tsbOzo3r16nz++edFvqeCILwcROIkCILRuHXrFoDWneXMzEx69uxJu3bt2LZtG1999RXp6em0bduWFStWMGHCBHbt2sUbb7zB//73P/r27Zun3e3bt/PTTz8xffp0Nm3aRPny5Xn99dfZtGlTofEcOnSIV155hfj4eBYtWsS2bduoV68eAwcOzHeOy9tvv42DgwPr1q1jypQprFmzhpEjR9KtWzfq1q3Lpk2bGDp0KD/88AM///xzocd+8OABd+/epXbt2nz++ed4eHhgampKzZo1+eOPP/Ls36ZNG2QyWaFt5rR06VIsLCwYPHgww4cPRyaTsXTpUq190tLSaN++PX/99RczZ85k48aNeHp6MnDgwDzt3b17l2rVqjF37lz27dvH999/T0REBI0bNyYmJkav7xVoLu67du3K0qVL2bt3Lx988AEbNmygR48eefZVKBT07NmToKAgtm3bxvDhw/nxxx/5/vvvn+u1F+Srr76iZs2aDBs2jE8//ZSbN2+yYsUKLCwsOHz4MO3atSMhIYGlS5eydu1a7Ozs6NGjB+vXr89u459//qFDhw6kpqbyxx9/sGjRIs6dO8eMGTOKHUdhvLy82Lt3LwAjRozgxIkTnDhxgi+++AKAq1ev0rhxYy5fvswPP/zAzp076datG+PHj88eVpfTpEmTuH//PosWLWLHjh24u7sTFxcHwLRp09i1axfLly+nYsWKtGnTJns+U1Fx5GfmzJmMGDGCmjVrsnnzZubNm8c///xDYGAgN2/e1Nq3OP/m69atY+zYsbRu3ZotW7awdetWPvzwQ53PRRMEoRSTBEEQStjy5cslQDp58qSkUCikpKQkaefOnZKbm5tkZ2cnRUZGSpIkSUOHDpUAadmyZVrPX7RokQRIGzZs0Nr+/fffS4D0119/ZW8DJCsrq+w2JUmSlEqlVL16daly5crZ2w4dOiQB0qFDh7K3Va9eXapfv76kUCi0jtO9e3fJy8tLUqlUWq/nvffe09qvd+/eEiDNmTNHa3u9evWkBg0aFPoenThxQgIke3t7qUaNGtKGDRukffv2Sf3795cA6ddff9Xav127dpJcLi+0zSx3796VTExMpNdeey17W+vWrSUbGxspMTExe9vChQslQNq2bZvW80eOHCkB0vLlyws8hlKplJKTkyUbGxtp3rx52dv18V7lplarJYVCIR0+fFgCpIsXL2Y/lvWZyv3Z6dq1q1StWrXs31/kted04cIFydzcXAKkr7/+Ont7s2bNJHd3dykpKSl7m1KplGrVqiX5+PhIarVakiRJevXVVyUbGxspOjo6ez+VSiXVqFFDAqQ7d+5kb2/durXUunXrImPKvV90dLQESNOmTcuzb6dOnSQfHx8pISFBa/u7774rWVpaSnFxcZIkPf37adWqVZHHVyqVkkKhkIKCgqQ+ffoUK46sz03W6338+LFkZWUlde3aVWu/+/fvSxYWFtKgQYOytxX33/zdd9+VHB0di4xfEISXl+hxEgTBYJo1a4aZmRl2dnZ0794dT09P9uzZg4eHh9Z+/fr10/r94MGD2NjY0L9/f63tWcN6Dhw4oLU9KChIq025XM7AgQO5desWYWFh+cZ269Ytrl+/nj3XSqlUZv907dqViIgIbty4ofWc7t27a/0eEBAA5J2TFBAQkGcIW25qtRqA9PR0du/ezauvvkrHjh3ZsGEDDRo0yB7uleXAgQMolcpC28yyfPly1Go1w4cPz942fPhwUlJStHo7Dh06hJ2dHT179tR6/qBBg/K0mZyczKeffkrlypUxNTXF1NQUW1tbUlJS8gydAt2+VwD//fcfgwYNwtPTE7lcjpmZGa1btwbIc3yZTJanJ6pOnTpax3mW116YunXr0rdvX6ysrJg0aRIAKSkpnDp1iv79+2Nra5u9r1wu58033yQsLCz7s5XVM+Xq6pq9n4mJCQMGDHimOJ5Heno6Bw4coE+fPlhbW+f5G0hPT+fkyZNaz8n9t5pl0aJFNGjQAEtLS0xNTTEzM+PAgQP5fjaK48SJE6SlpWkNKQTw9fWlXbt2eb4DivNv3qRJE+Lj43n99dfZtm1bvj2lgiC83ETiJAiCwaxYsYLQ0FDOnz/Pw4cP+eeff/JMmre2tsbe3l5rW2xsLJ6ennmGprm7u2NqakpsbKzWdk9PzzzHztqWe98sjx49AmDixImYmZlp/YwdOxYgz4WVs7Oz1u9Z1dTy256enp7vcbO4uLgAUL16dcqXL5+9XSaT0alTJ8LCwoiKiiq0jfyo1Wp+//13ypUrR8OGDYmPjyc+Pp727dtjY2OjNVwvNjY2TxIL+b+fgwYN4pdffuHtt99m3759nD59mtDQUNzc3EhLS8uzvy7fq+TkZFq2bMmpU6f45ptvCA4OJjQ0lM2bNwPkOb61tTWWlpZa2ywsLLSO8yyvvSgWFhaYmJggl8sBePz4MZIk5Vshrly5ctnHLyyO/LbpWmxsLEqlkp9//jnP30DXrl2BvH8D+b2mOXPmMGbMGJo2bcqff/7JyZMnCQ0NpXPnzvl+NoobW0HHK1euXJ6/6+L8m7/55pssW7aMe/fu0a9fP9zd3WnatCn79+9/rhgFQSh7RFU9QRAMJiAgILuqXkHym7fj4uLCqVOnkCRJ6/GoqCiUSqXW3Xkg3wn9WduyEpTcstqYNGlSvvOmAKpVq1Zo7C+iUqVKBVYPlCQJ0PQ8PKu///47+y57fq/95MmTXL16lRo1auDi4sLp06fz7JP7/UxISGDnzp1MmzaNzz77LHt7RkZG9vwWfTp48CAPHz4kODg4u5cJID4+/rnbLO5rfx5OTk6YmJjku3bWw4cPgaefPxcXl+wkXtdxFMXJySm7F2zcuHH57lOhQgWt3/P7e121ahVt2rRh4cKFWtuTkpKeO7asz25B72Hu74Dieuutt3jrrbdISUnhyJEjTJs2je7du/Pvv/9q3cAQBOHlJHqcBEEodYKCgkhOTmbr1q1a21esWJH9eE4HDhzQuvhUqVSsX7+eSpUq4ePjk+8xqlWrRpUqVbh48SKNGjXK98fOzk63LywHU1NTevXqxbVr17QqiUmSxN69e6lUqdJzXRwuXboUExMTtm7dyqFDh7R+Vq5cCWgqlYGmRHxSUhLbt2/XamPNmjVav8tkMiRJwsLCQmv7kiVLUKlUzxzjs8q6WM99/MWLine truncated
|
||||
"text/plain": [
|
||||
"Figure(PyObject <Figure size 1000x400 with 2 Axes>)"
|
||||
]
|
||||
},
|
||||
"metadata": {},
|
||||
"output_type": "display_data"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"PyObject Text(0.5, 1.0, '|error| in x')"
|
||||
]
|
||||
},
|
||||
"execution_count": 75,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"Random.seed!(314159)\n",
|
||||
"\n",
|
||||
"x = [0.0,0.0]\n",
|
||||
"X_yogi = [x]\n",
|
||||
"ΔX_yogi = [norm(x - x̂₀)]\n",
|
||||
"m = zero(x) # initial 1st-moment (mean gradient)\n",
|
||||
"v = zero(x) # initial 2nd-moment (mean gradient.^2)\n",
|
||||
"for t = 1:iterations\n",
|
||||
" # generate a mini-batch of size M:\n",
|
||||
" ã = randn(M)\n",
|
||||
" b̃ = @. 1 + 2ã + randn()\n",
|
||||
" Ã = [one.(ã) ã]\n",
|
||||
"\n",
|
||||
" # approximate gradient from this mini-batch\n",
|
||||
" ∇f̃ = (2/M)*(Ã'*(Ã*x - b̃))\n",
|
||||
"\n",
|
||||
" # Adam update\n",
|
||||
" m = β₁*m + (1-β₁)*∇f̃\n",
|
||||
" v = @. β₂*v - (1-β₂) * sign(v - ∇f̃^2) * (∇f̃^2)\n",
|
||||
" m̂ = m / (1 - β₁^t) # de-bias: normalize by total weights\n",
|
||||
" v̂ = v / (1 - β₂^t) # de-bias: normalize by total weights\n",
|
||||
" x = @. x - α * m̂ / (√v̂ + ϵ)\n",
|
||||
"\n",
|
||||
" push!(X_yogi, x)\n",
|
||||
" push!(ΔX_yogi, norm(x - x̂₀))\n",
|
||||
"end\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"# plot progress\n",
|
||||
"\n",
|
||||
"figure(figsize=(10,4))\n",
|
||||
"plt.suptitle(\"Problem 6: Adam and Yogi iterations\")\n",
|
||||
"\n",
|
||||
"subplot(1,2,1)\n",
|
||||
"plt.subplots_adjust(wspace=0.4)\n",
|
||||
"contour(x1', x2, fvals)\n",
|
||||
"plot(x̂₀[1], x̂₀[2], \"r*\")\n",
|
||||
"plot(getindex.(X,1), getindex.(X,2), \"ko-\")\n",
|
||||
"plot(getindex.(X_yogi,1), getindex.(X_yogi,2), \"b*-\")\n",
|
||||
"legend([\"Adam\", \"Yogi\"])\n",
|
||||
"xlabel(L\"x_1\")\n",
|
||||
"ylabel(L\"x_2\")\n",
|
||||
"title(\"gradient-descent iterates x\")\n",
|
||||
"\n",
|
||||
"subplot(1,2,2)\n",
|
||||
"semilogy(ΔX, \"ko-\")\n",
|
||||
"semilogy(ΔX_yogi, \"b*-\")\n",
|
||||
"ylim(1e-2, 2)\n",
|
||||
"legend([\"Adam\", \"Yogi\"])\n",
|
||||
"xlabel(\"iteration\")\n",
|
||||
"ylabel(L\"\\Vert x - \\hat{x} \\Vert\")\n",
|
||||
"title(\"|error| in x\")"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "ee14d9bc",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"When we plot them like this, we can see that there is virtually no difference! That's **because this example problem is so well-conditioned** (the $f$ contours are nearly circles), so the momentum terms make almost no difference — momentum-based algorithms should be essentially equivalent to steepest descent for well-conditioned problems.\n",
|
||||
"\n",
|
||||
"**Optional**:\n",
|
||||
"\n",
|
||||
"To actually see any difference between the iterates, let's plot the difference $\\Vert x_\\text{Adam} - x_\\text{Yogi} \\Vert$:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 77,
|
||||
"id": "54da14e6",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAHFCAYAAADi7703AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB4GUlEQVR4nO3deVxUVf8H8M8wrIIgm4BAYGqamZpLpUaillla5Dhu2aO2l5SiZj09pEhqmppJifaYpS0uCIzaXlpglFZm2mPZooWmhBsYuCuX8/vD350YmBnuDLPcgc/79ZrXC+6cuffMncvcL2f5Ho0QQoCIiIiIzPJydwWIiIiI1IzBEhEREZEVDJaIiIiIrGCwRERERGQFgyUiIiIiKxgsEREREVnBYImIiIjICgZLRERERFYwWCIiIiKygsGSDVatWgWNRmN8eHt7Iy4uDvfddx9KSkoceiyNRoPHH3+83nKFhYXQaDQoLCx06PEb6tKlS1i0aBGuvfZaBAQEoEWLFujduze2bdvm0nqcPXsWM2fOVN35qUm+rg4cOGDclpycjOTkZJNyBw4cwODBgxEWFgaNRoO0tDQAwK5du9C3b1+EhIRAo9Fg8eLFLqu72sh/D3l5eQ7b5969ezFz5kyTz6cx69atGzQaDRYuXGjT68xds01Ffe/966+/hre3N6ZOnWr2+eeffx4ajQYff/yxQ+tl7rtFabk1a9ao4rtELfXwdncFPNHKlSvRoUMHnDt3Dl988QXmzp2LrVu3Ys+ePQgMDHR39dxOkiQMHToUX375JZ566in07t0bZ86cwc6dO3HmzBmX1uXs2bPIzMwEAI/6Il+6dGmdbZMnT8Y333yDN954A9HR0YiJiQEA3H///Thz5gzWrVuH0NBQJCYmuri2jdvevXuRmZmJ5OTkRn9ud+/ejV27dgEAXn/9dTz55JNurlHjcOONN+Lpp5/GvHnzMHToUNx0003G53788UdkZmbikUcewaBBgxx63MGDB2P79u3G7wpbyq1ZswY//vij8Z8yd1FLPRgs2aFTp07o0aMHAKBfv36QJAmzZs3Cxo0bMWbMGLOvOXv2LJo1a+bKarrNK6+8go8++ghfffUVbrzxRuP2wYMHu7FWyqjlc+rYsWOdbT/++COuv/563H333XW2P/TQQ7j99tsdcmxJklBVVQU/Pz+H7I88x4oVKwBc/lv94IMPsG3bNvTu3dvNtWocMjIy8MEHH2D8+PH43//+h2bNmqGqqgrjx49HXFyczS15SkRGRiIyMtJh5Rzh3LlzCAgIcMmxHIndcA4gBwQHDx4EAIwfPx5BQUHYs2cPBg4ciObNm2PAgAEAgPLyckyYMAGxsbHw9fXFlVdeifT0dFy4cMHsvv/73//iqquugp+fHzp27Ih169YpqtN3332Hu+66C2FhYfD398d1112H9evXm5SRm14///xzPPTQQwgPD0dwcDDGjh2LM2fO4MiRIxgxYgRatGiBmJgYPPnkk7h06VK9x87KysLNN99sEijZQgiBpUuXomvXrggICEBoaCj0ej3++OMPY5l169ZBo9FgyZIlJq/NyMiAVqvF5s2bceDAAeMXQGZmprH7dPz48QCAmTNnQqPR4Pvvv4der0doaCjatGljPH+jRo1CYmIiAgICkJiYiNGjRxs/Y3t9/fXX6NOnD/z9/dGqVSs888wzZs9pzWZ9uWtp//79+Oijj4zvQ/78qqqqsGzZMuN22ZEjR/DII48gLi4Ovr6+aN26NTIzM1FVVWUsc+DAAWg0GsyfPx+zZ89G69at4efnh4KCAuN5UHodFRQU4LHHHkNERATCw8Oh0+nw119/1Xlva9asQa9evRAUFISgoCB07doVr7/+ukmZLVu2YMCAAQgODkazZs3Qp08ffPbZZ4rP8/nz5zFlyhRER0cjICAAffv2NbaY1FTf+1u1ahWGDx8O4PI/RjXPfXZ2Nry8vHDs2DFj+RdffBEajQapqanGbdXV1QgNDTXpgrl48SJmz56NDh06wM/PD5GRkbjvvvtw/PjxOnXMyclBr169EBgYiKCgINx222113ov8nbN//37ccccdCAoKQnx8PKZOnWrxu8XcOVuzZg26d++Ol156CQDwxhtv1CknhMD8+fORkJAAf39/dOvWDR999JHZ/U2dOhVdu3ZFSEgIwsLC0KtXL2zatKlOWXnYwcqVK9G+fXsEBASgR48e+PrrryGEwIIFC9C6dWsEBQWhf//+2L9/f73vZ//+/bjvvvvQrl07NGvWDLGxsbjzzjuxZ88ek3Ly39fatWuRnp6OVq1aITg4GLfccgt+/fVXu967Ob6+vnjrrbdw6NAhPP300wCAuXPnYteuXVi1ahWCgoIU3x/+/vtvPPDAAwgLC0NQUBAGDx6MP/74AxqNBjNnzjSWs7cbLjk5GR988AEOHjxoMvREpvT6TUxMxJAhQ2AwGHDdddfB39/f2NKfnZ2Nm2++GS1btkRgYCCuvfZazJ8/3+T70FH1+Pzzz5GcnIzw8HAEBATgiiuuwLBhw3D27FnrH1pNghRbuXKlACB27Nhhsj0rK0sAEMuXLxdCCDFu3Djh4+MjEhMTxdy5c8Vnn30mPvnkE3Hu3DnRuXNnERgYKBYuXCg+/fRTMX36dOHt7S3uuOMOk30CEPHx8aJjx45i7dq14t133xWDBg0SAERubq6xXEFBgQAgCgoKjNs+//xz4evrK5KSkkROTo74+OOPxfjx4wUAsXLlyjrvp3Xr1mLq1Kni008/FS+88ILQarVi9OjRolu3bmL27Nli8+bN4umnnxYAxIsvvmj1HP35558CgHjiiSfEM888I1q2bCm0Wq3o2LGjWLVqlaLz/NBDDwkfHx8xdepU8fHHH4s1a9aIDh06iKioKHHkyBFjuUcffVT4+voaP4/PPvtMeHl5iWeffVYIIcT58+fFxx9/LACIBx54QGzfvl1s375d7N+/XwghREZGhgAgEhISxNNPPy02b94sNm7cKIQQIjc3V8yYMUNs2LBBbN26Vaxbt0707dtXREZGiuPHj5vUNyEhQSQkJNT7vn766SfRrFkz42e6adMmcdttt4krrrhCABDFxcXGsn379hV9+/YVQghRUVEhtm/fLqKjo0WfPn2M7+PIkSNi+/btAoDQ6/XG7UIIUVpaKuLj40VCQoL473//K7Zs2SJmzZol/Pz8xPjx443HKS4uFgBEbGys6Nevn8jLyxOffvqpKC4utvk6uvLKK8UTTzwhPvnkE7FixQoRGhoq+vXrZ3IOpk+fLgAInU4ncnNzxaeffioWLVokpk+fbizz9ttvC41GI+6++25hMBjEe++9J4YMGSK0Wq3YsmWL1XMs/z3Ex8eLlJQU8d5774l33nlHtG3bVgQHB4vff//dWFbJ+zt27Jh4/vnnBQCRnZ1tPMfHjh0Tv/zyiwAg1qxZY9znoEGDREBAgGjXrp1x2zfffCMAiA8//FAIIYQkSWLQoEEiMDBQZGZmis2bN4sVK1aI2NhY0bFjR3H27Fnja+fMmSM0Go24//77xfvvvy8MBoPo1auXCAwMFD/99JOx3Lhx44Svr6+4+uqrxcKFC8WWLVvEjBkzhEajEZmZmVbPmWz16tXG9ymEEDfddJMICgoSp06dMikn/9088MAD4qOPPhLLly8XsbGxIjo62njNCiHE33//LcaPHy/efvtt8fnnn4uPP/5YPPnkk8LLy0u8+eabJvuU/w579+4tDAaD2LBhg7jqqqtEWFiYmDx5skhJSRHvv/++WL16tYiKihKdO3cW1dXVVt/P1q1bxdSpU0VeXp7YunWr2LBhg7j77rtFQECA+OWXX4zl5GsmMTFRjBkzRnzwwQdi7dq14oorrhDt2rUTVVVVNr93a+bOnSs0Go1YvHix8PHxEU8++aQQQii+P0iSJG666Sbh7+8v5s2bJz799FORmZkp2rVrJwCIjIwMY1n5b7Pmd4s5tcv99NNPok+fPiI6Otp4zcvfLbZcvwkJCSImJkZceeWV4o033hAFBQXi22+/FUIIMXnyZLFs2TLx8ccfi88//1y89NJLIiIiQtx3333G1zuiHsXFxcLf31/ceuutYuPGjaKwsFCsXr1a/Otf/xInT55U9JkJIQSDJRvIF9TXX38tLl26JE6dOiXef/99ERkZKZo3b268kY8bN04AEG+88YbJ61999VUBQKxfv95k+wsvvCAAiE8//dS4DYAICAgwCQ6qqqpEhw4dRNu2bY3bzAVLHTp0ENddd524dOmSyXGGDBkiYmJihCRJJu/niSeeMCl39913CwBi0aJFJtu7du0qunXrZvUcyTfv4OBg0bFjR7F+/XrxySefCL1ebxJQ1vf62kHZoUOHREBAgHjqqaeM286fPy+uu+460bp1a7F3714RFRUl+vbta/Lldvz48TpfIDL5i2/GjBlW6yTE5XN/+vRpERgYKLKyskyea9OmjWjTpk29+xg5cqTFz9RasCRLSEgQgwcPrrNfACI1NdVk2yOPPCKCgoLEwYMHTbYvXLhQADDeaOVgqU2bNuLixYsmZW29jiZMmGBSbv78+QKAKC0tFUII8ccffwitVivGjBlj6RSJM2fOiLCwMHHnnXeabJckSXTp0kVcf/31Fl8rxD9/D926dTO5mR44cED4+PiIBx980Ob3l5ubW+dvTBYXFyfuv/9+IYQQFy5cEIGBgcZ/LORzP2fOHOHj4yNOnz4thBBi7dq1AoDIz8832deOHTsEALF06VIhxOV/PLy9vev8fZ46dUpER0eLESNGGLfJ3zm1v1vuuOMO0b59e6vnTNa/f3/h7+9vvIHIn+vrr79uLHPy5Enh7+8vhg4davLar776SgCwGjBUVVWJS5cuiQceeEBcd911Js8BENHR0cZzJIQQGzduFABE165dTT7LxYsXCwDif//7n6L3VfP4Fy9eFO3atROTJ082bpevmdr/sK5fv14AMN6cG/LLine truncated
|
||||
"text/plain": [
|
||||
"Figure(PyObject <Figure size 640x480 with 1 Axes>)"
|
||||
]
|
||||
},
|
||||
"metadata": {},
|
||||
"output_type": "display_data"
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"PyObject Text(0.5, 1.0, 'Problem 6 extra: difference between Adam and Yogi iterates')"
|
||||
]
|
||||
},
|
||||
"execution_count": 77,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"semilogy(norm.(X .- X_yogi), \"ko-\")\n",
|
||||
"xlabel(\"iteration\")\n",
|
||||
"ylabel(L\"\\Vert x_\\mathrm{Adam} - x_\\mathrm{Yogi} \\Vert\")\n",
|
||||
"title(\"Problem 6 extra: difference between Adam and Yogi iterates\")"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "bea01a59",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The differences are bigger than roundoff errors, but they are still quite small!"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "9dd42d55",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": []
|
||||
}
|
||||
],
|
||||
"metadata": {
|
||||
"@webio": {
|
||||
"lastCommId": null,
|
||||
"lastKernelId": null
|
||||
},
|
||||
"kernelspec": {
|
||||
"display_name": "Julia 1.8.3",
|
||||
"language": "julia",
|
||||
"name": "julia-1.8"
|
||||
},
|
||||
"language_info": {
|
||||
"file_extension": ".jl",
|
||||
"mimetype": "application/julia",
|
||||
"name": "julia",
|
||||
"version": "1.8.5"
|
||||
}
|
||||
},
|
||||
"nbformat": 4,
|
||||
"nbformat_minor": 5
|
||||
}
|
||||
Reference in new issue
Block a user