SPOT (version 2.5.18)

buildKriging: Build Kriging Model

Description

This function builds a Kriging model based on code by Forrester et al.. By default exponents (p) are fixed at a value of two, and a nugget (or regularization constant) is used. To correct the uncertainty estimates in case of nugget, re-interpolation is also by default turned on.

Usage

buildKriging(x, y, control = list())

Arguments

x

design matrix (sample locations)

y

vector of observations at x

control

(list), with the options for the model building procedure: types a character vector giving the data type of each variable. All but "factor" will be handled as numeric, "factor" (categorical) variables will be subject to the hamming distance. thetaLower lower boundary for theta, default is 1e-4 thetaUpper upper boundary for theta, default is 1e2 algTheta algorithm used to find theta, default is optimDE. budgetAlgTheta budget for the above mentioned algorithm, default is 200. The value will be multiplied with the length of the model parameter vector to be optimized. optimizeP boolean that specifies whether the exponents (p) should be optimized. Else they will be set to two. Default is FALSE useLambda whether or not to use the regularization constant lambda (nugget effect). Default is TRUE. lambdaLower lower boundary for log10lambda, default is -6 lambdaUpper upper boundary for log10lambda, default is 0 startTheta optional start value for theta optimization, default is NULL reinterpolate whether (TRUE,default) or not (FALSE) reinterpolation should be performed target target values of the prediction, a vector of strings. Each string specifies a value to be predicted, e.g., "y" for mean, "s" for standard deviation, "ei" for expected improvement. See also predict.kriging. This can also be changed after the model has been built, by manipulating the respective object$target value.

Value

an object of class kriging. Basically a list, with the options and found parameters for the model which has to be passed to the predictor function: x sample locations (scaled to values between 0 and 1) y observations at sample locations (see parameters) thetaLower lower boundary for theta (see parameters) thetaUpper upper boundary for theta (see parameters) algTheta algorithm to find theta (see parameters) budgetAlgTheta budget for the above mentioned algorithm (see parameters) optimizeP boolean that specifies whether the exponents (p) were optimized (see parameters) normalizeymin minimum in normalized space normalizeymax maximum in normalized space normalizexmin minimum in input space normalizexmax maximum in input space dmodeltheta vector of activity parameters Theta log_10 vector of activity parameters (i.e. log10(dmodeltheta)) dmodellambda regularization constant (nugget) Lambda log_10 of regularization constant (nugget) (i.e. log10(dmodellambda)) yonemu Ay-ones*mu ssq sigma square mu mean mu Psi matrix large Psi Psinv inverse of Psi nevals number of Likelihood evaluations during MLE

Details

The model uses a Gaussian kernel: k(x,z)=exp(-sum(theta_i * |x_i-z_i|^p_i)). By default, p_i = 2. Note that if dimension x_i is a factor variable (see parameter types), Hamming distance will be used instead of |x_i-z_i|.

References

Forrester, Alexander I.J.; Sobester, Andras; Keane, Andy J. (2008). Engineering Design via Surrogate Modelling - A Practical Guide. John Wiley & Sons.

See Also

predict.kriging

Examples

# NOT RUN {
## Create design points
set.seed(1)
x <- cbind(runif(20)*15-5,runif(20)*15)
## Compute observations at design points (for Branin function)
# y <- as.matrix(apply(x,1,braninFunction))
y <- funBranin(x)
## Create model with default settings
fit <- buildKriging(x,y,control = list(algTheta=optimLHD))
## Print model parameters
print(fit)
## Predict at new location
predict(fit,cbind(1,2))
## True value at location
funBranin(matrix(c(1,2), 1))
##
## Next Example: Handling factor variables
# }
# NOT RUN {
## create a test function:
braninFunctionFactor <- function (x) {
y <- (x[2]  - 5.1/(4 * pi^2) * (x[1] ^2) + 5/pi * x[1]  - 6)^2 +
		10 * (1 - 1/(8 * pi)) * cos(x[1] ) + 10
if(x[3]==1)
		y <- y +1
else if(x[3]==2)
		y <- y -1
y
}
## create training data
set.seed(1)
x <- cbind(runif(50)*15-5,runif(50)*15,sample(1:3,50,replace=TRUE))
y <- as.matrix(apply(x,1,braninFunctionFactor))
## fit the model (default: assume all variables are numeric)
fitDefault <- buildKriging(x,y,control = list(algTheta=optimDE))
## fit the model (give information about the factor variable)
fitFactor <- buildKriging(x,y,control =
	list(algTheta=optimDE,types=c("numeric","numeric","factor")))
## create test data
xtest <- cbind(runif(200)*15-5,runif(200)*15,sample(1:3,200,replace=TRUE))
ytest <- as.matrix(apply(xtest,1,braninFunctionFactor))
## Predict test data with both models, and compute error
ypredDef <- predict(fitDefault,xtest)$y
ypredFact <- predict(fitFactor,xtest)$y
mean((ypredDef-ytest)^2)
mean((ypredFact-ytest)^2)
# }