Higher-order multivariate derivatives

Question

Is there a packaged way to compute higher-order multivariate derivatives (using finite differences, not symbolic calculations) in Python?

For example, if f computes the function cos(x)*y from R^2 to R, i.e., f takes numpy arrays of shape 2 and returns floats (or arrays of shape ()), is there a function partial such that partial([2,1])(f) computes the function (d^2/dx^2)(d/dy)f = -cos(x)*1 from R^2 to R, e.g.

np.isclose(partial([2,1])(f)(0),1.0)

There is a host of finite difference tools in the default libraries (so much for ''one way to do it''):

https://docs.scipy.org/doc/scipy/reference/generated/scipy.misc.derivative.html https://docs.scipy.org/doc/numpy/reference/generated/numpy.gradient.html https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.optimize.approx_fprime.html https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.misc.central_diff_weights.html https://github.com/statsmodels/statsmodels/blob/master/statsmodels/tools/numdiff.py

However, none of them can handle both multivariate functions and higher-order derivatives, i.e., none of them can handle the job above.

(There is also https://pypi.org/project/numdifftools/ but it seems to fall just short of doing what I need. The authors haven't responded to my question.)

It's easy to write the tool myself. However, it seems to be surprisingly difficult to do so well, i.e., in an accurate and stable fashion. For example, a straightforward recursive implementation is unstable for mesh-widths smaller than 1e-3, even for simple functions like the one above and for mixed derivatives of only second order.

PS: I am not asking for finite differences of an array (https://docs.scipy.org/doc/numpy/reference/generated/numpy.diff.html and https://github.com/maroba/findiff do that). I need to evaluate derivatives at arbitrary points without computing the values of the function on a full cartesian grid.

score 0 · Answer 1 · answered Jul 15 '19 at 12:41

0

There's scipy.optimize._approx_derivative which does it. This is however not a public function, so if you end up using it, you're on your own.

answered Jul 15 '19 at 12:41

ev-br

24,968
9
65
78

Shouldn't it be `scipy.optimize._numdiff.approx_derivative`? – joni Jul 15 '19 at 18:04
How do you compute mixed/cross derivatives with that? – Bananach Aug 08 '19 at 19:09

score 0 · Answer 2 · answered Jul 15 '19 at 17:07

0

You can write a wrapper and use the scipy.misc.derivative function.

For a simple x, y derivative, use this answer

For a vector derivative, you can define g(t, (x, y), vector) = f((x,y) + t * vector), which yields g'(0, args=((x,y), vector)) = directional derivative in the chosen vector direction

from scipy.misc import derivative

f = lambda x: x[0] * np.cos(x[1])

def vector_derivative(f, x0, vector, delta=1):
    def wrapper(x, x0, vector):
        return f(np.asarray(x0) + x * np.asarray(vector))
    return derivative(wrapper, 0, args=(x0, vector), dx=delta)

vector_derivative(f, [1, np.pi/2], [0, 1], delta=0.01)
>> -0.9999833334166673

answered Jul 15 '19 at 17:07

Tarifazo

4,118
1
9
22

This can't handle mixed derivatives, can it? – Bananach Jul 15 '19 at 19:22
The approach here is to transform the directional derivative (x, y, or generic direction) into a scalar function. You can always vectorize the function to work with different directions. – Tarifazo Jul 15 '19 at 20:51
I'm assuming by mixed you mean cross derivatives. In that case you need to wrap it again (but it is easier to do it with the approach in the link provided, since this approach is aimed to being generic in terms of direction) – Tarifazo Jul 15 '19 at 20:54

Higher-order multivariate derivatives

2 Answers2