[R] how lambda is computed in smoot.spline given _df_
Martin Maechler
maechler at stat.math.ethz.ch
Mon Nov 8 14:51:39 CET 2004
>>>>> "Eryk" == Eryk Wolski <wolski at molgen.mpg.de>
>>>>> on Mon, 08 Nov 2004 10:39:26 +0100 writes:
Eryk> Hi,
Eryk> I posted some days ago a question concerning the
Eryk> computation of lambda in the smooth.spline function
Eryk> (which I repreat at the bottom of the mail) given _df_.
I had thought one should be able to figure it out.
The help page cannot give all the theory behind smoothing splines!
The 'References' of help(smooth.spline) do explain how "df" is
defined, and even the help page mentions that
df = trace( <smoother matrix> ).
Eryk> Unfortunately the documentation is not clear to me.
Eryk> Maybee someone can help to answer in my view the
Eryk> basic question:
Eryk> If the penalized log likelihood is
Eryk> L = (y - f)' W (y - f) + lambda c' Sigma c
Eryk> how the _lambda_ in the above equation is computed if _df_ is
Eryk> given and _spar_ not?
Well, we have a computer for solving such problems:
spar (or lambda, equivalently) is determined such that the
resulting df matches the desired df by a traditional zero-finder
(Brent). Basically the same algorithm as used by uniroot().
Eryk> And, Is there a way to define lambda directly?
Almost {note: this is high school math}:
Since the help page gives the linear relation
between log(lambda) and spar, you call smooth.spline() once with
e.g. lambda = 0.5; and then know the formula how to compute spar
from lambda:
Using the 'cars' example from help(smooth.spline):
assume you want to set
my.lambda <- 1e-4
data(cars)# needed in R versions before 2.0.0
attach(cars)
## Call it with an arbitrary fixed 'spar':
(sspl <- smooth.spline(speed, dist, spar= .5))
## now solve the linear relationship for the only unknown s0 :
s0 <- with(sspl, spar - 0.0601*log(lambda))
ssp. <- smooth.spline(speed, dist, spar= s0 + 0.0601*log(my.lambda))
str(ssp.)## confirms that lambda *is* 1e-4 = my.lambda
---------
Hoping this helped:
Martin Maechler Seminar fuer Statistik, ETH Zurich
More information about the R-help
mailing list