interpolation - Is there an R library that estimates a multivariate natural cubic spline (or similar) function?

Question

note: originally posted on Cross Validated (stats SE) on 07-26-2011, with no correct answers to date.

Background

I have a model, f, where Y=f(X)

X is an n x m matrix of samples from m parameters and Y is the n x 1 vector of model outputs.

f is computationally intensive, so I would like to approximate f using a multivariate cubic spline through (X,Y) points, so that I can evaluate Y at a larger number of points.

Question

Is there an R function that will calculate an arbitrary relationship between X and Y?

Specifically, I am looking for a multivariate version of the splinefun function, which generates a spline function for the univariate case.

e.g. this is how splinefun works for the univariate case

x <- 1:100
y <- runif(100)
foo <- splinefun(x,y, method = "monoH.FC")
foo(x) #returns y, as example

The test that the function interpolates exactly through the points is successful:

all(y == foo(1:100))
## TRUE

What I have tried

I have reviewed the mda package, and it seems that the following should work:

library(mda)
x   <- data.frame(a = 1:100, b = 1:100/2, c = 1:100*2)
y   <- runif(100)
foo <- mars(x,y)
predict(foo, x) #all the same value

however the function does not interpolate exactly through the design points:

all(y == predict(foo,x))
## FALSE

I also could not find a way to implement a cubic-spline in either the gam, marss, or earth packages.

Answer 1

Actually several packages can do it. The one I use is the "rms" package which has rcs, but the survival package also has pspline and the splines package has the ns function {}. "Natural splines" (constructed with ns) are also cubic splines. You will need to form multivariate fitting function with the '*' operator in the multivariate formula creating "crossed" spline terms. that the example you offered was not sufficiently rich.

I guess I am confused that you want exact fits. R is a statistical package. Approximate estimation is the goal. Generally exact fits are more of a problem because they lead to multicollinearity.