G'day Rolf,
>>>>> "RT" == Rolf Turner <rolf at math.unb.ca>
writes:
RT> Is there a way of calculating the derivative of a function
RT> returned by splinefun()?
Yes, of course. :)
As somebody else once said: "this is R, anything can be done, the
question is just how easy it is to do so".
It may depend on what kind of spline you are fitting with
splinefun()...
RT> Such a function is a cubic spline, whence it has a calculable
RT> derivative, but is there a (simple) way of getting at it?
Depends on your definition of simple. ;-))
RT> One workaround that I have thought of is to take a fine grid
RT> of points, evaluate the function returned by splinefun() at
RT> these points, put an interpolating spline through these points
RT> using smooth.spline() with spar=0, and then calculate the
RT> derivative with predict(...,deriv=1).
Why easy if one can also make it complicated? There is the function
interpSpline() in the package splines for which, according to its help
page, the predict method also works. So using smooth.spline with
spar=0 does indeed look a bit kludgy. :)
RT> Is there a better way?
splinefun() essentially returns a function that contains one object in
its environment, namely z, which contains all the necessary
information for calculating the spline and whose body is just a call
to an internal C function.
The components y, b, c and d of the object z contain the piecewise
polynomial representation between knots of the spline if you choose
the methods "fmm" or "natural". If you use the method
"periodic" then
the component e of the object z is used too, but I don't know whether
it is just used to hold intermediate results while calculating the
spline or whether it is actually necessary for calculating the
spline. (An analysis of the C code that is called to calculate the
spline and to evaluate the spline should answer that question but
until now I never bothered of doing this since I didn't have use for
periodic splines yet.)
Thus, since the spline is stored in a piecewise polynomial fashion, it
shouldn't be hard to define a function that calculates the derivative
of the function returned by splinefun(). If my calculus is correct,
the following code should construct such a function:
> par(mfrow=c(1,2))
> n <- 9
> x <- 1:n
> y <- rnorm(n)
> f <- splinefun(x, y)
> ls(envir = environment(f))
[1] "z"> splinecoef <- eval(expression(z), envir = environment(f))
> curve(f(x), 1, 10, col = "green", lwd = 1.5)
> points(splinecoef, col = "purple", cex = 2)
##
## construct a function that calculates the derivative of f()
##> g <- f
> environment(g) <- new.env(parent=baseenv())
> z <- get("z", envir=environment(f))
> z$y <- z$b
> z$b <- 2*z$c
> z$c <- 3*z$d
> z$d <- rep(0, z$n)
> assign("z", z, envir=environment(g))
> curve(g(x), 1, 10, col = "green", lwd = 1.5)
As indicated above, I would only vouch for this method if splinefun()
is called with method="fmm" or method="natural". And the
second
caveat is, of course, that you are assuming that I get my calculus
correct on a Saturday morning which might be a bit much to ask. ;-))
Hope that helps.
Cheers,
Berwin