qsheppInt: Shepard's (modified) quadratic interpolation method

Description Usage Arguments Details Value Note Author(s) References Examples

Description

Shepard's (modified) quadratic interpolation method.

Usage

1
2
3
4
5
6
qsheppInt(X, y, XNew = NULL, nQ = 15L, nW = 33L, nR,
          checkX = TRUE)
qsheppInt2d(X, y, XNew = NULL, nQ = 13L, nW = 19L, nR,
            deriv = 0L, checkX = TRUE)
qsheppInt3d(X, y, XNew = NULL, nQ = 17L, nW = 32L, nR,
            deriv = 0L, checkX = TRUE)

Arguments

X

Design matrix. A matrix with n rows and d columns where n is the number of nodes and d is the dimension.

y

Numeric vector of length n containing the function values at design points.

XNew

Numeric matrix of new design points where interpolation is computed. It must have d columns.

nQ

Number of points used in the local polynomial fit. This is the parameter NQ of the Fortran routine qshepmd and the parameter N_p of the referenced article.

nW

Number of nodes within (and defining) the radii of influence, which enter into the weights. This is parameter NW of the Fortran routine qshepmd and N_w of the reference article.

nR

Number of divisions in each dimension for the cell grid defined in the subroutine storem. A hyperbox containing the nodes is partitioned into cells in order to increase search efficiency. The recommended value (n/3)^(1/d) is used as default, where n is the number of nodes and d is the dimension.

checkX

If TRUE, the dimensions and colnames of X and XNew will be checked.

deriv

Logical or integer value in c(0, 1). When the (coerced) integer value is 1, the derivatives are computed. This is is only possible for the dimensions 2 and 3, and only when the specific interpolation functions are used. The result is then a matrix with one column for the function and one column by derivative.

Details

Shepard's modified interpolation method works for scattered data in arbitrary dimension. It relies on a local polynomial fit and the quadratic version uses a polynomial of degree 2 (a quadratic form) as local approximation of the function.

Value

A list with several objects related to the computation method. The vector of interpolated value is in the list element named yNew.

Note

This function is an R interface to the qshepmd routine in the SHEPPACK Fortran package available on netlib http://www.netlib.org as algorithm 905A.

The qshepInt function is an interface for the QSHEPMD Fortran routine, while qshepInt2d and qshepInt3d are interfaces to the QSHEPM2D and QSHEPM3D Fortran routines. The general interpolation of qshepInt can be used also for the dimensions 2 and 3. However, this function does not allow the computation of the derivatives as qshepInt2d and qshepInt3d do.

Author(s)

Fortran code by William I. Thacker; Jingwei Zhang; Layne T. Watson; Jeffrey B. Birch; Manjula A. Iyer; Michael W. Berry. See References below. The Fortran code is a translation of M.W Berry's C++ code.

Code adaptation and R interface by Yves Deville.

References

W.I. Thacker, J. Zhang, L.T. Watson, J.B. Birch, M.A. Iyer and M.W. Berry (2010). Algorithm 905: SHEPPACK: Modified Shepard Algorithm for Interpolation of Scattered Multivariate Data ACM Trans. on Math. Software (TOMS) Vol. 37, n. 3. link

M.W. Berry and K.S. Minser (1999). Algorithm 798: High-dimensional interpolation using the modified Shepard method. ACM Trans. Math. Software (TOMS) Vol. 25, n. 3, pp. 353-366. link

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
n <- 1500; nNew <- 100; d <- 4
fTest <- function(x)((x[1] + 2 * x[2]   + 3 * x[3] + 4 * x[4]) / 12)^2
set.seed(12345)
X <- matrix(runif(n*d), nrow = n, ncol = d)
y <- apply(X, 1, FUN = fTest)
XNew <- matrix(runif(nNew * d), nrow = nNew, ncol = d)
yNew <- apply(XNew, 1, FUN = fTest)
system.time(res <- qsheppInt(X = X, XNew = XNew, y = y, nQ = 40,
                             checkX = FALSE))
## check errors
max(abs(res$yNew - yNew))

##=========================================================================
## Use SHEPPACK test functions see Thacker et al. section 7 'PERFORMANCE'
##=========================================================================
## Not run: 
   set.seed(1234)
   d <- 3
   k <- 0:4; n0 <- 100 * 2^k; n1 <- 4
   GD <- Grid(nlevels = rep(n1, d))
   XNew <- as.matrix(GD)
   RMSE <- array(NA, dim = c(5, length(k)),
                 dimnames = list(fun = 1:5, k = k))
   for (iFun in 1:5) {
      yNew <- apply(XNew, 1, ShepFuns[[iFun]])
      for (iN in 1:length(n0)) {
         X <- matrix(runif(n0[iN] * d), ncol = d)
         y <- apply(X, 1, ShepFuns[[iFun]])
         res <- qsheppInt(X = X, XNew = XNew, y = y, nQ = 40, checkX = FALSE)
         RMSE[iFun, iN] <- mean((res$yNew - yNew)^2)
      }
   }
   cols <- c("black", "SteelBlue2", "orangered", "SpringGreen3", "purple")
   pchs <-  c(16, 21, 22, 23, 24)
   matplot(k, t(RMSE), type = "o", lwd = 2, lty = 1,
           col = cols, xaxt = "n", pch = pchs, cex = 1.4,
           bg = "white",
           main = sprintf("dim = %d SHEPPACK test functions", d),
           xlab = "number of nodes", ylab = "RMSE")
   axis(side = 1, at = k, labels = n0)
   legend("topright", legend = paste("shepFun", 1:5),
           col = cols, pch = pchs, lwd = 2, pt.lwd = 2, pt.bg = "white")

## End(Not run)

smint documentation built on April 14, 2017, 1:49 p.m.