This module provides cubic smoothing spline for univariate/multivariate/gridded data approximation. The smoothing parameter can be calculated automatically or it can be set manually.
The smoothing parameter should be in range [0, 1]
where bounds are:
- 0: The smoothing spline is the least-squares straight line fit to the data
- 1: The natural cubic spline interpolant
The calculation of the smoothing spline requires the solution of a linear system
whose coefficient matrix has the form p*A + (1 - p)*B
, with the matrices A
and B
depending on the data sites X
. The automatically computed smoothing parameter makes
p*trace(A) equal (1 - p)*trace(B)
.
Python 3.5 or above is supported.
pip install csaps
The module depends only on NumPy and SciPy.
On Windows we highly recommend to use unofficial builds NumPy+MKL and SciPy from Christoph Gohlke.
> mkdir depends
download numpy (https://www.lfd.uci.edu/~gohlke/pythonlibs/#numpy) to "depends" directory
download scipy (https://www.lfd.uci.edu/~gohlke/pythonlibs/#scipy) to "depends" directory
> pip install --find-links depends/ csaps
Univariate data are two vectors: X and Y with the same size. X is data sites, Y is data values.
For univariate case X-values must satisfy the condition: x1 < x2 < ... < xN
.
You can use UnivariateCubicSmoothingSpline
class for uivariate data smoothing.
import numpy as np
import matplotlib.pyplot as plt
import csaps
np.random.seed(1234)
x = np.linspace(-5., 5., 25)
y = np.exp(-(x/2.5)**2) + (np.random.rand(25) - 0.2) * 0.3
sp = csaps.UnivariateCubicSmoothingSpline(x, y, smooth=0.85)
xs = np.linspace(x[0], x[-1], 150)
ys = sp(xs)
plt.plot(x, y, 'o', xs, ys, '-')
plt.show()
The algorithm supports weighting. You can set weights vector that will determine weights for all data values:
import csaps
x = [1., 2., 4., 6.]
y = [2., 4., 5., 7.]
w = [0.5, 1, 0.7, 1.2]
sp = csaps.UnivariateCubicSmoothingSpline(x, y, w)
...
The algorithm supports vectorization. You can compute smoothing splines for
X
, Y
data where X
is data site vector and Y
is ND-array of data value vectors.
The shape of Y
array must be: (d0, d1, ..., dN)
where dN
must equal of X
vector size.
In this case the smoothing splines will be computed for all Y
data vectors at a time.
For example:
import numpy as np
import csaps
# data sites
x = [1, 2, 3, 4]
# two data vectors
y = np.array([(2, 4, 6, 8),
(1, 3, 5, 7)])
sp = csaps.UnivariateCubicSmoothingSpline(x, y)
xi = np.linspace(1, 4, 10)
yi = sp(xi)
print(yi.shape) # (2, 10)
assert yi.shape[:-1] == y.shape[:-1]
assert yi.shape[-1] == xi.size
Important: The same weights vector and the same smoothing parameter will be used for all Y data.
We can easily smooth multivariate data using univariate smoothing spline, vectorization and parametrization.
We should define parametric data sites vector t
with condition t1 < t2 < ... < tN
and make spline for each dimension X(t), Y(t), ..., M(t)
.
The module provides MultivariateCubicSmoothingSpline
class for this case.
This class is a simple wrapper and it just automatically computes t
vector by default and calls UnivariateCubicSmoothingSpline
:
# Construct multivariate spline from t and X, Y, Z, ..., M
sx = UnivariateCubicSmoothingSpline(t, data[x])
sy = UnivariateCubicSmoothingSpline(t, data[y])
sz = UnivariateCubicSmoothingSpline(t, data[z])
...
sm = UnivariateCubicSmoothingSpline(t, data[m])
# Or the same with using vectorization
sxyz_m = UnivariateCubicSmoothingSpline(t, data)
In this example we make parametric 3D curve from theta
parameter and compute spline from the same theta
parameter.
In fact we could use UnivariateCubicSmoothingSpline
with vectorization instead of MultivariateCubicSmoothingSpline
in this case.
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import csaps
n = 100
theta = np.linspace(-4 * np.pi, 4 * np.pi, n)
z = np.linspace(-2, 2, n)
r = z ** 2 + 1
np.random.seed(1234)
x = r * np.sin(theta) + np.random.randn(n) * 0.3
np.random.seed(5678)
y = r * np.cos(theta) + np.random.randn(n) * 0.3
data = np.vstack((x, y, z))
sp_theta = csaps.MultivariateCubicSmoothingSpline(data, theta, smooth=0.95)
# or the same
# sp_theta = csaps.UnivariateCubicSmoothingSpline(theta, data, smooth=0.95)
theta_i = np.linspace(-4 * np.pi, 4 * np.pi, 250)
data_i = sp_theta(theta_i)
xi = data_i[0, :]
yi = data_i[1, :]
zi = data_i[2, :]
fig = plt.figure()
ax = fig.gca(projection='3d')
plt.show()
The algorithm can make smoothing splines for ND-gridded data approximation. In this case we use coordinatewise smoothing (tensor-product of univariate splines coefficients).
You can use NdGridCubicSmoothingSpline
class for ND-gridded data smoothing.
You also can set weights and smoothing parameters for each dimension.
Currently the implementation does not support vectorization for ND-gridded data.
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import csaps
xdata = [np.linspace(-3, 3, 61), np.linspace(-3.5, 3.5, 51)]
i, j = np.meshgrid(*xdata, indexing='ij')
ydata = (3 * (1 - j)**2. * np.exp(-(j**2) - (i + 1)**2)
- 10 * (j / 5 - j**3 - i**5) * np.exp(-j**2 - i**2)
- 1 / 3 * np.exp(-(j + 1)**2 - i**2))
np.random.seed(12345)
noisy = ydata + (np.random.randn(*ydata.shape) * 0.75)
sp = csaps.NdGridCubicSmoothingSpline(xdata, noisy, smooth=0.988)
ysmth = sp(xdata)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_wireframe(j, i, noisy, linewidths=0.5, color='r')
ax.scatter(j, i, noisy, s=5, c='r')
ax.plot_surface(j, i, ysmth, linewidth=0, alpha=1.0)
plt.show()
Please look through csaps.ipynb file for more examples.
We use pytest and tox (on Travis CI) for testing. Please see test_csaps.py file.
csaps
is a Python modified port of MATLAB CSAPS function that is an implementation of
Fortran routine SMOOTH from PGS (originally written by Carl de Boor).
csaps-cpp C++11 Eigen based implementation of the algorithm.
C. de Boor, A Practical Guide to Splines, Springer-Verlag, 1978.