Search code examples
pythonpython-2.7sympygradient-descent

multi-variable Gradient descent using Sympy


I'm practicing on Gradient descent algorithm implementation for two variables in Sympy library in Python 2.7.

My goal is to find minimum of two variable function using vector of derivatives according to following steps:

  1. For function f(a,b) of two varibale define the Matrix of first partial differentials - M.
  2. Then, I pass the the statring point of a,b (for instance V0 = (1.0,1.0)) into M and multiply this by step - this gives M0 * 𝜼 matrix.
  3. Next, substract the result calculated above from vector of starting values - V0. This gives a new vector of variables a,b - V1.
  4. Finally, values of V1 are put in M again. If results of matrix M are less than epsilon = > continue the iterations.

The code and steps description is attached. I suppose the problem lies in while loop, since it gives the worng values and algorithm takes the first i iteration only. enter image description here

Could you please advise on improvements?

import sympy as sp
from sympy import *
from sympy import lambdify

a=Symbol('a')
b=Symbol('b')
alpha=Symbol('alpha')

def test_f(a, b):
    # define a test function of a,b
    test_f=4*a**2 + 25*b**2 
    return test_f

def deriv_matrix(f, a, b):
    # define the matrix of derivatives
    d1=diff(f(a,b), a, 1)
    d2=diff(f(a,b), b, 1)
    M=Matrix([d1,d2])
    return M

epsilon=0.02
alpha=0.1
i=1 # strat the iteration from 1
vector1=Matrix([1,1]) # starting point 

while (i<100) & (f(vector1[0], vector1[1]).all()> epsilon):
    f = lambdify((a,b), deriv_matrix(test_f, a, b))
    vector=vector1
    N=-(alpha/i)*f(vector[0],vector[1])
    vector1=vector1+N
    i+=i
print vector1

Solution

  • The code has too many problems to list. As a small sample: f is undefined, but apparently its return value has method .all() (which in NumPy returns a bool), which is then compared to epsilon? Makes no sense. Generally, when using SymPy, one does not need to encode symbolic functions as Python functions. They are represented as SymPy expressions, which can be lambdified when speed is desirable, or evaluated directly with subs when speed is not an issue. For example:

    from sympy import *    
    a, b = symbols('a b')
    f = 4*a**2 + 25*b**2     # this is a SymPy expression, not a function
    grad_f = lambdify((a, b), derive_by_array(f, (a, b)))   # lambda from an expression for graduent
    
    epsilon = 0.02
    alpha = 0.01      # the value of 0.1 is too large for convergence here
    vector = Matrix([1, 1])  
    
    for i in range(100):    #  to avoid infinite loop
        grad_value = Matrix(grad_f(vector[0], vector[1]))
        vector += -alpha*grad_value
        if grad_value.norm() < epsilon:   #  gradient is small, we are done 
            break          
    
    if grad_value.norm() < epsilon:
        print('Converged to {}'.format(vector))    # only print vector if successful
    else:
        print('Failed to converge')