Second derivatives are very noisy near endpoints

[hersle]

The second derivative is very noisy near the spline endpoints:

using DataInterpolations, Plots
t = 0.0:0.01:1.0
y = @. t^2 # -> y′ = 2*t, y″ = 2
yspl = CubicSpline(y, t)
plot(t, yspl(t); label = "y")
plot!(t, t -> DataInterpolations.derivative(yspl, t, 1); label = "y′")
plot!(t, t -> DataInterpolations.derivative(yspl, t, 2); label = "y″")

The first derivative looks much better, but also looks more kinky than it should near the ends (?) 🙂

[hersle]

For comparison, it does not happen with cubic splines from e.g. Dierckx:

using Dierckx, Plots
t = 0.0:0.01:1.0
y = @. t^2 # -> y′ = 2*t, y″ = 2
yspl = Spline1D(t, y; k=3)
plot(t, yspl(t); label = "y")
plot!(t, t -> Dierckx.derivative(yspl, t; nu=1); label = "y′")
plot!(t, t -> Dierckx.derivative(yspl, t, nu=2); label = "y″")

Now that I have thought more about this: is this just an effect of the choice of boundary conditions on the splines, as in this thread? ?DataInterpolations.CubicSpline does say it is using natural boundary conditions, which means the second derivative is forced to 0 at the ends. I chose a function whose derivative is exactly 2 everywhere, so that must create some "tension" at the ends?

[SouthEndMusic]

You only plot the interpolation in the data points, which is misleading. Have a look at this:

using DataInterpolations, Plots
t = 0.0:0.1:1.0
y = @. t^2 # -> y′ = 2*t, y″ = 2
yspl = CubicSpline(y, t)

t_eval = 0.0:0.001:1.0
plot(yspl; label = "yspl") # Uses plotting recipe
plot!(t_eval, DataInterpolations.derivative.(Ref(yspl), t_eval, 1), label = "yspl′")
plot!(t_eval, DataInterpolations.derivative.(Ref(yspl), t_eval, 2), label = "yspl″")
plot!(t_eval, t_eval.^2, label = "y", ls = :dash)
plot!(t_eval, 2*t_eval, label = "y′", ls = :dash)
plot!(t_eval, fill(2, length(t_eval)), label = "y″", ls = :dash)

Because of how CubicSpline is initialized, the interpolation is close to, but not an exact representation of the function $y = x^2$. As you noticed, it couldn't be because of the boundary conditions. If you want more control, you could use CubicHermiteSpline (or QuinticHermiteSpline) which does give you an exact match:

using DataInterpolations, Plots
t = 0.0:0.1:1.0
y = t.^2
y′ = 2*t

A = CubicHermiteSpline(y′, y, t)

t_eval = 0.0:0.001:1.0
plot(A; label = "CubicHermiteSpline") # Uses plotting recipe
plot!(t_eval, DataInterpolations.derivative.(Ref(A), t_eval, 1); label = "CubicHermiteSpline′")
plot!(t_eval, DataInterpolations.derivative.(Ref(A), t_eval, 2); label = "CubicHermiteSpline″")
plot!(t_eval, t_eval.^2; ls = :dash, label = "y")
plot!(t_eval, 2 * t_eval; ls = :dash, label = "y′")
plot!(t_eval, fill(2, length(t_eval)); ls = :dash, label = "y″")

[hersle]

Thank you very much. I'm integrating ODEs, so I know the derivative, and will try to adapt to hermite splines.