# High-order derivatives of sampled data

**URL:** <https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045>\
**Category:** General Usage\
**Tags:** finitediff, splines\
**Created:** [August 26, 2021, 3:52pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045 "2021-08-26T15:52:34Z")\
**Posts on this page:** 10\
**Page:** 1

<div class="post-metadata">

**Author:** ![DanielIngraham](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielingraham/32/13055_2.png) [@DanielIngraham](https://discourse.julialang.org/u/DanielIngraham)\
**Post date:** [August 26, 2021, 3:52pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/1 "2021-08-26T15:52:34Z")

</div>

Hi,

I have discrete data sampled at a constant rate. I need the first, second, and third derivatives with respect to time. Is there a Julia package that could help me do that? Things that I’ve considered so far:

- Use finite differences from either [FiniteDiff.jl](https://github.com/JuliaDiff/FiniteDiff.jl) or [FiniteDifferences.jl](https://github.com/JuliaDiff/FiniteDifferences.jl), but both appear to differentiate Julia functions, not discrete data.
- Recursively use spline interpolation and then differentiate the spline. That led me to a discontinuous second derivative, which is no good.
- Use [ApproxFun.jl](https://github.com/JuliaApproximation/ApproxFun.jl), but that again seems to work only for continuous functions, not discrete data.

Anyone have any better ideas?

Thanks!

Daniel

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [August 26, 2021, 3:59pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/2 "2021-08-26T15:59:50Z")

</div>

> [@DanielIngraham](#):
>
> I have discrete data sampled at a constant rate. I need the first, second, and third derivatives with respect to time. Is there a Julia package that could help me do that? Things that I’ve considered so far:

Is O(\Delta t^2) accuracy sufficient? If so, you could just use low-order finite-difference formulas, e.g. see [here](https://www.mech.kth.se/~ardeshir/courses/literature/fd.pdf).

---

<div class="post-metadata">

**Author:** ![DanielIngraham](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielingraham/32/13055_2.png) [@DanielIngraham](https://discourse.julialang.org/u/DanielIngraham)\
**Post date:** [August 26, 2021, 4:43pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/3 "2021-08-26T16:43:14Z")

</div>

Yes, I think that finite differences would work fine—I was just surprised that I couldn’t find an existing implementation in the Julia ecosystem.

Thanks!

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [August 26, 2021, 6:47pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/4 "2021-08-26T18:47:06Z")

</div>

@DanielIngraham, have you tried [Dierckx.jl](https://github.com/kbarbary/Dierckx.jl/blob/master/README.md) or [BSplineKit.jl](https://jipolanco.github.io/BSplineKit.jl/stable/) for this?

A very simple example is provided below using _Dierckx_, but would recommend _BSplineKit_ as it is native Julia and greatly supported by the author:

```julia
using Dierckx
f(x) = exp(-2x)
x = 0:0.1:2; y = f.(x)
spl = Spline1D(x, y; k=5)
D1f = derivative(spl, x; nu=1)
D2f = derivative(spl, x; nu=2)
D3f = derivative(spl, x; nu=3)

```

 ![Dierckx_high_order_derivatives](https://global.discourse-cdn.com/julialang/original/3X/9/b/9bf5171508d5d37b10b9e839eb0447074f61d8cc.png)

---

<div class="post-metadata">

**Author:** ![arunsm](https://avatars.discourse-cdn.com/v4/letter/a/91b2a8/32.png) [@arunsm](https://discourse.julialang.org/u/arunsm)\
**Post date:** [August 26, 2021, 8:28pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/5 "2021-08-26T20:28:18Z")

</div>

You could take derivatives using a Gaussian Process regression, if you select a differentiable kernel. You could the Gaussian process julia package to do the regression.

[https://stats.stackexchange.com/questions/373446/computing-gradients-via-gaussian-process-regression](https://stats.stackexchange.com/questions/373446/computing-gradients-via-gaussian-process-regression)

---

<div class="post-metadata">

**Author:** ![rafael.guerra](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/rafael.guerra/32/216610_2.png) [@rafael.guerra](https://discourse.julialang.org/u/rafael.guerra)\
**Post date:** [August 26, 2021, 8:41pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/6 "2021-08-26T20:41:14Z")

</div>

> [@arunsm](#):
>
> Gaussian Process regression

Sounds like using a massive sledgehammer to crack a nut, but might be wrong.

---

<div class="post-metadata">

**Author:** ![jipolanco](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jipolanco/32/12129_2.png) [@jipolanco](https://discourse.julialang.org/u/jipolanco)\
**Post date:** [August 27, 2021, 6:57am UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/7 "2021-08-27T06:57:35Z")

</div>

Below is @rafael.guerra’s example using BSplineKit.

```julia
using BSplineKit
f(x) = exp(-2x)
x = 0:0.1:2
y = f.(x)
spl = interpolate(x, y, BSplineOrder(6))
D1f = diff(spl, Derivative(1))
D2f = diff(spl, Derivative(2))
D3f = diff(spl, Derivative(3))

using CairoMakie
D3f_exact(x) = -8f(x)
fig = Figure(resolution = (400, 300))
ax = Axis(fig[1, 1])
plot!(ax, 0..2, D3f_exact; label = "Exact")
scatter!(ax, x, D3f.(x); color = :red, label = "BSplineKit f‴(x)")
axislegend(position = :rb)

```

![interp](https://global.discourse-cdn.com/julialang/original/3X/9/2/92c472430554b2da350e74a6ddf3f4f69d5100e1.png)

Here I used B-splines of order k = 6 (piecewise polynomials of degree 5), but the order can be arbitrarily increased or decreased.

---

<div class="post-metadata">

**Author:** ![arunsm](https://avatars.discourse-cdn.com/v4/letter/a/91b2a8/32.png) [@arunsm](https://discourse.julialang.org/u/arunsm)\
**Post date:** [August 27, 2021, 1:25pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/8 "2021-08-27T13:25:31Z")

</div>

A simpler basis would be polynomial, which we could get using ApproxFun:

```julia
using ApproxFun

function vandermonde(S,n,x::AbstractVector)
    V=Array{Float64, 2}(undef, length(x), n)
    for k=1:n
        V[:,k]=Fun(S,[zeros(k-1);1]).(x)
    end
    V
end

x=collect(0:0.1:2)
dom = domain(Fun(identity,minimum(x)..maximum(x)))
v=exp.(-2x)

S=Taylor()
k = 9 # polynomial order
V=vandermonde(S,k,x)
f=Fun(S,V\v)

order = 3 # derivative order
g = Derivative(dom,order) * f

plot(x,-8 .* v, label="Exact")
scatter!(x, g.(x), label="ApproxFun", legend=:bottomright)

```

This is the least-squares fit of the data used in this thread using the Taylor polynomial. ApproxFun has a bunch of other bases. ApproxFun also includes piecewise domains for combining multiple bases, which you may be able to use if the data is piecewise continuous and you apply a change-detection algorithm.

---

<div class="post-metadata">

**Author:** ![arunsm](https://avatars.discourse-cdn.com/v4/letter/a/91b2a8/32.png) [@arunsm](https://discourse.julialang.org/u/arunsm)\
**Post date:** [August 28, 2021, 2:41pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/9 "2021-08-28T14:41:14Z")

</div>

As an aside, regarding piecewise continuous data, e.g., jump discontinuities, I just found a really nice julia package: [NoiseRobustDifferentiation](https://github.com/adrhill/NoiseRobustDifferentiation.jl). It’s based off of work from Rick Chartrand (a Los Alamos researcher):

- Rick Chartrand, “Numerical differentiation of noisy, nonsmooth data,” ISRN Applied Mathematics, Vol. 2011, Article ID 164564, 2011.

A while back I was interested in real-time numerical differentiation for nonlinear systems using embedded controllers. At the time, I developed a variant of the numerical differentiator described in:

- M. Mboup, C. Join, and M. Fliess, “Numerical differentiation with annihilators in noisy environment.” Numerical Algorithms, 50 (4), 439–467, 2009.

There’s a nice python/Matlab BSD-licensed implementation of these algebraic differentiators [here](https://github.com/aothmane-control/Algebraic-differentiators). I didn’t see an implementation in julia. I may end up translating it because it allows for explicitly setting the low-pass characteristics of the differentiator at higher frequencies. This is useful if you want as output a causal LTI filter with coefficients that could implement for real-time use.

Numerical differentiation is a very interesting rabbit hole with so many ways of attacking the problem, each with its own advantages / disadvantages.

---

<div class="post-metadata">

**Author:** ![DanielIngraham](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielingraham/32/13055_2.png) [@DanielIngraham](https://discourse.julialang.org/u/DanielIngraham)\
**Post date:** [September 1, 2021, 5:38pm UTC](https://discourse.julialang.org/t/high-order-derivatives-of-sampled-data/67045/10 "2021-09-01T17:38:53Z")

</div>

Excellent, thanks everyone! I ended up going with BSplineKit. This community is great.