Differentiation of functions in Python and R

Calculating value of derivative function is a part of many applied tasks such
– gradient descent/ascent for finding minimum/maximum of a function
– Newton’s method for finding approximations to the roots (or zeroes) of a real-valued function

There are three main ways to calculate derivatives
1. Finite differences i.e. calculate approximation using (f(x+Δx)-f(x))/Δx when Δx is close to zero.
2. Symbolic differentiation – generates exact formula for derivative function.
3. Automatic differentiation – generates evaluations (and not formulas) of the derivatives. All intermediate expressions are evaluated as soon as possible.
A brilliant explanation of automatic differentiation may be found here in a paper by Warwick Tucker “One-dimensional, first-order, and Taylor-mode automatic differentiation with programming examples”.

So let’s have a closer look on some technics for derivation in Python/R.

  • Python
  • Symbolic differentiation (computation) is implemented in a number of libraries for python. I’ll use sympy for demonstration.

    from sympy.core.function import diff
    from sympy.core.symbol import Symbol
    f_str = '(x - 1)*(x - 10)'
    x = Symbol('x')
    exec 'f0 = lambda x: ' + f_str
    exec 'f = ' + f_str
    print '*** functions ***'
    print f0
    print f
    print f.diff(x) # or diff(f) or diff(f_str)
    print '*** types ***'
    print type(f0)
    print type(f)
    print type(diff(f_str))
    print '*** calculated values ***'
    exec 'f1 = lambda x: ' + str(f.diff(x))
    f0_res = [f0(i) for i in range(1, 5+1)]
    f1_res = [f1(i) for i in range(1, 5+1)]
    print f0_res
    print f1_res


    *** functions ***
    <function <lambda> at 0x0000000008C4D278>
    (x - 10)*(x - 1)
    2*x - 11
    *** types ***
    <type 'function'>
    <class 'sympy.core.mul.Mul'>
    <class 'sympy.core.add.Add'>
    *** calculated values ***
    [0, -8, -14, -18, -20]
    [-9, -7, -5, -3, -1]

    To be able to get derivative we need to define variable of a function (using Symbol function).
    After this it’s possible to call diff method of ‘sympy.core.mul.Mul’ object to get symbolic representation of a derivative. On the other hand it can also be calculated from string.
    Crucial point is that f and it’s derivative are not callable objects of sympy classes Mul/Add.
    In order to work with original function and its derivative like with other mathematical functions we need to define them as lambda function using correspondent strings.
    In such case correspondent values can be calculated simply like f0(i) and f1(1).
    Another approach to calculate values would be to replace variable in a string and call eval function.

  • R
  • In case of R both symbolic and automatic (algorithmic) derivatives are built into engine.

    f_str <- '(x - 1)*(x - 10)'
    f_expr <- parse(text = f_str)
    print ('symbolic differentiation')
    print (D(f_expr, "x"))
    print ('automatic differentiation')
    f <- deriv(f_expr, "x", func=TRUE)
    print (f)
    print (f(0))
    print ('calculated values')
    f0_res = c()
    f1_res = c()
    for(i in 1:5) {f0_res[i] <- f(i)[1]}
    for(i in 1:5) {f1_res[i] <- attributes(f(i))$gradient}
    print (f0_res)
    print (f1_res)


    [1] "symbolic differentiation"
    (x - 10) + (x - 1)
    [1] "automatic differentiation"
    function (x)
        .expr1 <- x - 1
        .expr2 <- x - 10
        .value <- .expr1 * .expr2
        .grad <- array(0, c(length(.value), 1L), list(NULL, c("x")))
        .grad[, "x"] <- .expr2 + .expr1
        attr(.value, "gradient") <- .grad
    [1] 10
    [1,] -11
    [1] "calculated values"
    [1]   0  -8 -14 -18 -20
    [1] -9 -7 -5 -3 -1

    function D returns symbolic expression for derivative. As you can see it’s written in a different form compare to the symbolic result by sympy.
    deriv returns function (func=TRUE) which can be used for computing the original expression and its derivatives, simultaneously.
    R language does not support such concept as lambda functions so in order to calculate derivative for given point using symbolic differentiation one may need to replace variable in a derived symbolic string with specific value call eval function.

  • Quick test AD vs SD. Python
  • from ad import adnumber
    from ad.admath import *
    def AD(iter):
        res = []
        for i in range(1, iter+1):
            x = adnumber(i)
            y = sin(x**3+10*x) + exp(cos(x)*cos(x)+5)
        return res    
    from sympy.core.function import diff
    from sympy.core.symbol import Symbol
    def SD(iter):
        res = []
        x = Symbol('x')
        exec 'fd = lambda x: ' + str(diff('sin(x**3+10*x) + exp(cos(x)*cos(x)+5)', x))
        for i in range(1, iter+1):
        return res
    %timeit AD(1)
    print AD(1)
    %timeit SD(1)
    print SD(1)
    %timeit AD(1000)
    %timeit SD(1000)
    1000 loops, best of 3: 482 µs per loop
    100 loops, best of 3: 3.22 ms per loop
    1 loop, best of 3: 478 ms per loop
    10 loops, best of 3: 32.4 ms per loop

    Using lambda functions allows calculate derivative symbolic computations only once and than evaluate it for all required values.
    So single evaluation took almost 10 times longer in case of SD while 1000 evaluations completed 10 times faster.
    I used ad library for automatic differentiation, it’s implemented purely on Python.

    Leave a Reply

    Fill in your details below or click an icon to log in:

    WordPress.com Logo

    You are commenting using your WordPress.com account. Log Out /  Change )

    Google photo

    You are commenting using your Google account. Log Out /  Change )

    Twitter picture

    You are commenting using your Twitter account. Log Out /  Change )

    Facebook photo

    You are commenting using your Facebook account. Log Out /  Change )

    Connecting to %s