{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Intuition for Second Order Partial Derivatives and the Hessian Matrix\n", "> But what is a Hessian, really? \n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Julia Version 1.5.4\n", "Commit 69fcb5745b (2021-03-11 19:13 UTC)\n", "Platform Info:\n", " OS: macOS (x86_64-apple-darwin18.7.0)\n", " CPU: Apple M1 Max\n", " WORD_SIZE: 64\n", " LIBM: libopenlibm\n", " LLVM: libLLVM-9.0.1 (ORCJIT, westmere)\n", "Active project: /Users/alex/Documents/back_of_my_envelope/_notebooks/Project.toml\n", "Sys.BINDIR: /Applications/Julia-1.5.app/Contents/Resources/julia/bin\n", "Base.julia_cmd(): `\u001b[4m/Applications/Julia-1.5.app/Contents/Resources/julia/bin/julia\u001b[24m \u001b[4m-Cnative\u001b[24m \u001b[4m-J/Applications/Julia-1.5.app/Contents/Resources/julia/lib/julia/sys.dylib\u001b[24m \u001b[4m-g1\u001b[24m \u001b[4m--color=yes\u001b[24m`\n", "\u001b[32m\u001b[1mStatus\u001b[22m\u001b[39m `~/Documents/back_of_my_envelope/_notebooks/Project.toml`\n", " \u001b[90m [a03496cd] \u001b[39m\u001b[37mPlotlyBase v0.7.3\u001b[39m\n" ] } ], "source": [ "#hide\n", "import Pkg, InteractiveUtils\n", "\n", "InteractiveUtils.versioninfo()\n", "\n", "println(\"Active project: \", Base.active_project())\n", "println(\"Sys.BINDIR: \", Sys.BINDIR)\n", "println(\"Base.julia_cmd(): \", Base.julia_cmd())\n", "\n", "Pkg.status() # must show PlotlyBase in THIS env" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "render_plotly_fastpages (generic function with 1 method)" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#collapse\n", "# imports\n", "using PlotlyBase, IJulia\n", "\n", "# Define once at the top of the notebook\n", "function render_plotly_fastpages(p::PlotlyBase.Plot; include_plotlyjs::Union{Bool,String}=\"cdn\", full_html::Bool=false)\n", " html = sprint(io -> PlotlyBase.to_html(io, p; include_plotlyjs=include_plotlyjs, full_html=full_html))\n", "\n", " # Replace the literal plotly div with a script that inserts it at runtime\n", " if (m = match(r\"]*class=\\\"plotly-graph-div\\\"[^>]*>\\s*\"s, html)) !== nothing\n", " div_str = String(m.match)\n", " div_js = replace(replace(replace(div_str, \"\\\\\" => \"\\\\\\\\\"), \"\\\"\" => \"\\\\\\\"\"), \"\\n\" => \"\")\n", " inject = \"\"\n", " html = replace(html, m.match => inject; count=1)\n", " end\n", "\n", " display(\"text/html\", html)\n", " nothing\n", "end" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The Hessian matrix appears in the optimization literature, but the intuition for how the Hessian and its inverse transform vectors is opaque to me. Let's review second order partial derivatives, and then try to build intuition for the Hessian matrix. \n", "\n", "For the purpose of this intuition-building exercise, we'll work with functions $\\Reals^2 \\mapsto \\Reals^1$. I'll also use partial derivative notations $\\frac{\\partial}{\\partial y} f(x, y) = \\frac{\\partial f}{\\partial y} = f_y$ interchangeably. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 1. Partial Derivatives\n", "\n", "Take the $\\Reals^2 \\mapsto \\Reals^1$ function $f(x, y) = x^2 + 2y^2$.\n", "\n", "A partial derivative is the change in an \"output\" variable (in our case, $f$) with respect to infinitesimal changes in an \"input\" variable (in our case, $x$ or $y$). For example, $\\frac{\\partial}{\\partial y} f(x, y) = 4y$, which is to say, for any point in the domain, moving infinitsimally in the y direction changes f propotional to 4 times the y coordinate of the starting point. " ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "f_y (generic function with 1 method)" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f(x, y) = x^2 + 2y^2\n", "\n", "x = 6\n", "\n", "xlim=[-10, x]\n", "ylim=[-10, 10]\n", "\n", "xs = LinRange(xlim..., 101)\n", "ys = LinRange(ylim..., 101)\n", "zs = [f(x, y) for x in xs, y in ys]\n", "\n", "y = 4\n", "dy = 4\n", "f_y(y) = 4y" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# Build traces (PlotlyBase.* instead of PlotlyJS.*)\n", "traces = GenericTrace[]\n", "\n", "# main surface\n", "push!(traces, surface(x=xs, y=ys, z=zs;\n", " showscale=false, opacity=0.8))\n", "\n", "# thin \"plane\" at x\n", "ylim = [minimum(ys), maximum(ys)]\n", "push!(traces, surface(x=[x, x+0.001], y=ylim,\n", " z=[[maximum(zs) minimum(zs)];\n", " [maximum(zs) minimum(zs)]];\n", " showscale=false, colorscale=\"Greys\", opacity=0.2))\n", "\n", "# red curve at x\n", "push!(traces, scatter3d(x=fill(x, length(ys)), y=ys,\n", " z=[x^2 + 2y^2 for y in ys];\n", " showlegend=false, mode=\"lines\", line=attr(color=\"red\", width=2)))\n", "\n", "# orange tangent segments (one per slider step)\n", "step_ys = ys[1:5:end]\n", "for y in step_ys\n", " push!(traces, scatter3d(\n", " x=fill(x, 2),\n", " y=[y-dy, y+dy],\n", " z=[f(x,y) - f_y(y)*dy, f(x,y) + f_y(y)*dy],\n", " visible=false, showlegend=false, mode=\"lines\",\n", " line=attr(color=\"orange\", width=5)\n", " ))\n", "end\n", "\n", "# slider: toggle exactly one orange segment visible per step\n", "nsteps = length(step_ys)\n", "slider_steps = [\n", " attr(\n", " label = round(y, digits=2),\n", " method = \"update\",\n", " args = [attr(\n", " visible = vcat(fill(true, 3), [j == i for j in 1:nsteps])\n", " )]\n", " )\n", " for (i, y) in enumerate(step_ys)\n", "]\n", "\n", "scene = attr(\n", " xaxis = attr(range=[-10, 10]),\n", " yaxis = attr(range=[-10, 10]),\n", " zaxis = attr(range=[-50, 300]),\n", " aspectratio = attr(x=1, y=1, z=1)\n", ")\n", "\n", "layout = Layout(\n", " sliders = [attr(\n", " steps = slider_steps,\n", " active = 0, # start at first step\n", " currentvalue_prefix = \"x = $(x), y = \"\n", " )],\n", " scene = scene\n", ")\n", "\n", "p = Plot(traces, layout)\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# interactive plot\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=xs, y=ys, z=zs,\n", " showscale=false, opacity=0.8))\n", "push!(traces, surface(x=[x, x+0.001], y=ylim, z=[[maximum(zs), minimum(zs)] [maximum(zs), minimum(zs)]],\n", " showscale=false, colorscale=\"Greys\", opacity=0.2))\n", "push!(traces, scatter3d(x=fill(x, length(ys)), y=ys, z=[x^2 + 2y^2 for y in ys],\n", " showlegend=false, mode=\"lines\", line=attr(color=\"red\", width=2)))\n", "\n", "for y in ys[1:5:end]\n", " push!(traces, scatter3d(x=fill(x, 2), y=[y-dy, y+dy], z=[f(x,y)-f_y(y)*dy, f(x,y)+f_y(y)*dy],\n", " visible=false, showlegend=false, mode=\"lines\", line=attr(color=\"orange\", width=5)))\n", "end\n", "\n", "scene = attr(\n", " xaxis = attr(range=[-10,10]),\n", " yaxis = attr(range=[-10,10]),\n", " zaxis = attr(range=[-50,300]),\n", " aspectratio = attr(x=1, y=1, z=1)\n", " )\n", "\n", "layout = Layout(\n", " sliders=[attr(\n", " steps=[\n", " attr(\n", " label=round(y, digits=2),\n", " method=\"update\",\n", " args=[attr(visible=[fill(true, 3); fill(false, i-1); true; fill(false, 101-i)])]\n", " )\n", " for (i, y) in enumerate(ys[1:5:end])\n", " ],\n", " active = 0,\n", " currentvalue = attr(prefix=\"x = 6, y = \"),\n", " # pad_t=40\n", " )],\n", " scene = scene,\n", ")\n", "\n", "p = Plot(traces, layout)\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can plot the function $f_y$ for every starting point:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# plot partial derivative of f with respect to y, f_y\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=xs, y=ys, z=zs,\n", " showscale=false, opacity=0.8))\n", "push!(traces, surface(x=ylim, y=ylim, z=[[0, 0] [0, 0]],\n", " showscale=false, colorscale=\"Greys\", opacity=0.3))\n", "push!(traces, surface(x=xs, y=ys, z=[f_y(y) for x in xs, y in ys],\n", " showscale=false))\n", "\n", "p = Plot(traces, Layout(scene=scene))\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can do the exact same exercise with $f_x$:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "f_x (generic function with 1 method)" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "f_x(x) = 2x" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=xs, y=ys, z=zs,\n", " showscale=false, opacity=0.8))\n", "push!(traces, surface(x=ylim, y=ylim, z=[[0, 0] [0, 0]],\n", " showscale=false, colorscale=\"Greys\", opacity=0.3))\n", "push!(traces, surface(x=xs, y=ys, z=[f_x(x) for x in xs, y in ys],\n", " showscale=false))\n", "\n", "p = Plot(traces, Layout(scene=scene))\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "So the way the second order partial derivative is defined is as a composition, e.g. $f_{xx} = \\frac{\\partial}{\\partial x} \\left( \\frac{\\partial}{\\partial x}\\left( f(x, y) \\right) \\right) $. As second derivatives do, it captures the [change in the [change in the [output variable]]] with respect to infinitesimal changes in the input variable. This notion coincides with the curvature of the function: a positive second derivative at a particular point indicates that the output variable is concave up at a that point, and a negative second derivative indicates the output variable is concave down at that point. \n", "\n", "In the case of the function we've chosen, $f_{xx} = 2$ and $f_{yy} = 4$ which informs us that $f$ is concave up everywhere in the domain, which makes sense from looking at the plot. \n", "\n", "However, we've omitted the \"mixed\" partial derivatives here: $f_{xy}$ and $f_{yx}$. We can compute them to both be zero for this particular function. What does that tell us? " ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# plot partial derivative of f with respect to y, f_y\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=ylim, y=ylim, z=[[0, 0] [0, 0]],\n", " showscale=false, colorscale=\"Greys\", opacity=0.3))\n", "push!(traces, surface(x=xs, y=ys, z=[f_y(y) for x in xs, y in ys],\n", " showscale=false))\n", "push!(traces, scatter3d(x=[-7, 4], y=[5, 5], z=[f_y(5), f_y(5)],\n", " showlegend=false, mode=\"lines\", line=attr(color=\"orange\", width=5)))\n", "\n", "p = Plot(traces, Layout(scene=scene))\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It's clear from the graph above that infinitesimal changes in $x$ do not influence the value of $f_y$. Interpreting $f_{yx} = \\frac{\\partial}{\\partial x} \\left( \\frac{\\partial}{\\partial y}\\left( f(x, y) \\right) \\right)$ as $\\frac{\\partial}{\\partial x} f_y $, it's clear $f_{yx} = 0$. But how can we interpret that in terms of the original function f? " ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# interactive plot\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=xs, y=ys, z=zs,\n", " showscale=false, opacity=0.8))\n", "\n", "init_i = cld(length(xs), 2) # pick a middle index\n", "x0 = xs[init_i]\n", "\n", "for x in xs\n", " vis = (x == x0)\n", " push!(traces, surface(x=[x, x+0.001], y=ylim,\n", " z=[[maximum(zs), minimum(zs)] [maximum(zs), minimum(zs)]],\n", " visible=vis, showscale=false, colorscale=\"Greys\", opacity=0.2))\n", "\n", " push!(traces, scatter3d(x=fill(x, length(ys)), y=ys,\n", " z=[x^2 + 2y^2 for y in ys],\n", " visible=vis, showlegend=false, mode=\"lines\",\n", " line=attr(color=\"red\", width=2)))\n", "\n", " push!(traces, scatter3d(x=fill(x, 2), y=[y-dy, y+dy],\n", " z=[f(x,y)-f_y(y)*dy, f(x,y)+f_y(y)*dy],\n", " visible=vis, showlegend=false, mode=\"lines\",\n", " line=attr(color=\"orange\", width=5)))\n", "end\n", "\n", "layout = Layout(\n", " sliders = [attr(\n", " steps = [\n", " attr(\n", " label = round(x, digits=2),\n", " method = \"update\",\n", " args = [attr(visible = vcat(\n", " fill(true, 1), # main surface\n", " fill(false, 3*(i-1)), # before i\n", " fill(true, 3), # at i\n", " fill(false, 3*(length(xs)-i)) # after i\n", " ))]\n", " )\n", " for (i, x) in enumerate(xs)\n", " ],\n", " active = init_i - 1, # 0-based index for initial slider position\n", " currentvalue_prefix = \"x = \"\n", " )],\n", " scene = scene,\n", ")\n", "\n", "p = Plot(traces, layout)\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here is where we find the intuition. The mixed second order partial derivatives tell us how the slope along one coordinate axis changes as we move infinitesimally along an orthogonal coordinate axis. In the function we've chosen, the slice of the graph at any x coordinate is the same parabola (just vertically offset by $x^2$) and thus has the same slope for any y. The mixed second order partial derivatives could thus be said to relay information about the \"pucker\" of the graph of the function (I made that name up) which is the concavity of the graph with respect to two coordinate axes. \n", "\n", "In order to see this more clearly, let's add a mixed term to our function, to produce non-zero mixed second order partial derivatives:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "g_y (generic function with 1 method)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g(x, y) = x^2 + 2y^2 - x*y\n", "\n", "zs = [g(x, y) for x in xs, y in ys]\n", "\n", "g_x(x, y) = 2x - y\n", "g_y(x, y) = 4y - x" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", " \n", " \n", "\n", "\n", " \n", " \n", "
\n" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#collapse\n", "# built interactive plot\n", "traces = GenericTrace[]\n", "\n", "push!(traces, surface(x=xs, y=ys, z=zs, showscale=false, opacity=0.8))\n", "\n", "init_i = cld(length(xs), 2) # middle index\n", "x0 = xs[init_i]\n", "\n", "for x in xs\n", " vis = (x == x0)\n", " push!(traces, surface(x=[x, x+0.001], y=ylim,\n", " z=[[maximum(zs), minimum(zs)] [maximum(zs), minimum(zs)]],\n", " visible=vis, showscale=false, colorscale=\"Greys\", opacity=0.2))\n", "\n", " push!(traces, scatter3d(x=fill(x, length(ys)), y=ys,\n", " z=[g(x, y) for y in ys],\n", " visible=vis, showlegend=false, mode=\"lines\",\n", " line=attr(color=\"red\", width=2)))\n", "\n", " push!(traces, scatter3d(x=fill(x, 2), y=[y-dy, y+dy],\n", " z=[g(x,y)-g_y(x,y)*dy, g(x,y)+g_y(x,y)*dy],\n", " visible=vis, showlegend=false, mode=\"lines\",\n", " line=attr(color=\"orange\", width=5)))\n", "end\n", "\n", "layout = Layout(\n", " sliders = [attr(\n", " steps = [\n", " attr(\n", " label = round(x, digits=2),\n", " method = \"update\",\n", " args = [attr(visible = vcat(\n", " fill(true, 1),\n", " fill(false, 3*(i-1)),\n", " fill(true, 3),\n", " fill(false, 3*(length(xs)-i))\n", " ))]\n", " ) for (i, x) in enumerate(xs)\n", " ],\n", " active = init_i - 1, # <<< set initial position\n", " currentvalue_prefix = \"x = \"\n", " )],\n", " scene = scene\n", ")\n", "\n", "p = Plot(traces, layout)\n", "render_plotly_fastpages(p)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The slope along the y coordinate changes as we vary x. This function is \"puckered\" up. \n", "\n", "- Note: For twice-continuously differentiable functions, $f_{xy}$ = $f_{yx}$. This is called the [symmetry property of second derivatives](https://en.wikipedia.org/wiki/Symmetry_of_second_derivatives). And leads to a symmetric Hessian matrix. We'll assume that property of our functions throughout. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 2. Hessian Matrix" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Having now built intuition for each of the entries in the Hessian Matrix, we can begin to build intuition for the matrix itself: \n", "\n", "$$ \\bold{H}_f = \\begin{bmatrix}\n", "\\frac{\\partial^2 f}{\\partial x^2} & \\frac{\\partial^2 f}{\\partial x \\partial y} \\\\\n", "\\\\\n", "\\frac{\\partial^2 f}{\\partial y \\partial x} & \\frac{\\partial^2 f}{\\partial y^2}\n", "\\end{bmatrix} = \\begin{bmatrix}\n", "f_{xx} & f_{yx} \\\\\n", "\\\\\n", "f_{xy} & f_{yy}\n", "\\end{bmatrix}$$\n", "\n", "In some sense, the Hessian is a matrix-valued function we can evaluate for any point $\\bold{v_0} = \\begin{bmatrix} x_0 \\\\ y_0 \\end{bmatrix}$ in the domain of $f$:\n", "\n", "$$ \\bold{H}_f(\\bold{v_0}) = \\begin{bmatrix}\n", "f_{xx}(\\bold{v_0}) & f_{yx}(\\bold{v_0}) \\\\\n", "\\\\\n", "f_{xy}(\\bold{v_0}) & f_{yy}(\\bold{v_0})\n", "\\end{bmatrix}$$\n", "\n", "The Hessian shows up in the quadratic approximation (second order Taylor expansion) of multivariate functions around a particular point $\\bold{v_0}$:\n", "\n", "$$ Q_f(\\bold{v_0}) = \\color{green} f(\\bold{v_0}) \\color{black} + \\color{blue} \\nabla_f(\\bold{v_0}) \\cdot (\\bold{v} - \\bold{v_0}) \\color{black} + \\color{indigo} {1 \\over 2}(\\bold{v} - \\bold{v_0})^\\intercal \\left[ \\bold{H}_f(\\bold{v_0}) \\right] (\\bold{v} - \\bold{v_0})$$\n", "\n", "And we might ask: how does this matrix scale $(\\bold{v} - \\bold{v_0})$? Let's try to build some intuition: \n", "\n", "Suppose we multiply some vector against that Hessian Matrix: \n", "\n", "$$\\begin{bmatrix}\n", "f_{xx}(\\bold{v_0}) & f_{yx}(\\bold{v_0}) \\\\\n", "\\\\\n", "f_{xy}(\\bold{v_0}) & f_{yy}(\\bold{v_0})\n", "\\end{bmatrix} \\cdot \\begin{bmatrix} x \\\\ y \\end{bmatrix} = \n", "\\begin{bmatrix} \n", "f_{xx} \\cdot x + f_{xy} \\cdot y \\\\ \n", "f_{yy} \\cdot y + f_{yx} \\cdot x \\end{bmatrix} = \n", "$$\n", "\n", "$$\n", "\\begin{bmatrix} \n", "( \\text{the rate of change of the slope in the x direction as you move in the x direction} ) ( \\text{a distance in the x direction} ) + \\\\ ( \\text{the rate of change of the slope in the x direction as you move in the y direction} ) ( \\text{a distance in the y direction} ) \\\\ \\\\\n", "( \\text{the rate of change of the slope in the y direction as you move in the y direction} ) ( \\text{a distance in the y direction} ) + \\\\( \\text{the rate of change of the slope in the y direction as you move in the x direction} ) ( \\text{a distance in the x direction} )\n", "\\end{bmatrix} = $$\n", "\n", "$$ \n", "\\begin{bmatrix} \n", "\\text{the approximate slope in the x direction at the distance x from }x_0 \\\\\n", "\\text{the approximate slope in the y direction at the distance y from }y_0\n", "\\end{bmatrix} = \n", "$$\n", "\n", "$$ \\text{what a second-order approximation of } f \\text{ suggests the gradient is at } \\bold{v_0} + \\begin{bmatrix} x \\\\ y \\end{bmatrix} $$\n", "\n", "However, let's note that the quadratic approximation of $f$ includes as [quadratic form](https://en.wikipedia.org/wiki/Quadratic_form) $(\\bold{v} - \\bold{v_0})^\\intercal \\left[ \\bold{H}_f(\\bold{v_0}) \\right] (\\bold{v} - \\bold{v_0})$, not just a matrix-vector product. Applying this second product operation results in a scalar, which is the change in the value of the function as a result of extrapolating the (second-order) curvature of the function at $\\bold{v_0}$ out towards $\\bold{v} - \\bold{v_0}$. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3. Newton Method" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Newton's method for optimization is an iteration of \n", "- forming the quadratic approximation $Q_f(\\bold{v_0})$ of a function around the current point. \n", "- jumping to the min/max of the paraboloid approximation for the next iteration. \n", "\n", "This requires us to identify the min/max of the second-order approximation of $f$. Recall, the quadratic approximation (second order taylor expansion) of $f$ at a particular point $\\bold{v_0}$ is \n", "\n", "$$ Q_f(\\bold{v_0}) = \\color{green} f(\\bold{v_0}) \\color{black} + \\color{blue} \\nabla_f(\\bold{v_0}) \\cdot (\\bold{v} - \\bold{v_0}) \\color{black} + \\color{indigo} {1 \\over 2}(\\bold{v} - \\bold{v_0})^\\intercal \\left[ \\bold{H}_f(\\bold{v_0}) \\right] (\\bold{v} - \\bold{v_0})$$\n", "\n", "Searching for the extrema of this approximation, we can take $\\nabla \\left[ Q_f(v_0) \\right] = 0$. \n", "\n", "- $ \\nabla \\left[ \\color{green} f(\\bold{v_0}) \\color{black} \\right] = 0$\n", "- $ \\nabla \\left[ \\color{blue} \\nabla_f(\\bold{v_0}) \\cdot (\\bold{v} - \\bold{v_0}) \\color{black} \\right] = \\nabla_f(\\bold{v_0})$\n", "- $ \\nabla \\left[ \\color{indigo} {1 \\over 2}(\\bold{v} - \\bold{v_0})^\\intercal \\left[ \\bold{H}_f(\\bold{v_0}) \\right] (\\bold{v}- \\bold{v_0}) \\color{black} \\right] = {1 \\over 2} \\left( \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v}+ \\left[ \\bold{H}_f(\\bold{v_0}) \\right]^\\intercal \\right) \\bold{v}$ \n", "\n", "We can simplify ${1 \\over 2} \\left( \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v}+ \\left[ \\bold{H}_f(\\bold{v_0}) \\right]^\\intercal \\right) \\bold{v}$ further, because $\\left[ \\bold{H}_f(\\bold{v_0}) \\color{black} \\right] = \\left[ \\bold{H}_f(\\bold{v_0}) \\right]^\\intercal$, because of the symmetry property of second derivatives, so this expression simplifies to ${1 \\over 2} \\left( 2 \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\right) \\bold{v} = \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v}$. \n", "\n", "So all together, $\\nabla \\left[ Q_f(\\bold{v_0}) \\right] = \\color{blue} \\nabla_f(\\bold{v_0}) \\color{black} + \\color{indigo} \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v}$. \n", "\n", "Well, $\\nabla_f(\\bold{v_0}) + \\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v} = 0$ when:\n", "\n", "$$\n", "\\begin{aligned}\n", "\\left[ \\bold{H}_f(\\bold{v_0}) \\right] \\bold{v} &= -\\nabla_f(\\bold{v_0}) \\\\\n", "\\bold{v} &= - \\left[ \\bold{H}_f(\\bold{v_0}) \\right]^{-1} \\cdot \\nabla_f(\\bold{v_0})\n", "\\end{aligned}\n", "$$\n", "\n", "So $\\bold{v}$ is the vector which denotes the extremum of the quadratic approximation of $f$ around $\\bold{v_0}$. The familiar Newton method expression pops out of our above derivation: \n", "\n", "$$\\bold{v_{k+1}} = \\bold{v_k} - \\left[ \\bold{H}_f(\\bold{v_k}) \\right]^{-1} \\nabla_f(\\bold{v_k})$$\n", "\n", "And now we know that this expression which seemed to come out of nowhere $\\left[ \\bold{H}_f(\\bold{v_k}) \\right]^{-1} \\nabla_f(\\bold{v_k})$ including an inverse Hessian (!!) is really just the $\\bold{v}$ which solves $\\nabla \\left[ Q_f(v_0) \\right] = 0$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "- [Khan Academy: Second-Order Partial Derivatives](https://www.khanacademy.org/math/multivariable-calculus/multivariable-derivatives/partial-derivative-and-gradient-articles/a/second-partial-derivatives)\n", "- [Walter Schreiner's course notes](http://facstaff.cbu.edu/wschrein/media/M232%20Notes/M232L73.pdf): intuition for second order partial derivatives\n", "- [Khan Academy: The Hessian](https://www.khanacademy.org/math/multivariable-calculus/applications-of-multivariable-derivatives/quadratic-approximations/a/the-hessian)\n", "- [Khan Academy: Quadratic Approximations](https://www.khanacademy.org/math/multivariable-calculus/applications-of-multivariable-derivatives#quadratic-approximations)\n", "- [Kris Hauser's course notes](http://people.duke.edu/~kh269/teaching/b553/newtons_method.pdf): derivation of Newton's Method (especially page 4)\n" ] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.5.4", "language": "julia", "name": "julia-1.5" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.5.4" } }, "nbformat": 4, "nbformat_minor": 4 }