rest of pset 4

This commit is contained in:
Steven G. Johnson committed 2023-03-26 15:24:09 -04:00
1 parent d1edba8ddc
commit a96fdc7659
1 file changed
+138 -17
+138 -17
View File
@@ -66,21 +66,10 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 4, "execution_count": null,
"id": "ac419822", "id": "ac419822",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [],
{
"data": {
"text/plain": [
"-25.41770831865196"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [ "source": [
"using LinearAlgebra\n", "using LinearAlgebra\n",
"\n", "\n",
@@ -136,22 +125,154 @@
}, },
{ {
"cell_type": "markdown", "cell_type": "markdown",
"id": "9d19b287", "id": "6e9798d4",
"metadata": {}, "metadata": {},
"source": [ "source": [
"## Problems 5, 6, … coming soon" "## Problem 5\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.)"
]
},
{
"cell_type": "markdown",
"id": "c7f174a6",
"metadata": {},
"source": [
"## Problem 6"
]
},
{
"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": 1,
"id": "2181180f",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA04AAAGICAYAAAB2jrABAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd1QT2RcH8G8IEHqv0kFRsSt2Uey9Ifa6tnXtsuqu4tpdd22L7lrWvnZFUXftFcWOir2jiCC9d0gyvz9C+ElLJjCTSeB9zuEciZN5N5Qwd9599/EoiqJAEARBEARBEARBlEuD6wAIgiAIgiAIgiBUHUmcCIIgCIIgCIIg5CCJE0EQBEEQBEEQhBwkcSIIgiAIgiAIgpCDJE4EQRAEQRAEQRBykMSJIAiCIAiCIAhCDpI4EQRBEARBEARByEESJ4IgCIIgCIIgCDlI4kQQBEEQBEEQBCEHSZyISnN2dsa4ceOKPg8ODgaPx0NwcDCr427ZsgV79+6t9HkiIiLA4/EYOZcqOXToEAICAmgf7+3tDW9v76LPs7OzsXTpUta/j/KoShwEQRDqyNnZGUuXLuU6DIXRjXvp0qXg8XjsB0QQADS5DoCoepo2bYq7d+/Cw8OD1XG2bNkCCwuLYkkb8X+HDh3CixcvMHv2bFrHb9mypdjn2dnZWLZsGQAUS6iUTVXiIAiCIFTPxIkT0aNHD67DIKoJkjhVUxRFITc3F7q6uoyf28jICK1atWL8vAS72E50pQoKCsDj8aCpSd5+CIIg1E12djb09PRKPc7EdUVOTg50dHQUmkGyt7eHvb19hcckCEWQUj01d/r0aTRs2BACgQCurq7YuHFjmdPWPB4P06dPx7Zt21C3bl0IBAL8888/AIBly5ahZcuWMDMzg5GREZo2bYpdu3aBoqhi5ygoKMD8+fNhY2MDPT09tGvXDg8ePCgVU3mleg8fPkS/fv1gZmYGHR0dNGnSBMeOHSt2zN69e8Hj8XD9+nX88MMPsLCwgLm5OXx8fPD169ei45ydnfHy5UvcuHEDPB4PPB4Pzs7Ocr9eX79+xZAhQ2BoaAhjY2MMHToUsbGxZR5LJ97s7GzMnTsXLi4u0NHRgZmZGTw9PXH48OFix92/fx99+/aFubk5dHR04ObmVmom6P379xgxYgSsrKwgEAhQt25dbN68ucyv7eHDh+Hv748aNWrAyMgIXbp0wdu3b4uO8/b2xtmzZ/H58+eir4+8P0TflupFRETA0tISgOTnQ/r8b2f3FIl3//79+PHHH2FnZweBQIAPHz4gISEBU6dOhYeHBwwMDGBlZYVOnTohJCSk6PlMxSEWi7Fy5UrUrl0burq6MDExQcOGDbFx40aZX5MpU6ZAR0cHjx49Knauzp07w9raGjExMTKfTxAEoaqOHj2K1q1bQ19fHwYGBujevTvCwsKKHTNu3DgYGBjg+fPn6NatGwwNDdG5c2cAsq8rbt26hc6dO8PQ0BB6enpo06YNzp49W+zc0r/3ly5dwvjx42FpaQk9PT3k5eUp9DrKuuZxdnZGnz59cOHCBTRt2hS6urqoU6cOdu/eLfd8v/32GzQ0NPDff/+V+lro6enh+fPnCsVHVC3klq8au3DhAnx8fNC+fXscPXoUQqEQ69atQ1xcXJnHnzp1CiEhIVi8eDFsbGxgZWUFQHJx+v3338PR0REAcO/ePcyYMQPR0dFYvHhx0fMnTZqEffv2Ye7cuejatStevHgBHx8fZGRkyI31+vXr6NGjB1q2bIlt27bB2NgYR44cwdChQ5GdnV2q3G7ixIno3bs3Dh06hC9fvmDevHkYNWoUrl27BgA4efIkfH19YWxsXFRiJhAIZMaQk5ODLl264OvXr1i9ejXc3d1x9uxZDB06tMLx+vn5Yf/+/Vi5ciWaNGmCrKwsvHjxAklJSUXnunjxIvr27Yu6detiw4YNcHR0REREBC5dulR0zKtXr9CmTRs4Ojpi/fr1sLGxwcWLFzFz5kwkJiZiyZIlxeJbuHAh2rZti507dyI9PR0//fQT+vbti9evX4PP52PLli2YPHkywsPDcfLkSbnfn5JsbW1x4cIF9OjRAxMmTMDEiRMBoCiJUTTeBQsWoHXr1ti2bRs0NDRgZWWFhIQEAMCSJUtgY2ODzMxMnDx5Et7e3rh69Sq8vb0Zi2PNmjVYunQpFi1ahPbt26OgoABv3rxBamqqzK9DQEAA7t+/jyFDhuDRo0cwMTHBsmXLEBwcjAsXLsDW1lbhry1BEATXfv31VyxatAjfffcdFi1ahPz8fKxduxZeXl548OBBsQqE/Px89OvXD99//z1+/vlnCIXCov8r67rixo0b6Nq1Kxo2bIhdu3ZBIBBgy5Yt6Nu3Lw4fPlzqb+748ePRu3dv7N+/H1lZWdDS0mLkNT59+hQ//vgjfv75Z1hbW2Pnzp2YMGECatasifbt25f7vJ9++gkhISEYO3YswsLC4OTkhD179uCff/7Bzp070aBBA0biI9QURait5s2bUw4ODlReXl7RYxkZGZS5uTlV8lsLgDI2NqaSk5NlnlMkElEFBQXU8uXLKXNzc0osFlMURVGvX7+mAFBz5swpdvzBgwcpANTYsWOLHrt+/ToFgLp+/XrRY3Xq1KGaNGlCFRQUFHt+nz59KFtbW0okElEURVF79uyhAFBTp04tdtyaNWsoAFRMTEzRY/Xq1aM6dOgg8/V8a+vWrRQA6vTp08UenzRpEgWA2rNnj8Lx1q9fnxowYIDMcd3c3Cg3NzcqJyen3GO6d+9O2dvbU2lpacUenz59OqWjo1P0fZN+bXv16lXsuGPHjlEAqLt37xY91rt3b8rJyUlmbN/q0KFDsa9nQkICBYBasmRJpeNt37693PGFQiFVUFBAde7cmRo4cCCjcfTp04dq3Lix3BjK8v79e8rIyIgaMGAAdeXKFUpDQ4NatGhRhc5FEAShbE5OTsXePyMjIylNTU1qxowZxY7LyMigbGxsqCFDhhQ9NnbsWAoAtXv37lLnLe+6olWrVpSVlRWVkZFR9JhQKKTq169P2dvbF11XSP/ejxkzhlbc5VmyZEmpax4nJydKR0eH+vz5c9FjOTk5lJmZGfX999/LPWdiYiJlb29PtWjRgnr8+DGlp6dHjRo1Su7ziKqPlOqpqaysLDx8+BADBgyAtrZ20eMGBgbo27dvmc/p1KkTTE1NSz1+7do1dOnSBcbGxuDz+dDS0sLixYuRlJSE+Ph4AJIZGAAYOXJksecOGTJE7lqVDx8+4M2bN0XPFQqFRR+9evVCTExMsTIzAOjXr1+xzxs2bAgA+Pz5s8yxAEAkEhUbQywWF70GQ0PDUuceMWJEheNt0aIFzp8/j59//hnBwcHIyckpdq53794hPDwcEyZMgI6OTpnx5ubm4urVqxg4cCD09PRKjZebm4t79+4x9vWprIrEO2jQoDLPtW3bNjRt2hQ6OjrQ1NSElpYWrl69itevXzMaR4sWLfD06VNMnToVFy9eRHp6Ou3XW7NmTezYsQOnTp1Cnz594OXlpZYdqgiCIABJFYRQKMSYMWOKvW/q6OigQ4cOZXYwLe89vOR1RVZWFu7fvw9fX18YGBgUPc7n8zF69GhERUWV+ntf3rkrq3HjxkWVNACgo6MDd3d3Wn8nzc3NcfToUTx+/LioqmHbtm2sxEmoF5I4qamUlBRQFAVra+tS/1fWYwDKLCt68OABunXrBgDYsWMHbt++jdDQUPj7+wNAUSIgLT2zsbEp9nxNTU2Ym5vLjFVaOjh37lxoaWkV+5g6dSoAIDExsdhzSp5TWoZXMjEpi5ubW7Exli9fXvQayvralHxNisS7adMm/PTTTzh16hQ6duwIMzMzDBgwAO/fvweAonI0WQtXk5KSIBQK8eeff5Yar1evXox/fSqrIvGW9bO3YcMG/PDDD2jZsiVOnDiBe/fuITQ0FD169KD1OhSJY8GCBVi3bh3u3buHnj17wtzcHJ07d8bDhw9pvebevXvD2toaubm58PPzA5/Pp/U8giAIVSP9G9e8efNS751Hjx4t9f6tp6cHIyOjMs9V8r1dem1S1nt+jRo1AKBYKXtZ52BKWdcmAoGA9t/Jli1bol69esjNzcUPP/wAfX19pkMk1BBZ46SmTE1NwePxylzPVF6zg7KaAxw5cgRaWlo4c+ZMsRmRU6dOFTtO+gYUGxsLOzu7oseFQmGpN8GSLCwsAEguXn18fMo8pnbt2jLPoYj//vuv2OJS6Zu1ubl5mc0sSn69FIlXX18fy5Ytw7JlyxAXF1c0+9S3b1+8efOmaC1OVFRUufGampoW3Y2bNm1amce4uLiU+3xlq0i8Zf3sHThwAN7e3ti6dWuxx+msmVM0Dk1NTfj5+cHPzw+pqam4cuUKFi5ciO7du+PLly9ldoj61pQpU5CRkYF69eph5syZ8PLyKnP2liAIQtVJ/8YdP34cTk5Oco+X1Vio5P+ZmppCQ0OjzMY50gZP0vHpnJ9LS5YswfPnz9GsWTMsXrwYffr0gaurK9dhERwjiZOa0tfXh6enJ06dOoV169YVletlZmbizJkztM8jbQv97R30nJwc7N+/v9hx0m5rBw8eRLNmzYoeP3bsWLGFomWpXbs2atWqhadPn+LXX3+lHZs85d05Km/hZseOHXHs2DH8+++/xUrdDh06xEi81tbWGDduHJ4+fYqAgABkZ2fD3d0dbm5u2L17N/z8/MpsYKGnp4eOHTsiLCwMDRs2LFZ6WRmK3Fkr7/lA6VkspuLl8Xilvh7Pnj3D3bt34eDgwFocJiYm8PX1RXR0NGbPno2IiAiZrdh37tyJAwcOYPfu3ejQoQOaNm2K7777rtTNBYIgCHXQvXt3aGpqIjw8nPEyOX19fbRs2RJBQUFYt25dUWtysViMAwcOwN7eHu7u7oyOyYbLly9j9erVWLRoEWbPno3GjRtj6NChuH37NmN/own1RBInNbZ8+XLLine 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": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using LinearAlgebra, PyPlot\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",
"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": "code", "cell_type": "code",
"execution_count": null, "execution_count": null,
"id": "6ac54afe", "id": "167a6e22",
"metadata": {}, "metadata": {},
"outputs": [], "outputs": [],
"source": [] "source": []
} }
], ],
"metadata": { "metadata": {
"@webio": {
"lastCommId": null,
"lastKernelId": null
},
"kernelspec": { "kernelspec": {
"display_name": "Julia 1.8.3", "display_name": "Julia 1.8.3",
"language": "julia", "language": "julia",
@@ -161,7 +282,7 @@
"file_extension": ".jl", "file_extension": ".jl",
"mimetype": "application/julia", "mimetype": "application/julia",
"name": "julia", "name": "julia",
"version": "1.8.3" "version": "1.8.5"
} }
}, },
"nbformat": 4, "nbformat": 4,