[R] smooth.spline gives different results from sreg ?

Duncan Murdoch murdoch.duncan at gmail.com
Sun Jul 28 01:10:51 CEST 2013


On 13-07-27 2:50 PM, Jean-Luc Dupouey wrote:
> Dear R-helpers,
>
> I compared various programs for cubic spline smoothing, and it appeared
> that smooth.spline ( stats version 3.0.1) seems to behave surprisingly.
> For enough long series and low values of lambda (or spar), the results
> of smooth.spline seem to be different from those of sreg ( package
> fields version 6.8), Octave (=MATLAB) or SAS. These three last softwares
> always gave the same results.

Did you read the ?smooth.spline help page, in particular the Details and 
Note sections?  They indicate that the default computation makes some 
efficiency simplifications.

Duncan Murdoch

>
> Here is a script which shows the problem:
>
> #generate a random series of 2000 values
>
> set.seed(1)
> MyData=data.frame(Time=1:2000,Val=runif(1000))
>
> #calculate the sreg cubic smoothing spline with a given lambda parameter
> (0.006 here)
>
> library(fields)
>
> SplineFields=sreg(MyData$Time,MyData$Val,lambda=0.006)
>
> #keep the minimim fitted value (or any other from a long list of
> possible values)
>
> ValMin=min(SplineFields$fitted.values)
> TimeValMin=which.min(SplineFields$fitted.values)
>
> #calculations of all possible fitted values at the TimeValMin point with
> smooth.spline,
> #varying the spar parameter in the range of all its possible values
>
> SplineRValMin=sapply(seq(-0.5,2.5,0.1),
>     function(Ispar) {
>       SplineR=smooth.spline(MyData$Time,MyData$Val,spar=Ispar)
>       SplineR$y[TimeValMin]})
>
> #None of the smooth.spline fitted values reach the one calculated with
> sreg !
>
> Lim=range(ValMin,SplineRValMin)
>
> #smooth.spline values
> plot(seq(-0.5,2.5,0.1),SplineRValMin,type="l",ylim=Lim)
>
> #sreg value
> abline(h=ValMin)
>
> I hope there is no real problem here, but only some misunderstanding
> from my side, because cubic splines are very often used. Best regards,
>
> Jean-Luc Dupouey
>



More information about the R-help mailing list