gridInt: Grid interpolation in arbitrary dimension

In IRSN/smint: Smooth Multivariate Interpolation for Gridded and Scattered Data

Grid interpolation in arbitrary dimension

Description

***************************************************************************** Grid interpolation in arbitrary dimension.

Usage

gridInt(
  X,
  Y,
  Xout,
  interpFun = function(x, y, xout) approx(x = x, y = y, xout = xout)$y,
  intOrder = NULL,
  useC = TRUE,
  trace = 1L,
  out_of_bounds = stop,
  ...
)

Arguments

X	An object that can be coerced into Grid. This can be a data frame or a matrix in Scattered Data style: the number of columns is then equal to the spatial dimension d, and the number of rows is then equal to the number of nodes n. But it can also be a Grid object previously created. A data frame or matrix X will first be coerced into Grid by using the the S3 method as.Grid.
Y	Response to be interpolated. It must be a vector with length n equal to the number of nodes. When X is a matrix or a data frame, the elements in Y match rows of X, so Y[i] is the value of the interpolated function at X[i, ]. If instead X is an object with call "Grid" the order elements of Y must be given in the order associated with the object X.
Xout	Interpolation locations. Can be a vector or a matrix. In the first case, the length of the vector must be equal to the spatial dimension d as given by the number of columns of X if this is a matrix, or by dim(X) if X is a Grid object. In the second case, each row will be considered as a response to be interpolated, and the number of columns of Xout must be equal to the spatial dimension.
interpFun	The function to interpolate. This function must have as its first 3 formals 'x', 'y', and 'xout', as does approx. It must also return the vector of interpolated values as DOES NOT approx. In most cases, a simple wrapper will be enough to use an existing method of interpolation, see Examples.
intOrder	Order of the one-dimensional interpolations. Must be a permutation of 1:d where d is the spatial dimension. By default, the interpolation order is d, d-1, ..., 1, corresponding to intOrder = d:1L. Note that for the first element of intOrder, a vector of length nrow(Xout) is passed as xout formal to interpFun, while the subsequent interFun calls use a xout of length 1. So the choice of the first element can have an impact on the computation time.
useC	Logical. If TRUE the computation is done using a C program via the .Call. Otherwise, the computation is done entirely in R by repeated use of the apply function. Normally useC = TRUE should be faster, but it is not always the case.
trace	Level of verbosity.
out_of_bounds	Function to handle Xout outside x (default is stop). Then Xout will be bounded by x range.
...	Further arguments to be passed to interpFun.

Details

The grid interpolation is performed by looping over dimensions. For each dimension, a one-dimensional interpolation is carried out, leading to a collection of interpolation problems each with a dimension reduced by one.

Value

A single interpolated value when Xout is either a vector or a row matrix. If Xout is a matrix with several rows, the result is a vector of interpolated values, in the order of the rows of Xout.

Note

A future multivariate version to come will allow the simultaneous interpolation of several responses. Most probably, this possibility will be used by using a different rule for the object Y (matrix or list).

When X is a Grid object and Y is a vector, the user must ensure that the order of nodes is the same for inputs and output. If so, the order of the dimensions in X can be changed using the generalised transposition aperm. This will change the order of the univariate interpolations with possible effects on the computation time and on the result when the univariate interpolation method is not linear (w.r.t. the response).

When X is a data frame or a matrix and Y is a vector, the values in Y are matched to the rows of X is the same as the order. This situation typically arises when X and Y are columns extracted from a same data frame.

Author(s)

Yves Deville

Examples

set.seed(12345)
##========================================================================
## Select Interpolation Function. This function must have its first 3
## formals 'x', 'y', and 'xout', as does 'approx'. It must also return
## the vector of interpolated values as DOES NOT 'approx'. So a wrapper
## must be xritten.
##=======================================================================
myInterp <- function(x, y, xout) approx(x = x, y = y, xout = xout)$y

##=======================================================================
## ONE interpolation, d = 2. 'Xout' is a vector.
##=======================================================================
myFun1 <- function(x) exp(-x[1]^2 - 3 * x[2]^2)
myGD1 <- Grid(nlevels = c("X" = 8, "Y" = 12))
Y1 <- apply_Grid(myGD1, myFun1)
Xout1 <- runif(2)
GI1 <- gridInt(X = myGD1,  Y = Y1, Xout = Xout1, interpFun = myInterp)
c(true = myFun1(Xout1), interp = GI1)

##=======================================================================
## ONE interpolation, d = 7. 'Xout' is a vector.
##=======================================================================
d <- 7; a <- runif(d); myFun2 <- function(x) exp(-crossprod(a, x^2))
myGD2 <- Grid(nlevels = rep(4L, time = d))
Y2 <- apply_Grid(myGD2, myFun2)
Xout2 <- runif(d)
GI2 <- gridInt(X = myGD2,  Y = Y2, Xout = Xout2, interpFun = myInterp)
c(true = myFun2(Xout2), interp = GI2)

##=======================================================================
## n interpolations, d = 7. 'Xout' is a matrix. Same grid data and
## response as before
##=======================================================================
n <- 30
Xout3 <- matrix(runif(n * d), ncol = d)
GI3 <- gridInt(X = myGD2,  Y = Y2, Xout = Xout3, interpFun = myInterp)
cbind(true = apply(Xout3, 1, myFun2), interp = GI3)

##======================================================================
## n interpolation, d = 5. 'Xout' is a matrix. Test the effect of the
## order of interpolation.
##=======================================================================
d <- 5; a <- runif(d); myFun4 <- function(x) exp(-crossprod(a, x^2))
myGD4 <- Grid(nlevels = c(3, 4, 5, 2, 6))
Y4 <- apply_Grid(myGD4, myFun4)
n <- 100
Xout4 <- matrix(runif(n * d), ncol = d)
t4a <- system.time(GI4a <- gridInt(X = myGD4,  Y = Y4, Xout = Xout4,
                                   interpFun = myInterp))
t4b <- system.time(GI4b <- gridInt(X = myGD4,  Y = Y4, Xout = Xout4,
                                   interpFun = myInterp,
                                   intOrder = 1L:5L))
cbind(true = apply(Xout4, 1, myFun4), inta = GI4a, intb = GI4b)

IRSN/smint documentation built on Dec. 9, 2023, 9:53 p.m.