Approximating derivatives

In this workbook, we'll see a few ideas for doing derivative computations. You'll see more on the homework.

Automatic differentiation

The first type of derivative computation is a source code transform. This is called automatic differentiation, and consists of taking in a program and outputing a program. We won't talk too much about this idea, but it's used frequently in optimization. Purdue Prof. Alex Pothen is currently doing research on these ideas.

Julia has a package to do this (actually, a few). There are some bugs with the current one that mean it isn't quite as useful.

In [31]:
Pkg.add("ReverseDiffSource")
using ReverseDiffSource
using Plots
using Interact
gadfly();
INFO: Nothing to be done
Out[31]:
Plots.GadflyPackage()
INFO: METADATA is out-of-date — you may not have the latest version of ReverseDiffSource
INFO: Use `Pkg.update()` to get the latest versions of your packages
In [34]:
rdiff( :((1 - x[1])^2 + 100(x[2] - x[1]^2)^2) , x=ones(2))
Out[34]:
quote 
    _tmp1 = zeros(size(x))
    _tmp2 = 1 - x[1]
    _tmp3 = x[2] - x[1] ^ 2
    _tmp4 = 2 * (_tmp3 * 100.0)
    _tmp1[1] = 2 * (x[1] * -_tmp4) + -(2_tmp2)
    _tmp1[2] = _tmp1[2] + _tmp4
    (_tmp2 ^ 2 + 100 * _tmp3 ^ 2,_tmp1)
end
In [35]:
rdiff( :(exp(x^2)*(sin(100x)^2)) , x=1.)
Out[35]:
quote 
    _tmp1 = 100x
    _tmp2 = exp(x ^ 2)
    _tmp3 = sin(_tmp1)
    _tmp4 = _tmp2 * _tmp3 ^ 2
    (_tmp4,100 * (cos(_tmp1) * (2 * (_tmp3 * _tmp2))) + 2 * (x * _tmp4))
end

You should be able to directly evaluate these functions, but I couldn't get that working before class.

In [ ]:
In [36]:
# This should work, but as of 2016-02-23, it doesn't :(
rosenbrock(x) = (1 - x[1])^2 + 100(x[2] - x[1]^2)^2
rdiff( rosenbrock , x=ones(2))  # 'x=2.' indicates the type of x to rdiff
LoadError: [tograph] unmanaged type rosenbrock (Function)
while loading In[36], in expression starting on line 3

 in error at /Applications/Julia-0.4.3.app/Contents/Resources/julia/lib/julia/sys.dylib
 in explore at /Users/dgleich/.julia/ReverseDiffSource/src/tograph.jl:45
 in tograph at /Users/dgleich/.julia/ReverseDiffSource/src/tograph.jl:257
 in rdiff at /Users/dgleich/.julia/ReverseDiffSource/src/rdiff.jl:25

Finite differences

We can also use a finite difference approximation. The idea is that the derivative is formally $$ f'(x) = \lim_{h\to 0} \frac{f(x+h) - f(x)}{h} $$ and so we should be able to get a good approximation by using $$ f'(x) \approx \frac{f(x+h) - f(x)}{h} $$ with $h$ very small

In [37]:
h = 0.01
expderiv(x) = (exp(x+h)-exp(x))./h
xx = linspace(-2,2,1000)
plot(xx,exp(xx))
plot!(xx,expderiv(xx))
title!(string(norm(exp(xx) - expderiv(xx),Inf)))
-2 -1 0 1 2 y1 y2 0 2 4 6 8 0.037068739923602934
Out[37]:
In [38]:
@manipulate for h=logspace(-16,-1,16)
    expderiv(x) = (exp(x+h)-exp(x))./h
    xx = linspace(-2,2,1000)
    plot(xx,exp(xx))
    plot!(xx,expderiv(xx))
    title!(string(norm(exp(xx) - expderiv(xx),Inf))).o
end
Out[38]:
-7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 -6.0 -5.8 -5.6 -5.4 -5.2 -5.0 -4.8 -4.6 -4.4 -4.2 -4.0 -3.8 -3.6 -3.4 -3.2 -3.0 -2.8 -2.6 -2.4 -2.2 -2.0 -1.8 -1.6 -1.4 -1.2 -1.0 -0.8 -0.6 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2.4 2.6 2.8 3.0 3.2 3.4 3.6 3.8 4.0 4.2 4.4 4.6 4.8 5.0 5.2 5.4 5.6 5.8 6.0 -6 -3 0 3 6 -6.0 -5.5 -5.0 -4.5 -4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 y1 y2 -10 -8 -6 -4 -2 0 2 4 6 8 10 12 14 16 18 -8.0 -7.5 -7.0 -6.5 -6.0 -5.5 -5.0 -4.5 -4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 6.5 7.0 7.5 8.0 8.5 9.0 9.5 10.0 10.5 11.0 11.5 12.0 12.5 13.0 13.5 14.0 14.5 15.0 15.5 16.0 -10 0 10 20 -8.0 -7.5 -7.0 -6.5 -6.0 -5.5 -5.0 -4.5 -4.0 -3.5 -3.0 -2.5 -2.0 -1.5 -1.0 -0.5 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 6.5 7.0 7.5 8.0 8.5 9.0 9.5 10.0 10.5 11.0 11.5 12.0 12.5 13.0 13.5 14.0 14.5 15.0 15.5 16.0 7.38905609893065