6

Currently I have two numpy arrays: x and y of the same size.

I would like to write a function (possibly calling numpy/scipy... functions if they exist):

def derivative(x, y, n = 1):
    # something
    return result

where result is a numpy array of the same size of x and containing the value of the n-th derivative of y regarding to x (I would like the derivative to be evaluated using several values of y in order to avoid non-smooth results).

Vincent
  • 57,703
  • 61
  • 205
  • 388
  • 2
    Possible solution: http://stackoverflow.com/a/18993405/190597 – unutbu Nov 18 '13 at 09:44
  • How to use gaussian filter when the y data are not distributed regularly (dx between two consecutive points is not constant)? – Vincent Nov 18 '13 at 10:24
  • [This article](http://en.wikipedia.org/wiki/Finite_difference_coefficients) gives the coefficients you would use to compute finite differences for regularly spaced grids. That article links to [a paper by Bengt Fornberg](http://www.ams.org/journals/mcom/1988-51-184/S0025-5718-1988-0935077-0/home.html) which shows how to compute the coefficients for **arbitrarily spaced grids**. Here is a link to [the PDF](http://www.ams.org/journals/mcom/1988-51-184/S0025-5718-1988-0935077-0/S0025-5718-1988-0935077-0.pdf). – unutbu Nov 18 '13 at 11:01
  • If you know that your data can be modeled by some function, then you can do better by [fitting your data to that function](http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html). If the function is differentiable, then you can compute the derivative analytically. If it is not differentiable, you could sample the function on a regularly spaced grid, compute the finite differences, and interpolate to find an approximation to the derivative at your irregular grid points. – unutbu Nov 18 '13 at 13:41

3 Answers3

6

This is not a simple problem, but there are a lot of methods that have been devised to handle it. One simple solution is to use finite difference methods. The command numpy.diff() uses finite differencing where you can specify the order of the derivative.

Wikipedia also has a page that lists the needed finite differencing coefficients for different derivatives of different accuracies. If the numpy function doesn't do what you want.

Depending on your application you can also use scipy.fftpack.diff which uses a completely different technique to do the same thing. Though your function needs a well defined Fourier transform.

There are lots and lots and lots of variants (e.g. summation by parts, finite differencing operators, or operators designed to preserve known evolution constants in your system of equations) on both of the two ideas above. What you should do will depend a great deal on what the problem is that you are trying to solve.

The good thing is that there is a lot of work has been done in this field. The Wikipedia page for Numerical Differentiation has some resources (though it is focused on finite differencing techniques).

Stephen Rauch
  • 47,830
  • 31
  • 106
  • 135
Ben Whale
  • 493
  • 2
  • 9
  • `scipy.fftpack.diff` seems extremely useful. However, if I have function given as array `ampl` on a time grid `tgrid`, I found (through trial and error) that I have to do `T = tgrid[-1] - tgrid[0]; deriv = scipy.fftpack.diff(ampl) * (2.0*pi / T)`. I was thinking that the same effect could probably be achieved by passing the correct `period` to the `diff` function, but I couldn't figure out the correct value. Is there an example for this use case somewhere? – Michael Goerz Aug 12 '14 at 00:40
1

The findiff project is a Python package that can do derivatives of arrays of any dimension with any desired accuracy order (of course depending on your hardware restrictions). It can handle arrays on uniform as well as non-uniform grids and also create generalizations of derivatives, i.e. general linear combinations of partial derivatives with constant and variable coefficients.

maroba
  • 45
  • 6
0

Would something like this solve your problem?

def get_inflection_points(arr, n=1):
    """
    returns inflextion points from array
        arr: array
        n: n-th discrete difference
    """
    inflections = []
    dx = 0
    for i, x in enumerate(np.diff(arr, n)):
        if x >= dx and i > 0:
            inflections.append(i*n)
        dx = x
    return inflections

Yaakov Bressler
  • 9,056
  • 2
  • 45
  • 69