Describe the workflow you want to enable

As B-splines naturally have a sparse structure, I'd like to have the option that SplineTransformer returns a sparse matrix instead of always an ndarray.

import numpy as np
from sklearn.preprocessing import SplineTransformer

X = np.arange(6).reshape(6, 1)
spline = SplineTransformer(degree=2, n_knots=3, sparse=True)
spline.fit_transform(X)

Describe your proposed solution

With scipy >= 1.8 (yet to be released), we can use design_matrix from scipy/scipy#14344.

Describe alternatives you've considered, if relevant

No response

Additional context

No response