diff --git a/doc/pub/week39/html/week39-bs.html b/doc/pub/week39/html/week39-bs.html index 4d8240a34..011e4abdf 100644 --- a/doc/pub/week39/html/week39-bs.html +++ b/doc/pub/week39/html/week39-bs.html @@ -248,10 +248,11 @@ doconce format html week39.do.txt --html_style=bootstrap --pygments_html_style=d 2, None, 'same-code-but-now-with-momentum-gradient-descent'), - ('Same problem but now with AdaGrad', + ('Similar (second order function now) problem but now with ' + 'AdaGrad', 2, None, - 'same-problem-but-now-with-adagrad'), + 'similar-second-order-function-now-problem-but-now-with-adagrad'), ('And Logistic Regression', 2, None, 'and-logistic-regression'), ('Introducing "JAX":"https://jax.readthedocs.io/en/latest/"', 2, @@ -373,7 +374,7 @@ MathJax.Hub.Config({
  • But noen of these can compete with Newton's method
  • Including Stochastic Gradient Descent with Autograd
  • Same code but now with momentum gradient descent
  • -
  • Same problem but now with AdaGrad
  • +
  • Similar (second order function now) problem but now with AdaGrad
  • And Logistic Regression
  • Introducing "JAX":"https://jax.readthedocs.io/en/latest/"
  • Weekend challenge
  • @@ -406,7 +407,7 @@ MathJax.Hub.Config({
    -

    Sep 30, 2022

    +

    Oct 4, 2022


    diff --git a/doc/pub/week39/html/week39-reveal.html b/doc/pub/week39/html/week39-reveal.html index c5d92e266..4c46ed246 100644 --- a/doc/pub/week39/html/week39-reveal.html +++ b/doc/pub/week39/html/week39-reveal.html @@ -184,7 +184,7 @@ MathJax.Hub.Config({
    -

    Sep 30, 2022

    +

    Oct 4, 2022


    @@ -3486,7 +3486,7 @@ delta_momentum = 0.3
    -

    Same problem but now with AdaGrad

    +

    Similar (second order function now) problem but now with AdaGrad

    @@ -3494,8 +3494,8 @@ delta_momentum = 0.3
    -
    # Using Autograd to calculate gradients using SGD
    -# OLS example with AdaGrad
    +  
    # Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent
    +# OLS example
     from random import random, seed
     import numpy as np
     import autograd.numpy as np
    @@ -3506,53 +3506,45 @@ delta_momentum = 0.3
     def CostOLS(y,X,theta):
         return np.sum((y-X @ theta)**2)
     
    -n = 100
    -x = 2*np.random.rand(n,1)
    -y = 4+3*x+np.random.randn(n,1)
    +n = 10000
    +x = np.random.rand(n,1)
    +y = 2.0+3*x +4*x*x# +np.random.randn(n,1)
     
    -X = np.c_[np.ones((n,1)), x]
    +X = np.c_[np.ones((n,1)), x, x*x]
     XT_X = X.T @ X
     theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y)
     print("Own inversion")
     print(theta_linreg)
    -# Hessian matrix
    -H = (2.0/n)* XT_X
    -EigValues, EigVectors = np.linalg.eig(H)
    -print(f"Eigenvalues of Hessian Matrix:{EigValues}")
     
    -theta = np.random.randn(2,1)
    -eta = 1.0/np.max(EigValues)
    -Niterations = 100
     
     # Note that we request the derivative wrt third argument (theta, 2 here)
     training_gradient = grad(CostOLS,2)
    -
    -for iter in range(Niterations):
    -    gradients = (1.0/n)*training_gradient(y, X, theta)
    -    theta -= eta*gradients
    -print("theta from own gd")
    -print(theta)
    -print(np.size(gradients))
    -
    +# Define parameters for Stochastic Gradient Descent
     n_epochs = 50
     M = 5   #size of each minibatch
     m = int(n/M) #number of minibatches
    -theta = np.random.randn(2,1)
    +# Guess for unknown parameters theta
    +theta = np.random.randn(3,1)
     
    -# Including AdaGrad
    -delta  = 0.000001
    -r = [0.0 for _ in range(gradients.shape[0])]
    +# Value for learning rate
    +eta = 0.01
    +# Including AdaGrad parameter to avoid possible division by zero
    +delta  = 1e-8
     for epoch in range(n_epochs):
    +    # The outer product is calculated from scratch for each epoch
    +    Giter = np.zeros(shape=(3,3))
         for i in range(m):
             random_index = M*np.random.randint(m)
             xi = X[random_index:random_index+M]
             yi = y[random_index:random_index+M]
             gradients = (1.0/M)*training_gradient(yi, xi, theta)
    -        # calculate squared gradient by Hadamard multiplication
    -        r -= gradients*gradients
    +	# Calculate the outer product of the gradients
    +        Giter +=gradients @ gradients.T
    +	# Simpler algorithm with only diagonal elements
    +        Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))]
             # compute update
    -        update = 1.0/(delta+np.sqrt(r))*gradients
    -        theta -= eta*update
    +        update = np.multiply(Ginverse,gradients)
    +        theta -= update
     print("theta from own AdaGrad")
     print(theta)
     
    @@ -3569,6 +3561,8 @@ r = [0.0
    -

    Sep 30, 2022

    +

    Oct 4, 2022


    @@ -3413,7 +3414,7 @@ delta_momentum = 0.3









    -

    Same problem but now with AdaGrad

    +

    Similar (second order function now) problem but now with AdaGrad

    @@ -3421,8 +3422,8 @@ delta_momentum = 0.3
    -
    # Using Autograd to calculate gradients using SGD
    -# OLS example with AdaGrad
    +  
    # Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent
    +# OLS example
     from random import random, seed
     import numpy as np
     import autograd.numpy as np
    @@ -3433,53 +3434,45 @@ delta_momentum = 0.3
     def CostOLS(y,X,theta):
         return np.sum((y-X @ theta)**2)
     
    -n = 100
    -x = 2*np.random.rand(n,1)
    -y = 4+3*x+np.random.randn(n,1)
    +n = 10000
    +x = np.random.rand(n,1)
    +y = 2.0+3*x +4*x*x# +np.random.randn(n,1)
     
    -X = np.c_[np.ones((n,1)), x]
    +X = np.c_[np.ones((n,1)), x, x*x]
     XT_X = X.T @ X
     theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y)
     print("Own inversion")
     print(theta_linreg)
    -# Hessian matrix
    -H = (2.0/n)* XT_X
    -EigValues, EigVectors = np.linalg.eig(H)
    -print(f"Eigenvalues of Hessian Matrix:{EigValues}")
     
    -theta = np.random.randn(2,1)
    -eta = 1.0/np.max(EigValues)
    -Niterations = 100
     
     # Note that we request the derivative wrt third argument (theta, 2 here)
     training_gradient = grad(CostOLS,2)
    -
    -for iter in range(Niterations):
    -    gradients = (1.0/n)*training_gradient(y, X, theta)
    -    theta -= eta*gradients
    -print("theta from own gd")
    -print(theta)
    -print(np.size(gradients))
    -
    +# Define parameters for Stochastic Gradient Descent
     n_epochs = 50
     M = 5   #size of each minibatch
     m = int(n/M) #number of minibatches
    -theta = np.random.randn(2,1)
    +# Guess for unknown parameters theta
    +theta = np.random.randn(3,1)
     
    -# Including AdaGrad
    -delta  = 0.000001
    -r = [0.0 for _ in range(gradients.shape[0])]
    +# Value for learning rate
    +eta = 0.01
    +# Including AdaGrad parameter to avoid possible division by zero
    +delta  = 1e-8
     for epoch in range(n_epochs):
    +    # The outer product is calculated from scratch for each epoch
    +    Giter = np.zeros(shape=(3,3))
         for i in range(m):
             random_index = M*np.random.randint(m)
             xi = X[random_index:random_index+M]
             yi = y[random_index:random_index+M]
             gradients = (1.0/M)*training_gradient(yi, xi, theta)
    -        # calculate squared gradient by Hadamard multiplication
    -        r -= gradients*gradients
    +	# Calculate the outer product of the gradients
    +        Giter +=gradients @ gradients.T
    +	# Simpler algorithm with only diagonal elements
    +        Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))]
             # compute update
    -        update = 1.0/(delta+np.sqrt(r))*gradients
    -        theta -= eta*update
    +        update = np.multiply(Ginverse,gradients)
    +        theta -= update
     print("theta from own AdaGrad")
     print(theta)
     
    @@ -3497,6 +3490,7 @@ r = [0.0 And Logistic Regression diff --git a/doc/pub/week39/html/week39.html b/doc/pub/week39/html/week39.html index 2d7932742..511b312c4 100644 --- a/doc/pub/week39/html/week39.html +++ b/doc/pub/week39/html/week39.html @@ -352,10 +352,11 @@ div.toc p,a { 2, None, 'same-code-but-now-with-momentum-gradient-descent'), - ('Same problem but now with AdaGrad', + ('Similar (second order function now) problem but now with ' + 'AdaGrad', 2, None, - 'same-problem-but-now-with-adagrad'), + 'similar-second-order-function-now-problem-but-now-with-adagrad'), ('And Logistic Regression', 2, None, 'and-logistic-regression'), ('Introducing "JAX":"https://jax.readthedocs.io/en/latest/"', 2, @@ -399,7 +400,7 @@ MathJax.Hub.Config({
    -

    Sep 30, 2022

    +

    Oct 4, 2022


    @@ -3490,7 +3491,7 @@ delta_momentum = Same problem but now with AdaGrad +

    Similar (second order function now) problem but now with AdaGrad

    @@ -3498,8 +3499,8 @@ delta_momentum =
    -
    # Using Autograd to calculate gradients using SGD
    -# OLS example with AdaGrad
    +  
    # Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent
    +# OLS example
     from random import random, seed
     import numpy as np
     import autograd.numpy as np
    @@ -3510,53 +3511,45 @@ delta_momentum = def CostOLS(y,X,theta):
         return np.sum((y-X @ theta)**2)
     
    -n = 100
    -x = 2*np.random.rand(n,1)
    -y = 4+3*x+np.random.randn(n,1)
    +n = 10000
    +x = np.random.rand(n,1)
    +y = 2.0+3*x +4*x*x# +np.random.randn(n,1)
     
    -X = np.c_[np.ones((n,1)), x]
    +X = np.c_[np.ones((n,1)), x, x*x]
     XT_X = X.T @ X
     theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y)
     print("Own inversion")
     print(theta_linreg)
    -# Hessian matrix
    -H = (2.0/n)* XT_X
    -EigValues, EigVectors = np.linalg.eig(H)
    -print(f"Eigenvalues of Hessian Matrix:{EigValues}")
     
    -theta = np.random.randn(2,1)
    -eta = 1.0/np.max(EigValues)
    -Niterations = 100
     
     # Note that we request the derivative wrt third argument (theta, 2 here)
     training_gradient = grad(CostOLS,2)
    -
    -for iter in range(Niterations):
    -    gradients = (1.0/n)*training_gradient(y, X, theta)
    -    theta -= eta*gradients
    -print("theta from own gd")
    -print(theta)
    -print(np.size(gradients))
    -
    +# Define parameters for Stochastic Gradient Descent
     n_epochs = 50
     M = 5   #size of each minibatch
     m = int(n/M) #number of minibatches
    -theta = np.random.randn(2,1)
    +# Guess for unknown parameters theta
    +theta = np.random.randn(3,1)
     
    -# Including AdaGrad
    -delta  = 0.000001
    -r = [0.0 for _ in range(gradients.shape[0])]
    +# Value for learning rate
    +eta = 0.01
    +# Including AdaGrad parameter to avoid possible division by zero
    +delta  = 1e-8
     for epoch in range(n_epochs):
    +    # The outer product is calculated from scratch for each epoch
    +    Giter = np.zeros(shape=(3,3))
         for i in range(m):
             random_index = M*np.random.randint(m)
             xi = X[random_index:random_index+M]
             yi = y[random_index:random_index+M]
             gradients = (1.0/M)*training_gradient(yi, xi, theta)
    -        # calculate squared gradient by Hadamard multiplication
    -        r -= gradients*gradients
    +	# Calculate the outer product of the gradients
    +        Giter +=gradients @ gradients.T
    +	# Simpler algorithm with only diagonal elements
    +        Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))]
             # compute update
    -        update = 1.0/(delta+np.sqrt(r))*gradients
    -        theta -= eta*update
    +        update = np.multiply(Ginverse,gradients)
    +        theta -= update
     print("theta from own AdaGrad")
     print(theta)
     
    @@ -3574,6 +3567,7 @@ r = [0.0
    +

    Running this code we note an almost perfect agreement with the results from matrix inversion.











    And Logistic Regression

    diff --git a/doc/pub/week39/ipynb/ipynb-week39-src.tar.gz b/doc/pub/week39/ipynb/ipynb-week39-src.tar.gz index a6ec70f87..09221a7db 100644 Binary files a/doc/pub/week39/ipynb/ipynb-week39-src.tar.gz and b/doc/pub/week39/ipynb/ipynb-week39-src.tar.gz differ diff --git a/doc/pub/week39/ipynb/week39.ipynb b/doc/pub/week39/ipynb/week39.ipynb index 552f83332..6120e7b51 100644 --- a/doc/pub/week39/ipynb/week39.ipynb +++ b/doc/pub/week39/ipynb/week39.ipynb @@ -2,8 +2,10 @@ "cells": [ { "cell_type": "markdown", - "id": "70c42524", - "metadata": {}, + "id": "e24bca62", + "metadata": { + "editable": true + }, "source": [ "\n", @@ -12,21 +14,25 @@ }, { "cell_type": "markdown", - "id": "cffad5f6", - "metadata": {}, + "id": "e9d046e5", + "metadata": { + "editable": true + }, "source": [ "# Week 39: Optimization and Gradient Methods\n", "**Morten Hjorth-Jensen**, Department of Physics, University of Oslo and Department of Physics and Astronomy and Facility for Rare Isotope Beams, Michigan State University\n", "\n", - "Date: **Sep 30, 2022**\n", + "Date: **Oct 4, 2022**\n", "\n", "Copyright 1999-2022, Morten Hjorth-Jensen. Released under CC Attribution-NonCommercial 4.0 license" ] }, { "cell_type": "markdown", - "id": "34c8fbcd", - "metadata": {}, + "id": "a42fc301", + "metadata": { + "editable": true + }, "source": [ "## Plan for week 39\n", "\n", @@ -49,8 +55,10 @@ }, { "cell_type": "markdown", - "id": "9e3863c0", - "metadata": {}, + "id": "86a444d3", + "metadata": { + "editable": true + }, "source": [ "## Optimization, the central part of any Machine Learning algortithm\n", "\n", @@ -68,8 +76,10 @@ }, { "cell_type": "markdown", - "id": "aac1849d", - "metadata": {}, + "id": "18c198a7", + "metadata": { + "editable": true + }, "source": [ "## Revisiting our Logistic Regression case\n", "\n", @@ -83,8 +93,10 @@ }, { "cell_type": "markdown", - "id": "eda84591", - "metadata": {}, + "id": "a78a7971", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\begin{align*}\n", @@ -96,16 +108,20 @@ }, { "cell_type": "markdown", - "id": "e5e6a3fc", - "metadata": {}, + "id": "5f7ba8f7", + "metadata": { + "editable": true + }, "source": [ "where $\\boldsymbol{\\beta}$ are the weights we wish to extract from data, in our case $\\beta_0$ and $\\beta_1$." ] }, { "cell_type": "markdown", - "id": "f5a2093d", - "metadata": {}, + "id": "7c688cd4", + "metadata": { + "editable": true + }, "source": [ "## The equations to solve\n", "\n", @@ -118,8 +134,10 @@ }, { "cell_type": "markdown", - "id": "ababb255", - "metadata": {}, + "id": "f8d89f53", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\frac{\\partial \\mathcal{C}(\\boldsymbol{\\beta})}{\\partial \\boldsymbol{\\beta}} = -\\boldsymbol{X}^T\\left(\\boldsymbol{y}-\\boldsymbol{p}\\right).\n", @@ -128,8 +146,10 @@ }, { "cell_type": "markdown", - "id": "81ee053a", - "metadata": {}, + "id": "6520ed2c", + "metadata": { + "editable": true + }, "source": [ "If we in addition define a diagonal matrix $\\boldsymbol{W}$ with elements \n", "$p(y_i\\vert x_i,\\boldsymbol{\\beta})(1-p(y_i\\vert x_i,\\boldsymbol{\\beta})$, we can obtain a compact expression of the second derivative as" @@ -137,8 +157,10 @@ }, { "cell_type": "markdown", - "id": "22359c61", - "metadata": {}, + "id": "d4eceb89", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\frac{\\partial^2 \\mathcal{C}(\\boldsymbol{\\beta})}{\\partial \\boldsymbol{\\beta}\\partial \\boldsymbol{\\beta}^T} = \\boldsymbol{X}^T\\boldsymbol{W}\\boldsymbol{X}.\n", @@ -147,16 +169,20 @@ }, { "cell_type": "markdown", - "id": "f41b1042", - "metadata": {}, + "id": "d329b027", + "metadata": { + "editable": true + }, "source": [ "This defines what is called the Hessian matrix." ] }, { "cell_type": "markdown", - "id": "e7ff9f7e", - "metadata": {}, + "id": "3e676963", + "metadata": { + "editable": true + }, "source": [ "## Solving using Newton-Raphson's method\n", "\n", @@ -167,8 +193,10 @@ }, { "cell_type": "markdown", - "id": "488f1674", - "metadata": {}, + "id": "c1c25a40", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\beta}^{\\mathrm{new}} = \\boldsymbol{\\beta}^{\\mathrm{old}}-\\left(\\frac{\\partial^2 \\mathcal{C}(\\boldsymbol{\\beta})}{\\partial \\boldsymbol{\\beta}\\partial \\boldsymbol{\\beta}^T}\\right)^{-1}_{\\boldsymbol{\\beta}^{\\mathrm{old}}}\\times \\left(\\frac{\\partial \\mathcal{C}(\\boldsymbol{\\beta})}{\\partial \\boldsymbol{\\beta}}\\right)_{\\boldsymbol{\\beta}^{\\mathrm{old}}},\n", @@ -177,16 +205,20 @@ }, { "cell_type": "markdown", - "id": "fc58bce3", - "metadata": {}, + "id": "49a9b988", + "metadata": { + "editable": true + }, "source": [ "or in matrix form as" ] }, { "cell_type": "markdown", - "id": "9063bdf3", - "metadata": {}, + "id": "9ff443da", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\beta}^{\\mathrm{new}} = \\boldsymbol{\\beta}^{\\mathrm{old}}-\\left(\\boldsymbol{X}^T\\boldsymbol{W}\\boldsymbol{X} \\right)^{-1}\\times \\left(-\\boldsymbol{X}^T(\\boldsymbol{y}-\\boldsymbol{p}) \\right)_{\\boldsymbol{\\beta}^{\\mathrm{old}}}.\n", @@ -195,8 +227,10 @@ }, { "cell_type": "markdown", - "id": "b60c1a1e", - "metadata": {}, + "id": "478c51a9", + "metadata": { + "editable": true + }, "source": [ "The right-hand side is computed with the old values of $\\beta$. \n", "\n", @@ -205,8 +239,10 @@ }, { "cell_type": "markdown", - "id": "1836b157", - "metadata": {}, + "id": "4584f036", + "metadata": { + "editable": true + }, "source": [ "## Brief reminder on Newton-Raphson's method\n", "\n", @@ -223,8 +259,10 @@ }, { "cell_type": "markdown", - "id": "e1ecf421", - "metadata": {}, + "id": "0713a4c6", + "metadata": { + "editable": true + }, "source": [ "## The equations\n", "\n", @@ -237,8 +275,10 @@ }, { "cell_type": "markdown", - "id": "d1696b80", - "metadata": {}, + "id": "564a1130", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -251,8 +291,10 @@ }, { "cell_type": "markdown", - "id": "a9018f84", - "metadata": {}, + "id": "523c58c6", + "metadata": { + "editable": true + }, "source": [ "For small enough values of the function and for well-behaved\n", "functions, the terms beyond linear are unimportant, hence we obtain" @@ -260,8 +302,10 @@ }, { "cell_type": "markdown", - "id": "d8f2841f", - "metadata": {}, + "id": "89e62fe1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "f(x)+(s-x)f'(x)\\approx 0,\n", @@ -270,16 +314,20 @@ }, { "cell_type": "markdown", - "id": "9623e5a1", - "metadata": {}, + "id": "83c6b223", + "metadata": { + "editable": true + }, "source": [ "yielding" ] }, { "cell_type": "markdown", - "id": "624339a9", - "metadata": {}, + "id": "596fc91e", + "metadata": { + "editable": true + }, "source": [ "$$\n", "s\\approx x-\\frac{f(x)}{f'(x)}.\n", @@ -288,16 +336,20 @@ }, { "cell_type": "markdown", - "id": "642946f7", - "metadata": {}, + "id": "8aea8c22", + "metadata": { + "editable": true + }, "source": [ "Having in mind an iterative procedure, it is natural to start iterating with" ] }, { "cell_type": "markdown", - "id": "6c05a90c", - "metadata": {}, + "id": "ebd738c5", + "metadata": { + "editable": true + }, "source": [ "$$\n", "x_{n+1}=x_n-\\frac{f(x_n)}{f'(x_n)}.\n", @@ -306,8 +358,10 @@ }, { "cell_type": "markdown", - "id": "edcd1d95", - "metadata": {}, + "id": "cbade8a1", + "metadata": { + "editable": true + }, "source": [ "## Simple geometric interpretation\n", "\n", @@ -326,8 +380,10 @@ }, { "cell_type": "markdown", - "id": "2928dbd1", - "metadata": {}, + "id": "54284f78", + "metadata": { + "editable": true + }, "source": [ "## Extending to more than one variable\n", "\n", @@ -337,8 +393,10 @@ }, { "cell_type": "markdown", - "id": "14edc1f1", - "metadata": {}, + "id": "05e8bd6c", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\begin{array}{cc} f_1(x_1,x_2) &=0\\\\\n", @@ -348,16 +406,20 @@ }, { "cell_type": "markdown", - "id": "c6fcded8", - "metadata": {}, + "id": "75eb8466", + "metadata": { + "editable": true + }, "source": [ "which we Taylor expand to obtain" ] }, { "cell_type": "markdown", - "id": "5d46ef00", - "metadata": {}, + "id": "8c3a8a36", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\begin{array}{cc} 0=f_1(x_1+h_1,x_2+h_2)=&f_1(x_1,x_2)+h_1\n", @@ -372,16 +434,20 @@ }, { "cell_type": "markdown", - "id": "29230b7d", - "metadata": {}, + "id": "0d7ef225", + "metadata": { + "editable": true + }, "source": [ "Defining the Jacobian matrix ${\\bf \\boldsymbol{J}}$ we have" ] }, { "cell_type": "markdown", - "id": "4050e8df", - "metadata": {}, + "id": "314d130a", + "metadata": { + "editable": true + }, "source": [ "$$\n", "{\\bf \\boldsymbol{J}}=\\left( \\begin{array}{cc}\n", @@ -393,16 +459,20 @@ }, { "cell_type": "markdown", - "id": "b9acd139", - "metadata": {}, + "id": "d16e653c", + "metadata": { + "editable": true + }, "source": [ "we can rephrase Newton's method as" ] }, { "cell_type": "markdown", - "id": "378b9253", - "metadata": {}, + "id": "8e0d1824", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\left(\\begin{array}{c} x_1^{n+1} \\\\ x_2^{n+1} \\end{array} \\right)=\n", @@ -413,16 +483,20 @@ }, { "cell_type": "markdown", - "id": "0f5c3c0d", - "metadata": {}, + "id": "060c713d", + "metadata": { + "editable": true + }, "source": [ "where we have defined" ] }, { "cell_type": "markdown", - "id": "bb12d0e0", - "metadata": {}, + "id": "c6af6eb4", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\left(\\begin{array}{c} h_1^{n} \\\\ h_2^{n} \\end{array} \\right)=\n", @@ -433,8 +507,10 @@ }, { "cell_type": "markdown", - "id": "93274086", - "metadata": {}, + "id": "78b1b1a2", + "metadata": { + "editable": true + }, "source": [ "We need thus to compute the inverse of the Jacobian matrix and it\n", "is to understand that difficulties may\n", @@ -446,8 +522,10 @@ }, { "cell_type": "markdown", - "id": "122ee26a", - "metadata": {}, + "id": "0cdc4f21", + "metadata": { + "editable": true + }, "source": [ "## Steepest descent\n", "\n", @@ -461,8 +539,10 @@ }, { "cell_type": "markdown", - "id": "4235ab9d", - "metadata": {}, + "id": "f1abab41", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{x}_{k+1} = \\mathbf{x}_k - \\gamma_k \\nabla F(\\mathbf{x}_k),\n", @@ -471,8 +551,10 @@ }, { "cell_type": "markdown", - "id": "4fc1a870", - "metadata": {}, + "id": "9edac1f7", + "metadata": { + "editable": true + }, "source": [ "with $\\gamma_k > 0$.\n", "\n", @@ -483,8 +565,10 @@ }, { "cell_type": "markdown", - "id": "57d01841", - "metadata": {}, + "id": "d70a9f7a", + "metadata": { + "editable": true + }, "source": [ "## More on Steepest descent\n", "\n", @@ -496,8 +580,10 @@ }, { "cell_type": "markdown", - "id": "1d84085f", - "metadata": {}, + "id": "68089b45", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{x}_{k+1} = \\mathbf{x}_k - \\gamma_k \\nabla F(\\mathbf{x}_k), \\ \\ k \\geq 0.\n", @@ -506,8 +592,10 @@ }, { "cell_type": "markdown", - "id": "a2a841e7", - "metadata": {}, + "id": "2a3ac951", + "metadata": { + "editable": true + }, "source": [ "The parameter $\\gamma_k$ is often referred to as the step length or\n", "the learning rate within the context of Machine Learning." @@ -515,8 +603,10 @@ }, { "cell_type": "markdown", - "id": "cfe42a8d", - "metadata": {}, + "id": "5fa04a8c", + "metadata": { + "editable": true + }, "source": [ "## The ideal\n", "\n", @@ -541,8 +631,10 @@ }, { "cell_type": "markdown", - "id": "37e43c57", - "metadata": {}, + "id": "a0b7d4f5", + "metadata": { + "editable": true + }, "source": [ "## The sensitiveness of the gradient descent\n", "\n", @@ -561,8 +653,10 @@ }, { "cell_type": "markdown", - "id": "d73f9ec1", - "metadata": {}, + "id": "09df91da", + "metadata": { + "editable": true + }, "source": [ "## Convex functions\n", "\n", @@ -581,8 +675,10 @@ }, { "cell_type": "markdown", - "id": "9ef05d17", - "metadata": {}, + "id": "ffdbe770", + "metadata": { + "editable": true + }, "source": [ "## Convex function\n", "\n", @@ -591,8 +687,10 @@ }, { "cell_type": "markdown", - "id": "84ad8b8f", - "metadata": {}, + "id": "5e8f7025", + "metadata": { + "editable": true + }, "source": [ "## Conditions on convex functions\n", "\n", @@ -626,8 +724,10 @@ }, { "cell_type": "markdown", - "id": "b73096ec", - "metadata": {}, + "id": "034bba01", + "metadata": { + "editable": true + }, "source": [ "## More on convex functions\n", "\n", @@ -652,8 +752,10 @@ }, { "cell_type": "markdown", - "id": "01363dd6", - "metadata": {}, + "id": "2c6e6d8e", + "metadata": { + "editable": true + }, "source": [ "## Some simple problems\n", "\n", @@ -680,8 +782,10 @@ }, { "cell_type": "markdown", - "id": "3035712a", - "metadata": {}, + "id": "df7481e8", + "metadata": { + "editable": true + }, "source": [ "## Standard steepest descent\n", "\n", @@ -698,8 +802,10 @@ }, { "cell_type": "markdown", - "id": "5499b54c", - "metadata": {}, + "id": "8ef11a23", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{A}\\boldsymbol{x} = \\boldsymbol{b}.\n", @@ -708,16 +814,20 @@ }, { "cell_type": "markdown", - "id": "c564beea", - "metadata": {}, + "id": "de4dd468", + "metadata": { + "editable": true + }, "source": [ "In the iterative process we end up with a problem like" ] }, { "cell_type": "markdown", - "id": "2d78636b", - "metadata": {}, + "id": "a2024839", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{r}= \\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x},\n", @@ -726,8 +836,10 @@ }, { "cell_type": "markdown", - "id": "ae1fa6c3", - "metadata": {}, + "id": "eb954cfc", + "metadata": { + "editable": true + }, "source": [ "where $\\boldsymbol{r}$ is the so-called residual or error in the iterative process.\n", "\n", @@ -736,8 +848,10 @@ }, { "cell_type": "markdown", - "id": "a4f97f65", - "metadata": {}, + "id": "f2873926", + "metadata": { + "editable": true + }, "source": [ "## Gradient method\n", "\n", @@ -746,8 +860,10 @@ }, { "cell_type": "markdown", - "id": "a3990eaa", - "metadata": {}, + "id": "408639fb", + "metadata": { + "editable": true + }, "source": [ "$$\n", "P(\\boldsymbol{x})=\\frac{1}{2}\\boldsymbol{x}^T\\boldsymbol{A}\\boldsymbol{x} - \\boldsymbol{x}^T\\boldsymbol{b},\n", @@ -756,8 +872,10 @@ }, { "cell_type": "markdown", - "id": "c0d95c8a", - "metadata": {}, + "id": "ad2e03f0", + "metadata": { + "editable": true + }, "source": [ "with the constraint that the matrix $\\boldsymbol{A}$ is positive definite and\n", "symmetric. This defines also the Hessian and we want it to be positive definite." @@ -765,8 +883,10 @@ }, { "cell_type": "markdown", - "id": "da0edbb9", - "metadata": {}, + "id": "c4f86219", + "metadata": { + "editable": true + }, "source": [ "## Steepest descent method\n", "\n", @@ -776,8 +896,10 @@ }, { "cell_type": "markdown", - "id": "6129f857", - "metadata": {}, + "id": "4e634352", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x}_0=0,\n", @@ -786,16 +908,20 @@ }, { "cell_type": "markdown", - "id": "7ddc7a65", - "metadata": {}, + "id": "fce2a84a", + "metadata": { + "editable": true + }, "source": [ "or consider the system" ] }, { "cell_type": "markdown", - "id": "bad0a0c3", - "metadata": {}, + "id": "bb1c0818", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{A}\\boldsymbol{z} = \\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_0,\n", @@ -804,16 +930,20 @@ }, { "cell_type": "markdown", - "id": "ac85206a", - "metadata": {}, + "id": "28191eb6", + "metadata": { + "editable": true + }, "source": [ "instead." ] }, { "cell_type": "markdown", - "id": "52646ec0", - "metadata": {}, + "id": "efb36ffb", + "metadata": { + "editable": true + }, "source": [ "## Steepest descent method\n", "One can show that the solution $\\boldsymbol{x}$ is also the unique minimizer of the quadratic form" @@ -821,8 +951,10 @@ }, { "cell_type": "markdown", - "id": "c88c7f41", - "metadata": {}, + "id": "cbb54980", + "metadata": { + "editable": true + }, "source": [ "$$\n", "f(\\boldsymbol{x}) = \\frac{1}{2}\\boldsymbol{x}^T\\boldsymbol{A}\\boldsymbol{x} - \\boldsymbol{x}^T \\boldsymbol{x} , \\quad \\boldsymbol{x}\\in\\mathbf{R}^n.\n", @@ -831,8 +963,10 @@ }, { "cell_type": "markdown", - "id": "b817b6a2", - "metadata": {}, + "id": "02f3a580", + "metadata": { + "editable": true + }, "source": [ "This suggests taking the first basis vector $\\boldsymbol{r}_1$ (see below for definition) \n", "to be the gradient of $f$ at $\\boldsymbol{x}=\\boldsymbol{x}_0$, \n", @@ -841,8 +975,10 @@ }, { "cell_type": "markdown", - "id": "999b8908", - "metadata": {}, + "id": "e74d016a", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{A}\\boldsymbol{x}_0-\\boldsymbol{b},\n", @@ -851,8 +987,10 @@ }, { "cell_type": "markdown", - "id": "a86f550e", - "metadata": {}, + "id": "b221f8cd", + "metadata": { + "editable": true + }, "source": [ "and \n", "$\\boldsymbol{x}_0=0$ it is equal $-\\boldsymbol{b}$." @@ -860,8 +998,10 @@ }, { "cell_type": "markdown", - "id": "c6fad45a", - "metadata": {}, + "id": "d1fcda0e", + "metadata": { + "editable": true + }, "source": [ "## Final expressions\n", "We can compute the residual iteratively as" @@ -869,8 +1009,10 @@ }, { "cell_type": "markdown", - "id": "d90ca11a", - "metadata": {}, + "id": "5e580524", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{r}_{k+1}=\\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_{k+1},\n", @@ -879,16 +1021,20 @@ }, { "cell_type": "markdown", - "id": "0771db2b", - "metadata": {}, + "id": "b18d9941", + "metadata": { + "editable": true + }, "source": [ "which equals" ] }, { "cell_type": "markdown", - "id": "6c2d2fb6", - "metadata": {}, + "id": "ddfb86f1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{b}-\\boldsymbol{A}(\\boldsymbol{x}_k+\\alpha_k\\boldsymbol{r}_k),\n", @@ -897,16 +1043,20 @@ }, { "cell_type": "markdown", - "id": "528ce8f1", - "metadata": {}, + "id": "c15c0c1e", + "metadata": { + "editable": true + }, "source": [ "or" ] }, { "cell_type": "markdown", - "id": "03f1a1a7", - "metadata": {}, + "id": "1e70c362", + "metadata": { + "editable": true + }, "source": [ "$$\n", "(\\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_k)-\\alpha_k\\boldsymbol{A}\\boldsymbol{r}_k,\n", @@ -915,16 +1065,20 @@ }, { "cell_type": "markdown", - "id": "1c963f61", - "metadata": {}, + "id": "eda47934", + "metadata": { + "editable": true + }, "source": [ "which gives" ] }, { "cell_type": "markdown", - "id": "3b341cb6", - "metadata": {}, + "id": "a8807ca8", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\alpha_k = \\frac{\\boldsymbol{r}_k^T\\boldsymbol{r}_k}{\\boldsymbol{r}_k^T\\boldsymbol{A}\\boldsymbol{r}_k}\n", @@ -933,16 +1087,20 @@ }, { "cell_type": "markdown", - "id": "24a25a5b", - "metadata": {}, + "id": "797d70a7", + "metadata": { + "editable": true + }, "source": [ "leading to the iterative scheme" ] }, { "cell_type": "markdown", - "id": "e1f575b2", - "metadata": {}, + "id": "91253590", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x}_{k+1}=\\boldsymbol{x}_k+\\alpha_k\\boldsymbol{r}_{k},\n", @@ -951,8 +1109,10 @@ }, { "cell_type": "markdown", - "id": "733353ab", - "metadata": {}, + "id": "bb586698", + "metadata": { + "editable": true + }, "source": [ "## Steepest descent example" ] @@ -960,8 +1120,11 @@ { "cell_type": "code", "execution_count": 1, - "id": "29e79c22", - "metadata": {}, + "id": "c99ffffe", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "%matplotlib inline\n", @@ -990,8 +1153,10 @@ }, { "cell_type": "markdown", - "id": "9f0b7abf", - "metadata": {}, + "id": "23383c7d", + "metadata": { + "editable": true + }, "source": [ "And then as countor plot" ] @@ -999,8 +1164,11 @@ { "cell_type": "code", "execution_count": 2, - "id": "d9b290ba", - "metadata": {}, + "id": "c0bf99b8", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "pt.axis(\"equal\")\n", @@ -1010,8 +1178,10 @@ }, { "cell_type": "markdown", - "id": "16c54be4", - "metadata": {}, + "id": "2073d43c", + "metadata": { + "editable": true + }, "source": [ "Find guesses" ] @@ -1019,8 +1189,11 @@ { "cell_type": "code", "execution_count": 3, - "id": "3afefbb0", - "metadata": {}, + "id": "098e1651", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "x = guesses[-1]\n", @@ -1029,8 +1202,10 @@ }, { "cell_type": "markdown", - "id": "623773d6", - "metadata": {}, + "id": "8f3c24d6", + "metadata": { + "editable": true + }, "source": [ "Run it!" ] @@ -1038,8 +1213,11 @@ { "cell_type": "code", "execution_count": 4, - "id": "14404bfb", - "metadata": {}, + "id": "b53b03db", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "def f1d(alpha):\n", @@ -1053,8 +1231,10 @@ }, { "cell_type": "markdown", - "id": "938e1614", - "metadata": {}, + "id": "65e9016f", + "metadata": { + "editable": true + }, "source": [ "What happened?" ] @@ -1062,8 +1242,11 @@ { "cell_type": "code", "execution_count": 5, - "id": "03906b76", - "metadata": {}, + "id": "90a32f75", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "pt.axis(\"equal\")\n", @@ -1074,16 +1257,20 @@ }, { "cell_type": "markdown", - "id": "bef38425", - "metadata": {}, + "id": "6ecc0eaf", + "metadata": { + "editable": true + }, "source": [ "Note that we did only one iteration here. We can easily add more using our previous guesses." ] }, { "cell_type": "markdown", - "id": "974eb8f3", - "metadata": {}, + "id": "391ea54c", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "In the CG method we define so-called conjugate directions and two vectors \n", @@ -1094,8 +1281,10 @@ }, { "cell_type": "markdown", - "id": "93e05093", - "metadata": {}, + "id": "0f827f59", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{s}^T\\boldsymbol{A}\\boldsymbol{t}= 0.\n", @@ -1104,8 +1293,10 @@ }, { "cell_type": "markdown", - "id": "a7b2c6eb", - "metadata": {}, + "id": "8d85252f", + "metadata": { + "editable": true + }, "source": [ "The philosophy of the CG method is to perform searches in various conjugate directions\n", "of our vectors $\\boldsymbol{x}_i$ obeying the above criterion, namely" @@ -1113,8 +1304,10 @@ }, { "cell_type": "markdown", - "id": "1be3929c", - "metadata": {}, + "id": "4b962c6a", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x}_i^T\\boldsymbol{A}\\boldsymbol{x}_j= 0.\n", @@ -1123,8 +1316,10 @@ }, { "cell_type": "markdown", - "id": "9436c590", - "metadata": {}, + "id": "836212de", + "metadata": { + "editable": true + }, "source": [ "Two vectors are conjugate if they are orthogonal with respect to \n", "this inner product. Being conjugate is a symmetric relation: if $\\boldsymbol{s}$ is conjugate to $\\boldsymbol{t}$, then $\\boldsymbol{t}$ is conjugate to $\\boldsymbol{s}$." @@ -1132,8 +1327,10 @@ }, { "cell_type": "markdown", - "id": "c41bb077", - "metadata": {}, + "id": "c86350b5", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "An example is given by the eigenvectors of the matrix" @@ -1141,8 +1338,10 @@ }, { "cell_type": "markdown", - "id": "3aa89379", - "metadata": {}, + "id": "5303838b", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{v}_i^T\\boldsymbol{A}\\boldsymbol{v}_j= \\lambda\\boldsymbol{v}_i^T\\boldsymbol{v}_j,\n", @@ -1151,16 +1350,20 @@ }, { "cell_type": "markdown", - "id": "f930e2c1", - "metadata": {}, + "id": "a61193cf", + "metadata": { + "editable": true + }, "source": [ "which is zero unless $i=j$." ] }, { "cell_type": "markdown", - "id": "9fd91284", - "metadata": {}, + "id": "9acdaf47", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "Assume now that we have a symmetric positive-definite matrix $\\boldsymbol{A}$ of size\n", @@ -1169,8 +1372,10 @@ }, { "cell_type": "markdown", - "id": "9d942a7d", - "metadata": {}, + "id": "86d52ca2", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x}_{i+1}=\\boldsymbol{x}_{i}+\\alpha_i\\boldsymbol{p}_{i}.\n", @@ -1179,8 +1384,10 @@ }, { "cell_type": "markdown", - "id": "302a3c2c", - "metadata": {}, + "id": "8818bb65", + "metadata": { + "editable": true + }, "source": [ "We assume that $\\boldsymbol{p}_{i}$ is a sequence of $n$ mutually conjugate directions. \n", "Then the $\\boldsymbol{p}_{i}$ form a basis of $R^n$ and we can expand the solution \n", @@ -1189,8 +1396,10 @@ }, { "cell_type": "markdown", - "id": "b3fbd072", - "metadata": {}, + "id": "974957b1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x} = \\sum^{n}_{i=1} \\alpha_i \\boldsymbol{p}_i.\n", @@ -1199,8 +1408,10 @@ }, { "cell_type": "markdown", - "id": "be71d25b", - "metadata": {}, + "id": "0974b834", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "The coefficients are given by" @@ -1208,8 +1419,10 @@ }, { "cell_type": "markdown", - "id": "c4fd737f", - "metadata": {}, + "id": "9f4062b2", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{A}\\mathbf{x} = \\sum^{n}_{i=1} \\alpha_i \\mathbf{A} \\mathbf{p}_i = \\mathbf{b}.\n", @@ -1218,16 +1431,20 @@ }, { "cell_type": "markdown", - "id": "8d649170", - "metadata": {}, + "id": "d3ba22b1", + "metadata": { + "editable": true + }, "source": [ "Multiplying with $\\boldsymbol{p}_k^T$ from the left gives" ] }, { "cell_type": "markdown", - "id": "6b329c10", - "metadata": {}, + "id": "17f98074", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{p}_k^T \\boldsymbol{A}\\boldsymbol{x} = \\sum^{n}_{i=1} \\alpha_i\\boldsymbol{p}_k^T \\boldsymbol{A}\\boldsymbol{p}_i= \\boldsymbol{p}_k^T \\boldsymbol{b},\n", @@ -1236,16 +1453,20 @@ }, { "cell_type": "markdown", - "id": "bc1420b8", - "metadata": {}, + "id": "02eb7c67", + "metadata": { + "editable": true + }, "source": [ "and we can define the coefficients $\\alpha_k$ as" ] }, { "cell_type": "markdown", - "id": "3d65e7db", - "metadata": {}, + "id": "6b6e2535", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\alpha_k = \\frac{\\boldsymbol{p}_k^T \\boldsymbol{b}}{\\boldsymbol{p}_k^T \\boldsymbol{A} \\boldsymbol{p}_k}\n", @@ -1254,8 +1475,10 @@ }, { "cell_type": "markdown", - "id": "bfb356ea", - "metadata": {}, + "id": "581e48c0", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method and iterations\n", "\n", @@ -1272,8 +1495,10 @@ }, { "cell_type": "markdown", - "id": "0123ff7d", - "metadata": {}, + "id": "8d80d39f", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{x}_0=0,\n", @@ -1282,16 +1507,20 @@ }, { "cell_type": "markdown", - "id": "1b4d83fe", - "metadata": {}, + "id": "c85bf481", + "metadata": { + "editable": true + }, "source": [ "or consider the system" ] }, { "cell_type": "markdown", - "id": "6e66acbd", - "metadata": {}, + "id": "4d752584", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{A}\\boldsymbol{z} = \\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_0,\n", @@ -1300,16 +1529,20 @@ }, { "cell_type": "markdown", - "id": "9786cafd", - "metadata": {}, + "id": "d0c84272", + "metadata": { + "editable": true + }, "source": [ "instead." ] }, { "cell_type": "markdown", - "id": "47cba065", - "metadata": {}, + "id": "faff592c", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "One can show that the solution $\\boldsymbol{x}$ is also the unique minimizer of the quadratic form" @@ -1317,8 +1550,10 @@ }, { "cell_type": "markdown", - "id": "a2ac921c", - "metadata": {}, + "id": "fa4879d1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "f(\\boldsymbol{x}) = \\frac{1}{2}\\boldsymbol{x}^T\\boldsymbol{A}\\boldsymbol{x} - \\boldsymbol{x}^T \\boldsymbol{x} , \\quad \\boldsymbol{x}\\in\\mathbf{R}^n.\n", @@ -1327,8 +1562,10 @@ }, { "cell_type": "markdown", - "id": "6aad75d4", - "metadata": {}, + "id": "4b823ab8", + "metadata": { + "editable": true + }, "source": [ "This suggests taking the first basis vector $\\boldsymbol{p}_1$ \n", "to be the gradient of $f$ at $\\boldsymbol{x}=\\boldsymbol{x}_0$, \n", @@ -1337,8 +1574,10 @@ }, { "cell_type": "markdown", - "id": "040738cc", - "metadata": {}, + "id": "ed317a38", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{A}\\boldsymbol{x}_0-\\boldsymbol{b},\n", @@ -1347,8 +1586,10 @@ }, { "cell_type": "markdown", - "id": "956aa919", - "metadata": {}, + "id": "d1ef8688", + "metadata": { + "editable": true + }, "source": [ "and \n", "$\\boldsymbol{x}_0=0$ it is equal $-\\boldsymbol{b}$.\n", @@ -1358,8 +1599,10 @@ }, { "cell_type": "markdown", - "id": "ca75dc60", - "metadata": {}, + "id": "968986eb", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "Let $\\boldsymbol{r}_k$ be the residual at the $k$-th step:" @@ -1367,8 +1610,10 @@ }, { "cell_type": "markdown", - "id": "dfc34c66", - "metadata": {}, + "id": "e984713c", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{r}_k=\\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_k.\n", @@ -1377,8 +1622,10 @@ }, { "cell_type": "markdown", - "id": "4403d672", - "metadata": {}, + "id": "ffe69df7", + "metadata": { + "editable": true + }, "source": [ "Note that $\\boldsymbol{r}_k$ is the negative gradient of $f$ at \n", "$\\boldsymbol{x}=\\boldsymbol{x}_k$, \n", @@ -1391,8 +1638,10 @@ }, { "cell_type": "markdown", - "id": "4b9213fd", - "metadata": {}, + "id": "dd86af7f", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{p}_{k+1}=\\boldsymbol{r}_k-\\frac{\\boldsymbol{p}_k^T \\boldsymbol{A}\\boldsymbol{r}_k}{\\boldsymbol{p}_k^T\\boldsymbol{A}\\boldsymbol{p}_k} \\boldsymbol{p}_k.\n", @@ -1401,8 +1650,10 @@ }, { "cell_type": "markdown", - "id": "e05d442d", - "metadata": {}, + "id": "8ba70f4e", + "metadata": { + "editable": true + }, "source": [ "## Conjugate gradient method\n", "We can also compute the residual iteratively as" @@ -1410,8 +1661,10 @@ }, { "cell_type": "markdown", - "id": "3b5b19e7", - "metadata": {}, + "id": "8031d92b", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{r}_{k+1}=\\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_{k+1},\n", @@ -1420,16 +1673,20 @@ }, { "cell_type": "markdown", - "id": "41c98480", - "metadata": {}, + "id": "b2ec51f3", + "metadata": { + "editable": true + }, "source": [ "which equals" ] }, { "cell_type": "markdown", - "id": "500ef35c", - "metadata": {}, + "id": "85c0b7b4", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{b}-\\boldsymbol{A}(\\boldsymbol{x}_k+\\alpha_k\\boldsymbol{p}_k),\n", @@ -1438,16 +1695,20 @@ }, { "cell_type": "markdown", - "id": "019056be", - "metadata": {}, + "id": "38aad522", + "metadata": { + "editable": true + }, "source": [ "or" ] }, { "cell_type": "markdown", - "id": "4eeba33e", - "metadata": {}, + "id": "ecf1a239", + "metadata": { + "editable": true + }, "source": [ "$$\n", "(\\boldsymbol{b}-\\boldsymbol{A}\\boldsymbol{x}_k)-\\alpha_k\\boldsymbol{A}\\boldsymbol{p}_k,\n", @@ -1456,16 +1717,20 @@ }, { "cell_type": "markdown", - "id": "f558af91", - "metadata": {}, + "id": "7b32860c", + "metadata": { + "editable": true + }, "source": [ "which gives" ] }, { "cell_type": "markdown", - "id": "1b906f66", - "metadata": {}, + "id": "f024d4eb", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{r}_{k+1}=\\boldsymbol{r}_k-\\boldsymbol{A}\\boldsymbol{p}_{k},\n", @@ -1474,8 +1739,10 @@ }, { "cell_type": "markdown", - "id": "c802c90e", - "metadata": {}, + "id": "8a8bc39c", + "metadata": { + "editable": true + }, "source": [ "## Revisiting our first homework\n", "\n", @@ -1496,8 +1763,11 @@ { "cell_type": "code", "execution_count": 6, - "id": "f3f4306e", - "metadata": {}, + "id": "500755fd", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "x = 2*np.random.rand(m,1)\n", @@ -1506,8 +1776,10 @@ }, { "cell_type": "markdown", - "id": "d617011d", - "metadata": {}, + "id": "181dc8d0", + "metadata": { + "editable": true + }, "source": [ "with $x_i \\in [0,1] $ is chosen randomly using a uniform distribution. Additionally we have a stochastic noise chosen according to a normal distribution $\\cal {N}(0,1)$. \n", "The linear regression model is given by" @@ -1515,8 +1787,10 @@ }, { "cell_type": "markdown", - "id": "bf98c067", - "metadata": {}, + "id": "33fdef1d", + "metadata": { + "editable": true + }, "source": [ "$$\n", "h_\\beta(x) = \\boldsymbol{y} = \\beta_0 + \\beta_1 x,\n", @@ -1525,16 +1799,20 @@ }, { "cell_type": "markdown", - "id": "28cc86a9", - "metadata": {}, + "id": "1488d4ed", + "metadata": { + "editable": true + }, "source": [ "such that" ] }, { "cell_type": "markdown", - "id": "c015617b", - "metadata": {}, + "id": "b60cb10a", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{y}_i = \\beta_0 + \\beta_1 x_i.\n", @@ -1543,8 +1821,10 @@ }, { "cell_type": "markdown", - "id": "f2eca32c", - "metadata": {}, + "id": "10455f35", + "metadata": { + "editable": true + }, "source": [ "## Gradient descent example\n", "\n", @@ -1555,8 +1835,10 @@ }, { "cell_type": "markdown", - "id": "479ef1b2", - "metadata": {}, + "id": "0e7097e1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "X \\equiv \\begin{bmatrix}\n", @@ -1569,16 +1851,20 @@ }, { "cell_type": "markdown", - "id": "b562c596", - "metadata": {}, + "id": "c337500f", + "metadata": { + "editable": true + }, "source": [ "The cost/loss/risk function is given by (" ] }, { "cell_type": "markdown", - "id": "391bf862", - "metadata": {}, + "id": "539a67c8", + "metadata": { + "editable": true + }, "source": [ "$$\n", "C(\\beta) = \\frac{1}{n}||X\\beta-\\mathbf{y}||_{2}^{2} = \\frac{1}{n}\\sum_{i=1}^{100}\\left[ (\\beta_0 + \\beta_1 x_i)^2 - 2 y_i (\\beta_0 + \\beta_1 x_i) + y_i^2\\right]\n", @@ -1587,16 +1873,20 @@ }, { "cell_type": "markdown", - "id": "f4a16fba", - "metadata": {}, + "id": "bac485c6", + "metadata": { + "editable": true + }, "source": [ "and we want to find $\\beta$ such that $C(\\beta)$ is minimized." ] }, { "cell_type": "markdown", - "id": "f3add448", - "metadata": {}, + "id": "66ac2ee5", + "metadata": { + "editable": true + }, "source": [ "## The derivative of the cost/loss function\n", "\n", @@ -1605,8 +1895,10 @@ }, { "cell_type": "markdown", - "id": "002e739c", - "metadata": {}, + "id": "d0fa5500", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\nabla_{\\beta} C(\\beta) = \\frac{2}{n}\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\beta_0+\\beta_1x_i-y_i\\right) \\\\\n", @@ -1617,16 +1909,20 @@ }, { "cell_type": "markdown", - "id": "c5d0ba0a", - "metadata": {}, + "id": "a0f49a4d", + "metadata": { + "editable": true + }, "source": [ "where $X$ is the design matrix defined above." ] }, { "cell_type": "markdown", - "id": "0d3f481d", - "metadata": {}, + "id": "0d6b270c", + "metadata": { + "editable": true + }, "source": [ "## The Hessian matrix\n", "The Hessian matrix of $C(\\beta)$ is given by" @@ -1634,8 +1930,10 @@ }, { "cell_type": "markdown", - "id": "77136fc8", - "metadata": {}, + "id": "e296a37d", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{H} \\equiv \\begin{bmatrix}\n", @@ -1647,16 +1945,20 @@ }, { "cell_type": "markdown", - "id": "84130614", - "metadata": {}, + "id": "f182250d", + "metadata": { + "editable": true + }, "source": [ "This result implies that $C(\\beta)$ is a convex function since the matrix $X^T X$ always is positive semi-definite." ] }, { "cell_type": "markdown", - "id": "2cb8e687", - "metadata": {}, + "id": "fbedfc1a", + "metadata": { + "editable": true + }, "source": [ "## Simple program\n", "\n", @@ -1665,8 +1967,10 @@ }, { "cell_type": "markdown", - "id": "38e7d2b8", - "metadata": {}, + "id": "a903fc25", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\beta_{k+1} = \\beta_k - \\gamma \\nabla_\\beta C(\\beta_k), \\ k=0,1,\\cdots\n", @@ -1675,8 +1979,10 @@ }, { "cell_type": "markdown", - "id": "365c99d5", - "metadata": {}, + "id": "0457e928", + "metadata": { + "editable": true + }, "source": [ "We can use the expression we computed for the gradient and let use a\n", "$\\beta_0$ be chosen randomly and let $\\gamma = 0.001$. Stop iterating\n", @@ -1688,8 +1994,10 @@ }, { "cell_type": "markdown", - "id": "54fe7054", - "metadata": {}, + "id": "9bbfc263", + "metadata": { + "editable": true + }, "source": [ "## Gradient Descent Example\n", "\n", @@ -1699,8 +2007,11 @@ { "cell_type": "code", "execution_count": 7, - "id": "457e4af8", - "metadata": {}, + "id": "c39e3e8b", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "\n", @@ -1753,8 +2064,10 @@ }, { "cell_type": "markdown", - "id": "97c93fc7", - "metadata": {}, + "id": "72247c35", + "metadata": { + "editable": true + }, "source": [ "## And a corresponding example using **scikit-learn**" ] @@ -1762,8 +2075,11 @@ { "cell_type": "code", "execution_count": 8, - "id": "ef498ea5", - "metadata": {}, + "id": "61003ccc", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "# Importing various packages\n", @@ -1786,8 +2102,10 @@ }, { "cell_type": "markdown", - "id": "9275952a", - "metadata": {}, + "id": "8168d80a", + "metadata": { + "editable": true + }, "source": [ "## Gradient descent and Ridge\n", "\n", @@ -1796,8 +2114,10 @@ }, { "cell_type": "markdown", - "id": "5645bcc8", - "metadata": {}, + "id": "946a7392", + "metadata": { + "editable": true + }, "source": [ "$$\n", "C_{\\text{ridge}}(\\beta) = \\frac{1}{n}||X\\beta -\\mathbf{y}||^2 + \\lambda ||\\beta||^2, \\ \\lambda \\geq 0.\n", @@ -1806,16 +2126,20 @@ }, { "cell_type": "markdown", - "id": "8b2a9206", - "metadata": {}, + "id": "47969ccc", + "metadata": { + "editable": true + }, "source": [ "In order to minimize $C_{\\text{ridge}}(\\beta)$ using GD we adjust the gradient as follows" ] }, { "cell_type": "markdown", - "id": "f9078f32", - "metadata": {}, + "id": "efdf37a2", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\nabla_\\beta C_{\\text{ridge}}(\\beta) = \\frac{2}{n}\\begin{bmatrix} \\sum_{i=1}^{100} \\left(\\beta_0+\\beta_1x_i-y_i\\right) \\\\\n", @@ -1826,16 +2150,20 @@ }, { "cell_type": "markdown", - "id": "6e817a7b", - "metadata": {}, + "id": "480c4e6a", + "metadata": { + "editable": true + }, "source": [ "We can easily extend our program to minimize $C_{\\text{ridge}}(\\beta)$ using gradient descent and compare with the analytical solution given by" ] }, { "cell_type": "markdown", - "id": "ae8fee07", - "metadata": {}, + "id": "fc4e15ae", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\beta_{\\text{ridge}} = \\left(X^T X + n\\lambda I_{2 \\times 2} \\right)^{-1} X^T \\mathbf{y}.\n", @@ -1844,8 +2172,10 @@ }, { "cell_type": "markdown", - "id": "d9b0fe0f", - "metadata": {}, + "id": "1d7d09cd", + "metadata": { + "editable": true + }, "source": [ "## The Hessian matrix for Ridge Regression\n", "The Hessian matrix of Ridge Regression for our simple example is given by" @@ -1853,8 +2183,10 @@ }, { "cell_type": "markdown", - "id": "1ac722e2", - "metadata": {}, + "id": "5194a928", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{H} \\equiv \\begin{bmatrix}\n", @@ -1866,8 +2198,10 @@ }, { "cell_type": "markdown", - "id": "e25866ff", - "metadata": {}, + "id": "f50f65e6", + "metadata": { + "editable": true + }, "source": [ "This implies that the Hessian matrix is positive definite, hence the stationary point is a\n", "minimum.\n", @@ -1878,8 +2212,10 @@ }, { "cell_type": "markdown", - "id": "564ee0e6", - "metadata": {}, + "id": "d0903aad", + "metadata": { + "editable": true + }, "source": [ "## Program example for gradient descent with Ridge Regression" ] @@ -1887,8 +2223,11 @@ { "cell_type": "code", "execution_count": 9, - "id": "e19e3034", - "metadata": {}, + "id": "f77d8aec", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "from random import random, seed\n", @@ -1945,8 +2284,10 @@ }, { "cell_type": "markdown", - "id": "519c3099", - "metadata": {}, + "id": "eb135b12", + "metadata": { + "editable": true + }, "source": [ "## Using gradient descent methods, limitations\n", "\n", @@ -1965,8 +2306,10 @@ }, { "cell_type": "markdown", - "id": "c0224a76", - "metadata": {}, + "id": "76ccfd64", + "metadata": { + "editable": true + }, "source": [ "## Improving gradient descent with momentum\n", "\n", @@ -1976,8 +2319,11 @@ { "cell_type": "code", "execution_count": 10, - "id": "3ec96499", - "metadata": {}, + "id": "3106ae68", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "from numpy import asarray\n", @@ -2039,8 +2385,10 @@ }, { "cell_type": "markdown", - "id": "013518f9", - "metadata": {}, + "id": "4a09ca84", + "metadata": { + "editable": true + }, "source": [ "## Same code but now with momentum gradient descent" ] @@ -2048,8 +2396,11 @@ { "cell_type": "code", "execution_count": 11, - "id": "1177cebe", - "metadata": {}, + "id": "afdc854d", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "from numpy import asarray\n", @@ -2119,8 +2470,10 @@ }, { "cell_type": "markdown", - "id": "aa3e42e1", - "metadata": {}, + "id": "187eefc1", + "metadata": { + "editable": true + }, "source": [ "## Overview video on Stochastic Gradient Descent\n", "\n", @@ -2129,8 +2482,10 @@ }, { "cell_type": "markdown", - "id": "d2804e9d", - "metadata": {}, + "id": "6dd22999", + "metadata": { + "editable": true + }, "source": [ "## Batches and mini-batches\n", "\n", @@ -2148,8 +2503,10 @@ }, { "cell_type": "markdown", - "id": "0cb1f030", - "metadata": {}, + "id": "1c21e598", + "metadata": { + "editable": true + }, "source": [ "## Stochastic Gradient Descent (SGD)\n", "\n", @@ -2178,8 +2535,10 @@ }, { "cell_type": "markdown", - "id": "07c69389", - "metadata": {}, + "id": "55a20c5d", + "metadata": { + "editable": true + }, "source": [ "## Stochastic Gradient Descent\n", "\n", @@ -2193,8 +2552,10 @@ }, { "cell_type": "markdown", - "id": "ceb4ddeb", - "metadata": {}, + "id": "f5966077", + "metadata": { + "editable": true + }, "source": [ "$$\n", "C(\\mathbf{\\beta}) = \\sum_{i=1}^n c_i(\\mathbf{x}_i,\n", @@ -2204,8 +2565,10 @@ }, { "cell_type": "markdown", - "id": "b74d079f", - "metadata": {}, + "id": "e2ce9ab8", + "metadata": { + "editable": true + }, "source": [ "## Computation of gradients\n", "\n", @@ -2215,8 +2578,10 @@ }, { "cell_type": "markdown", - "id": "123b53cc", - "metadata": {}, + "id": "d9333453", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\nabla_\\beta C(\\mathbf{\\beta}) = \\sum_i^n \\nabla_\\beta c_i(\\mathbf{x}_i,\n", @@ -2226,8 +2591,10 @@ }, { "cell_type": "markdown", - "id": "3960f61f", - "metadata": {}, + "id": "ef1ee4d3", + "metadata": { + "editable": true + }, "source": [ "Stochasticity/randomness is introduced by only taking the\n", "gradient on a subset of the data called minibatches. If there are $n$\n", @@ -2238,8 +2605,10 @@ }, { "cell_type": "markdown", - "id": "0d84e11d", - "metadata": {}, + "id": "09390a19", + "metadata": { + "editable": true + }, "source": [ "## SGD example\n", "As an example, suppose we have $10$ data points $(\\mathbf{x}_1,\\cdots, \\mathbf{x}_{10})$ \n", @@ -2258,8 +2627,10 @@ }, { "cell_type": "markdown", - "id": "a23014d3", - "metadata": {}, + "id": "b75555e1", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\nabla_{\\beta}\n", @@ -2271,8 +2642,10 @@ }, { "cell_type": "markdown", - "id": "51163212", - "metadata": {}, + "id": "b29e4d92", + "metadata": { + "editable": true + }, "source": [ "## The gradient step\n", "\n", @@ -2281,8 +2654,10 @@ }, { "cell_type": "markdown", - "id": "1ac440c6", - "metadata": {}, + "id": "43c24131", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\beta_{j+1} = \\beta_j - \\gamma_j \\sum_{i \\in B_k}^n \\nabla_\\beta c_i(\\mathbf{x}_i,\n", @@ -2292,8 +2667,10 @@ }, { "cell_type": "markdown", - "id": "271fd2f1", - "metadata": {}, + "id": "f5bcccfe", + "metadata": { + "editable": true + }, "source": [ "where $k$ is picked at random with equal\n", "probability from $[1,n/M]$. An iteration over the number of\n", @@ -2304,8 +2681,10 @@ }, { "cell_type": "markdown", - "id": "1bdeb533", - "metadata": {}, + "id": "5c2a21d3", + "metadata": { + "editable": true + }, "source": [ "## Simple example code" ] @@ -2313,8 +2692,11 @@ { "cell_type": "code", "execution_count": 12, - "id": "4ee8002d", - "metadata": {}, + "id": "f6c39f8d", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import numpy as np \n", @@ -2335,8 +2717,10 @@ }, { "cell_type": "markdown", - "id": "41a7a26d", - "metadata": {}, + "id": "0dc5d595", + "metadata": { + "editable": true + }, "source": [ "Taking the gradient only on a subset of the data has two important\n", "benefits. First, it introduces randomness which decreases the chance\n", @@ -2349,8 +2733,10 @@ }, { "cell_type": "markdown", - "id": "ab92166f", - "metadata": {}, + "id": "b72fbc30", + "metadata": { + "editable": true + }, "source": [ "## When do we stop?\n", "\n", @@ -2368,8 +2754,10 @@ }, { "cell_type": "markdown", - "id": "a20c6717", - "metadata": {}, + "id": "35db732b", + "metadata": { + "editable": true + }, "source": [ "## Slightly different approach\n", "\n", @@ -2386,8 +2774,10 @@ }, { "cell_type": "markdown", - "id": "d25cf8bd", - "metadata": {}, + "id": "705c8d83", + "metadata": { + "editable": true + }, "source": [ "## Time decay rate\n", "\n", @@ -2403,8 +2793,11 @@ { "cell_type": "code", "execution_count": 13, - "id": "2dc23da5", - "metadata": {}, + "id": "d30150f7", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import numpy as np \n", @@ -2435,8 +2828,10 @@ }, { "cell_type": "markdown", - "id": "068ff42b", - "metadata": {}, + "id": "a5fab6e0", + "metadata": { + "editable": true + }, "source": [ "## Code with a Number of Minibatches which varies\n", "\n", @@ -2445,39 +2840,13 @@ }, { "cell_type": "code", - "execution_count": 4, - "id": "1d9d8857", - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Own inversion\n", - "[[3.83223285]\n", - " [3.16783558]]\n", - "Eigenvalues of Hessian Matrix:[0.2892 4.69411478]\n", - "theta from own gd\n", - "[[3.83223285]\n", - " [3.16783558]]\n", - "theta from own sdg\n", - "[[3.81685178]\n", - " [3.19056895]]\n" - ] - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEWCAYAAAB1xKBvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAsWElEQVR4nO3deZgU5bn38e89w7AMLsCwig5gRAQ0bqMmSKLGOBhjNMub44KJW4LLa1azk0QTDzlZfE8Wj4mSxGhkYkxM4olmYVAHRFwHBIVBEJBNdhCQHWbu94/uwWbo7qmeru6ubn6f6+qLnqrqqruqm7rrWeopc3dERETaU1boAEREpDgoYYiISCBKGCIiEogShoiIBKKEISIigShhiIhIIEoYUrLM7HYzm1ToOMJkZuea2cpCxyGHJiUMySszW2pmO81sm5mtMbP7zeywQsclIu1TwpBC+Ii7HwacApwKfLOw4RyazKxToWOQ4qKEIQXj7muAycQSBwBm9g0zW2xmb5tZk5l9LGHeNWb2jJndaWZvmdkbZvahhPlDzGxa/LNTgN6J2zOzS8xsnpltNrOpZjY8Yd5SM/uqmb1iZtvN7Ldm1s/M/hVf3xNm1jPZfrRWE5nZrWa2zsxWm9m1CfOnmtln2u5Hwt9uZjeb2evxbd1hZu8ys2fNbKuZ/cnMOrfZ5rfMbEM87rEJ07vEj89yM1trZveYWbc2cX7dzNYAvzOz3mb2ePyYbDKz6Wam84IkpR+GFIyZHQ18CFiUMHkx8D7gSOB7wCQzG5Aw/yxgAbFk8GPgt2Zm8Xl/AGbG590BXJ2wreOBh4AvAn2AfwKPtTkRfwK4ADge+AjwL+Bb8eXLgM+n2Z3+8ZgHAtcDd6dKMCmMAU4H3gN8DZgIXAUcA5wIXNFmW73j27oamGhmw+LzfhiP/xTguPgy323z2V7AIGAccCuwMr6P/eL7q/GCJCklDCmER83sbWAFsA64rXWGu//Z3Ve5e4u7Pwy8DpyZ8Nll7v5rd28GHgAGAP3MrBo4A/iOu+9296eBxxI+dxnwD3ef4u57gTuBbsCohGXucve17v4mMB14wd1fdvddwN+IVZ+lshf4vrvvdfd/AtuAYWmWb+vH7r7V3ecBc4F6d1/i7luIJa62227dz2nAP4D/iCfOccCX3H2Tu78N/AC4POFzLcBt8c/ujMc9ABgUj326a4A5SUEJQwrho+5+OHAucAIJVUdm9mkzmx2vItlM7Oo6sWppTesbd98Rf3sYcBTwlrtvT1h2WcL7oxL/dvcWYglrYMIyaxPe70zyd7rG+Y3uvi/h7x3tLN9WJttOtp9HESslVAIzE47fv+PTW62PJ8BWPyFWwqs3syVm9o0MYpZDjBKGFEz86vh+Ylf7mNkg4NfALUCVu/cgdrVtKVaRaDXQ08y6J0yrTni/ilg1DPFtGbHqnjc7vgeBbSd2Im/VP8v1JdvPVcAGYsllpLv3iL+OjHcwaHVA6cHd33b3W939WOAS4Mtmdn6W8UmJUsKQQvsZcIGZnQx0J3ZCWw8Qbzg+MchK3H0Z0Ah8z8w6m9loYu0Qrf4EfNjMzjezCmJ197uBZ8PakTRmAx83s0ozO45YG0e2WvfzfcDFwJ/jpaZfAz81s74AZjbQzMakWomZXWxmx8UT6BagmVi1lchBlDCkoNx9PfB74Lvu3gT8P+A5YlUyJwEzMljdlcQaxTcRaxf5fcJ2FhBrRL6L2JX4R4h1790Twm6056fAHmL79ABQl+X61gBvEStV1AE3uvtr8XlfJ1bF9LyZbQWeIH1bytD4MtuIHfdfuntDlvFJiTK1b4mISBAqYYiISCBKGCIiEogShoiIBKKEISIigRTV4GO9e/f2wYMHFzoMEZGiMnPmzA3u3qf9JdMrqoQxePBgGhsbCx2GiEhRMbNl7S/VPlVJiYhIIEoYIiISiBKGiIgEooQhIiKBKGGIiEggShgiIhKIEoaIiASihCEiIoEoYYiISCBKGCIiEkjOE4aZ3Wdm68xsbpJ5t5qZm1nvXMchIiLZyUcJ437gwrYTzewYoBZYnocYREQkSzlPGO7+NLFnLLf1U+BrgJ4RKyJSBArShmFmlwJvuvucAMuOM7NGM2tcv359HqITEZFk8p4wzKwS+Bbw3SDLu/tEd69x95o+fbIezl1ERDqoECWMdwFDgDlmthQ4GphlZv0LEIuIiASU9wcoufurQN/Wv+NJo8bdN+Q7FhERCS4f3WofAp4DhpnZSjO7PtfbFBGR8OW8hOHuV7Qzf3CuYxARkezpTm8REQlECUNERAJRwhARkUCUMEREJBAlDBERCUQJQ0REAlHCEBGRQJQwREQkECUMEREJRAlDREQCUcIQEZFAlDBERCQQJQwREQlECUNERAJRwhARkUCUMEREJBAlDBERCUQJQ0REAlHCEBGRQJQwREQkkJwnDDO7z8zWmdnchGk/MbPXzOwVM/ubmfXIdRwiIpKdfJQw7gcubDNtCnCiu78bWAh8Mw9xiIgUp7o6GDwYyspi/9bVFSSMnCcMd38a2NRmWr2774v/+TxwdK7jEBEpSnV1MG4cLFsG7rF/x40rSNKIQhvGdcC/Us00s3Fm1mhmjevXr89jWCIiETB+POzYceC0HTti0/OsoAnDzMYD+4CUqdLdJ7p7jbvX9OnTJ3/BiYhEwfLlmU3PoYIlDDO7BrgYGOvuXqg4REQirbo6s+k5VJCEYWYXAl8DLnH3He0tLyJyyJowASorD5xWWRmbnmf56Fb7EPAcMMzMVprZ9cD/AIcDU8xstpndk+s4RESK0tixMHEiDBoEZrF/J06MTc8zK6baoJqaGm9sbCx0GCIiRcXMZrp7TbbriUIvKRERKQJKGCIiEogShoiIBKKEISIigShhiIgUWkTGimpPp0IHICJySGsdK6p1+I/WsaKgIF1n01EJQ0SkkDIZK6rAJRElDBGRfEh1sg86VlQERq1VwhARybV0J/ugY0VFYNRaJQwRkVxLd7IPOlZUBEatVcIQEcm1dCf7oGNFRWDUWiUMEZFca+9kP3YsLF0KLS2xf5P1jorAqLVKGCIiuRbGyT4Co9YqYYhIcSmSm9wOENbJPkhJJIeUMESkeESga2mHhX2yL0DiVMIQkeIRga6lkVCgxKmEISLFI+yupcVYvQUFS5xKGCJSPMLsWhq16q1MkleB7slQwhCR4tFeb6NMTrpRqt7KNHkV6J4MJQwRKR7pehslO+leey307p08gUTgzun9Mk1eAbvpbl25lb+PfyG0MM3dQ1tZ0g2Y3QdcDKxz9xPj03oBDwODgaXAf7j7W+2tq6amxhsbG3MXrIgUr8GDY0kincrKdxJMquUHDYr1YsqnsrJYkmvLLNarKpm6ulhCWb48VrKYMIHmT17OzLrXqK9bT/1LPXlu6wj2UQHYTHevyTrMbFcQwP3AhW2mfQN40t2HAk/G/xYpvGJtBM21YjguQUoGiVftEbhzer+OVDHFu+mueG4lvz33QS775hD6dt3CWdeN5DtPnsuOfRV89b0zaPjp7PDidPecv4iVJOYm/L0AGBB/PwBYEGQ9p59+uovkzKRJ7pWV7rFrvdirsjI2/VBWLMdl0KADY0z1MnvnM5MmxT5nFvu3UPuUwTHetnab/+P2F/0Lp0z1Ezov2r/4UWWr/JrjnvY/3DLD1zWtP+AzQKOHcS4PYyXtbuTghLE54b0l/p3ks+OARqCxurq6I1+FSDCpTjiDBhU6ssIqluOS7KQbtbjTJagU85r3NvusP8z3H17Y4Of1mOWf4n5/g0HejPmasv7+j9PG+9xHX/eW5paUmy2ZhBH/+60g61EJQ3LKLPkJJvGK9FBUTMcl8aRbVeXeufOBMReyZJRBKWLVy2v8gc9O9ysHPeN9bN3+xb/R6Se+u6xLxvtU7AlDVVISPcVyJd0R2VS9FPNxiUqVk3va47hj4w6f/INGv/X0Bj+p64L9s/raOr9qyHT//Q3TfdXLazr8XRR7wvgJ8I34+28APw6yHiUMyaliqavPVLb7VarHpVWYSSXdulKU1Jox78oOB/fO7PLze870H32owV/+42vevLf5wPVnWtqLx3M6uBdDwgAeAlYDe4GVwPVAFbHeUa8DTwC9gqxLCUNyLkpXpGEJo4RQ7MclVfxhJsN21rVvwNFJv4eVHOVfPHWq//N7L/r29dvTbyPVd1lefnDMCfEUTcII86WEIZJCB65sI9kGkahte0RVVccSVqrG8NZ1hlXdluJkvqlzPz+1W5NfwSTfxoFxNHftFs6+JEt0CfEoYYhITHtXycXYBtFej6dMSgFBu9tmmUxb0lQ5vf/Il33CBQ2++IYfeUt1dXYltUmTYiWK9r7ThHiUMEQkpr2EUIxtEEFO8kETXqoSVgjr3rhokz/8xRl+/fFP+woGJk8YA48Jv0ovSKlRJQwlDJGDBDl5FFsbRJCTfNBSQKYljDTJdM/2Pf70XbP926Mb/Mzur7rR7OB+JJv9pz1u8z3lSbq83nRT+Ak7SKlRbRhKGCIHKcYqp/aEWcJor3qrqiptMn39iaV+92VT/dL+z/vhbHFwL2Ofjzp8jt9+boM/e+8rvnfn3ne21XZdufh+gpYai62XVJgvJQyRJIqxyqk9YbZhtK4vWQN3kvVsXrbZ//q15/zGEdP82E5L9y86uNNyv2H4NP/LV5/zt5ZuDr7tXHU6yKDUWFT3YYT1UsIQSaHYqpyCCKuXVKp1xtezd+def+7Xr/r3zmvwUYfP8XL2Orgfxla/pP/zfvdlU/31J5amHXojrXSlpTx9V2EljJwPbx4mDW8uUgKSDMvN2LF5DWHZjJXU37uEyU9V8OSq4Wz2Hhgt1HSfz5hT11N7WU/ec90IKiorst9Y63M62j7volXikOs5YhbO8OadwghGRCSQtifP1ifLQU5PmNvWbKPh7ibq/3cn9a9Vs3DvEOBoji5fxSeGvkrtReWcf/MJVA0dmd2GUiXDGTPgnnti5Yq2Wodcz3PS7AiVMEQkf/L00KKWfS3M+sNr1E9aR/2LPXh2ywj20plKtnNunyZqR++g9rqjOeGiY7EyC2ejyUoSraWH8ePTP9wp3YOSQhBWCUMJQ0TypyNPlgto5UurmfKrRdQ/Wc6UFcPY6FUAnNptPrXvXkvtJ4/k7M+OoMsRXbLaTkrpkuHy5cn3O3GZHD7lT1VSIlJ8qquTn1TTPVkuhR0bdvD0PU3U/20b9fMGMm/3UGAA/cvW8uEh86kdY1zwf4+n78jhwPCsQ29XumeEp9pvKNxT/jogH49oFRGJyeKxqN7izPnTAn580VQ+2GsWPfuU86Hv1PDLWe/hqO5bufPiqbzyyEJW7e3LA4tHM/aXZ9N3ZJ8c7UgS6R6zmmy/AaqqgjV4R+URuWF0tcrXS91qRUpABl2AV89Z67+/YbpfNWS69ytbu7836oldFvqtpzf45B80+o6NO/IWelrt3Q/T0a7PIdxng7rVikjB5Khr7K7Nu3hmYhP1j2yl/tUBzNk1DIDetoHa6gXUXuBccNNxHHVa/6y3lRO5OC4hdBRQo7eIFEa63kAZnhy9xWl6bDH1v1tJ/YzuTNswkp1UUsEeRveYR+1ZW6j9VD9OuWwYZZ0O0Rr0EDoKKGGISGFkecW7YcFGnrh7AfX/bqZ+8XG82TIAgBM6L6Z2+EpqP1rJOTcO57D+h4Ubd7GKUAnjEE3ZUlKi0iAYRDHFmkq63kBJ9m/Ptj1M+/lsxp89lZruTfQ9oSdX3DWKRxedyKijlvKbq6ez7Nk3mb/7Xfx89jl8+PYzDu1k0fYYXnRRhzsKhC6MhpB8vdToLQcppoH3iinWdFKNjVRVddD+7aSrX8NvHdzL2eujj5jtd5zf4C/cN9f37d5X6D2JnlS/kZtuymqsMDT4oIgX19DeUYs1Xa+d9ua1Oam1dOvmO7sekXT/NlT080e/+bxvWbElv/tXjHL0GwkrYagNQ4pbDu8cDl2UYk3XcA3tNmo33/8ge2/9Ol02rWGN9ecr/iMe5GrKiMj+Fasc/UZKog3DzL5kZvPMbK6ZPWRmXQsZjxShdDdLRU2UYh0//uDRU1sHwUsxb++tX+eeK5/m40c9T69rL6HbplV0Yi8f7z6Z495/DHt6pejqGsXvIpkotC9F6TeSTBjFlI68gIHAG0C3+N9/Aq5J9xlVSclBiqldIEqxpnuoT4p5zZiDe3X5Cv/sCdP8z19+1jcteSua+5epqMSeozgo9jaMeMJYAfQiNqbV40Btus8oYUhSxfTwoKjEmqKufFfVAN/WKXlbxPaKI3zBv5fEHiSUaj+isn+ZilL7Ug6OYdEnjNg+8AVgG7AeqEuxzDigEWisrq7O+sDJIaRYT175kORKdjvd/Aom+TqSPMq0tRdUis8WTUkilVw9RjUiwkoYBWvDMLOewKXAEOAooLuZXdV2OXef6O417l7Tp08eBxKT4tbaqLtsWey/fuuDeorxvocQbV+3nX9+7yW+8JOj+cru/2Qpg2jBWMlAHuz3FS753BB626bkH94Un56u/aNYRb3tICIK2ej9QeANd1/v7nuBvwKjChiPREFYDY+leFLrgJZ9Lcyqm88PL5zKB3q+TM9+FXz49jP49ZwzmNdzNH+79H7mP7qQgc0ruGHN97n8F6Ow9k6e6W7cK1ZZjKJ7SAmjmNKRF3AWMA+oBAx4APhcus+oDaPEhVnVUeJVDOm8OXO13/+Z6X7loGe8j63bv+snd33Nv3pGg0/50Uzf+dbO1Cto73uIUn1/mHJVhRmBqlFKpA3je8BrwFzgQaBLuuWVMEpcmCeiUj2pJbFj4w6f/INGv/X0Bj+xy8L9u9rX1vlVQ6b772+Y7qvnrM1spRneuBd6G0YETrKhiEh7T0kkjExfShglLsxSQUT+o+ZCS3OLv/LIAr/z4ga/oFejd2Gng3tndvn5PWf6jz7U4LMffs2b9zbnLohcntBL6buLyIVL3hIGMAU4OYyNZftSwihxYf/nKpWrVHdfO3ed1938jH/62Ok+oGz1/kMzosvr/qXTGvxfd7zk29dvz3zFUTxGETnJhiIiVaP5TBinAQ3A74ABYWy0oy8ljBJXDFeWeTrB7tqyy5+8c5Z//awGP7Vb0zs9W22DX149w++79mlf8eKq7DYSpeOdeFyTnWCLtf0pIskv71VSwCeAOcBtxO/OzvdLCeMQMGlSrL9/Yt//fJ7A8lF3n2QbLc0t3vTYIv/Zx6b6RX1e9Eq2Obh3Yo+fc+TLPuGCBn/pgXnhjvAakZNZ0uMahbjCEJGknNeEEe/FdCJwI7ABWAl8KowAMnkpYRwC8v0fLPHkXVXlXlGRetthnGCT7N8u6+K32F37Jx1fscRvOWmq//3bL/jWN7fmYKfjIlJdkvK4RqHkE4YIVPvls0pqBrAq3pZxB3AxcBxwFzAxjCCCvpQwDgH5vOrN9Mo23TIB7Nm+x3f2HJD082vL+/m9Y6f5G9NXdHxfMj0pBTnW+TjZtVcNFZW2lSKWz4QxkvijXJPMmx9GEEFfShhFIpuTTD6veoNc2SZuu7w8+fzy8qSrb2lu8defWOp3XzbVL+n/vB/OFm8mw/0Lciw7Wipr73P5Ku1FpWqshEWiWy1wbBhBBH0pYRSBbE8y+Tx5pLuy7WAJ462lm/2vX3vObxg+zYd0WrZ/9pBOy/zGEdN8+xH9g+9f0GOZzTFLl5Dy9V1EpJ6/lEUiYeT7pYRRBLI9yeTz5JFJ3fmkSSlLGDt7DfDbz23wUYfP8XL2Orgfzha/tP/zfvdlU/31J5bGRnjNdP+CHstMSmWZlP7yWdqLQD1/KVPCkGgK4ySTr5NHspN3586xxu/Ebadp62gd4dVo9jO7v+rfHt3gT9812/ds35P9/gU9lkETS6bJWFVFJUMJQ6Kp2E4yAU7ezQOPSbpPeyn3if2/4w9/cYZvWLgx/NjCTgSZfjcdKe2ppBBJShgSTSVQH928t9lfemCeT7igwc858uWUDdUtue5+msmxDHKi7kjpL5MEUALffalSwpDoKuRVZqYnuPiye/sP9Ib3f8cvO2aGV9mG/ee7U7s1+abO/ZKfaPNRagrzWOa69FdspctDiBKGSFvprnDbnHj3XPtZ31vR9YBlt1HpN3G3X/2u6V538zO+du669tdbTHK9H1G5EVAOooQhuZOPEkIutpHqCreqylvanChTVjOlegxwqdTN53I/VMKILCUMyY1cX4W2HSsqzG2kuMJtSXYSS/XS1XDHlUpJrASFlTAK+YhWiaJcPtq0rg6uuw42bjx4XgjbaBl4dFafB/QM52yMHQsTJ8KgQWAW+3fixNh0KQkWSz7FoaamxhsbGwsdRmkzSz0v299K797Jk0XitltaAq/OW5ymxxZT/7uV1M/oTr8Nc7mbW+jOOwmvuXNXyg6rxDZtSr69xH2qrNQJTkqSmc1095ps16MSRimrq4PBg6GsLPZvXV37nykvz2x6JtIlCwh0db9hwUb++PlnuXbodI6uWMOJHz2OL//vuSzd2osepxzLwo9/g5ajj9l/hVt+32+wX/wilgwSVVbCjTfqalgkA50KHYDkSF0djBv3TvXSsmWxvyH9SbG5ObPpYamshAkTDpq8Z9senvttE5Mf3kz9nL7M2nECzih62ltccPR8aj+wiAtuOJbq974LeBdwDvCd5NsYPx6WL48lpgkTlBxEMhVGQ0i+Xmr0zkBHe6zksqdLssZucC8r298w2tLc4gv+vcR/8YmpfnHfF7w7bzvEHiT0viNm+x3nN/gL980N9iChUunZJMUjor85SqGXFNADeAR4DZgPvDfd8koYGehon/hc9nSZNOngBxRVVPjbP53of/7ys/7ZE6b5oPIV+2cdV/GG33ziVH/0m8/7lhVbMt9Wrnt7RfDEIAUU4V5ipZIwHgA+E3/fGeiRbnkljAzkasjrbE2a5C3HVHsL5m916eff7fJfXsY+B/cj2OwfG/Cc/+qKab64YVl228llSUnJSJKJ8H0oYSWMgvWSMrMjgdnEnqkRKAj1kspA2zYMKGgvoCVTl1M/cSmTp3bmqdXD2cqRlNHMmYc1UXvaRsZcWcWZVw+nU9eQmtXKypL36sqwJ1ZSgwfH2oTaGjQIli7Nbt0R+94kA7n8zWUprF5ShUwYpwATgSbgZGAm8AV3395muXHAOIDq6urTlyX7jyrJ1dUVrKF368qtNPxyPpP/vpv6BYNYvG8QANXlKxkzdAljLq7gAzcPp+eQHrkJIJcn9WJNRpJbEf7uSiFh1ADPA2e7+wtm9nNgq7un6OKiEkaUNe9pZmbda0yetJ76xp48t3UkzXSiO9s4r28TY96/k9rPVDP0gsFYWZp7PcKSyyv1Yk1GklsRLh2Wwn0YK4GV7v5C/O9HgNMKGE80dOTeiQJZ/tyb/Obq6fzHMc/Rp+tWzrpuJLc99X527avg66OeYerPZrPp7c48tvZMbvnzORw/Zkh+kgW8c9dxVdU707p1C2fdEyYkv68jSbfgjKW6F0V3oEffoXCnexgNIR19AdOBYfH3twM/Sbd8yTd6R7iXhbv726vf9sdve9E/f/JUP6Hzov0hDixb5dcOfdof+twMX//ahsxXnKtG3lz3+MpXzK093tQALh1EifSSOgVoBF4BHgV6pls+sgkjrJNHxHpZNO9t9pmTmvy/xjT4eT1meQW7Hdy7sd0v7P2i//elDT7vf19/53nVHZHLk3rEjmdgrb+nxGQRwQsIKR5hJQyNJZWtMOstI1B/vWrWGqb8ahGTpxhTlg9jg/cG4OSuC6g9aTVjPnkEZ392BF17dA1ng2oPSC3CjahSXMJqw9DQINlKNbrr1VfH3meSNKqrk58gclh/vXPTTqbf28Tkv7xN/dyBzN09FOhPX1vPhYNfY8wY44M3DaX/u4cBw8IPYPnyzKZnogDHM1S5PDYiHaDBB7OV6j9vc3Os5JFJo3UuG1PjvMV55ZGF3HnxVGqrZtKzyhjzrdO5e+Z76N99Kz++aCqzH17A6gcm82DLVVx17/vof8mZuWt8z2Ujbx6OZ06pAVyiJox6rXy9ItmGkaqevKP15TloTF07d51PuukZ//Sx071/2Zr9oY3sstC/dFqD/+uOl3z7+u0HxpBNu0LiPlRVxV6p9kd3TacW8U4QUjwohUbvTF+RTBjJ/lNnMnZTDuzassufvHOWf/2sBj+1W9P+UKpsg19ePcPvu/ZpX/HiqtQryHZYkXTHI9kJr5hP6rmmYyMhCCthqNE7DHV1sTaLZEOA56GB0luc1/65hPr7VlL/TCVT149gB93pxF7OPnIetWduZsyn+3Hq5cMo6xSgFjJoY3GyO8nHj0/ebpBIjbYieaVG7yhpbdhO1lsqR/XlG1/fxMLP/4KhT9xDr33r6EY1LzCBRRWjuO6kRsZ8tBvn3Dicw486JfOVp2os7tUr1nNn+fLY+7ffhj17YvNan7fRtgNAMmq0FSlOYRRT8vXqcJVUvor1OdzOnu17/Om7Zvv4sxv8jO5z/Uoe9G0cWPXT3KVbeDelta1Wqqhw79w5dVVT66u8vP1lon4fhEiJQVVSAUV4fJd0vMVZ9OQy6n+9jPrpXXlqzQi2cTjl7OM9hzfx+J5aeuxee/AHw6ruaVvdtG1b+49YbVVZmbqkUQTHXqTUFP3ggx3RoYRRRDc/bV62hafunk/943uof30wb+yLdZ8c0mk5Y45fSu3FnTnv5uH0GHRk/m9KS7W9tgYNeqcto7XqCmDTJj0aVaRASmHwwfzoyM1PeRoAcN+ufTw38VW+d95Uzj7iFXoP7s4nfvIe/jD/FN7dezV3XzaN159YxpK91fzqWyv42MOX02NIz1hMrSfitnLVRz/IelvbbMaOjSXjlhbYsCH2ammJTVOyECleYdRr5evVoTaMVF1Ey8uTtzXkuO/7G9NX+L1jp/knBj7rPeytWM9bmv3M7q/6t0c3+PS75/ie7XsO/FCymDp3Pvhxp7nso5+qXSPdPRYiEgnoPoyA2rsvoO2JNuQB67a+udX//u0X/JaTpvrQiiX7V3dM+Ur/zLBp/qcvPesbF21Kv5JUMVVV5bePvu4JEClKYSWM0m/DgAMbcMvK0t8vkWXbQPOeZl7+4wImP7iO+pd68OyWkeyjgkq2c17fedSO3knt9ccw7MIMng1R7IPoiUhBqdG7o9o7+XagkXzlS6up/+Ui6p8s54mVw9josYf2nNZtPrUnr6X2k0cy6jMj6HJEl47FXEQN9yISPbpxr6PaG8F0woR2b8Dbvm47T9/TRP2jO6hvGkjT7uOAAQwoW8PFx86ndozxwZuPp+/I4cDw7GMOEJOISK6Vfi+pttobwTTJYxZb7rmX2Z1q+PFFU/lgr1n06teJi247g3tePpOjD9vMnRdP5ZVHFvLm3n7cv2g0V45aSt8PnxFeL6tD4dGPIhJ5h16VFCQfA6nNyXfNK+uY8suF1NfDlGXHs7alLwAndV1I7chVjPk/hzN63Ai69ep28LqL8EZBESldasMIKkByANi1eRfPTGxi8p+3Uj93AK/sij0sqI+t54LqhYypdT5443EcdVr/9NtTe4OIRIzaMNpKlhjgwKv91gHyAL/iSpoeW8zk+1ZSP+Mwpm0cyS5OozO7Gd1zHj88dyq1n+7PyZ88nrJOZwePQ09JE5ESVRoljFTVQN26JR3/aEOnfpzc8jKrWgYAMLzzYmpHrKT2o5Wcc9MIuvft3vEgVcIQkYgpmaFBzKzczF42s8c7vJIUz9X2FIPl9dq3jtED3+A3V09n2bNv0rT7Xfzs5XO46LYzsksWUPyPBe2IsIdSydPQLCKSoTDu/svmBXwZ+APweHvLnn766cnvNjZLeid0S4o7u1uqq7O6a7Jdh9Id0WEPpaLHkoqEjlK409vMjgYeACYAX3b3i9MtX3Pssd64du0BpQnv1o3dLRV03b31oOW3lR9BV9tDp3273pmoHkvhCrsKTlV6IqErlSqpnwFfA1KOb2Fm48ys0cwaW5YvP6jqyXbu5O3dFeygTffWykoOe+CXdLr/N7p/IZfCbuRXpwGRyCpYwjCzi4F17j4z3XLuPtHda9y9pizZGFBAbzbR5Xf3Jk8MiUNtR2F47VKrn0817HlHh1kPe30iEppCljDOBi4xs6XAH4EPmNmkdB/YQ0XS6TaomvJrPhWtxJBMa2+uZctitfOt3XyLOWmE3ch/KHYaECkWYTSEZPsCziVIo/eQIcXdIBry0OmREXYj/6HUaUAkDyiFRu9WZnYu8BVvr9G7psYbv/SlQHduR5KGKReRAtDQIMVIPYBEpABKpZdUdvLdgJzt9gpVP19qDe0iUhhh1Gvl63XAI1rzfYNXWNvLd/28boQTOeRRSm0YQR1QJZVJ9U7AEWvTSrW9qirYsCGzdeWTqsFEDnmqkgp6g1dYXVlTbW/jxmhX8ehGOBEJSfEmjKA3eKUYmJDx48PZXus2oko3wolISIo3YQRtQA7rCjtdw3SUr9Z1I5yIhKR4E0bQ51yHdYU9dmysvSKMdeWTngcuIiEp3oQBwcaJCvMK++c/L86r9aiNpyUiRam4E0YQYV5h62pdRA5hxdutVkREAlG3WhERySslDBERCUQJQ0REAlHCEBGRQJQwREQkECUMEREJRAlDREQCUcIQEZFAlDBERCQQJQwREQmkYAnDzI4xswYzazKzeWb2hULFIiIi7etUwG3vA25191lmdjgw08ymuHtTAWMSEZEUClbCcPfV7j4r/v5tYD4wsFDxiIhIepFowzCzwcCpwAtJ5o0zs0Yza1y/fn3eYxMRkZiCJwwzOwz4C/BFd9/adr67T3T3Gnev6dOnT/4DFBERoMAJw8wqiCWLOnf/ayFjERGR9ArZS8qA3wLz3f2/CxWHiIgEU8gSxtnAp4APmNns+OuiAsYjIiJpFKxbrbs/A1ihti8iIpkpeKO3iIgUByUMEREJRAlDREQCUcIQEZFAlDBERCQQJQwREQlECUNERAJRwhARkUCUMEREJBAlDBERCUQJQ0REAlHCEBGRQJQwREQkECUMEREJRAlDREQCUcIQEZFAlDBERCQQJQwREQlECUNERAJRwhARkUAKmjDM7EIzW2Bmi8zsG4WMRURE0itYwjCzcuBu4EPACOAKMxtRqHhERCS9QpYwzgQWufsSd98D/BG4tIDxiIhIGp0KuO2BwIqEv1cCZ7VdyMzGAePif+42s7l5iC1bvYENhQ4iAMUZnmKIERRn2IolzmFhrKSQCSMQd58ITAQws0Z3rylwSO1SnOEqhjiLIUZQnGErpjjDWE8hq6TeBI5J+Pvo+DQREYmgQiaMl4ChZjbEzDoDlwN/L2A8IiKSRsGqpNx9n5ndAkwGyoH73H1eOx+bmPvIQqE4w1UMcRZDjKA4w3ZIxWnuHsZ6RESkxOlObxERCUQJQ0REAolMwmhvmBAz62JmD8fnv2BmgxPmfTM+fYGZjSlgjF82syYze8XMnjSzQQnzms1sdvyV08b9AHFeY2brE+L5TMK8q83s9fjr6gLH+dOEGBea2eaEeXk5nmZ2n5mtS3X/j8X8Ir4Pr5jZaQnz8nks24tzbDy+V83sWTM7OWHe0vj02WF1v8wiznPNbEvCd/vdhHl5G0ooQJxfTYhxbvz32Cs+Ly/H08yOMbOG+Dlnnpl9Icky4f4+3b3gL2KN3ouBY4HOwBxgRJtlbgbuib+/HHg4/n5EfPkuwJD4esoLFON5QGX8/U2tMcb/3hahY3kN8D9JPtsLWBL/t2f8fc9Cxdlm+c8R6xiR7+P5fuA0YG6K+RcB/wIMeA/wQr6PZcA4R7Vun9hwPC8kzFsK9I7I8TwXeDzb30uu42yz7EeAp/J9PIEBwGnx94cDC5P8Xw/19xmVEkaQYUIuBR6Iv38EON/MLD79j+6+293fABbF15f3GN29wd13xP98nti9JfmWzZArY4Ap7r7J3d8CpgAXRiTOK4CHchRLSu7+NLApzSKXAr/3mOeBHmY2gPwey3bjdPdn43FA4X6bQY5nKnkdSijDOAv121zt7rPi798G5hMbQSNRqL/PqCSMZMOEtN3x/cu4+z5gC1AV8LP5ijHR9cQye6uuZtZoZs+b2UdzEF+roHF+Il5EfcTMWm+gzNexzGhb8aq9IcBTCZPzdTzbk2o/8nksM9X2t+lAvZnNtNhQPIX2XjObY2b/MrOR8WmRPJ5mVknsRPuXhMl5P54Wq6I/FXihzaxQf5+RHxqkGJnZVUANcE7C5EHu/qaZHQs8ZWavuvviwkTIY8BD7r7bzG4gVnL7QIFiCeJy4BF3b06YFqXjWTTM7DxiCWN0wuTR8WPZF5hiZq/Fr7ALYRax73abmV0EPAoMLVAsQXwEmOHuiaWRvB5PMzuMWML6ortvzdV2IDoljCDDhOxfxsw6AUcCGwN+Nl8xYmYfBMYDl7j77tbp7v5m/N8lwFRiVwO50G6c7r4xIbbfAKcH/Ww+40xwOW2K/Hk8nu1JtR+RG/rGzN5N7Pu+1N03tk5POJbrgL+RmyrdQNx9q7tvi7//J1BhZr2J4PGMS/fbzPnxNLMKYsmizt3/mmSRcH+fuW6YCdh404lYo8sQ3mnQGtlmmf/LgY3ef4q/H8mBjd5LyE2jd5AYTyXWMDe0zfSeQJf4+97A6+SowS5gnAMS3n8MeN7faQh7Ix5vz/j7XoWKM77cCcQaEa0QxzO+jcGkbqT9MAc2Kr6Y72MZMM5qYu17o9pM7w4cnvD+WeDCAsbZv/W7JnaiXR4/toF+L/mKMz7/SGLtHN0LcTzjx+X3wM/SLBPq7zNnB7sDO38RsVb+xcD4+LTvE7tSB+gK/Dn+o38RODbhs+Pjn1sAfKiAMT4BrAVmx19/j08fBbwa/5G/Clxf4GP5X8C8eDwNwAkJn70ufowXAdcWMs7437cDP2zzubwdT2JXj6uBvcTqea8HbgRujM83Yg8CWxyPpaZAx7K9OH8DvJXw22yMTz82fhznxH8T4wsc5y0Jv83nSUhwyX4vhYozvsw1xDrcJH4ub8eTWLWiA68kfK8X5fL3qaFBREQkkKi0YYiISMQpYYiISCBKGCIiEogShoiIBKKEISIigShhiIhIIEoYIiISiBKGSBbizyO4IP7+P83srkLHJJIrGnxQJDu3Ad+PDzR3KnBJgeMRyRnd6S2SJTObBhwGnOux5xKIlCRVSYlkwcxOIvbksz1KFlLqlDBEOij+5LI6Yk8122ZmOXuinkgUKGGIdED8SWt/BW519/nAHcTaM0RKltowREQkEJUwREQkECUMEREJRAlDREQCUcIQEZFAlDBERCQQJQwREQlECUNERAL5/88Z6O+sYRe6AAAAAElFTkSuQmCC\n", - "text/plain": [ - "
    " - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], + "execution_count": 14, + "id": "a00269ad", + "metadata": { + "collapsed": false, + "editable": true + }, + "outputs": [], "source": [ "# Importing various packages\n", "from math import exp, sqrt\n", @@ -2549,8 +2918,10 @@ }, { "cell_type": "markdown", - "id": "9e64e8f4", - "metadata": {}, + "id": "1c87e8f8", + "metadata": { + "editable": true + }, "source": [ "## Replace or not\n", "\n", @@ -2562,8 +2933,10 @@ }, { "cell_type": "markdown", - "id": "4a9ac2cb", - "metadata": {}, + "id": "0a059042", + "metadata": { + "editable": true + }, "source": [ "## Momentum based GD\n", "\n", @@ -2575,8 +2948,10 @@ }, { "cell_type": "markdown", - "id": "7f63c4a3", - "metadata": {}, + "id": "6390a129", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{v}_{t}=\\gamma \\mathbf{v}_{t-1}+\\eta_{t}\\nabla_\\theta E(\\boldsymbol{\\theta}_t) \\nonumber\n", @@ -2585,8 +2960,10 @@ }, { "cell_type": "markdown", - "id": "6b126246", - "metadata": {}, + "id": "819f92b8", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -2601,8 +2978,10 @@ }, { "cell_type": "markdown", - "id": "98cdbb32", - "metadata": {}, + "id": "f00de63c", + "metadata": { + "editable": true + }, "source": [ "where we have introduced a momentum parameter $\\gamma$, with\n", "$0\\le\\gamma\\le 1$, and for brevity we dropped the explicit notation to\n", @@ -2618,8 +2997,10 @@ }, { "cell_type": "markdown", - "id": "9da9b8d8", - "metadata": {}, + "id": "2b7f4f9b", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\Delta \\boldsymbol{\\theta}_{t+1} = \\gamma \\Delta \\boldsymbol{\\theta}_t -\\ \\eta_{t}\\nabla_\\theta E(\\boldsymbol{\\theta}_t),\n", @@ -2628,16 +3009,20 @@ }, { "cell_type": "markdown", - "id": "9926b3c4", - "metadata": {}, + "id": "9204976e", + "metadata": { + "editable": true + }, "source": [ "where we have defined $\\Delta \\boldsymbol{\\theta}_{t}= \\boldsymbol{\\theta}_t-\\boldsymbol{\\theta}_{t-1}$." ] }, { "cell_type": "markdown", - "id": "57779254", - "metadata": {}, + "id": "29e64510", + "metadata": { + "editable": true + }, "source": [ "## More on momentum based approaches\n", "\n", @@ -2650,8 +3035,10 @@ }, { "cell_type": "markdown", - "id": "348a6e95", - "metadata": {}, + "id": "a4ddeb6b", + "metadata": { + "editable": true + }, "source": [ "$$\n", "m {d^2 \\mathbf{w} \\over dt^2} + \\mu {d \\mathbf{w} \\over dt }= -\\nabla_w E(\\mathbf{w}).\n", @@ -2660,16 +3047,20 @@ }, { "cell_type": "markdown", - "id": "87aae119", - "metadata": {}, + "id": "72f692f0", + "metadata": { + "editable": true + }, "source": [ "We can discretize this equation in the usual way to get" ] }, { "cell_type": "markdown", - "id": "a4669289", - "metadata": {}, + "id": "8f41eeba", + "metadata": { + "editable": true + }, "source": [ "$$\n", "m { \\mathbf{w}_{t+\\Delta t}-2 \\mathbf{w}_{t} +\\mathbf{w}_{t-\\Delta t} \\over (\\Delta t)^2}+\\mu {\\mathbf{w}_{t+\\Delta t}- \\mathbf{w}_{t} \\over \\Delta t} = -\\nabla_w E(\\mathbf{w}).\n", @@ -2678,16 +3069,20 @@ }, { "cell_type": "markdown", - "id": "8554a357", - "metadata": {}, + "id": "bc96de6b", + "metadata": { + "editable": true + }, "source": [ "Rearranging this equation, we can rewrite this as" ] }, { "cell_type": "markdown", - "id": "1fa72611", - "metadata": {}, + "id": "c9cf2f01", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\Delta \\mathbf{w}_{t +\\Delta t}= - { (\\Delta t)^2 \\over m +\\mu \\Delta t} \\nabla_w E(\\mathbf{w})+ {m \\over m +\\mu \\Delta t} \\Delta \\mathbf{w}_t.\n", @@ -2696,8 +3091,10 @@ }, { "cell_type": "markdown", - "id": "8a825fb9", - "metadata": {}, + "id": "5fc6491d", + "metadata": { + "editable": true + }, "source": [ "## Momentum parameter\n", "\n", @@ -2710,8 +3107,10 @@ }, { "cell_type": "markdown", - "id": "36b5f3de", - "metadata": {}, + "id": "fc566a20", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\gamma= {m \\over m +\\mu \\Delta t }, \\qquad \\eta = {(\\Delta t)^2 \\over m +\\mu \\Delta t}.\n", @@ -2720,8 +3119,10 @@ }, { "cell_type": "markdown", - "id": "49c47929", - "metadata": {}, + "id": "890fa973", + "metadata": { + "editable": true + }, "source": [ "Thus, as the name suggests, the momentum parameter is proportional to\n", "the mass of the particle and effectively provides inertia.\n", @@ -2751,8 +3152,10 @@ }, { "cell_type": "markdown", - "id": "cfd6472d", - "metadata": {}, + "id": "e958e4c2", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{v}_{t}=\\gamma \\mathbf{v}_{t-1}+\\eta_{t}\\nabla_\\theta E(\\boldsymbol{\\theta}_t +\\gamma \\mathbf{v}_{t-1}) \\nonumber\n", @@ -2761,8 +3164,10 @@ }, { "cell_type": "markdown", - "id": "ccb3b789", - "metadata": {}, + "id": "277711eb", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -2777,16 +3182,20 @@ }, { "cell_type": "markdown", - "id": "52a50fee", - "metadata": {}, + "id": "efc50ab5", + "metadata": { + "editable": true + }, "source": [ "One of the major advantages of NAG is that it allows for the use of a larger learning rate than GDM for the same choice of $\\gamma$." ] }, { "cell_type": "markdown", - "id": "60b98b21", - "metadata": {}, + "id": "a5abcf0b", + "metadata": { + "editable": true + }, "source": [ "## Second moment of the gradient\n", "\n", @@ -2814,8 +3223,10 @@ }, { "cell_type": "markdown", - "id": "a3f65369", - "metadata": {}, + "id": "972eeee6", + "metadata": { + "editable": true + }, "source": [ "## RMS prop\n", "\n", @@ -2827,8 +3238,10 @@ }, { "cell_type": "markdown", - "id": "36e483ff", - "metadata": {}, + "id": "f3e5aef6", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -2843,8 +3256,10 @@ }, { "cell_type": "markdown", - "id": "c565ac54", - "metadata": {}, + "id": "0fba3e6a", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{s}_t =\\beta \\mathbf{s}_{t-1} +(1-\\beta)\\mathbf{g}_t^2 \\nonumber\n", @@ -2853,8 +3268,10 @@ }, { "cell_type": "markdown", - "id": "da98e855", - "metadata": {}, + "id": "53a91982", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\theta}_{t+1}=\\boldsymbol{\\theta}_t - \\eta_t { \\mathbf{g}_t \\over \\sqrt{\\mathbf{s}_t +\\epsilon}}, \\nonumber\n", @@ -2863,8 +3280,10 @@ }, { "cell_type": "markdown", - "id": "3f14f027", - "metadata": {}, + "id": "50547f38", + "metadata": { + "editable": true + }, "source": [ "where $\\beta$ controls the averaging time of the second moment and is\n", "typically taken to be about $\\beta=0.9$, $\\eta_t$ is a learning rate\n", @@ -2879,8 +3298,10 @@ }, { "cell_type": "markdown", - "id": "2717fec0", - "metadata": {}, + "id": "86fe9ed3", + "metadata": { + "editable": true + }, "source": [ "## [ADAM optimizer](https://arxiv.org/abs/1412.6980)\n", "\n", @@ -2900,8 +3321,10 @@ }, { "cell_type": "markdown", - "id": "170fb635", - "metadata": {}, + "id": "7cd27183", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -2916,8 +3339,10 @@ }, { "cell_type": "markdown", - "id": "faa2bccc", - "metadata": {}, + "id": "ab967150", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{m}_t = \\beta_1 \\mathbf{m}_{t-1} + (1-\\beta_1) \\mathbf{g}_t \\nonumber\n", @@ -2926,8 +3351,10 @@ }, { "cell_type": "markdown", - "id": "9ea53bf7", - "metadata": {}, + "id": "5c70db98", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\mathbf{s}_t =\\beta_2 \\mathbf{s}_{t-1} +(1-\\beta_2)\\mathbf{g}_t^2 \\nonumber\n", @@ -2936,8 +3363,10 @@ }, { "cell_type": "markdown", - "id": "7338ddf3", - "metadata": {}, + "id": "0503bdfd", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\mathbf{m}}_t={\\mathbf{m}_t \\over 1-\\beta_1^t} \\nonumber\n", @@ -2946,8 +3375,10 @@ }, { "cell_type": "markdown", - "id": "b585e12c", - "metadata": {}, + "id": "88165a6f", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\mathbf{s}}_t ={\\mathbf{s}_t \\over1-\\beta_2^t} \\nonumber\n", @@ -2956,8 +3387,10 @@ }, { "cell_type": "markdown", - "id": "d1fe61e6", - "metadata": {}, + "id": "d295a6ac", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\boldsymbol{\\theta}_{t+1}=\\boldsymbol{\\theta}_t - \\eta_t { \\boldsymbol{\\mathbf{m}}_t \\over \\sqrt{\\boldsymbol{\\mathbf{s}}_t} +\\epsilon}, \\nonumber\n", @@ -2966,8 +3399,10 @@ }, { "cell_type": "markdown", - "id": "4b64f46a", - "metadata": {}, + "id": "405ccb0a", + "metadata": { + "editable": true + }, "source": [ "\n", "
    \n", @@ -2981,8 +3416,10 @@ }, { "cell_type": "markdown", - "id": "69212d81", - "metadata": {}, + "id": "15512517", + "metadata": { + "editable": true + }, "source": [ "where $\\beta_1$ and $\\beta_2$ set the memory lifetime of the first and\n", "second moment and are typically taken to be $0.9$ and $0.99$\n", @@ -2998,8 +3435,10 @@ }, { "cell_type": "markdown", - "id": "e4ec2ca8", - "metadata": {}, + "id": "9f3fe2ce", + "metadata": { + "editable": true + }, "source": [ "$$\n", "\\Delta \\theta_{t+1}= -\\eta_t { \\boldsymbol{m}_t \\over \\sqrt{\\sigma_t^2 + m_t^2 }+\\epsilon}.\n", @@ -3008,8 +3447,10 @@ }, { "cell_type": "markdown", - "id": "13ebaeaa", - "metadata": {}, + "id": "fe96771f", + "metadata": { + "editable": true + }, "source": [ "## Practical tips\n", "\n", @@ -3026,8 +3467,10 @@ }, { "cell_type": "markdown", - "id": "cc132324", - "metadata": {}, + "id": "380ee721", + "metadata": { + "editable": true + }, "source": [ "## Automatic differentiation\n", "\n", @@ -3062,8 +3505,10 @@ }, { "cell_type": "markdown", - "id": "47ad5167", - "metadata": {}, + "id": "f6648761", + "metadata": { + "editable": true + }, "source": [ "$$\n", "f(x) = \\sin\\left(2\\pi x + x^2\\right)\n", @@ -3072,16 +3517,20 @@ }, { "cell_type": "markdown", - "id": "9ee6225f", - "metadata": {}, + "id": "df6b487a", + "metadata": { + "editable": true + }, "source": [ "which has the following derivative" ] }, { "cell_type": "markdown", - "id": "cc49c1e1", - "metadata": {}, + "id": "cd355c9d", + "metadata": { + "editable": true + }, "source": [ "$$\n", "f'(x) = \\cos\\left(2\\pi x + x^2\\right)\\left(2\\pi + 2x\\right)\n", @@ -3090,38 +3539,23 @@ }, { "cell_type": "markdown", - "id": "ea0c79f1", - "metadata": {}, + "id": "ffb3aed8", + "metadata": { + "editable": true + }, "source": [ "Using **autograd** we have" ] }, { "cell_type": "code", - "execution_count": 5, - "id": "aba1be09", - "metadata": {}, - "outputs": [ - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdwAAAEWCAYAAADM0CYnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABDoklEQVR4nO3dd3gU5drH8e+dAiEkhEBC70gnECAI6BEsWLCAYsMGVizH+loPNjzHip5jV+wFu6CIoiKoWAFJIITeWwglBAKhhLTn/WMmuMZ0kn223J/rypXdndnZ3zw7M/e0nRFjDEoppZSqXSG2AyillFLBQAuuUkop5QVacJVSSikv0IKrlFJKeYEWXKWUUsoLtOAqpZRSXmCl4IrIcSKyspaGPU5EXq+NYQc6EXlbRB4up/vDIrJTRLZ5M1cgEpHjRSTddg6bRMSIyFFHOIylInJ8Od1ni8jVRzD88SLyXnXf7y1H0pYi8o2IjDnCz79cRH6tQv8bRGRoNT+rjYjsE5HQ6ry/gmFfIiLf1fRwi1VYcN2GOSgiOSKSLSK/i8h1IlLtYm2M+cUY06W67/fI9reFljHmUWNMtWcwf1YTC7Byht0GuB3oboxpVhufUckcl7vjeWEV31fuyoTyT8aYHsaY2XDkxTEYVoJKayNjzDBjzDu2MlWVMWaTMSbKGFN4JMMRkXbusiTMY9jvG2NOOfKUpats0TzLGBMNtAUeB+4G3qjOB3qOnPIrbYAsY8yO0jp68XsdA+wCRnvp8yqk07S2gfIOv5/OjDHl/gEbgKElXjsaKAJ6us/rAk8Bm4DtwESgntvteCAdp0hvAyYVv+Z2vxuYXGL4zwLPuY+vAJYDOcA64Fr39frAQTfHPvevBTAeeM/t5xvgxhLDXgSMdB93BWbiLMBXAheU0w6NgLeADGA3MNWj2zXAGnc404AWHt0McAOw2h2H/wAdgd+BvcAnQJ0SbTUO2Om2/SUew5oNXO3x/HLgV/fxz+5n7Xfb4kL39TOBVCDb/cxeHu/vAyxwc30MfAQ8XMq4Dy3R1m8D7dzPu8r93n/GWYG7D9gI7ADeBWLcYRT3fwWw2W3D64D+QJqb74UKpsW2boZzgQKgWWltUaLtjwLGAvlAnpv/S7d7N7dNs4GlwHCP9zYGvnS/o/nAw57Dd4f9T/d7Xe8x3W5235MCHOfRfz233XYDy4A7ceeBMsa1B39Om9uBcR7z2jM402GG+7huiennLrf9twJnA6cDq9xhjfP4jPHAZPe7z3Gnhd4e3e8B1rrdlgHnlGjv34CngSy3fcpcDrjvudPNlAFcWfz9lDLuJwCLPZ7PBOZ7PP8FONtz+QSc5n6/+e53vMhjnvmPmzUH+A6IK+Uzy1uefIIzLee400mSx/taAFOATGA9cHM53+kZwEJ3+tgMjPfo1s5tjzFu++0E7i2xzJ2DM61uBV7AXW6UmNb7u20f6tFtJM5yr7w28lyuXMOfy9xlQN9KTg+/ljPul+EsF7KAe/GoKzjLjeJhZ7nt3ahEu3guZ4pfCwMuBJJLfNZtwLRKtPkmdzjF3/cg/rpMfRl4qsSwvwD+r6rf/eH3V9hDKQXXI+z17uOncQpNIyAaZ0H1mMdCoAB4AmeGrMdfC25b4AAQ7T4PdSeogR4N1hEQYIjbb1+PYaeXyDWePwvuaOA3j27dcSbYujgz2GacAhCGU3x24uwyLa0dpuMsmGKBcGCI+/qJ7vv6usN9Hvi5xIzwBdAAZyF6CPge6ADE4Ey4Y0q01f/cYQ3BKaBdypgxDk8cnjOdx/M+OAveAW67jnG/z7pAHZwZ4DZ3fM7DmRH/VnBLa2v+nOjfdduyHs5CdI07blHAZ8CkEv1PBCKAU4BcYCrQBGjpZh1SzrR4P/CH+3gxcHt5M7xne+AUu4c9uoW7Wce5bXEizoKkuK0/cv8icaabzaW09Uycab545fJSnEIdhrP7fRsQ4XZ7HKdQNAJaA0soo+DizENb3WFEuM8HuN3+Dcx12yweZyXqPyWmnwfc8bsGZ2HwgTuMHjhFpb3HvJLvfvfhwB04C45wt/v5OAuVEJwF236guUd7FwA3ueNbj/KXA6fhFIKeONPLB5RdcOu500acm2s7sMUdZj13HBqXXD7hMe97DGs2zoK8s/ve2cDjlZnGPYaZi7PSEgo8Bsz1KBQpbnvXwZnu1wGnljP8BPd9vdzxKl5xaOe2x2tuzt44y4pubvd+wEC3rdvhFMRby5jWlwHDPLp9jjuvlNNGV3t851twCrfgFPG2lZweSi24OPPPPmAwzrLnfzjTTvH3dgvONN3K7f4K8GE5y5ni18Jw5s8coJPH580HRlWhzcNKW464eTcD4j6PxZn2itug0t/94eGX17HkBF3i9bk4ayriNnxHj26D+HOt/3ictaqIchbevwKj3ccnA2vLyTMVuKWCGaS44Ea72YonmEeAN93HFwK/lHjvK8CDpXxmc5w139hSur0BTPB4HoWzEGvnMSMc69E9Bbjb4/l/gWdKLDDre3T/BLi/5IxR2kTO3wvuy7gLY4/XVuIU8sE4Wxri0e13ql5wO3i89j1wg8fzLm5bhHn039Kjexbulrj7fAoeC5FSMqwu7g78C3cNvawZnvIL7nE4BTHE47UP3ekn1M3dxaNbaVu4J1Yw7+zG3WLEmRlP8+g2lrIL7kXAwjK6rQVO93h+KrDB4zs6iLt1gzP9G9xi7TH9FS9wxuMWD/d5CE6hP66Mz04FRni09yaPbhUtB97Eo9DhFMBSC67b/RecLbOBOFuln+AU7ROANI/+NlBxwb3P4/kNwLeVmcY9hjnL43l34KD7eIBnG3hMl2+VN1149PsM8HSJ+amVR/c/cAtHKe+9Ffi8jGn9buB993EjnI2U5h7jU17BnYG7fK1E/pLTQ1kF9wHgI4/n9XFqQvH3thw4yaN7c/6+3PBczhS/FuY+fw94wH3cCacAR1ahzcsquIKzYTnYfX4N8MORfPdHcpZyS5xdVPE4axkp7klV2cC37uvFMo0xueUM6wOchQzAxe5zAERkmIjMFZFd7rBPx1nzrZAxJgdny3SU+9JFwPvu47bAgOLM7rAvAUo7Iag1sMsYs7uUbi1wthSLP3MfTiFp6dHPdo/HB0t5HuXxfLcxZr/H843uZ1RHW+D2EuPY2h1eC2CLcacUj8+qqs0ej//SFu7jMKCpx2tVaYvDRORYoD3OVic400iCiCRWI3Nx1s3GmKISeVviTLth/HXcPB+X+pqI3CEiy0Vkj9vWMfw5rbYo0X95bd0ap7CWlbtkG3tOH1nmz5NJDrr/y2vjw5nctkgvHp6IjBaRVI9ppyd/nfc8x6ei5UBVxh/gJ5wCONh9PBtnRXGI+7wqPM+qP0AZ01gV3h/hHktsC7QoMX+N46/T+2EiMkBEfhSRTBHZg3NIpeSyrNSsItJZRL4SkW0ishd4tJT3FnsPOEtE6gMX4GxYbK3kuJY57VVieijLX757d/mW5dG9LfC5x3CXA4X8tR1Lm/+KlawfU40xB9zMlWnzUrnLxo9KDNuzflT6uy9WrYIrIv1xFky/4uxOPQj0MMY0dP9ijDGeE7UpbTgePgWOF5FWwDm4BVdE6uJs9TwFNDXGNAS+xlnzqMxwwdlquUhEBuHsnvvRfX0z8JNH5obGOfPt+lKGsRloJCINS+mWgdP4uJnr4+xW3FKJbKWJdYdRrI37GeBsQUR6dKvobOHNwCMlxjHSGPMhzpZMSxERj/7bVCOv53fwl7Zwh1fAXxf41TUG53tPdX+WNM/jdSjRNiJSsm1KTisZQOsSZ9u3wfneMt3crTy6tS4l0+FhishxOMdOL8DZE9IQ2MOf0+rWEsMor6034+yiKk1pbZxRRr+VcTiT2xatgAwRaYuze/NGnN23DXF2g3tOL55tWtFyoCrjD38vuD9RccGtzPKgPFV9/2acLXjP+SvaGHN6Gf1/gLPLvbUxJgbn8IqU0W9JLwMrcHadNsBZuJf6XmPMFpzjvSNxjp1O8uxciXHqWPLFSk4PZfnLdy8ikTjLSM/PHFaiHSPc8ahM7plAvLvyfREeG2yU3+aVrR/nueM/AKceFWeuyncPVLHgikgDETkTp+q/Z4xZ7K4VvwY8LSJN3P5aisiplR2uMSYTZw32LXcklrud6uDs088ECkRkGM6xv2LbgcYiElPO4L/GWUD9G/jYY4vmK6CziFwmIuHuX38R6VZKvq04J2C9JCKxbr+D3c4fAleISKK7gvAoMM8Ys6Gy41+Kh0SkjrsQPxNnhQScXTgjRSTS/fnPVSXet52/LqhfA65z1/JEROqLyBkiEo0zQxYAN7vjMxLnxIwj8SFwm4i0F5EonLb42BhTcCQDFZEInEI2Fkj0+LsJuNjd2lgE9HC/hwicXWeeSrbNPJwtiLvc8T8eOAtn11chzvHn8W5bd6Xis6KjcdozEwgTkQdwjtsX+wT4lzv9tHKzl+UroLmI3CoidUUkWkQGuN0+BO4TkXgRicPZXXckvxPtJyIj3Ta8Fee44Vyc3X7GHR9E5AqcLZpSVWI58AlwuYh0dxe4D1aQ63ecQxJH4xy3X4q7VwrnxJnSbAfaSfV/sliZ5YmnP4AcEblbROqJSKiI9HQ3SEoTjbOnLFdEjsbZYqqsaJwTf/a502NpGwae3sVZAUzAmZaLVdRGrwN3iEg/d5lxlFtsqjQ9lDAZOFNE/iEidXCWxZ6fPxF4xP0c3Gl7RCWHjTEmH2cZ+STOLvSZHp3La/NMnEOFZa3cYoxZiLMy+TowwxiT7Xaq6ncPVL7gfikiOThV/V6cg95XeHS/G+cElLni7O6YhTOzVMUHOGcbHl47cXcJ34wzs+7GaaxpHt1X4CyA1rmb9X/b9WqMOYQzwZU27FNwdjdn4OzKKT6xqzSX4RxXWIFzcs+t7nBm4ZzMMwVnTa4jf+7Cro5tOOOagbP74jp3PME5KSUPZ6Z5hz93bxQbD7zjtsUFxphknOMOL7jDXINzjAJjTB7OGvDlOIcGLuSvM2Z1vImzNv0zzsk3uZRfWCrrbJytp3eNMduK/9zPC8M5NroKZ0aehXOst+SP8N8AurttM9Ud/7OAYTgz1Es45xEUt/WNOLuEi8+s/xCnGJVlBs4u1FU4u0tz+etusIfc19fjHJOcVHIAxdxp82Q33zZ3fE5wOz8MJOOc2b0Y58ziI/l98Rc43/1unGl8pDEm3xizDOf8gjk401sCzpm+5SlzOWCM+Qbn+NkPbj8/lDcgd7fjAmCp+13hZtloyvhpGn+umGaJyIIKspb2mRUuT0r0X4izQpyI870WL5jLKtg3AP92l6UP4CzXKusOnOVfDs6KzccV9P857q7a4t2rrnLbyBjzKc65Lh+4nzUV54zh6kwPxcNcinNG/wc4y8jdOIcuij2Ls1z/zm2buTgrVlVRXD8+LbGCX2abu+3yCPCb+30PrGDYnvWjqt898OfZV8oHuFtZ7xljWlXQq/IyEXkC52dIYyrs2U+IyHicE20utZ1F1TwRWYvzM8pZtrMoh15LWalSiEhXEenl7lY7Gmf3/ee2cylVGSJyLs4u4HL3JCjv8u+rdihVe6Jxdi+2wNmF9l+c3a9K+TQRmY3z86XLSpyFryzTXcpKKaWUF+guZaWUUsoLdJdyOeLi4ky7du1sx1BKKb+SkpKy0xgTX3GfwUULbjnatWtHcnKy7RhKKeVXRKQ6V60LeLpLWSmllPICLbhKKaWUF2jBVUoppbxAC65SSinlBVpwlVJKKS/QgquUUkp5QVAVXBG5TUSWisgSEfnQvY2bUkopVeuC5ne4ItIS51Z/3Y0xB0XkE5zb6L1tNZhSKugU5OexccUCdq1bQGHODkzeAaRuFGENmxPfsR+tjupFaFjQLJ6DRrB9o2FAPRHJByJx7jmrlFK1zhQVsfinz8hb+CFd9/xGRzlIx9J6nA9ZxLCm8fHEHnslnfse7+WkqrYETcE1xmwRkaeATTg3M//OGPNdyf5EZCwwFqBNmzbeDamUCjimqIgF375Do+Sn6VW0kWyiWNboREI6DCG+8wBim7Ymsn40+/ftJWvLWnau/oPQtbNI2PktkdO+YOm3CYSe8hBdk06yPSrqCAXN3YJEJBaYAlwIZAOfApONMe+V9Z6kpCSjl3ZUSlXXplWp7J1yCz0PpbIhpA2Zva+n92lXUqduxaeP7Nu7myVfPkfH1W8Rz27+aHg6XcY8T0xsnBeSHxkRSTHGJNnO4WuC6aSpocB6Y0ymMSYf+Aw4xnImpVQAMkVF/DHlaZq8P5Q2h1Yzr9s4Wo9bSP8RN1Sq2AJENYhl4CUPEnl7KnOaX0qf3TM48OxAVsyfVcvpVW0JpoK7CRgoIpEiIsBJwHLLmZRSASb34H6Snx3F0YvHszqiJ3nXzWXAhXdX+ySo+tENGXTti6w7azJGhI5fXcD8qS/WcGrlDUFTcI0x84DJwAJgMc64v2o1lFIqoOzZlcm6p0+l/54ZzGlzLd3vnEVcs5o5F6RL0onUv3kOKyMS6J86jjlv3oUpKqqRYSvvCJpjuNWhx3CVUpW1c9tm9r16Oi0Kt5B29BMknXFNrXxO3qFcFr00hv57vmVOi9EMvPpZJMS3tp30GG7pfOtbUkopP7Q7cys5r55Bk8LtrBr6Vq0VW4A6dSNIuuVD5jUewaCMd5n79j219lmqZmnBVUqpI7Bn906yJp5Bi8IM1g19jZ7Hjaj1z5SQEPrf8BbzGw5j0KZXmPfx47X+merIacFVSqlqyjuUS/rL59CmYAMrjn/ZK8W2WEhoKH1vfI+FkceQtOxx0n6c7LXPVtWjBVcpparBFBWR+vIV9MhLI63fo/Q+4XyvZwgNC6Pz9R+yIaw97WffyMblKV7PoCpPC65SSlXDvA8e4ujsr5nb6iqShl9nLUf96IZEXTGZQ1IXPh3D/pxsa1lU+bTgKqVUFS35dRpHr36WBVFDOPqKJ23HoWmrjmwb+jytC9NZ9vo1+nMhH6UFVymlqmDnts00m3UTm0Nb0uXadwkJDbUdCYCe/xjOvLbX0H/Pd8z/4gXbcVQptOAqpVQlFRUWsvWt0USZ/RSd+yb1oxvajvQXR49+jKV1Euie+ijbNq22HUeVoAVXKaUq6Y9J95NwaAFpCeNo32OA7Th/ExoWRsNRrxFCEZnv665lX6MFVymlKmFt2u/0Wz+RlOgT6D/yVttxytSyQzcW97iDhEML+WPK/2zHUR604CqlVAXyDuXCFzewR6LpOGaiz11KsaT+597OkrqJ9FzyJNvT19qOo1y+PdUopZQPSHnvPjoWrmfzMY/QMK6Z7TgVCgkNJfbCiYRSSPpH/2c7jnJpwVVKqXKsTfudpE1vktxgKH1OudR2nEpr2aEbC9teSb99s1n88xe24yi04CqlVJkKCwoonHYLeySaTmNesh2nyvpc9CDp0oyGP97DodwDtuMEPS24SilVhuTPn6FzwSo29P0XMY2b2o5TZRH16rNr8CO0Nhks/OQx23GCXpjtAIFo3qf/JWr1VMKLcgkvOkSBhJMfWo9D4TEcimqNNGpPw479addzIHUjIm3HVUqVYteOLXRd+j+W1ulFvzPH2o5Tbb1OOI/UP16j+7rXyd75T784Bh2otODWhsI8QkwhuWEN2B9SlxCTT1jhQRod3EjTffOI2J4PyyHvyzCW1u1GTtuTaT3oPFp26GE7uVLKteaDO+hjcok65xmfPyu5IrHDH6X+Rycx/5P7GXjDa7bjBC0xxtjO4DUi0hB4HegJGOBKY8ycsvpPSkoyycnJNZrBFBWxfcs6ti77jUPr59F0x2+0L9oAwLLwnuzvcTEJp4whIjKqRj9XKVV5K/6YSdevz2NO80sZdO2LtuPUiD+eu5TErK/JHPNLra/ci0iKMSapVj/EDwVbwX0H+MUY87qI1AEijTHZZfVfGwW3NBnrV7Dxl/dptX4yrU0GO2nImqOuIGHErT536TilAl1RYSFrHhtITEEWUbcvCJh5cGfGRiJf6c+K6IH0vWNarX6WFtzS+fd+kioQkRhgMPAGgDEmr7xi600t2ndl0Oj/0Or+pSwZOoltddsxcM3T5P03gXkfP05Bfp7tiEoFjQXfvEHnglVs6nNHwBRbgLgWbVnUZjR99/3EqgU/2Y4TlIJmC1dEEoFXgWVAbyAFuMUYs79Ef2OBsQBt2rTpt3HjRi8ndayYP4uCWf+h56FU1oe048DQx+hxzOlWsigVLHIP7mf3E73ZH9qADuPm+8ydgGpKzp5dFD6dwMZ63el998xa+xzdwi1d0Gzh4pwg1hd42RjTB9gP3FOyJ2PMq8aYJGNMUnx8vLczHta1/1B63P0jCwc9R0TRfnp8dxHzXriCA/v2WMukVKBL/fRxmpNJ7gnjA67YAkTHNGJ5+zH0PvgHK5N/sB0n6ARTwU0H0o0x89znk3EKsM+SkBD6nDqG2LtSmdv0Ivpnfs6u/w5gxfxZtqMpFXB2Z26lx9rXWFRvAD3/Mdx2nFqTcM6d7CaaQ7MetR0l6ARNwTXGbAM2i0gX96WTcHYv+7yIyCgGXj+R5ad+QAiFdPzqAuZ++IjeekupGrTqk/uIJJeGwwP7AhFRDWJZ0eEKeuXO15V3Lwuaguu6CXhfRNKARMCvVvF6HHM6UbfMZUn9AQxcOYGUp8/TXcxK1YCMDSvps+NzkhufRdtu/WzHqXW9zrmd3TQg73u/WgT6vaAquMaYVPf4bC9jzNnGmN22M1VVg4aN6X37V8xt90/67v2BjKdPYOe2TbZjKeXX0qeOxxBC+3P/bTuKV9SPbsiKDpfTKzeF1am/2I4TNIKq4AaKkNBQBl7+KIuHTKRFQTr5E09k44oFtmMp5Zc2r1lM393fsrDpOTRp2d52HK/pMfxWckw99s56ynaUoKEF14/1PnEUW86eTDj5xH50ph6PUaoatk8bTz5hHDXyfttRvKpBw8YsaXEuiTk/sWXdUttxgoIWXD/Xqc9g8i7/jr3SgNZfXcLyeTNsR1LKb2xcnkLfPd+zqPn5xDVrYzuO1x111p0UEkr69Am2owQFLbgBoEW7LtS5+huyQhvT9uvLWPrbdNuRlPILWdMf4iB16XJucG3dFotv0Y7URqfSe+d0sran244T8LTgBogmLdsTOfZbdoQ2ocN3l2vRVaoC65fOo+++n0hrdRGx8c1tx7Gm2bC7qEMBq6bpsdzapgU3gMQ1a0OD62awPbQpbb+7Ss8+VKocu759nP0mgu4j/2U7ilVtOieyKOpYemz5RH9mWMu04AaYRk1aEnnVl+RIFHFTL2bTqlTbkZTyOZvXLCZx74+ktTiPmMZNbcexLmLwzTRgP4u/0Xvl1iYtuAGoScv2FFz6OUUIdT44l+3pa21HUsqnbJ3+GPmE0WnE3y6nHpS69j+ZNaEdabr8bb2CXS3SghugWh+VQPbIj4gy+9n31nnsz8m2HUkpn7Bt8xoSd33LoviziGvW2nYcnyAhIezqeQXtijaz9LcvbccJWFpwA1jHXsew7vgXaFewnlUvX0RhQYHtSEpZt+HLJxCgzVm6deup12lXsosG5M+ZaDtKwNKCG+B6nXAeyd3ups+B35n/+s224yhl1a4dW+i9fSqpsafQvG2Xit8QRCLq1Wdlq/PovX8OW9Yttx0nIGnBDQIDRv2LeXHnMnDb+/wx5RnbcZSyZuUXE6hLPk1P163b0nQcdjOFhLB5xrO2owQkLbhBot+1E0mLSCIx7T+sWvCT7ThKed3e7Cx6pn9MavRg2nROtB3HJzVp2Z60BkPovv0LPe+jFmjBDRJh4XVoe80HZEksDaZdye7MrbYjKeVVy758lmg5SIOT77IdxafVH/xPGnCApbMm2Y4ScLTgBpGYxk3Zd/ZbxJo9pL9+sZ5EpYJGft4h2q99j6V1enNU73/YjuPTuvQ7kWWnfkTS8BtsRwk4WnCDTKfE41jU614SDi1g/tu6pq+Cw6Jv36IpWRQM1CJSEQkJofugYYSEhtqOEnCCruCKSKiILBSRr2xnsaX/ObfwR8PTGZj+Bot//tx2HKVqlSkqouGiV9kY0oqEIefbjqOCWNAVXOAWIKjPeZeQEBKueZWNIa1p9sNt7NqxxXYkpWrN0jnTOapwLdu7X6VbbcqqoCq4ItIKOAN43XYW2+rVj6bgnNeJMTlsevsqvZybCliFvz7PLhrQ64xrbUdRQS6oCi7wDHAXoNUF6JgwkAVdbiPxwBz++PRJ23GUqnEbVyyg98F5rGwzioh69W3HUUEuaAquiJwJ7DDGpFTQ31gRSRaR5MzMTC+ls2fAqHEsiuhP72VPsn7ZfNtxlKpR27/7H7kmnC5n3mo7ilLBU3CBY4HhIrIB+Ag4UUTeK9mTMeZVY0ySMSYpPj7e2xm9TkJCaHn5W+yXSMzkq8k7lGs7klI1Imt7Or2zvmVR42E0atLSdhylgqfgGmP+ZYxpZYxpB4wCfjDGXGo5lk+Ia9aazcc+ToeiDaS8d6/tOErViFXTn6Ou5NPs1NttR1EKCKKCq8qXePLFJDc4maRNb7Fm0W+24yh1RPIO5dJp08ekRfSnbZdE23GUAoK04BpjZhtjzrSdw9d0GvMieySakC9u0F3Lyq+lzZpEHNmYo8fajqLUYUFZcFXpYho3Jf0fzq7lBZPG2Y6jVLVFpb5JujQnYci5tqModZgWXPUXiUMvYn7MKSRtfos1i361HUepKluz6Fe65i8jvdMleqEL5VO04Kq/6TzmJXZLDDLtRgry82zHUapKdv/4AgdMXbqfrtdNVr5FC676m5hG8aQPeoiOhetJ/vhR23GUqrRdO7bQa/csFscNo0HDxrbjKPUXWnBVqRJPvoyFkcfQa/VLZGxYaTuOUpWy8psXnZ8CDb3ZdhSl/kYLriqVhITQfNRzAGR+fKNea1n5vIL8PNqv/5gldRNp262f7ThK/Y0WXFWmZm06kdb5Rnof/IOFM96xHUepci3+4UOasZP8ftfYjqJUqbTgqnIlXXAPa0I70nreQ+zNzrIdR6ky1Ul5nQxpQq8TR9mOolSptOCqcoWF18Gc9SyNTDbL37vDdhylSrV+6Tx65KWxqcNFhIaF2Y6jVKm04KoKdUo8juT4kSRlfs7axXNtx1Hqb3Z8/wK5Jpyuw/SnQMp3acFVldL14ifYK9Ecmna7nkClfMqe3TtJyJpBWqNTaBjXzHYcpcqkBVdVSkyjeFb3vI3u+UtImf6a7ThKHbZixmtEyiFih1xvO4pS5dKCqyqt39k3szr0KNqkPM7+nGzbcZTCFBXRdNUHrArrTKfE42zHUapcWnBVpYWGhVE4bAJN2EXaB/fbjqMUK+fPol3RJrK7XWI7ilIV0oKrqqRr0knMjzmNfhnvs3n1IttxVJDb99ur5Jh69Dz1CttRlKqQFlxVZe0vepI86rBrip5ApezJ3rmNhD2zWRZ/OpFRMbbjKFUhLbiqyuKatWFJ5+vpnTuftNmf2I6jgtSKGa9QV/JpcqKeLKX8Q9AUXBFpLSI/isgyEVkqIrfYzuTP+p53N5ulBQ1/fVhv4ae8zhQV0XLNRywP70777v1tx1GqUoKm4AIFwO3GmO7AQOCfItLdcia/VaduBDsH3Uvbos2kfP6M7TgqyCz9/Stamwz29bzMdhSlKi1oCq4xZqsxZoH7OAdYDrS0m8q/JQ69mKV1Eui87Hm9zrLyqkNz3yCbKBJOHm07ilKVFjQF15OItAP6APMsR/FrEhJCndMfI5a9LP14vO04Kkjs3LaZXjm/sKLpWURERtmOo1SlBV3BFZEoYApwqzFmbyndx4pIsogkZ2Zmej+gn+mUeBzzY06lX8YHeqN65RWrZ7xMuBTS4qTrbEdRqkqCquCKSDhOsX3fGPNZaf0YY141xiQZY5Li4+O9G9BPtTn/MQoJYeuUe2xHUQGuqLCQths+ZWmd3rTpnGg7jlJVEjQFV0QEeANYboz5n+08gaRpq46ktr6Ufjk/sCL5e9txVABb8vPntDA7yE0cYzuKUlUWNAUXOBa4DDhRRFLdv9NthwoUCRc8wE4awrf36sUwVK0pnP8GWcSQcJJeylH5n6ApuMaYX40xYozpZYxJdP++tp0rUEQ1iGVdz1voWrCc1JmTbMdRAWh7+loS9s9lVYsR1KkbYTuOUlUWNAVX1b6+I25kY0grGs97Qi+GoWrcuu8mEoKhzVC9ybzyT1pwVY0JC69D1oB7aFO0hQXTXrQdRwWQgvw8Om6azJJ6/WjZoZvtOEpVixZcVaP6nHwJK8O60m7xc+Qe2Gc7jgoQS2ZPpgm7KOyrdwVS/ksLrqpREhJCwYkP0oRdpE5+wnYcFSAk5U120IiEEy6wHUWpatOCq2pcj2NOZ1FEf7qve4M9u/TiIerIZGxYScLBZNa2HklYeB3bcZSqNi24qlZEnfEwUeYAyz79t+0oys9tnPkSBmh/it6GT/k3LbiqVnRMGMiCmJPok/EhO7astx1H+am8Q7l02jKVxfUH0qz1UbbjKHVEtOCqWtPinEcIoYgNU+63HUX5qcXff0gc2UjSlbajKHXEtOCqWtOifVcWNBlJ36zpbFqVajuO8kN1F73NVuLpOXik7ShKHTEtuKpWdTpvPIeow84vdCtXVc3m1YvoeSiVDe3OIzQszHYcpY6YFlxVqxo3bUVam8vou/9nVqf+YjuO8iNbZr1Mvgmlk54spQKEFlxV63qeN45sojgwQ89YVpWTe3A/XbZ/xeLoY4lr0dZ2HKVqhBZcVeuiYxqxosOV9D74Byv+mGk7jvIDS2ZOIpYcwo++2nYUpWqMFlzlFb1G3sFOGlI46z+2oyg/UH/xu6RLc3oce6btKErVGJ8suCJyk4jE2s6hak5kVAxruoylR94ilvw6zXYc5cM2LE+mW/5S0jtcQEhoqO04StUYnyy4QFNgvoh8IiKniYjYDqSOXOLZt7KdxoTNflRvUq/KtP2Hl8kzYXQ5TU+WUoHFJwuuMeY+oBPwBnA5sFpEHhWRjlaDqSMSUa8+G3veSNeC5aTN/sR2HOWDDu7PoVvm16TFHE9sfHPbcZSqUT5ZcAGMMQbY5v4VALHAZBGZUN1hulvLK0VkjYjcU0NRVRX0Gf5PtkhT6v82gaLCQttxlI9Z/N1bNOAAkcdcYzuKUjXOJwuuiNwiIinABOA3IMEYcz3QDzi3msMMBV4EhgHdgYtEpHsNRVaVFF6nLhmJt3JU4VpSZ06yHUf5mIZLJ7EhpDXdjj7FdhSlapxPFlygETDSGHOqMeZTY0w+gDGmCKjuaYtHA2uMMeuMMXnAR8CImomrqqLvGWPZGNKaxn88RWFBge04ykesWfQrnQtWsa3TRUiIry6alKo+n5yqjTEPGmM2ltFteTUH2xLY7PE83X3tL0RkrIgki0hyZqbey7U2hIaFkdX/dtoWbWbh16/ZjqN8RNZPr3DQ1KHbqWNtR1GqVvhkwbXJGPOqMSbJGJMUHx9vO07ASjxlNGtDO9Bs4TPk5x2yHUdZlrNnFwlZM1gcO5SYRjrfqcAUTAV3C9Da43kr9zVlQUhoKDnH3E0rs42F0160HUdZtmzGG0TKIRoOvtZ2FKVqTTAV3PlAJxFpLyJ1gFGAXoHBot4nXMDKsK60XfICuQf3246jLDFFRcStfJ+1oR3olDjYdhylak3QFFxjTAFwIzADWA58YoxZajdVcJOQEPKHjKMpWaROfcZ2HGXJygU/0rFwPTu7XqInS6mAFlRTtzHma2NMZ2NMR2PMI7bzKOhx7FksrdOLo1a+ysH9ObbjKAtyfnmV/SaCHqdeZTuKUrUqqAqu8j0SEkLISfcRRzaLPnvSdhzlZXt2ZZKQ/T1L4k4jqoFePl0FNi24yrpuA04lLSKJLmvfZN/e3bbjKC9a/u0rREg+cUP0ZCkV+LTgKp8QccoDxJLD4imP246ivMQUFdF8zYesDOtCx17H2I6jVK3Tgqt8Que+Q1gYeQw9Nr7Lnl16wZFgsGzut7QtSmdPj8tsR1HKK7TgKp8Rc/p4GnCAZZ89ajuK8oKDc15nL5EknHK57ShKeYUWXOUzOvQcQErU8fTa/AG7dug1SQJZ1vZ0eu2dzbL4M6hXP9p2HKW8Qguu8ilxZz5IBIdY9Zn+aiuQrZ7xCnWkkOYn6k3mVfDQgqt8StuufVnQ8BQSt37CzoxS71+h/FxhQQGt13/MsvCetO3Wz3YcpbxGC67yOS1GPEgYhaz9/N+2o6hasOTnybQ02zmYeIXtKEp5lRZc5XNadujBgsZn0GfHVLZtWm07jqppf7xOJrEkDL3UdhKlvEoLrvJJbc5+EIBNUx+ynETVpPQ1S+idO581bc6nTt0I23GU8iotuMonNWvTiYXxI+iT9TVb1uk9JgJF+swXyDehdDrtRttRlPI6LbjKZ3Uc+SCFhJDxhW7lBoKD+3Povn0aaQ0GE9eire04SnmdFlzls+JatCW1+QX0zf6OjSsW2I6jjtDib1+nAfupd+x1tqMoZYUWXOXTOo+8l1zqsvMr3cr1Z6aoiMbL3mFdSDu6HX2K7ThKWaEFV/m0Rk1aktbqIvrtm826JfNsx1HVtGL+TDoWriez22i9ybwKWjrlK5/X/dx72Uske74ebzuKqqYDv050rps87GrbUZSyJigKrog8KSIrRCRNRD4XkYa2M6nKi2kUz9K2o+lz4HdWLfjJdhxVRTszNtJr708sa3IWkVExtuMoZU1QFFxgJtDTGNMLWAX8y3IeVUUJ597DbqLJ/U6vPuVvVn/7IuFSSMtTbrIdRSmrgqLgGmO+M8YUuE/nAq1s5lFVF9UglpUdr6RXbjLL582wHUdVUn7eIY7a9AlpEf1pfVSC7ThKWRUUBbeEK4FvyuooImNFJFlEkjMz9UbovqT3yDvZSUOKvv8PpqjIdhxVCWkzJxHPbor667FbpQKm4IrILBFZUsrfCI9+7gUKgPfLGo4x5lVjTJIxJik+Pt4b0VUl1asfzZouY+mRt5ilv31pO46qhKiFr5Iuzel1/Pm2oyhlXcAUXGPMUGNMz1L+vgAQkcuBM4FLjDHGalhVbYln38p2GhP202O6levjVsyfRZeClWzpMoaQ0FDbcZSyLmAKbnlE5DTgLmC4MeaA7Tyq+iLq1WdDj3/StWA5abMn246jyrH/p+fZSyQ9z9CbzCsFQVJwgReAaGCmiKSKyETbgVT19R1xI1ukKfV/e1y3cn3Utk2r6Z3zM8uanUP96Ia24yjlE4Ki4BpjjjLGtDbGJLp/ejFXPxZepy4ZvW/mqMK1pM6cZDuOKsWGr58GoN3pt1lOopTvCIqCqwJPnzPGsimkJbHznqKwoKDiNyiv2Z+TTfdtn7MoejDN2nSyHUcpn6EFV/mlsPA6ZPb7P9oVbWLht2/ajqM8LJn+Mg04QP3jb7YdRSmfogVX+a0+p13B+pB2NE15moL8PNtxFFBUWEjLle+wMqwLXZNOsh1HKZ+iBVf5rZDQULIH3klrk8GCaS/ZjqOAtB8+ppXZyr4+Y21HUcrnaMFVfi1x6MWsDOtCu8XPcnB/ju04QS98/stsI47ep4y2HUUpn6MFV/k1CQmh8KSHaMIuUic/bjtOUFu14Cd65KWx4ahLCQuvYzuOUj5HC67ye90HDSO13kB6rH+T7J3bbMcJWvu+f4q9RNLjrFtsR1HKJ2nBVQGh4fBHqG8OsuLTB21HCUqbVy8icd8vLG1xPtExjWzHUconacFVAaFdtyRSYofRd9tkMjastB0n6Gz9egL5hNFp+J22oyjls7TgqoDR9ryHKULI+Gyc7ShBZWfGRhJ3fUtq3BnENWttO45SPksLrgoYTVt1ZGGLi0jaO4s1i36zHSdorP5yAqEU0uqMu2xHUcqnacFVAaX7BQ+STRQHvrnfdpSgsDc7i54ZU0iNHkLLDj1sx1HKp2nBVQElJjaOFZ3G0is3hcU/f2E7TsBbOu0ZouUgMSfrsVulKqIFVwWcPufeyVbiqffTQxQVFtqOE7AO7s+h07p3WVy3D0f1/oftOEr5PC24KuDUjYhkS9/bOapwLSlf6q2Pa8uiL54ljmxCj7/bdhSl/IIWXBWQ+p4xllVhnWmb+hT7c7Jtxwk4uQf303HV6yytk0D3QcNsx1HKL2jBVQEpJDSUolMfowm7SPv4IdtxAs6iac8Tz27MEN26VaqygqrgisjtImJEJM52FlX7uvYfSnL0SfTZPImtG/ViGDXlUO4B2i1/leXhPegx6AzbcZTyG0FTcEWkNXAKsMl2FuU9rS6YgEHImHyP7SgBI3XaizQli4Lj7kRCgmYRotQRC6a55WngLsDYDqK8p1nro0htM5p+OT+wYt53tuP4vbxDubRd9gorwrrR8x8jbMdRyq8ERcEVkRHAFmPMokr0O1ZEkkUkOTMz0wvpVG3rfeED7KARoTPH6c+EjtDCaS/SjEzyjr1Dt26VqqKAmWNEZJaILCnlbwQwDnigMsMxxrxqjEkyxiTFx8fXbmjlFZFRMWzqexedClaT8uXLtuP4rdwD+2i/9AVWhnUlYchI23GU8jsBU3CNMUONMT1L/gHrgPbAIhHZALQCFohIM5t5lXf1PWMsK8O60j71SfZmZ9mO45dSP3uKJuyi4IQHdOtWqWoI+LnGGLPYGNPEGNPOGNMOSAf6GmP0TuVBJCQ0lJAznyLW7GH5+3qR/arK2bOLLmteJy2iHz2O1TOTlaqOgC+4ShXrlHgcyfHnkLRjit5NqIqWTH6EWHKod5r+plmp6gq6gutu6e60nUPZ0fWSJ8mWBhR+eZueQFVJu3Zsodem91gQNZhOicfZjqOU3wq6gquCW0xsHOv63EOXgpUkf/6s7Th+YdXkh4jgEI3P+o/tKEr5NS24KugknXUdy+ok0HnJf9mdudV2HJ+WvmYJfbdPJqXR6bTtkmg7jlJ+TQuuCjoSEkK9s5+mvjnI6g9utx3Hp2V+dhcFhNHh/MdsR1HK72nBVUGpfff+pDQfxdG7p7P0t+m24/ikJb99SZ8Dv5HW/mriWrS1HUcpv6cFVwWt3pc9wRZpSsys/+Pg/hzbcXxKYUEB9X54gG3Ek3jBONtxlAoIWnBV0KpXP5rdQ/9LK7ONRZP0NnOeFkx7kY6F60jvfzcRkVG24ygVELTgqqDW89izmNdoOP23fsCqBT/ZjuMT9u3dTfu0/7EirBv9hl1lO45SAUMLrgp63UY/Q5bEEv7VTeQdyrUdx7ol791NI7OHkNMf10s4KlWDwmwH8Df5+fmkp6eTm6sL5uqIiIigVatWhIeH245yWIOGjVl33GMk/nItc967l0FX/dd2JGvWLPqN/ts/YX7cCAb0Pd52HKUCihbcKkpPTyc6Opp27dohIrbj+BVjDFlZWaSnp9O+fXvbcf4i8aRRzE+bQv9Nb7Iy+Qy6JJ1oO5LXFRYUUPTlrWRLA7pe8pTtOEoFHN1fVEW5ubk0btxYi201iAiNGzf22b0DXa6YyE5pTP3pN7A/J9t2HK9L/uxpOhesYn3fccQ00ltTKlXTtOBWgxbb6vPltmvQsDFZpzxHi6JtLHnrJttxvGrntk10W/Y/ltRNpN+ZY23HUSogacFVykOPY05nXotLGLBrGqmzPrQdxytMURGb372WuiafBuc+qydKKVVLdM4KAlOnTmXZsmVe/9zx48fz1FP+dyyw75gnWRvanja/3kVmxgbbcWpdylev0ufA7yzsdCNtOifajqNUwNKCGwRqsuAWFBTUyHB8Wd2ISMLOf5MIc4jMty6hID/PdqRaszNjI50W/JsV4d3pP+o+23GUCmh6lvIReOjLpSzL2Fujw+zeogEPntWjwv7OPvtsNm/eTG5uLrfccgtjx44lKiqKffv2ATB58mS++uorxo4dy7Rp0/jpp594+OGHmTJlCjk5OVx33XUcOHCAjh078uabbxIbG8v8+fO56qqrCAkJ4eSTT+abb75hyZIlvP3223z22Wfs27ePwsJCpk+fzogRI9i9ezf5+fk8/PDDjBgxAoBHHnmEd955hyZNmtC6dWv69etXo+3jLW279iU56T8kpdzFnDdvY9C1L9qOVONMURHp711LV5NH/QsmEhqmiwOlapNu4fqpN998k5SUFJKTk3nuuefIysoqtb9jjjmG4cOH8+STT5KamkrHjh0ZPXo0TzzxBGlpaSQkJPDQQw8BcMUVV/DKK6+QmppKaGjoX4azYMECJk+ezE8//URERASff/45CxYs4Mcff+T222/HGENKSgofffQRqampfP3118yfP7/W26E2JZ11LfMan82gre+ROvMD23Fq3B+fTiDxwBxSu9xK6069bcdRKuAFzSqtiNwE/BMoBKYbY+460mFWZku0tjz33HN8/vnnAGzevJnVq1dX6n179uwhOzubIUOGADBmzBjOP/98srOzycnJYdCgQQBcfPHFfPXVV4ffd/LJJ9OoUSPA+T3tuHHj+PnnnwkJCWHLli1s376dX375hXPOOYfIyEgAhg8fXmPja0vvq19izVOL6fDbHWzp2IeWHbrZjlQj1qb9Tp9lT7IocgADRunNCZTyhqDYwhWRE4ARQG9jTA/A/87k8TB79mxmzZrFnDlzWLRoEX369CE3N/cvP7mp6d+61q9f//Dj999/n8zMTFJSUkhNTaVp06Y++9vaIxVRrz6Rl36AAfLfO5+92aXvSfAn+3OyCZ96NdnSgNZXvKVnJSvlJcEyp10PPG6MOQRgjNlhOc8R2bNnD7GxsURGRrJixQrmzp0LQNOmTVm+fDlFRUWHt34BoqOjyclxbj8XExNDbGwsv/zyCwCTJk1iyJAhNGzYkOjoaObNmwfARx99VO7nN2nShPDwcH788Uc2btwIwODBg5k6dSoHDx4kJyeHL7/8slbG39tatO/K5pNfoWVhBhsmXujXJ1GZoiKWv3Y1rQoz2HHyczRq0tJ2JKWCRrAU3M7AcSIyT0R+EpH+ZfUoImNFJFlEkjMzM70YsfJOO+00CgoK6NatG/fccw8DBw4E4PHHH+fMM8/kmGOOoXnz5of7HzVqFE8++SR9+vRh7dq1vPPOO9x555306tWL1NRUHnjgAQDeeOMNrrnmGhITE9m/fz8xMTGlfv4ll1xCcnIyCQkJvPvuu3Tt2hWAvn37cuGFF9K7d2+GDRtG//5lNrPf6XnsWSxMuI9eufNJfu2ftuNU27wP/k3S3pnMa3ctPY89y3YcpYKKGGNsZ6gRIjILaFZKp3uBR4AfgZuB/sDHQAdTwcgnJSWZ5OTkv7y2fPlyunULjON4Je3bt4+oKOfep48//jhbt27l2WefrfHP8ec2nPvSWAbu+Ji5Xe5i4EX32o5TJWmzp9Djx6tYFHUcif83lZASJ8YpVVNEJMUYk2Q7h68JmJOmjDFDy+omItcDn7kF9g8RKQLiAN/chLVk+vTpPPbYYxQUFNC2bVvefvtt25F8Tv+xL7Hwf+kMXDmB5GmNSBp+ve1IlbJ59SLazb6JjaFt6XLdJC22SlkQMAW3AlOBE4AfRaQzUAfYaTWRD7rwwgu58MILbcfwaaFhYXS76VOWPj2MxJRxpNaPJfGkUbZjlSszYwNh759LAaFEjP6E+tENbUdSKigFyzHcN4EOIrIE+AgYU9HuZKXKElGvPm1v/IL14R3p+vONpM2eYjtSmfbsymTf68NpYHLIOvsDWrTrYjuSUkErKAquMSbPGHOpMaanMaavMeYH25mUf4tqEEvctV+yJaw1XX8cy8Lv3rMd6W9y9uxi68vDaVm4hfVDX6VT4nG2IykV1IKi4CpVG2LjmxP3z+/YEN6Bnr/dTMr0121HOmzPrky2Pn8qHfNWsvSY/9HzuBG2IykV9LTgKnUEYhrF0/ymGayu041+829n7rv3Y4qKrGbK2p5O5oun0i5/HUv/8QJ9Th1jNY9SyqEFN4i9/fbb3HjjjRX2k5GRcfj51VdfbeVWf74sOqYRHW6bQUr0CQxc9xzzn7uEvEN2rry1ftl8Dk08gVYFm1h+/EQST77YSg6l1N9pwVXlKllwX3/9dbp3724xkW+KiIyiz61TmNPqSo7O/poNTx7HlnXLvZph0Q8fEf/xWYSbfDaNmELvE8736ucrpcoXLD8Lqh3f3APbFtfsMJslwLDHK+ytrNvz3XLLLXz11VfUq1ePL774gqZNm/Lll1/y8MMPk5eXR+PGjXn//fdp2rTp4WHl5OTQq1cvVq1aRXh4OHv37qV3795MmDCB5ORkLrnkEurVq8ecOXMYNmwYTz31FElJSXz77beMGzeOwsJC4uLi+P7772u2LfxMSGgog65+mgUz+nLUnLuRd04gud94+p05tlavV3wo9wAL37yVgTs+Zk1YR6Kv+JTOrTrW2ucppapHt3D9VGm359u/fz8DBw5k0aJFDB48mNdeew2Af/zjH8ydO5eFCxcyatQoJkyY8JdhRUdHc/zxxzN9+nTAuY7yyJEjOf/880lKSuL9998nNTWVevXqHX5PZmYm11xzDVOmTGHRokV8+umn3ht5H9f31MvYf8WPbKnTnqQFd7N4wlC2rFtaK5+1fN4MMiYMZOCOj5kXdy6t7viFplpslfJJuoV7JCqxJVpbSrs9X506dTjzzDMB6NevHzNnzgQgPT2dCy+8kK1bt5KXl0f79u3/Nryrr76aCRMmcPbZZ/PWW28dLtZlmTt3LoMHDz48rOJb9ylH87ZdaHL3L8yb/CQ9lj9L+DtDmNvsXDqNvJ/GTVsd8fC3rFtOxhcP0n/PDLYRR+pxrzDAxy/AoVSw0y1cP1TW7fnCw8MP36IvNDSUgoICAG666SZuvPFGFi9ezCuvvFLqrfSOPfZYNmzYwOzZsyksLKRnz55eHadAFBoWxoBR/+LA2DmkNTyJ/ts+pt5LfZn3/GjWpv1e5eGZoiJW/DGT+U9fQNN3jqFX9g/MaX4p0ben+PzVrpRSuoXrl8q6PV95/bds6dyG7Z133imzv9GjR3PxxRdz//33H37N89Z+ngYOHMgNN9zA+vXrad++Pbt27dKt3DI0admeJrd9zKZVqWz/+nF67/yaiM++YOPUVmxtMoSIToNp3qU/TVq0/8ux3qLCQnZkrGfbqmRyV/1Iyx0/09VkcMDUJbnp+XQ8exyDWrSzN2JKqSrRguuHTjvtNCZOnEi3bt3o0qXL4dvzlWX8+PGcf/75xMbGcuKJJ7J+/fpS+7vkkku47777uOiiiw6/dvnll3PdddcdPmmqWHx8PK+++iojR46kqKiIJk2aHN6FrUrXpnMibTp/xJ5dmSz67nUi18+g79aPqLPtffgFDplw9kg0eVKXOuYQ0WYfzSSPZjjdVkX0ZFuX6+g2dDQDG8TaHh2lVBUFzO35akOw3Z5v8uTJfPHFF0yaNKlWPyeQ27Cq9udks3nZH+xZn4LJ3kxI7m5CCg9RFFqXoroxSNxRRLXsTsc+Q4ioV992XKUqRW/PVzrdwlWAc5z3m2++4euvv7YdJajUj25I1wGnwIBTbEdRStUyLbgKgOeff952BKWUCmh6lnI16G746tO2U0oFKy24VRQREUFWVpYWjmowxpCVlUVERITtKEop5XW6S7mKWrVqRXp6OpmZmbaj+KWIiAhatTryCz8opZS/0YJbReHh4aVeqUkppZQqj+5SVkoppbxAC65SSinlBVpwlVJKKS/QK02VQ0QygY3VfHscsLMG4/gDHefgoOMcHI5knNsaY+JrMkwg0IJbS0QkOdgubabjHBx0nINDMI5zbdNdykoppZQXaMFVSimlvEALbu151XYAC3Scg4OOc3AIxnGuVXoMVymllPIC3cJVSimlvEALrlJKKeUFWnCPkIicJiIrRWSNiNxTSve6IvKx232eiLSzELNGVWKc/09ElolImoh8LyJtbeSsSRWNs0d/54qIERG//zlFZcZZRC5wv+ulIvKBtzPWpEpM121E5EcRWehO26fbyFmTRORNEdkhIkvK6C4i8pzbJmki0tfbGQOKMUb/qvkHhAJrgQ5AHWAR0L1EPzcAE93Ho4CPbef2wjifAES6j68PhnF2+4sGfgbmAkm2c3vhe+4ELARi3edNbOeu5fF9Fbjefdwd2GA7dw2M92CgL7CkjO6nA98AAgwE5tnO7M9/uoV7ZI4G1hhj1hlj8oCPgBEl+hkBvOM+ngycJCLixYw1rcJxNsb8aIw54D6dC/j7/fgq8z0D/Ad4Asj1ZrhaUplxvgZ40RizG8AYs8PLGWtSZcbXAA3cxzFAhhfz1QpjzM/ArnJ6GQG8axxzgYYi0tw76QKPFtwj0xLY7PE83X2t1H6MMQXAHqCxV9LVjsqMs6ercNaQ/VmF4+zuamttjJnuzWC1qDLfc2egs4j8JiJzReQ0r6WreZUZ3/HApSKSDnwN3OSdaFZVdX5X5dD74apaIyKXAknAENtZapOIhAD/Ay63HMXbwnB2Kx+PsxfjZxFJMMZk2wxViy4C3jbG/FdEBgGTRKSnMabIdjDlH3QL98hsAVp7PG/lvlZqPyIShrMrKssr6WpHZcYZERkK3AsMN8Yc8lK22lLROEcDPYHZIrIB51jXND8/caoy33M6MM0Yk2+MWQ+swinA/qgy43sV8AmAMWYOEIFzgf9AVqn5XVWOFtwjMx/oJCLtRaQOzklR00r0Mw0Y4z4+D/jBuGcj+KkKx1lE+gCv4BRbfz6uV6zccTbG7DHGxBlj2hlj2uEctx5ujEm2E7dGVGbanoqzdYuIxOHsYl7nxYw1qTLjuwk4CUBEuuEU3EyvpvS+acBo92zlgcAeY8xW26H8le5SPgLGmAIRuRGYgXOW45vGmKUi8m8g2RgzDXgDZ9fTGpyTE0bZS3zkKjnOTwJRwKfu+WGbjDHDrYU+QpUc54BSyXGeAZwiIsuAQuBOY4xf7r2p5PjeDrwmIrfhnEB1uZ+vPCMiH+KsNMW5x6YfBMIBjDETcY5Vnw6sAQ4AV9hJGhj00o5KKaWUF+guZaWUUsoLtOAqpZRSXqAFVymllPICLbhKKaWUF2jBVUoppbxAC65SSinlBVpwlVJKKS/QgquUHxGR/u59SSNEpL57H9qetnMppSqmF75Qys+IyMM4lxWsB6QbYx6zHEkpVQlacJXyM+61fufj3Hf3GGNMoeVISqlK0F3KSvmfxjjXqo7G2dJVSvkB3cJVys+IyDTgI6A90NwYc6PlSEqpStC7BSnlR0RkNJBvjPlAREKB30XkRGPMD7azKaXKp1u4SimllBfoMVyllFLKC7TgKqWUUl6gBVcppZTyAi24SimllBdowVVKKaW8QAuuUkop5QVacJVSSikv+H/fAuFxNTDSQAAAAABJRU5ErkJggg==\n", - "text/plain": [ - "
    " - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - }, - { - "name": "stdout", - "output_type": "stream", - "text": [ - "The max absolute difference is: 1.77636e-15\n" - ] - } - ], + "execution_count": 15, + "id": "c01acfbf", + "metadata": { + "collapsed": false, + "editable": true + }, + "outputs": [], "source": [ "import autograd.numpy as np\n", "\n", @@ -3161,8 +3595,10 @@ }, { "cell_type": "markdown", - "id": "b692e704", - "metadata": {}, + "id": "bdd7abde", + "metadata": { + "editable": true + }, "source": [ "## Using autograd\n", "\n", @@ -3176,8 +3612,11 @@ { "cell_type": "code", "execution_count": 16, - "id": "85444dae", - "metadata": {}, + "id": "898937c9", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3201,8 +3640,10 @@ }, { "cell_type": "markdown", - "id": "70e49fb6", - "metadata": {}, + "id": "3b343ba9", + "metadata": { + "editable": true + }, "source": [ "## Autograd with more complicated functions\n", "\n", @@ -3214,8 +3655,11 @@ { "cell_type": "code", "execution_count": 17, - "id": "94bfdf53", - "metadata": {}, + "id": "334ce966", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3255,16 +3699,20 @@ }, { "cell_type": "markdown", - "id": "b4686bc1", - "metadata": {}, + "id": "8b661c94", + "metadata": { + "editable": true + }, "source": [ "Note that the grad function will not produce the true gradient of the function. The true gradient of a function with two or more variables will produce a vector, where each element is the function differentiated w.r.t a variable." ] }, { "cell_type": "markdown", - "id": "de9bbb0b", - "metadata": {}, + "id": "f6778038", + "metadata": { + "editable": true + }, "source": [ "## More complicated functions using the elements of their arguments directly" ] @@ -3272,8 +3720,11 @@ { "cell_type": "code", "execution_count": 18, - "id": "3e932adf", - "metadata": {}, + "id": "13700fa8", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3297,8 +3748,10 @@ }, { "cell_type": "markdown", - "id": "53c2fedf", - "metadata": {}, + "id": "daa46b77", + "metadata": { + "editable": true + }, "source": [ "Note that in this case, when sending an array as input argument, the\n", "output from Autograd is another array. This is the true gradient of\n", @@ -3310,8 +3763,10 @@ }, { "cell_type": "markdown", - "id": "afa8f89f", - "metadata": {}, + "id": "9d975664", + "metadata": { + "editable": true + }, "source": [ "## Functions using mathematical functions from Numpy" ] @@ -3319,8 +3774,11 @@ { "cell_type": "code", "execution_count": 19, - "id": "bc6756d2", - "metadata": {}, + "id": "2256de07", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3344,8 +3802,10 @@ }, { "cell_type": "markdown", - "id": "80e3c79f", - "metadata": {}, + "id": "4a6deaaf", + "metadata": { + "editable": true + }, "source": [ "## More autograd" ] @@ -3353,8 +3813,11 @@ { "cell_type": "code", "execution_count": 20, - "id": "af4ac652", - "metadata": {}, + "id": "76ab0823", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3375,8 +3838,10 @@ }, { "cell_type": "markdown", - "id": "85d76920", - "metadata": {}, + "id": "c7e1ff75", + "metadata": { + "editable": true + }, "source": [ "## And with loops" ] @@ -3384,8 +3849,11 @@ { "cell_type": "code", "execution_count": 21, - "id": "aa71a6e5", - "metadata": {}, + "id": "9b7bc7ac", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3417,8 +3885,11 @@ { "cell_type": "code", "execution_count": 22, - "id": "00f338d9", - "metadata": {}, + "id": "f48ea612", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3434,8 +3905,10 @@ }, { "cell_type": "markdown", - "id": "514775dc", - "metadata": {}, + "id": "41d0c6f4", + "metadata": { + "editable": true + }, "source": [ "## Using recursion" ] @@ -3443,8 +3916,11 @@ { "cell_type": "code", "execution_count": 23, - "id": "56292f34", - "metadata": {}, + "id": "96932b44", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3478,16 +3954,20 @@ }, { "cell_type": "markdown", - "id": "0821aef0", - "metadata": {}, + "id": "9f7eae7b", + "metadata": { + "editable": true + }, "source": [ "Note that if n is equal to zero or one, Autograd will give an error message. This message appears when the output is independent on input." ] }, { "cell_type": "markdown", - "id": "07bccded", - "metadata": {}, + "id": "5441446f", + "metadata": { + "editable": true + }, "source": [ "## Unsupported functions\n", "Autograd supports many features. However, there are some functions that is not supported (yet) by Autograd.\n", @@ -3498,8 +3978,11 @@ { "cell_type": "code", "execution_count": 24, - "id": "7e1c00c8", - "metadata": {}, + "id": "810cd2a6", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3517,16 +4000,20 @@ }, { "cell_type": "markdown", - "id": "91b32000", - "metadata": {}, + "id": "3699391b", + "metadata": { + "editable": true + }, "source": [ "Here, Autograd tells us that an 'ArrayBox' does not support item assignment. The item assignment is done when the program tries to assign x[2] to the value 3. However, Autograd has implemented the computation of the derivative such that this assignment is not possible." ] }, { "cell_type": "markdown", - "id": "5bb6505b", - "metadata": {}, + "id": "2db124a4", + "metadata": { + "editable": true + }, "source": [ "## The syntax a.dot(b) when finding the dot product" ] @@ -3534,8 +4021,11 @@ { "cell_type": "code", "execution_count": 25, - "id": "b698d800", - "metadata": {}, + "id": "3c3b74ba", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3553,8 +4043,10 @@ }, { "cell_type": "markdown", - "id": "7daefc20", - "metadata": {}, + "id": "37817c30", + "metadata": { + "editable": true + }, "source": [ "Here we are told that the 'dot' function does not belong to Autograd's\n", "version of a Numpy array. To overcome this, an alternative syntax\n", @@ -3564,8 +4056,11 @@ { "cell_type": "code", "execution_count": 26, - "id": "f4db3939", - "metadata": {}, + "id": "d8e9ea75", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -3586,8 +4081,10 @@ }, { "cell_type": "markdown", - "id": "63911856", - "metadata": {}, + "id": "64351a6b", + "metadata": { + "editable": true + }, "source": [ "## Recommended to avoid\n", "The documentation recommends to avoid inplace operations such as" @@ -3596,8 +4093,11 @@ { "cell_type": "code", "execution_count": 27, - "id": "fc661b03", - "metadata": {}, + "id": "b61840f6", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "a += b\n", @@ -3608,8 +4108,10 @@ }, { "cell_type": "markdown", - "id": "a873728b", - "metadata": {}, + "id": "d87544e5", + "metadata": { + "editable": true + }, "source": [ "## Using Autograd with OLS\n", "\n", @@ -3620,36 +4122,13 @@ }, { "cell_type": "code", - "execution_count": 6, - "id": "68ec2659", - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Own inversion\n", - "[[3.90257879]\n", - " [3.03235887]]\n", - "Eigenvalues of Hessian Matrix:[0.34535072 4.41651835]\n", - "theta from own gd\n", - "[[3.90257879]\n", - " [3.03235887]]\n" - ] - }, - { - "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEWCAYAAAB1xKBvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAsMElEQVR4nO3deZgU5bXH8e8Z9gFlGZBFmcEV0cQNonFHB6MxRrPcJBrMdQ3GbMaYxESu0cSoebKYmMTE4BKNELMYk9xovFHHGXEDMigqq7gAggjIIiADzHLuH9UDzdDdUz1d3dU98/s8Tz90V9dyqqapU1XnrbfM3REREWlPWdwBiIhIaVDCEBGRUJQwREQkFCUMEREJRQlDRERCUcIQEZFQlDCk0zKz681satxxRMnMxpvZ8rjjkK5JCUMKysyWmFmDmW02s7fN7B4z6xd3XCLSPiUMicNH3b0fcARwJPCdeMPpmsyse9wxSGlRwpDYuPvbwL8JEgcAZvZtM3vNzDaZ2Xwz+3jSdxea2dNm9hMzW29mb5jZh5O+39fMnkxM+xgwOHl5Zna2mc0zsw1mVmdmY5K+W2Jm3zSzl8zsPTO7y8yGmtkjifk9bmYDU61H62UiM7vKzFab2Uozuyjp+zozu7TteiR9djP7opktTizrBjPb38yeNbONZvZnM+vZZpnXmNk7ibgnJg3vldg+y8xslZndbmZ92sR5tZm9DfzOzAab2UOJbbLOzJ4yM+0XJCX9MCQ2ZrYP8GHg1aTBrwEnAv2B7wFTzWx40vfHAIsIksGPgLvMzBLf/QGYnfjuBuCCpGUdBNwPfA0YAvwL+GebHfEngdOAg4CPAo8A1yTGLwO+mmF1hiVi3hu4BLgtXYJJ43RgLPBB4FvAFOB8YCTwPuC8NssanFjWBcAUMxud+O6HifiPAA5IjPPdNtMOAqqAScBVwPLEOg5NrK/6C5KUlDAkDn83s03Am8Bq4LrWL9z9L+7+lru3uPufgMXA0UnTLnX3O9y9GbgXGA4MNbNK4APAte6+zd2nA/9Mmu4zwMPu/pi7NwI/AfoAxyWN80t3X+XuK4CngJnu/oK7bwX+RnD5LJ1G4Pvu3uju/wI2A6MzjN/Wj9x9o7vPA+YCj7r76+7+LkHiarvs1vV8EngY+HQicU4CrnT3de6+CbgJODdpuhbgusS0DYm4hwNVidifcnUwJ2koYUgcPubuewDjgYNJunRkZv9tZnMSl0g2EBxdJ19aerv1jbtvSbztB4wA1rv7e0njLk16PyL5s7u3ECSsvZPGWZX0viHF50zF+bXu3pT0eUs747eVzbJTrecIgrOEcmB20vb7v8TwVmsSCbDVjwnO8B41s9fN7NtZxCxdjBKGxCZxdHwPwdE+ZlYF3AF8Gahw9wEER9uWZhbJVgIDzaxv0rDKpPdvEVyGIbEsI7jcs6LjaxDaewQ78lbDcpxfqvV8C3iHILkc6u4DEq/+iQYGrXY5e3D3Te5+lbvvB5wNfN3MqnOMTzopJQyJ28+B08zscKAvwQ5tDUCicPy+MDNx96VAPfA9M+tpZicQ1CFa/Rn4iJlVm1kPgmv324Bno1qRDOYAnzCzcjM7gKDGkavW9TwROAv4S+Ks6Q7gZ2a2F4CZ7W1mp6ebiZmdZWYHJBLou0AzwWUrkd0oYUis3H0N8Hvgu+4+H/gp8BzBJZn3A89kMbvPEhTF1xHURX6ftJxFBEXkXxIciX+UoHnv9ghWoz0/A7YTrNO9wLQc5/c2sJ7grGIa8AV3X5j47mqCS0wzzGwj8DiZaykHJsbZTLDdf+3utTnGJ52Uqb4lIiJh6AxDRERCUcIQEZFQlDBERCQUJQwREQmlpDofGzx4sI8aNSruMERESsrs2bPfcfch7Y+ZWUkljFGjRlFfXx93GCIiJcXMlrY/Vvt0SUpEREJRwhARkVCUMEREJBQlDBERCUUJQ0REQlHCEBGRUJQwREQkFCUMEREJRQlDRERCUcIQEZFQ8p4wzOxuM1ttZnNTfHeVmbmZDc53HCIikptCnGHcA5zRdqCZjQQ+BCwrQAwiIpKjvCcMd59O8Izltn4GfAvQM2JFREpALDUMMzsHWOHuL4YYd5KZ1ZtZ/Zo1awoQnYiIpFLwhGFm5cA1wHfDjO/uU9x9nLuPGzIk5+7cRUSkg+I4w9gf2Bd40cyWAPsAz5vZsBhiERGRkAr+ACV3fxnYq/VzImmMc/d3Ch2LiIiEV4hmtfcDzwGjzWy5mV2S72WKiEj08n6G4e7ntfP9qHzHICIiudOd3iIiEooShoiIhKKEISIioShhiIhIKEoYIiISihKGiIiEooQhIiKhKGGIiEgoShgiIhKKEoaIiISihCEiIqEoYYiISChKGCIiEooShoiIhKKEISIioShhiIhIKEoYIiISihKGiIiEooQhIiKhKGGIiEgoeU8YZna3ma02s7lJw35sZgvN7CUz+5uZDch3HCIikptCnGHcA5zRZthjwPvc/TDgFeA7BYhDRERykPeE4e7TgXVthj3q7k2JjzOAffIdh4iI5KYYahgXA4+k+9LMJplZvZnVr1mzpoBhiYhIslgThplNBpqAaenGcfcp7j7O3ccNGTKkcMGJiMguuse1YDO7EDgLqHZ3jysOEREJJ5aEYWZnAN8CTnb3LXHEICIi2SlEs9r7geeA0Wa23MwuAX4F7AE8ZmZzzOz2fMchItIpTJsGo0ZBWVnw77S0V/Qjl/czDHc/L8Xgu/K9XBGRTmfaNJg0CbYkLswsXRp8Bpg4Me+LL4ZWUiIiEsbkyTuTRastW4LhBaCEISJSKpYty254xJQwRERKRWVldsMjpoQhIlIqbrwRyst3HVZeHgwvACUMEZFSMXEiTJkCVVVgFvw7ZUpBCt4Q4417IiIS0rRpQWF72bLg8tONNxYsSSRTwhARKWYxN6VNpktSIiLFLOamtMmUMEREilnMTWmTKWGISNcQY5caOYm5KW0yJQwR6fxa6wBLl4L7zjpAKSSNmJvSJlPCEJHOr4jqAFmLuSltMiulR1GMGzfO6+vr4w5DREpNWVlwZtGWGbS0FD6eAjOz2e4+Ltf56AxDRDq/IqoDlDIlDBHp/IqoDpBWCRTllTBEpPMrojpASiVSlFfCEJGuYeJEWLIkqFksWVI8yQJyL8oX6OxEXYOIiMQtl5vzCth1iM4wRETilktRvoBNhpUwRETilktRvoBdh+Q9YZjZ3Wa22szmJg0bZGaPmdnixL8D8x2HiEjRyqUoX8Amw4U4w7gHOKPNsG8DNe5+IFCT+Cwi0nV1tChfwCbDeU8Y7j4dWNdm8DnAvYn39wIfy3ccItJFlcD9DTlJc3bi532W155YypTzp0e2qLhqGEPdfWXi/dvA0HQjmtkkM6s3s/o1a9YUJjoR6RxK5P6GXXQkwSXOTla9tIr7P/oHLvleJfv2XMEB1VVcNu2kyEIrSF9SZjYKeMjd35f4vMHdByR9v97d261jqC8pEcnKqFFBkmirqiq47FNs2jaRheDyUpp6xsblG3nytwupeaiBmgUjmLvtQAAG2AZOGbaQ6hO2UX3hSMZ8ZP9I+pKKK2EsAsa7+0ozGw7Uufvo9uajhCEiWSm1TgfbSXDbNm7jubsXUPPXDdTMqWDW5jE0053eNHDCoPlMOHoT1eftxZHnjqZbz247Ji/1zgf/F7gg8f4C4B8xxSEiHVEqdYF8tCDK57qnaQrrS5dx+uB6BvZv5pQrj+Cmp0+kxY1vH/80T/z0BdavNx5bO5arHxnPuP8+ZJdkESl3z+sLuB9YCTQCy4FLgAqC1lGLgceBQWHmNXbsWBeRmE2d6l5e7h4cuwev8vJgeLGJOtY8r3tLZeWu80683qDKD+m12L96eJ3/45oZvmHphswxVlW5mwX/Tp3qQL1HsT+PYiaFeilhiBSBqqqUOzWvqoo7stRS7EA7LIp1bxPPuhtv8/u+8LRfeMB0/5L90jeza0Jq7N7L1914W/h5p0hog+F1j2AfrAcoiUh2Sq0uEKVc133aNPzzn8caGnYMeo9yPs8U/m0f5tQRi7is8l+c/Nrv6L7mLayyMrifIuw9GWlqIIfD9hfde4WbSXrqGkREshN1XaBU6iHQoXXfumErNT9+nmuOq2Pl5765S7IA6MsWflfxTdZsH8Bflh/LhGdvoMeq5VhHetVNUwPpAT3DzyQ9JQwRyU6UdxaXyn0SrUlt6dLgbCJZm3Vv3t7MrN/N4+bT65gw6HkGDnQmfOsofvTcCQz1t1POvte6tynrHsHuOE3iaoTtuc9cCUNEshXlw4gK2NNqhyUnNQgSW2vSqKrCfzuFBf0/yK8+9SQfHzGDwb03cczFh3LNo+NZ3dCPy4+ayUPX/Yf1Kxooq8pzv09pkvlbsCKS+UdRCCnUS0VvkU7GLHUR2Sy6ZeRa9E5T6N7Ub6j/esi1/iZ7ezPmb1DlX7Vb/dLRT/r9X3nGV81dnTqWfLcwUyspJQzp4qJs6ROVKGJqr9VRrsuIYgedJqk1g79Hn+zn3ZF1ynE7KGGIdBXFeN9DVDFlmk8Uy+hgM9j31rznj95c71cfU+tv2fCU82gp69aheWctgu2ghCHSVcR130PrUS24d+u2c5nJw6OIKd3RcxTLCHnJq7Gh0Z+742W/obrWxw943nuy1cG9B9v8e71v8m1lvXbfYaeab9SX09wj2Q5RJQzdhyFS7OK47yFVJ3itystTD486pijWO819CT5oEPPvnsnj966g5tnePLnqYDbSH4Aj+iyk+pC3qT67LydOGkO/Yf2C7TF5ctBstfXeiMmTC9OxYQTbIaq+pLrnOgMRybPKytQ7pjw8UW2HVK2XWm3ZAt26QXNzfmOKYr1vvBEuuggaG3cZvH3dJm782EzuZyL7d1/KuQe/yIQP9+CUyw5i8OiDgYN3nc/EialbgaXqWTbqBxfF8fdPQ81qRYpdAZ+otkN7z4Nubs5/TKnW2yzYeYa4we+dRWv5y+z92NRSvtt3vWhkSt+rWPL0cl5trOK3C07iU7ccy+DRFeHji7J5cSZx/P3TieK6VqFeqmFIl1XoVlLprpsnXz8vREzJ9ZK29Yg2hd/Nqzb7Izf8x78xrtaP7DPfjWYH92YK0HQ334qklZRqGCJdXarr85C5hpGPI+lM0tQitg4azo8O/wOPzx7IjI1jaKQnPdnGcQPmUz32Xao/XcEHf3AW9maKM6ZifYhSHpT68zBEpBik65oDdl5ugaBmAfm77NKeNJfIeq57m+trT2JLUw+uPPpZ/n3TbNavaaZ2/ZH8z+PjOXbS+7Gbb0p9SefMM0unD6tU4uiDK4rTlEK9dElKJGJF3lX5608u8zv+e7qv7jY0ZZyb9xzqa19d1/6M2l7Sufzy4ru3JRtZ3ptBRJekdIYh0pWlK263V/ROJ8ej3tXz1vDHrz7L5w+ezn49lrHfySP5/O9P5J/NZ7DbxfPycvr++qcM2n9g+zOeODG4/HTffcHn3/ym+PuwgvTbM64+uKLIOoV66QxDuoyoC8rZ3hyXXNjOZhmpbmirqEg7n40rNvpD183yK4+q9cN6L9wxSX82+DnDZvgvPlnny6+6xVvaztcsOEvIdhtkuuGu9VUsMp1FZNkHF7rTW6STKuRjRdvbiWaz3EzJJzGfbZu2+ZO/mOPfPanWj9/jRe/Odj+Pqb6Ekd6M+YbeQ33x53/ojQ2N7c8328tm7bX8at3hFstlqUzrneU2UcIQ6ayiriuE7eAv05lGGOmOehOvt8uGeTmbHdzLaPIP9J3rfzxgsjf16J05SWWabzbaiS/n7Ry1TGcRqmGIdHHJD+lJpaN1hfbqFK3X99s+GCjb5bZz5/GQllVc/P56/vbtmaxdsplZmw/lM41T6da4ddcR216LTzdfs+xqJGHvjO7odo5apqf7FeqmwbaiyDodfQFXAvOAucD9QO9M4+sMowQVY7fcxSjM9fV8nWG0N15FRca/4coXV/m0Lz7ttw+9dvfuvtuLP8y1+EzX7LPZJmFrGMVyhhHhpUlK/ZIUsDfwBtAn8fnPwIWZplHCKDHF2C13sWrv+nq+ahjtjdezp3uPHrsMa+nTx+vPvt6/enidH9rrlR1fDbR1fsuA631Lzz29JWz8YZNZuu2S7d3ayQcwFRXB+hXq9xnDczBadZaE8SYwiKATxIeAD2WaRgmjxBR5G/+ikun6eq5nZlOnBjvH5DOGdPNru4NKni7p9QZV3of3/LRB9f7DM2q9/r753rStKf18Mi0vTDLL128pH2fAqeYZ88FTySeMYB24AtgMrAGmpRlnElAP1FdWVka6ESXPCvH4zc4inzvELHdUTdua/D/3zvMfnlGbth+mFvCt727NLbbkGDPttNsmvGI+W023vdMk3kIdPJV8wgAGAk8AQ4AewN+B8zNNozOMEqMzjPDydQQa4m/Q0tziC//1mt/2mTr/xIjnfKCt2zHaCkakT/qF2Fl34N6OWIVpuhvDwVNUCSPOVlITgDfcfY27NwIPAsfFGI9ErZi6ZS52+Wr1kqbFjy9bxn1feJoLD3iayp4rOfjM/fjSn06mftU+fPyAuUz74jOsfHE1I6b+KHXrKffC3BWd7rkc/foVvj+rMLJtYRXDMy1yEVtvtWZ2DHA38AGgAbiHIAv+Mt006q22BKXqCbUY/6N3Vmma6S6hin1ZQoWt5dS9X6H6pEaqL65i/1MqsbI2CSJdc9t8PvGvVRxPG8xFumbRFRXQ0LD7w5YK1JFjyfdW6+4zgQeA54GXE7FMiSseyZPWNv4tLcG/ShYF0bCugZofP88fu59PA713+W4rvZl7xPk8/4eFrN4+kD+/eSyXTTuJA6qrdk8WsLPH2rYKcXSc6V6EYpTurPrWW+O5byJqUVzXKtRLNQwpWoW63yTNcpq2NfnMu+f6jafV+qkDZ3svGhzcu7Pdr+t1k6/vNdRbMG8ZWZl9bB2tr0SxTUqxaXYR3ntEoYrewGPA4VEsLNeXEoYUpULt1FIsZ3u3Xv6TPa/3/mzYMfiw3gv9yqNq/eHrZ/nGFRujW3Y2O8Eot0kR7oBLTVQJo90ahpkdBfwUWAJc4+4r83vOk55qGFKU0l23jviJbk3D96H72yt2G76cvfne6D8w4fTunDLpQPY6dEhky+ywAm0TCSeqGkbooreZfRL4LkFrph+5e0OuC8+WEoYUpTwVZte9tp7a2xdS80gjNa+MZEHj/pTt/lSI4iwAl1qxupMraNHbzAxYBPwG+Aqw2Mw+l+vCRTqFiAqzW97ZwqM3z+bqY+oY13c+gw/oz3/95Fjum3ckBw5cw+Y+e0WynIIotWK1hNJuwjCzZ4AVwM8IuvO4EBgPHG1matVUiuJ4FnCp6Mi26eD9Jk1bm3huysv8YEId4wfMYeCQbpx+zVh+Nus4+vbYzvWnTOfpX7/Euvd689Cqo9nzjp+Wzn0tugenc2qvyAEcSuLSVYrvFkRRSAn7UtE7AqXY6qRQctk2IQqzLc0t/vKDr/jPP17nZ+010/fg3eBmX5r9yD7z/Rvjav2RG/7jm1dtzmk5RaOUYu3kKFTROxMz28/dX48kc4WgGkYEVIxMLw/bZukzy3n8jjeoqSvjiTcPZFVLcFnpgB5LqD5gGdVn9OCUyw5i8OiKjsctuekCN5dGVcPonsvEhUwWEpH2HqbTlUWwbd5ZtJYnbl9Ezf81UfNqFa81VQH7MLRsNdUjF1N9yiKqL92XquNHAaPCx9YFdmqxmDYNJk3aeQf20qXBZ9D2TSWK05RCvXRJKgLqEDC9DmybTSs3+b++N8uvGlvrR/RZsGOSPXjXPzp0ht/6iTqf+/fF3tLc0vG4CnkZsatdRuoi/x8o9d5qO/JSwohAPnY+nWUnE2LbbH9vuz9124t+/fhaP3HPOd6DbQ7uPdnq4wc87z+YUOvP3fGyNzY0RhdXoXZqXbG+1UW64FfCkI6Lcgff2XYybbZN8+/v8zl/Wug/PbvWPzxklvdl045C9bjyeX71MbX+6M31/t6a9/IXU6F2al3kaHsXXWSdlTCkOOTjP1zMZyyv1S71KZ970j898hkfbGsc3M9jqi9nhLdgvnnPYb7plt8WLqBC7dS6yNH2LjrbAU8aShhSHKLeycTwH3jV3NV+/1ee8UsOmu6jui/bsdgRZW/55/Z7yp8cf6039+4T306lUNukixxt76azXFLNQAlDikPUO5kC7LQ2rtjoD103y792ZJ2/v/eiHYvozwb/2PDn/Jf/VecLHn5tZ6G6GHakhdipdZGj7a5ICUN2ivMIKeqdTDZnLCHXe+u7W73u5y/4tSfW+nF7vOjd2e7g3osGrx4422/6UK3PumeeN21ryj2mUtcFjra7IiUMCRTDUWGUO5mwR/MZ1ru5sdlnT53vPzqz1k+v+I+Xs9nBvYwmP7rvy/6dY2u95ifPe8P6hnAxVVSkjqmiIvr1F8kDJQwJFOJySSF3iGETYJr1Xt1tqA+ytTsGjen5qn/lsDr/+3dm+PolGzoWU6aEUQwJW6QdShgSyPflkjh2iGH6ZUqz3s2YX3jAdL/vC0/7uh/cFk2iy7SNi6G+kQ2dDXVJShgSyPcOq0h2iBuWbvB/XDPDv3JYnR/Sa7G/Qeq4WiorgwmiTHSZtkEp1Td0NtRlRZUwQj0PQ4pYvruRjqnvqa0btlJ7ywtMPr6OD/aby6Cqfpxz0zHc+dIH2KffBl4+7Hyae/bedaLycuymm4L3kyfv7B+o1ZYtwfBsZdrGpfTchyi3iXRNUWSdjr6AAcADwEJgAXBspvE79RlGLpcK8nmZIaozjHZibNrW5LPumec3n17rEwbVe2+2OLh3o9GP7feS/88JtV738xd867tbw80zH/eHpFpWKR21l9LZkESKznBJCrgXuDTxvicwINP4nTZhFPNOJ4rYUsyjpbzcV3zzFv/Vp+r848Of8wG2fsfX7+v1il9xRJ3/89qZ/u6b73Ys7kJeSiuVukCRXF6Uwiv5hAH0B94gzcOZUr06bcIo9v/Iue4Q06zfG1QFq9ntTb/4wOn+hy8/42+/vDq7eCoqglcpH/kXirZJl9UZEsYRwCzgHuAF4E6gb4rxJgH1QH1la0Gzs0l3qQDijiwn615f7w9+6zlvIfX6tWD+as2S7Lr+TrXTS7cDLJUj/0LSNumSokoYOT1xLxdmNg6YARzv7jPN7FZgo7tfm26aTvvEvXRPejOD++4rmQe5NKxr4Jk7F1Dzt408/tJePL9lNC10YymVVPLm7hMkP8ku7AOC0m2rdPMVkcieuBdnK6nlwHJ3n5n4/ABwVEEjmDYt2AGVlQX/TptW0MXvcOONQXJoy72oW7A0bW1i5l1zuelDdZw68AUGVhinXX0UP5lxPL27N3LtyU8x/ZcvMuKuGzK35Gp96tnSpcE6tz71LNXfI0zrLD09UCQ/ojhN6egLeAoYnXh/PfDjTONHWsPo6PXcfJ3Sp7vEUkQtWFqaW3zePxb7rZ+o87OHzfA92bAjzMN7L/Svj631h6+f5ZtWbtp94kzbLZsaTrpxi7H2I1IkKPUaRrAOHEFQn3gJ+DswMNP4kSaMjhSa81k0LORT1bJIeEufXe53XzTdJ4562oeXrdwR1n7dl/jnD37S//jVZ3z1/DW5xZRth4Nhaxgi4u6dJGFk+4o0YXSkTXo+d+qFaMESYhnvvLLW//L1Z/2uEf/jyxnhzZi/QZVfxm/83Mpn/M4LpvsbT70ZXUzu2W/XMK2kCkmFZClyShi56sjOP9skE3ZH0joeuHfrtjOOyy+PdkeUZp239B/m3/xArR/VZ74bzX4eU/09dn1gUEs+j9xLublnKccuXYYSRq468h89myQTdv7pxrv88uh3RBk67OvBNj+5/wv+/VNrvWHQ8OyTaa5K9Si92O+hEXEljGhku5PKJsmE3ZGkG6/1TCPHHVFLc4u/9MAi/9nH6nxV2bCU89zSf5hvXrV550TqQiI8bSspAVEljK7d+eDEiUF7/ZaW4N/27neYOBGmTIGKip3D+vSBZ57ZvXlu2E770o3X3Bxu+hSWPL2cuy58ivOqnmVY93c47L8O4sq/n8wPuZpt1mvXkcvL6XPbT+i7V9+dw0qpQ724aVtJF9K1E0ZHNTTsfL92LfzmN7vfQzBoUOpp2+5I0u1YunULNz2wZsE7/PnKZ5k0Zjr791jKvifuw6X3nkjd8v05reoV7r7oKZY+u4Jbmr9Gr/vuCm5sMwv+nTJl90SZ7x5wOxNtK+lKojhNKdSrKPqSCnMfALiXlbn36NH+5asO1DA2rdzkD18/y78+ttYP771wx9d7ssHPHjbDb/1Enc/7x+Lsutxoq1RrCnHQtpIih2oYWYjyP3Smfp/avnr2DNfkM1PX2VVV3mLmDQOH+wMHT/YT9pzj3dkezJ6tfsqA5/0HE2p9xp0ve2NDY8fXS0Q6ragSRmx9SXVEh/qSau12IvnBMeXlqS/FhBGmL6NkHejXqKWphZf+upiaqSupmdmX6WvG8B79MFoYW76QCYetpvrje3L8pWPoM6hPVvOOVNj+n0QkVlH1JRX7WUM2rw6dYUTd7LG9O4070FqmpbnFX61Z4r+d+KR/euQzPtjW7Jj84J6v+pfeX+cPfus5X/f6+vQxFfqSiO4/ECkZ6AwjpLKyYHfWllnQOqoj2h5Zn3lmcMaSqmVTmjOMVXPX8MRvF/P4v5upeX1fljbvA8DeZSupHvUa1dVQPWl/9h43vP1YojyDCivdmZZ6ihUpOlGdYXT+hFGoHVs7O+6NyzcyfcpCHv9nAzULRjB324EADLANnDJsIdUnbKP6wpGMPmNfrCxFz7XpxLXjzkciFpG8iCphdI8imKJ2442pd+RRN3tsPZpPnHn4PiNZcPJl3H/73tR84WVmbR5DM0fTmwZOGDSfiePrmPDZvTjy3NF06/nBji837P0eUausTJ2odP+BSKfV+RNGmx15voqzzdubmeNHUXPIPdRs7sdTbx5Cw9RyymjmA30XcPVxTzPhk/059uIx9B4wNroFx7XjLlQiFpGi0fkTRjo5tvDxFmfxY0uouedNap7qQe1bo1nnY4AxHNLrVS497D9MOKs3J192MP0r35e/9Yhrx12gRCwixaPz1zBS1RZ69gyuvzc27hwWolC8cs4qan77KjWPt1Dzxn682bw3ACO7rWDCfq9TPaGMUycdwPAjhmYXY67UvFVEMlDRO6xs7ptoUyh+d9m71N2+kJqHt/L4wn1YsH1/AAbZOk4dsYjqkxqpvqiSA6qrsitUi4gUkIreYWVR/PWlS8HK2NBrL27xK7lp+zdo4RjKeY8TKxZw0TF1VH92KEd8ZjRl3Y/NY9AiIsWn83c+mEXx1wDDGbhtFddsv54/jr6Oup/PYd273fm/d8bxzYfHc9TEMZR170Sbbdq03XvaFRFJoRPt+RLa7gDPPHO33kQ98cqkD1v51NapnHzFEfTas1c7Y5eo1vpO2552lTREJIXOVcNIUeBu6d2HRfuezpDFz1HRtArHKEtKF63vUlYgOvtNaLpbW6RLiKqGEfsZhpl1M7MXzOyhnGc2efKuraGAsq0N9FnwAmOa5/JOt6G7JAtIXIbK4tkTnUpcN/2JSEmKPWEAVwALOjpxw7oGHvvhbL79wTpalqbe0VXZMlZtH8SQltWpZ9Lc3DUfghP2aXGqc4gIMScMM9sH+AhwZ9hpmrY2MePOudx4Wh2nDnyBARVlfOg7Y/npzONZbcNSL6eyMihUp9tBtj55rr0n0XU2YZ4WpzqHiLSKosvbjr6AB4CxwHjgoTTjTALqgfo9yg7yPdmwozftI/os8KvG1vq/vjfLN63c1H6X2+qSe3ftdY0edffwIlJwlHr35mZ2FnCmu3/RzMYD33D3szJN09sO8wsO/hXVZ3TnlEkHMWTM4N1Hau+uZ90VnR31SitS8kr+Tm8zuxn4HNAE9Ab2BB509/PTTdOhO70lN2pJJVLySr6VlLt/x933cfdRwLnAE5mSRVHJVxG4GIvLYeocItIlFEMrqezEvVPNVxG4WIvLEyd2zQYBIrK7KAohhXqN3Xff9ovarQXciorg1d5zrrN9HnauReB0y1NxWUTyhIiK3qXV+eCKFbB9+67DtmwJitiw613ea9fuHKf1aB12L4BffPHOeS5dGnxuO16yjt7sNm0aXHFF+rh0E52IFLnS6hrEzFOWvM3SP3kuWdtC7eDBu+7AW1VUwDvvpJ5HR4rAqZ7J0XZaUHFZRPKi5IveHdKzZ+rhlZXhjsTbjpMqWWQaDh0rAqfosmS3uFRcFpEiV1oJY++90+9Uw/T7FEXfUGGLwMnF+fbOfCorVVwWkeIXRSGkUK+xY8emLxqnuou7vTu6KypSj1tRkX1VKVl7sehOcxEpICIqepfWGQYER9xLlgR3GS9ZsvMIvO0RekVF8Mp0tH7rrdCjx67DevQIhueivUtQrSoqdBYhIiWjtFpJtWfixOx2vq3jRt1VSKZ6SmuBXl2SiEiJ6VwJoyOyTTJhpGuxpRZPIlLCSu+SVClQiycR6YSUMPJBLZ5EpBPSJal8ycelLhGRGOkMQ0REQlHCEBGRUJQwREQkFCUMEREJRQlDRERCUcIQEZFQlDBERCQUJQwREQlFCUNEREKJLWGY2UgzqzWz+WY2z8yuiCsWERFpX5xdgzQBV7n782a2BzDbzB5z9/kxxiQiImnEdobh7ivd/fnE+03AAmDvuOIREZHMiqKGYWajgCOBmSm+m2Rm9WZWv2bNmoLHJiIigdgThpn1A/4KfM3dN7b93t2nuPs4dx83ZMiQwgcoIiJAzAnDzHoQJItp7v5gnLGIiEhmcbaSMuAuYIG73xJXHCIiEk6cZxjHA58DTjWzOYnXmTHGIyIiGcTWrNbdnwYsruWLiEh2Yi96i4hIaVDCEBGRUJQwREQkFCUMEREJRQlDRERCUcIQEZFQlDBERCQUJQwREQlFCUNEREJRwhARkVCUMEREJBQlDBERCUUJQ0REQlHCEBGRUJQwREQkFCUMEREJRQlDRERCUcIQEZFQlDBERCQUJQwREQkl1oRhZmeY2SIze9XMvh1nLCIikllsCcPMugG3AR8GDgHOM7ND4opHREQyi/MM42jgVXd/3d23A38EzokxHhERyaB7jMveG3gz6fNy4Ji2I5nZJGBS4uM2M5tbgNhyNRh4J+4gQlCc0SmFGEFxRq1U4hwdxUziTBihuPsUYAqAmdW7+7iYQ2qX4oxWKcRZCjGC4oxaKcUZxXzivCS1AhiZ9HmfxDARESlCcSaM/wAHmtm+ZtYTOBf43xjjERGRDGK7JOXuTWb2ZeDfQDfgbnef185kU/IfWSQUZ7RKIc5SiBEUZ9S6VJzm7lHMR0REOjnd6S0iIqEoYYiISChFkzDa6ybEzHqZ2Z8S3880s1FJ330nMXyRmZ0eY4xfN7P5ZvaSmdWYWVXSd81mNifxymtxP0ScF5rZmqR4Lk367gIzW5x4XRBznD9LivEVM9uQ9F1BtqeZ3W1mq9Pd/2OBXyTW4SUzOyrpu0Juy/binJiI72Uze9bMDk/6bkli+Jyoml/mEOd4M3s36W/73aTvCtaVUIg4v5kU49zE73FQ4ruCbE8zG2lmtYl9zjwzuyLFONH+Pt099hdB0fs1YD+gJ/AicEibcb4I3J54fy7wp8T7QxLj9wL2TcynW0wxngKUJ95f3hpj4vPmItqWFwK/SjHtIOD1xL8DE+8HxhVnm/G/QtAwotDb8yTgKGBumu/PBB4BDPggMLPQ2zJknMe1Lp+gO56ZSd8tAQYXyfYcDzyU6+8l33G2GfejwBOF3p7AcOCoxPs9gFdS/F+P9PdZLGcYYboJOQe4N/H+AaDazCwx/I/uvs3d3wBeTcyv4DG6e627b0l8nEFwb0mh5dLlyunAY+6+zt3XA48BZxRJnOcB9+cplrTcfTqwLsMo5wC/98AMYICZDaew27LdON392UQcEN9vM8z2TKegXQllGWdcv82V7v584v0mYAFBDxrJIv19FkvCSNVNSNsV3zGOuzcB7wIVIactVIzJLiHI7K16m1m9mc0ws4/lIb5WYeP8ZOIU9QEza72BslDbMqtlJS7t7Qs8kTS4UNuzPenWo5DbMlttf5sOPGpmsy3oiidux5rZi2b2iJkdmhhWlNvTzMoJdrR/TRpc8O1pwSX6I4GZbb6K9PdZ9F2DlCIzOx8YB5ycNLjK3VeY2X7AE2b2sru/Fk+E/BO43923mdllBGdup8YUSxjnAg+4e3PSsGLaniXDzE4hSBgnJA0+IbEt9wIeM7OFiSPsODxP8LfdbGZnAn8HDowpljA+Cjzj7slnIwXdnmbWjyBhfc3dN+ZrOVA8ZxhhugnZMY6ZdQf6A2tDTluoGDGzCcBk4Gx339Y63N1XJP59HagjOBrIh3bjdPe1SbHdCYwNO20h40xyLm1O+Qu4PduTbj2KrusbMzuM4O99jruvbR2etC1XA38jP5d0Q3H3je6+OfH+X0APMxtMEW7PhEy/zbxvTzPrQZAsprn7gylGifb3me/CTMjiTXeCosu+7CxoHdpmnC+xa9H7z4n3h7Jr0ft18lP0DhPjkQSFuQPbDB8I9Eq8HwwsJk8Fu5BxDk96/3Fghu8shL2RiHdg4v2guOJMjHcwQRHR4tieiWWMIn2R9iPsWlScVehtGTLOSoL63nFthvcF9kh6/yxwRoxxDmv9WxPsaJcltm2o30uh4kx835+gztE3ju2Z2C6/B36eYZxIf59529gdWPkzCar8rwGTE8O+T3CkDtAb+EviRz8L2C9p2smJ6RYBH44xxseBVcCcxOt/E8OPA15O/MhfBi6JeVveDMxLxFMLHJw07cWJbfwqcFGccSY+Xw/8sM10BdueBEePK4FGguu8lwBfAL6Q+N4IHgT2WiKWcTFty/bivBNYn/TbrE8M3y+xHV9M/CYmxxznl5N+mzNISnCpfi9xxZkY50KCBjfJ0xVsexJcVnTgpaS/65n5/H2qaxAREQmlWGoYIiJS5JQwREQkFCUMEREJRQlDRERCUcIQEZFQlDBERCQUJQwREQlFCUMkB4nnEZyWeP8DM/tl3DGJ5Is6HxTJzXXA9xMdzR0JnB1zPCJ5ozu9RXJkZk8C/YDxHjyXQKRT0iUpkRyY2fsJnny2XclCOjslDJEOSjy5bBrBU802m1nenqgnUgyUMEQ6IPGktQeBq9x9AXADQT1DpNNSDUNERELRGYaIiISihCEiIqEoYYiISChKGCIiEooShoiIhKKEISIioShhiIhIKP8PCU4B9LndZL4AAAAASUVORK5CYII=\n", - "text/plain": [ - "
    " - ] - }, - "metadata": { - "needs_background": "light" - }, - "output_type": "display_data" - } - ], + "execution_count": 28, + "id": "0b8b13c8", + "metadata": { + "collapsed": false, + "editable": true + }, + "outputs": [], "source": [ "# Using Autograd to calculate gradients for OLS\n", "from random import random, seed\n", @@ -3704,95 +4183,23 @@ }, { "cell_type": "markdown", - "id": "7a2e9992", - "metadata": {}, + "id": "7e0cfb4d", + "metadata": { + "editable": true + }, "source": [ "## Same code but now with momentum gradient descent" ] }, { "cell_type": "code", - "execution_count": 7, - "id": "b03e6d68", - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Own inversion\n", - "[[4.]\n", - " [3.]]\n", - "Eigenvalues of Hessian Matrix:[0.260336 4.68650355]\n", - "0 [-9.36327285] [-11.55399216]\n", - "1 [-0.03764415] [0.03029259]\n", - "2 [-0.03555301] [0.02860984]\n", - "3 [-0.03357804] [0.02702056]\n", - "4 [-0.03171277] [0.02551956]\n", - "5 [-0.02995112] [0.02410194]\n", - "6 [-0.02828733] [0.02276308]\n", - "7 [-0.02671597] [0.02149858]\n", - "8 [-0.02523189] [0.02030434]\n", - "9 [-0.02383026] [0.01917643]\n", - "10 [-0.02250648] [0.01811117]\n", - "11 [-0.02125624] [0.01710509]\n", - "12 [-0.02007545] [0.0161549]\n", - "13 [-0.01896026] [0.0152575]\n", - "14 [-0.01790701] [0.01440994]\n", - "15 [-0.01691228] [0.01360947]\n", - "16 [-0.0159728] [0.01285346]\n", - "17 [-0.01508551] [0.01213945]\n", - "18 [-0.0142475] [0.0114651]\n", - "19 [-0.01345605] [0.01082821]\n", - "20 [-0.01270857] [0.0102267]\n", - "21 [-0.0120026] [0.00965861]\n", - "22 [-0.01133586] [0.00912207]\n", - "23 [-0.01070615] [0.00861534]\n", - "24 [-0.01011142] [0.00813675]\n", - "25 [-0.00954973] [0.00768476]\n", - "26 [-0.00901924] [0.00725787]\n", - "27 [-0.00851822] [0.00685469]\n", - "28 [-0.00804503] [0.00647391]\n", - "29 [-0.00759813] [0.00611428]\n", - "theta from own gd\n", - "[[3.97243542]\n", - " [3.02218147]]\n", - "0 [-0.00717605] [0.00577464]\n", - "1 [-0.00677742] [0.00545385]\n", - "2 [-0.00628135] [0.00505466]\n", - "3 [-0.00578359] [0.00465411]\n", - "4 [-0.00531299] [0.00427541]\n", - "5 [-0.00487667] [0.0039243]\n", - "6 [-0.00447487] [0.00360097]\n", - "7 [-0.00410575] [0.00330394]\n", - "8 [-0.00376694] [0.00303129]\n", - "9 [-0.00345605] [0.00278111]\n", - "10 [-0.00317079] [0.00255157]\n", - "11 [-0.00290908] [0.00234096]\n", - "12 [-0.00266897] [0.00214774]\n", - "13 [-0.00244867] [0.00197047]\n", - "14 [-0.00224656] [0.00180783]\n", - "15 [-0.00206113] [0.00165861]\n", - "16 [-0.001891] [0.00152171]\n", - "17 [-0.00173492] [0.0013961]\n", - "18 [-0.00159172] [0.00128087]\n", - "19 [-0.00146034] [0.00117515]\n", - "20 [-0.0013398] [0.00107815]\n", - "21 [-0.00122922] [0.00098916]\n", - "22 [-0.00112776] [0.00090752]\n", - "23 [-0.00103467] [0.00083261]\n", - "24 [-0.00094927] [0.00076389]\n", - "25 [-0.00087092] [0.00070084]\n", - "26 [-0.00079903] [0.00064299]\n", - "27 [-0.00073308] [0.00058992]\n", - "28 [-0.00067257] [0.00054122]\n", - "29 [-0.00061706] [0.00049655]\n", - "theta from own gd wth momentum\n", - "[[3.9978254 ]\n", - " [3.00174992]]\n" - ] - } - ], + "execution_count": 29, + "id": "2ba6e4bb", + "metadata": { + "collapsed": false, + "editable": true + }, + "outputs": [], "source": [ "# Using Autograd to calculate gradients for OLS\n", "from random import random, seed\n", @@ -3851,37 +4258,23 @@ }, { "cell_type": "markdown", - "id": "d518ed04", - "metadata": {}, + "id": "671d74e3", + "metadata": { + "editable": true + }, "source": [ "## But noen of these can compete with Newton's method" ] }, { "cell_type": "code", - "execution_count": 3, - "id": "4f365ef3", - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Own inversion\n", - "[[4.13397208]\n", - " [2.83470568]]\n", - "Eigenvalues of Hessian Matrix:[0.28959697 4.92315107]\n", - "0 [-20.96904851] [-27.15233045]\n", - "1 [2.55464053e-14] [4.39327394e-14]\n", - "2 [5.68121938e-17] [-1.46584134e-16]\n", - "3 [9.86623977e-16] [1.19956128e-15]\n", - "4 [9.86623977e-16] [1.19956128e-15]\n", - "beta from own Newton code\n", - "[[4.13397208]\n", - " [2.83470568]]\n" - ] - } - ], + "execution_count": 30, + "id": "d282ba32", + "metadata": { + "collapsed": false, + "editable": true + }, + "outputs": [], "source": [ "# Using Newton's method\n", "from random import random, seed\n", @@ -3925,8 +4318,10 @@ }, { "cell_type": "markdown", - "id": "bb7a2b23", - "metadata": {}, + "id": "5504a30b", + "metadata": { + "editable": true + }, "source": [ "## Including Stochastic Gradient Descent with Autograd\n", "In this code we include the stochastic gradient descent approach discussed above. Note here that we specify which argument we are taking the derivative with respect to when using **autograd**." @@ -3935,8 +4330,11 @@ { "cell_type": "code", "execution_count": 31, - "id": "7be1bad2", - "metadata": {}, + "id": "3e36d821", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "# Using Autograd to calculate gradients using SGD\n", @@ -4016,8 +4414,10 @@ }, { "cell_type": "markdown", - "id": "998c490d", - "metadata": {}, + "id": "8115b858", + "metadata": { + "editable": true + }, "source": [ "## Same code but now with momentum gradient descent" ] @@ -4025,8 +4425,11 @@ { "cell_type": "code", "execution_count": 32, - "id": "8c8323d7", - "metadata": {}, + "id": "b4346837", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "# Using Autograd to calculate gradients using SGD\n", @@ -4100,21 +4503,26 @@ }, { "cell_type": "markdown", - "id": "e8e20118", - "metadata": {}, + "id": "944be76d", + "metadata": { + "editable": true + }, "source": [ - "## Same problem but now with AdaGrad" + "## Similar (second order function now) problem but now with AdaGrad" ] }, { "cell_type": "code", "execution_count": 33, - "id": "38b88401", - "metadata": {}, + "id": "2bd8357a", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ - "# Using Autograd to calculate gradients using SGD\n", - "# OLS example with AdaGrad\n", + "# Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent\n", + "# OLS example\n", "from random import random, seed\n", "import numpy as np\n", "import autograd.numpy as np\n", @@ -4125,61 +4533,65 @@ "def CostOLS(y,X,theta):\n", " return np.sum((y-X @ theta)**2)\n", "\n", - "n = 100\n", - "x = 2*np.random.rand(n,1)\n", - "y = 4+3*x+np.random.randn(n,1)\n", + "n = 10000\n", + "x = np.random.rand(n,1)\n", + "y = 2.0+3*x +4*x*x# +np.random.randn(n,1)\n", "\n", - "X = np.c_[np.ones((n,1)), x]\n", + "X = np.c_[np.ones((n,1)), x, x*x]\n", "XT_X = X.T @ X\n", "theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y)\n", "print(\"Own inversion\")\n", "print(theta_linreg)\n", - "# Hessian matrix\n", - "H = (2.0/n)* XT_X\n", - "EigValues, EigVectors = np.linalg.eig(H)\n", - "print(f\"Eigenvalues of Hessian Matrix:{EigValues}\")\n", "\n", - "theta = np.random.randn(2,1)\n", - "eta = 1.0/np.max(EigValues)\n", - "Niterations = 100\n", "\n", "# Note that we request the derivative wrt third argument (theta, 2 here)\n", "training_gradient = grad(CostOLS,2)\n", - "\n", - "for iter in range(Niterations):\n", - " gradients = (1.0/n)*training_gradient(y, X, theta)\n", - " theta -= eta*gradients\n", - "print(\"theta from own gd\")\n", - "print(theta)\n", - "print(np.size(gradients))\n", - "\n", + "# Define parameters for Stochastic Gradient Descent\n", "n_epochs = 50\n", "M = 5 #size of each minibatch\n", "m = int(n/M) #number of minibatches\n", - "theta = np.random.randn(2,1)\n", + "# Guess for unknown parameters theta\n", + "theta = np.random.randn(3,1)\n", "\n", - "# Including AdaGrad\n", - "delta = 0.000001\n", - "r = [0.0 for _ in range(gradients.shape[0])]\n", + "# Value for learning rate\n", + "eta = 0.01\n", + "# Including AdaGrad parameter to avoid possible division by zero\n", + "delta = 1e-8\n", "for epoch in range(n_epochs):\n", + " # The outer product is calculated from scratch for each epoch\n", + " Giter = np.zeros(shape=(3,3))\n", " for i in range(m):\n", " random_index = M*np.random.randint(m)\n", " xi = X[random_index:random_index+M]\n", " yi = y[random_index:random_index+M]\n", " gradients = (1.0/M)*training_gradient(yi, xi, theta)\n", - " # calculate squared gradient by Hadamard multiplication\n", - " r -= gradients*gradients\n", + "\t# Calculate the outer product of the gradients\n", + " Giter +=gradients @ gradients.T\n", + "\t# Simpler algorithm with only diagonal elements\n", + " Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))]\n", " # compute update\n", - " update = 1.0/(delta+np.sqrt(r))*gradients\n", - " theta -= eta*update\n", + " update = np.multiply(Ginverse,gradients)\n", + " theta -= update\n", "print(\"theta from own AdaGrad\")\n", "print(theta)" ] }, { "cell_type": "markdown", - "id": "0989e18a", - "metadata": {}, + "id": "8069cc8c", + "metadata": { + "editable": true + }, + "source": [ + "Running this code we note an almost perfect agreement with the results from matrix inversion." + ] + }, + { + "cell_type": "markdown", + "id": "990df1f5", + "metadata": { + "editable": true + }, "source": [ "## And Logistic Regression" ] @@ -4187,8 +4599,11 @@ { "cell_type": "code", "execution_count": 34, - "id": "c6622d90", - "metadata": {}, + "id": "ce42a59a", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import autograd.numpy as np\n", @@ -4228,8 +4643,10 @@ }, { "cell_type": "markdown", - "id": "442a319f", - "metadata": {}, + "id": "d4e18db5", + "metadata": { + "editable": true + }, "source": [ "## Introducing [JAX](https://jax.readthedocs.io/en/latest/)\n", "\n", @@ -4245,8 +4662,11 @@ { "cell_type": "code", "execution_count": 35, - "id": "8c6259ef", - "metadata": {}, + "id": "2fe09a96", + "metadata": { + "collapsed": false, + "editable": true + }, "outputs": [], "source": [ "import jax.numpy as jnp\n", @@ -4262,8 +4682,10 @@ }, { "cell_type": "markdown", - "id": "4495aa58", - "metadata": {}, + "id": "c1960c06", + "metadata": { + "editable": true + }, "source": [ "## Weekend challenge\n", "\n", @@ -4275,25 +4697,7 @@ ] } ], - "metadata": { - "kernelspec": { - "display_name": "Python 3 (ipykernel)", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.9.14" - } - }, + "metadata": {}, "nbformat": 4, "nbformat_minor": 5 } diff --git a/doc/src/week39/adagrad.py b/doc/src/week39/adagrad.py index a847267a7..e072281a3 100644 --- a/doc/src/week39/adagrad.py +++ b/doc/src/week39/adagrad.py @@ -4,53 +4,28 @@ from random import random, seed import numpy as np import autograd.numpy as np import matplotlib.pyplot as plt -from autograd import grad - -# Note change from previous example -def CostOLS(y,X,theta): - return np.sum((y-X @ theta)**2) n = 10000 x = np.random.rand(n,1) -y = 2.0+3*x +4*x*x# +np.random.randn(n,1) - +y = 4*x+3*x*x +# Setting up Design matrix X = np.c_[np.ones((n,1)), x, x*x] -XT_X = X.T @ X -theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y) +XTX = X.T @ X +XTy = X.T @ y +theta_linreg = np.linalg.pinv(XTX) @ (XTy) print("Own inversion") print(theta_linreg) -# Note that we request the derivative wrt third argument (theta, 2 here) -training_gradient = grad(CostOLS,2) -# Define parameters for Stochastic Gradient Descent -n_epochs = 50 -M = 5 #size of each minibatch -m = int(n/M) #number of minibatches -# Guess for unknown parameters theta -theta = np.random.randn(3,1) -gradients = np.zeros(theta.shape) -r = gradients*gradients -print(r.shape) -print(gradients.shape) -# Value for learning rate +beta = np.random.randn(3,1) eta = 0.01 -# Including AdaGrad parameter to avoid possible division by zero -delta = 1e-8 -for epoch in range(n_epochs): - for i in range(m): - random_index = M*np.random.randint(m) - xi = X[random_index:random_index+M] - yi = y[random_index:random_index+M] - gradients = (1.0/M)*training_gradient(yi, xi, theta) - # calculate squared gradient by Hadamard multiplication -# r += (gradients*gradients) - r = np.sum(gradients*gradients) - # compute update - update = 1.0/(delta+np.sqrt(r))*gradients - theta -= eta*update -print("theta from own AdaGrad") -print(theta) - - +delta = 1e-8 +Niterations = 10000 +Giter = np.zeros(shape=(3,3)) +for iter in range(Niterations): + gradient = (2.0/n)*(XTX @ beta - XTy) + Giter +=gradient @ gradient.T + Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))] + beta -= np.multiply(Ginverse,gradient) +print("Optimal parameters with AdaGrad",beta) diff --git a/doc/src/week39/adagradSGD.py b/doc/src/week39/adagradSGD.py new file mode 100644 index 000000000..157638fc8 --- /dev/null +++ b/doc/src/week39/adagradSGD.py @@ -0,0 +1,55 @@ +# Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent +# OLS example +from random import random, seed +import numpy as np +import autograd.numpy as np +import matplotlib.pyplot as plt +from autograd import grad + +# Note change from previous example +def CostOLS(y,X,theta): + return np.sum((y-X @ theta)**2) + +n = 10000 +x = np.random.rand(n,1) +y = 2.0+3*x +4*x*x# +np.random.randn(n,1) + +X = np.c_[np.ones((n,1)), x, x*x] +XT_X = X.T @ X +theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y) +print("Own inversion") +print(theta_linreg) + + +# Note that we request the derivative wrt third argument (theta, 2 here) +training_gradient = grad(CostOLS,2) +# Define parameters for Stochastic Gradient Descent +n_epochs = 50 +M = 5 #size of each minibatch +m = int(n/M) #number of minibatches +# Guess for unknown parameters theta +theta = np.random.randn(3,1) + +# Value for learning rate +eta = 0.01 +# Including AdaGrad parameter to avoid possible division by zero +delta = 1e-8 +for epoch in range(n_epochs): + Giter = np.zeros(shape=(3,3)) + for i in range(m): + random_index = M*np.random.randint(m) + xi = X[random_index:random_index+M] + yi = y[random_index:random_index+M] + gradients = (1.0/M)*training_gradient(yi, xi, theta) + Giter +=gradients @ gradients.T + Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))] + + # calculate squared gradient by Hadamard multiplication +# r += (gradients*gradients) +# r = np.sum(gradients*gradients) + # compute update + update = np.multiply(Ginverse,gradients) + theta -= update +print("theta from own AdaGrad") +print(theta) + diff --git a/doc/src/week39/week39.do.txt b/doc/src/week39/week39.do.txt index 33d3d3cf3..786fa02c9 100644 --- a/doc/src/week39/week39.do.txt +++ b/doc/src/week39/week39.do.txt @@ -2432,10 +2432,10 @@ print(theta) !split -===== Same problem but now with AdaGrad ===== +===== Similar (second order function now) problem but now with AdaGrad ===== !bc pycod -# Using Autograd to calculate gradients using SGD -# OLS example with AdaGrad +# Using Autograd to calculate gradients using AdaGrad and Stochastic Gradient descent +# OLS example from random import random, seed import numpy as np import autograd.numpy as np @@ -2446,57 +2446,50 @@ from autograd import grad def CostOLS(y,X,theta): return np.sum((y-X @ theta)**2) -n = 100 -x = 2*np.random.rand(n,1) -y = 4+3*x+np.random.randn(n,1) +n = 10000 +x = np.random.rand(n,1) +y = 2.0+3*x +4*x*x# +np.random.randn(n,1) -X = np.c_[np.ones((n,1)), x] +X = np.c_[np.ones((n,1)), x, x*x] XT_X = X.T @ X theta_linreg = np.linalg.pinv(XT_X) @ (X.T @ y) print("Own inversion") print(theta_linreg) -# Hessian matrix -H = (2.0/n)* XT_X -EigValues, EigVectors = np.linalg.eig(H) -print(f"Eigenvalues of Hessian Matrix:{EigValues}") -theta = np.random.randn(2,1) -eta = 1.0/np.max(EigValues) -Niterations = 100 # Note that we request the derivative wrt third argument (theta, 2 here) training_gradient = grad(CostOLS,2) - -for iter in range(Niterations): - gradients = (1.0/n)*training_gradient(y, X, theta) - theta -= eta*gradients -print("theta from own gd") -print(theta) -print(np.size(gradients)) - +# Define parameters for Stochastic Gradient Descent n_epochs = 50 M = 5 #size of each minibatch m = int(n/M) #number of minibatches -theta = np.random.randn(2,1) +# Guess for unknown parameters theta +theta = np.random.randn(3,1) -# Including AdaGrad -delta = 0.000001 -r = [0.0 for _ in range(gradients.shape[0])] +# Value for learning rate +eta = 0.01 +# Including AdaGrad parameter to avoid possible division by zero +delta = 1e-8 for epoch in range(n_epochs): + # The outer product is calculated from scratch for each epoch + Giter = np.zeros(shape=(3,3)) for i in range(m): random_index = M*np.random.randint(m) xi = X[random_index:random_index+M] yi = y[random_index:random_index+M] gradients = (1.0/M)*training_gradient(yi, xi, theta) - # calculate squared gradient by Hadamard multiplication - r -= gradients*gradients + # Calculate the outer product of the gradients + Giter +=gradients @ gradients.T + # Simpler algorithm with only diagonal elements + Ginverse = np.c_[eta/(delta+np.sqrt(np.diagonal(Giter)))] # compute update - update = 1.0/(delta+np.sqrt(r))*gradients - theta -= eta*update + update = np.multiply(Ginverse,gradients) + theta -= update print("theta from own AdaGrad") print(theta) !ec +Running this code we note an almost perfect agreement with the results from matrix inversion. !split ===== And Logistic Regression =====