espdev / csaps

Cubic spline approximation (smoothing)
https://csaps.readthedocs.io
MIT License
165 stars 28 forks source link

Is there an easy way to access the number of parameters and degrees of freedom, so one could calculate the GCV? #64

Closed ZeroCool2u closed 8 months ago

ZeroCool2u commented 1 year ago

I'm interested in optimizing the smoothing parameter using the GCV (page 244 7.52). Oftentimes, this requires the number of effective number of parameters and degrees of freedom for the generated spline. I've been looking at the source code in _sspumv.py and just wanted to confirm there's no easy/built-in way to get this from a spline that csaps has created.

espdev commented 1 year ago

Hello @ZeroCool2u,

Under the hood csaps uses scipy.interpolate.PPoly class as the base class for the spline representation. You can get the spline coefficients and the degree from the spline instance.

s = CubicSmoothingSpline(x, y, smooth=0.5)

s.spline
SplinePPForm
  breaks: [1. 2. 3. 4. 5.]
  coeffs shape: (4, 4)
  data shape: (5,)
  axis: 0
  pieces: 4
  order: 4
  ndim: 1

s.spline.coeffs
array([[-0.02813299,  0.08695652, -0.14578005,  0.08695652],
       [ 0.        , -0.08439898,  0.17647059, -0.26086957],
       [ 0.16879795,  0.08439898,  0.17647059,  0.09207161],
       [ 2.16879795,  2.30946292,  2.39641944,  2.60358056]])

s.spline.order
4

s.spline.pieces
4

That's what you need?

ZeroCool2u commented 1 year ago

Thanks for the reply! I did see the information here that the SplinePPForm provides, but it doesn't seem sufficient to calculate the GCV. I've modified my code to calculate it manually based on _sspumv.py and that does seem to work. If I have time I'll try to submit a PR that stores the variables needed for someone to calculate the GCV.

Calculating the GCV directly would require the package to generate a spline at each step along the valid range of the smoothing factor (0, 1), so maybe I can just provide an example of how to do it that folks can copy and paste if they'd like.

espdev commented 1 year ago

Thanks @ZeroCool2u

PRs are always welcome. :)

ZeroCool2u commented 1 year ago

Just submitted my PR!