Scipy b spline basis functions

Question

Is there a more efficient way in Scipy to generate the b-spline basis functions similar to the recursive code given in the example of the b-spline function:

scipy.interpolate.BSpline

def B(x, k, i, t):
   if k == 0:
      return 1.0 if t[i] <= x < t[i+1] else 0.0
   if t[i+k] == t[i]:
      c1 = 0.0
   else:
      c1 = (x - t[i])/(t[i+k] - t[i]) * B(x, k-1, i, t)
   if t[i+k+1] == t[i+1]:
      c2 = 0.0
   else:
      c2 = (t[i+k+1] - x)/(t[i+k+1] - t[i+1]) * B(x, k-1, i+1, t)
   return c1 + c2

I tried to use scipy.interpolate.BSpline.basis_element but was unable to generate the same results as the function "B".

user_5603657 · Answer

The way I have previously done this is by the following:

import scipy.interpolate as si
import numpy as np
import matplotlib.pyplot as plt

nt = 3 # the number of knots
x = np.linspace(0,24,10)
y = np.ones((x.shape)) # the function we wish to interpolate (very trivial in this case)
t = np.linspace(x[0],x[-1],nt)
t = t[1:-1]

tt = np.linspace(0.0, 24, 100)
y_rep = si.splrep(x, y, t=t, k=2)

y_i = si.splev(tt, y_rep)

spl = np.zeros((100,))

plt.figure()
for i in range(nt+1):
    vec = np.zeros(nt+1)
    vec[i] = 1.0
    y_list = list(y_rep)
    y_list[1] = vec.tolist()
    y_i = si.splev(tt, y_list) # your basis spline function
    spl = spl + y_i*y_rep[1][i] # the interpolated function
    plt.plot(tt, y_i)

You can have a look at the two functions scipy.interpolate.splrep and scipy.interpolate.splev. They respectively make a B-spline representation and evaluate the B-splines. The above code generates the (in this case) four basis spline functions:

I hope this is what you wanted.

Scipy b spline basis functions

1 Answers