Simple, concrete example in docs

[OliverEvans96]

Hello,

There's a whole page of benchmarks front-and-center in the README, but I've been searching for 10 minutes now to figure out how to do a 1D linear interpolation of a set of points similar to scipy's interp1d.

This page discusses (but doesn't give an example of) the first argument to interpolate but not the required options.

This page discusses the options but not the data (until the section on Parametric Splines, which is clear but not what I'm trying to do).

Also, The tests are abstracted such that they do not contain a simple, concrete example either.

The equivalent scipy documentation is abundantly clear and easy to find.

from scipy.interpolate import interp1d

x = np.linspace(0, 10, num=11, endpoint=True)
y = np.cos(-x**2/9.0)
f = interp1d(x, y)
f2 = interp1d(x, y, kind='cubic')
xnew = np.linspace(0, 10, num=41, endpoint=True)

import matplotlib.pyplot as plt
plt.plot(x, y, 'o', xnew, f(xnew), '-', xnew, f2(xnew), '--')
plt.legend(['data', 'linear', 'cubic'], loc='best')
plt.show()

Sorry to complain, but hopefully it's helpful complaining :) Oliver

[OliverEvans96]

The first Dierckx README is also great (but I'd rather avoid a fortran wrapper).

Example Usage

using Dierckx

Fit a 1-d spline to some input data (points can be unevenly spaced):

x = [0., 1., 2., 3., 4.]
y = [-1., 0., 7., 26., 63.]  # x.^3 - 1.
spl = Spline1D(x, y)

Evaluate the spline at some new points:

spl([1.5, 2.5])  # result = [2.375, 14.625]
spl(1.5)  # result = 2.375

[OliverEvans96]

This post provides a good but outdated example:

using Interpolations

function interp1(xpt, ypt, x; method="linear", extrapvalue=nothing)

    if extrapvalue == nothing
        y = zeros(x)
        idx = trues(x)
    else
        y = extrapvalue*ones(x)
        idx = (x .>= xpt[1]) .& (x .<= xpt[end])
    end
    
    if method == "linear"
        intf = interpolate((xpt,), ypt, Gridded(Linear()))
        y[idx] = intf[x[idx]]

    elseif method == "cubic"
        itp = interpolate(ypt, BSpline(Cubic(Natural())), OnGrid())
        intf = scale(itp, xpt)
        y[idx] = [intf[xi] for xi in x[idx]]
    end
    
    return y
end

x = 0:pi/4:2*pi
v = sin.(x)
xq = 0:pi/16:1.5*2*pi
vq = interp1(x, v, xq, method="cubic")

[mkitti]

Something like this maybe? http://juliamath.github.io/Interpolations.jl/latest/convenience-construction/#Convenience-notation-1

f(x) = log(x)
xs = 1:0.2:5
A = [f(x) for x in xs]

# linear interpolation
interp_linear = LinearInterpolation(xs, A)
interp_linear(3) # exactly log(3)
interp_linear(3.1) # approximately log(3.1)

# cubic spline interpolation
interp_cubic = CubicSplineInterpolation(xs, A)
interp_cubic(3) # exactly log(3)
interp_cubic(3.1) # approximately log(3.1)

[OliverEvans96]

Yes, exactly like that :facepalm: Thank you very much. Apparently "Convenience Constructors" didn't pop out to me, although I now see that it's also linked at the bottom of the General Usage page.

[vini-fda]

I'd also love if there was a plotting example. It took hours for me to discover that if I create an interpolation object itp, I can't just treat it as a function input to plots, e.g. plot(itp, x) like in the basic examples of Plots.jl. So, I created an example(directly mirroring the scipy docs mentioned by @OliverEvans96), which I think helps with visual learners like me:

using Interpolations, Plots

# Lower and higher bound of interval
a = 1.0
b = 10.0
# Interval definition
x = a:1.0:b
# Function application by broadcasting
y = @. cos(x^2 / 9.0)
# Interpolations
itp_linear = LinearInterpolation(x, y)
itp_cubic = CubicSplineInterpolation(x, y)
# Interpolation functions
f_linear(x) = itp_linear(x)
f_cubic(x) = itp_cubic(x)
# Plots
width, height = 1500, 800 # not strictly necessary
x_new = a:0.1:b # smoother interval, necessary for cubic spline

scatter(x, y, markersize=10,label="Data points")
plot!(f_linear, x_new, w=3,label="Linear interpolation")
plot!(f_cubic, x_new, linestyle=:dash, w=3, label="Cubic Spline interpolation")
plot!(size = (width, height))
plot!(legend = :bottomleft)

And the generated plot is:

[mkitti]

Could you submit it as a pull request?

[timholy]

Arguably Plots should support plot(itp, x). What's broken about it? Maybe they just put ::Function in places where it might be better not to constrain the type.

[vini-fda]

Done! @mkitti this is my first PR here, hope I haven't done anything wonky, :sweat_smile:. Feel free to modify as you please!