[R] smooth.spline

rkevinburton at charter.net rkevinburton at charter.net
Mon Jul 21 02:46:50 CEST 2008


Actually that was my next question. From the books that I have I see a "natural spline" and a clamped spline. I am assuming that "natural" (Umerical Analysis, Burden, et. all.) cooresponds to 'R''s "natural" method. I am not clear on what a clamped spline cooresponds to (fmm or perodic). Or what the difference between fmm and periodic.

Thank you.

Kevin
---- Spencer Graves <spencer.graves at pdf.com> wrote: 
>       Are you aware that there are many different kinds of splines?  
> With "spline" and "splinefun", you can use method = "fmm" (Forsyth, 
> Malcolm and Moler), "natural", or "periodic".  I'm not familiar with 
> "fmm", but it seems to be adequately explained by the "Manual spline 
> evaluation" you quoted from the documentation. 
> 
>       Natural splines are perhaps the simplest:  I(x-x0)*(x-x0)^j, where 
> x0 is a knot, and I(z) = 1 if z>0 and 0 otherwise. 
> 
>       However, computations using natural splines are numerically 
> unstable.  The standard solution to this problem is to use B-splines, 
> which are 0 outside a finite interval. 
> 
>       Let's look at your example: 
> 
> n <- 9
> x <- 1:n
> y <- rnorm(n)
> plot(x, y, main = paste("spline[fun](.) through", n, "points"))
> spl <- smooth.spline(x,y)
> lines(spl)
> 
>       The 'smooth.spline' function uses B-splines.  To see what they 
> look like, let's do the following: 
> 
> library(fda)
> Bspl.basis <- create.bspline.basis(unique(spl$fit$knot))
> 
> # Check to make sure: 
> all.equal(knots(Bspl.basis, interior=FALSE), spl$fit$knot)
> # TRUE
> 
> # What do B-splines look like? 
> plot(Bspl.basis)
> abline(v=knots(Bspl.basis), lty='dotted', col='red')
> #  7 interior knots, 2 end knots replicated 4 times each, for a spline 
> of order 4, degree 3 (cubic splines) 
> # total of 15 knots
> # Each spline uses 5 consecutive knots, which means there will be 11 
> basis functions. 
>      
> # NOTE:  'smooth.spline' rescaled the interval [1, 9] to [0, 1]. 
> # Evaluate the 11 B-splines at 'x'
> Bspl.basis.x <- eval.basis((x-1)/8, Bspl.basis)
> 
> round(Bspl.basis.x, 4)
> 
> # Now the manual computation: 
> y.spl <- Bspl.basis.x %*% spl$fit$coef
> 
> # Plot to confirm: 
> plot(x, y, main = paste("spline[fun](.) through", n, "points"))
> spl.xy <- spline(x, y)
> lines(spl.xy)
> points(x, y.spl, pch=2, col='red')
> 
>       Hope this helps. 
>       Spencer
> 
> rkevinburton at charter.net wrote:
> > Fair enough. FOr a spline interpolation I can do the following:
> >
> >   
> >> n <- 9
> >> x <- 1:n
> >> y <- rnorm(n)
> >> plot(x, y, main = paste("spline[fun](.) through", n, "points"))
> >> lines(spline(x, y))
> >>     
> >
> > Then look at the coefficients generated as:
> >
> >   
> >> f <- splinefun(x, y)
> >> ls(envir = environment(f))
> >>     
> > [1] "ties" "ux"   "z"   
> >   
> >> splinecoef <- get("z", envir = environment(f))
> >> slinecoef
> >>     
> > $method
> > [1] 3
> >
> > $n
> > [1] 9
> >
> > $x
> > [1] 1 2 3 4 5 6 7 8 9
> >
> > $y
> > [1]  0.93571604  0.44240485  0.45451903 -0.96207396 -1.13246522 -0.60032698
> > [7] -1.77506105 -0.09171419 -0.23262573
> >
> > $b
> > [1] -1.53673409  0.22775629 -0.81788209 -1.16966436  0.73558677 -0.68744178
> > [7]  0.08639287  1.86770869 -2.92992167
> >
> > $c
> > [1]  1.3657783  0.3987121 -1.4443504  1.0925682  0.8126830 -2.2357115  3.0095462
> > [8] -1.2282303 -3.5694000
> >
> > $d
> > [1] -0.32235542 -0.61435416  0.84563953 -0.09329507 -1.01613149  1.74841922
> > [7] -1.41259217 -0.78038989 -0.78038989
> >
> > WHen I look at ?spline there is even an example of "manually" using these coefficeients:
> >
> > ## Manual spline evaluation --- demo the coefficients :
> > .x <- get("ux", envir = environment(f))
> > u <- seq(3,6, by = 0.25)
> > (ii <- findInterval(u, .x))
> > dx <- u - .x[ii]
> > f.u <- with(splinecoef,
> >             y[ii] + dx*(b[ii] + dx*(c[ii] + dx* d[ii])))
> > stopifnot(all.equal(f(u), f.u))
> >
> >
> > For the smooth.spline as
> >
> > spl <- smooth.spline(x,y)
> >
> > I can also look at the coefficients:
> >
> > spl$fit
> > $knot
> >  [1] 0.000 0.000 0.000 0.000 0.125 0.250 0.375 0.500 0.625 0.750 0.875 1.000
> > [13] 1.000 1.000 1.000
> >
> > $nk
> > [1] 11
> >
> > $min
> > [1] 1
> >
> > $range
> > [1] 8
> >
> > $coef
> >  [1]  0.90345898  0.73823276  0.40777431 -0.08046715 -0.54625461 -0.85205147
> >  [7] -0.96233408 -0.91373830 -0.66529714 -0.47674774 -0.38246971
> >
> > attr(,"class")
> > [1] "smooth.spline.fit"
> >
> > But there isn't an example on how to "manual" use these coefficients. This is what I was asking about. Once I hae the coefficients how do I "manually" interpolate using the coefficients given and x.
> >
> > Thank you.
> >
> > Kevin
> >
> >
> > ---- Spencer Graves <spencer.graves at pdf.com> wrote: 
> >   
> >>       PLEASE do read the posting guide 
> >> http://www.R-project.org/posting-guide.html and provide commented, 
> >> minimal, self-contained, reproducible code.
> >>
> >>       I do NOT know how to do what you want, but with a self-contained 
> >> example, I suspect many people on this list -- probably including me -- 
> >> could easily solve the problem.  Without such an example, there is a 
> >> high probability that any answer might (a) not respond to your need, and 
> >> (b) take more time to develop, just because we don't know enough of what 
> >> you are asking. 
> >>
> >>       Spencer
> >>
> >> rkevinburton at charter.net wrote:
> >>     
> >>> Like I indicated. I understand the coefficients in a B-spline context. If I use the the 'spline' or 'splinefun' I can get the coefficients and they are grouped as 'a', 'b', 'c', and 'd' coefficients. But the coefficients for smooth.spline is just an array. I basically want to take these coefficients and outside of 'R' use them to form an interpolation. In other words I want 'R' to do the hard work and then export the results so they can be used else where.
> >>>
> >>> Thank you.
> >>>
> >>> Kevin
> >>>   
> >>>       
> >> Spencer Graves wrote:
> >>     
> >>>      I believe that a short answer to your question is that the 
> >>> "smooth" is a linear combination of B-spline basis functions, and the 
> >>> coefficients are the weights assigned to the different B-splines in 
> >>> that basis.
> >>>      Before offering a much longer answer, I would want to know what 
> >>> problem you are trying to solve and why you want to know.  For a brief 
> >>> description of B-splines, see 
> >>> "http://en.wikipedia.org/wiki/B-spline".  For a slightly longer 
> >>> commentary on them I suggest the "scripts\ch01.R" in the DierckxSpline 
> >>> package:  That script computes and displays some B-splines using 
> >>> "splineDesign", "spline.des" in the 'splines' package plus comparable 
> >>> functions in the 'fda' package.  For more info on this, I found the 
> >>> first chapter of Paul Dierckx (1993) Curve and Surface Fitting with 
> >>> Splines (Oxford U. Pr.).  Beyond that, I've learned a lot from the 
> >>> 'fda' package and the two companion volumes by Ramsay and Silverman 
> >>> (2006) Functional Data Analysis, 2nd ed. and (2002) Applied Functional 
> >>> Data Analysis (both Springer).
> >>>      If you'd like more help from this listserve, PLEASE do read the 
> >>> posting guide http://www.R-project.org/posting-guide.html and provide 
> >>> commented, minimal, self-contained, reproducible code.
> >>>         Hope this helps.      Spencer Graves
> >>>
> >>> rkevinburton at charter.net wrote:
> >>>       
> >>>> I like what smooth.spline does but I am unclear on the output. I can 
> >>>> see from the documentation that there are fit.coef but I am unclear 
> >>>> what those coeficients are applied to.With spline I understand the 
> >>>> "noraml" coefficients applied to a cubic polynomial. But these 
> >>>> coefficients I am not sure how to interpret. If I had a description 
> >>>> of the algorithm maybe I could figure it out but as it is I have this 
> >>>> question. Any help?
> >>>>
> >>>> Kevin
> >>>>
> >>>> ______________________________________________
> >>>> R-help at r-project.org mailing list
> >>>> https://stat.ethz.ch/mailman/listinfo/r-help
> >>>> PLEASE do read the posting guide 
> >>>> http://www.R-project.org/posting-guide.html
> >>>> and provide commented, minimal, self-contained, reproducible code.
> >>>>   
> >>>>         
> >
> > ______________________________________________
> > R-help at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
> >



More information about the R-help mailing list