On Wed, 2009-09-16 at 14:10 -0700, Chris Martin wrote:
> Dear list members,
> 
> I am a looking for a function that can calculate a "surface" from at least
> four predictor variables and one response variable. I would then like to
> calculate the first derivative of a specific point on this surface.
> 
> I have looked at many packages for nonparametric smoothing and kernel
> density estimation but have been unable to find any that fulfill both these
> criteria. For example, loess can handle multivariate data, but I do not how
> to extract the derivative from the resulting fit? Many smoothing splines
> offer predict functions to extract the derivative, but these functions can
> only handle univariate data (e.g. smooth.spline). Ideally, I would like to
> use local estimates of the surface (i.e. loess).

Depending upon how many observations you have (i.e. not very large
numbers) you could estimate this model using mgcv::gam which allows
multivariate smooths. You can then compute the derivatives for this
smooth using finite differences using the predict.gam method as a basis
from which to work (using type = "lpmatrix").

The trick is to get the bits from the model that you need. This code
snippet uses data from the semiPar package. I produced the script as I
was reading the text that semiPar is support for. Simon Wood, author of
the mgcv package, kindly helped me with the code to compute the
derivatives:

## Example from Ruppert, Wand and Carroll 2003 pp 156-8
## Semiparametric regression
## Cambridge University Press
##
## Data for the strontium example are in package SemiPar, support
## software for the book above.
require(SemiPar)
require(mgcv)
# '`
require(nlme)

## load the Strontium data
data(fossil)

mod1 <- gamm(strontium.ratio ~ s(age), data = fossil)

## Compute derivatives ###############################################
## where to evaluate derivatives
x.age <- with(fossil, seq(min(age), max(age), length = 200))
newd <- data.frame(age = x.age)
X0 <- predict(mod1$gam, newd, type = "lpmatrix")
eps <- 1e-7 ## finite difference interval
x.age <- x.age + eps ## shift the evaluation points
newd <- data.frame(age = x.age)
X1 <- predict(mod1$gam, newd, type = "lpmatrix")
Xp <- (X1 - X0) / eps ## maps coefficients to (fd approx.) derivatives
colnames(Xp)

Xi <- Xp*0
Xi[, 2:10] <- Xp[,2:10] ## Xi\%*\%coef(b) = smooth deriv i
## you'll need to change 2:10 in the above to match the cols for your 
## smoother
## ith smooth derivative
df <- Xi %*% coef(mod1$gam)
## cheap diag(Xi\%*\%b$Vp\%*\%t(Xi))^.5
df.sd <- rowSums(Xi %*% mod1$gam$Vp * Xi)^.5 ## SE of deriv

I haven't used this approach to compute the derivatives of a
multivariate smooth, but it should work as long as you can identify the
columns of X0 and X1 that arise from predict.gam.

HTH

G

> 
> I would appreciate any suitable functions or advice on where to look for
> functions that fulfill both these criteria.
> 
> Thank you very much for your time.
> 
> best wishes,
> Chris Martin
> Population Biology Graduate Group '12
> University of California, Davis
> 
>       [[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-ecology mailing list
> R-sig-ecology@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson             [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,          [f] +44 (0)20 7679 0565
 Pearson Building,             [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London          [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT.                 [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%

_______________________________________________
R-sig-ecology mailing list
R-sig-ecology@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-ecology

Reply via email to