Colocación de spline cúbico natural multivariante

17

nota: sin respuestas correctas después de un mes, he vuelto a publicar en SO

Antecedentes

Tengo un modelo, , donde Y = f ( X )FY=F(X)

es unamatriz de muestras n × m de m parámetros e Y es la n × 1Xn×mmYn×1 vector de las salidas del modelo.

es computacionalmente intensivo, por lo que me gustaría aproximar f usando una spline cúbica multivariada a través depuntos ( X , Y ) , para poder evaluar Y en un mayor número de puntos.ff(X,Y)Y

Pregunta

¿Existe una función R que calcule una relación arbitraria entre X e Y?

Específicamente, estoy buscando una versión multivariada de la splinefunfunción, que genera una función de spline para el caso univariante.

Por ejemplo, así es como splinefunfunciona el caso univariante

x <- 1:10
y <- runif(10)
foo <- splinefun(x,y)
foo(1:10) #returns y, as example
all(y == foo(1:10))
## TRUE

Lo que he intentado

He revisado el paquete mda , y parece que lo siguiente debería funcionar:

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

pero no pude encontrar ninguna manera de implementar una spline cúbica en mars

actualización desde la salida a la recompensa, I cambió el título - Si no hay una función R, aceptaría, en orden de preferencia: una función de R que da salida a una función gaussiana proceso, u otra función de interpolación multivariante que pasa a través de los puntos de diseño, preferiblemente en R, si no Matlab.

David LeBauer
fuente
tratar gam (función), permite a cualquier dimensión de splines cúbicos
user5563

Respuestas:

11

Este artículo presentado en UseR! 2009 parece abordar un problema similar

http://www.r-project.org/conferences/useR-2009/slides/Roustant+Ginsbourger+Deville.pdf

Sugiere el paquete DiceKriging http://cran.r-project.org/web/packages/DiceKriging/index.html

En particular, verifique las funciones km y prediga.

Aquí hay un ejemplo de interpolación tridimensional. Parece ser sencillo generalizar.

x <- c(0, 0.4, 0.6, 0.8, 1)
y <- c(0, 0.2, 0.3, 0.4, 0.5)
z <- c(0, 0.3, 0.4, 0.6, 0.8)

model <- function(param){
2*param[1] + 3*param[2] +4*param[3]
}


model.in <- expand.grid(x,y,z)
names(model.in) <- c('x','y','z')

model.out <- apply(model.in, 1, model)

# fit a kriging model 
m.1 <- km(design=model.in, response=model.out, covtype="matern5_2")

# estimate a response 
interp <- predict(m.1, newdata=data.frame(x=0.5, y=0.5, z=0.5), type="UK",    se.compute=FALSE)
# check against model output
interp$mean
# [1]  4.498902
model(c(0.5,0.5,0.5))
# [1] 4.5

# check we get back what we put in
interp <- predict(m.1, newdata=model.in, type="UK", se.compute=FALSE)
all.equal(model.out, interp$mean)
# TRUE
Tony Ladson
fuente
6

Necesita más datos para un ajuste de spline. mgcv de hecho es una buena opción. Para su solicitud específica, debe establecer la spline cúbica como la función básica bs = 'cr' y tampoco tenerla penalizada con fx = TRUE. Ambas opciones están configuradas para un término uniforme que se establece con s (). Predecir trabajos como se esperaba.

library(mgcv)
x <- data.frame(a = 1:100, b = 1:100/2, c = 1:100*2)
y <- runif(100)
foo <- gam(y~a+b+s(c,bs="cr",fx=TRUE),data=x)
plot(foo)
predict(foo,x)
Alex
fuente
Gracias por su ayuda, pero si esto fuera una ranura cúbica, ¿no debería esperar predict(foo,x)volver y?
David LeBauer
Lo sentimos, no me di cuenta de que quieres una aproximación perfecta. Entonces, aparentemente mgcv no es de mucha ayuda: stop ("Basis solo maneja suavizaciones 1D") (de svn.r-project.org/R-packages/trunk/mgcv/R/smooth.r )
Alex
0

No da detalles sobre la forma de la función F(X); Es posible que una función constante por partes sea una aproximación suficientemente buena, en cuyo caso es posible que desee ajustar un árbol de regresión (con el paquete, rpartpor ejemplo). De lo contrario, es posible que desee ver el paquete earth, además de lo que ya se ha sugerido.

F. Tusell
fuente
1
la forma de F(X)está más allá del alcance de este problema, pero es un modelo dinámico de vegetación global descrito en el apéndice de Medvigy et al 2009
David LeBauer