| Title: | Kriging Models using the 'libKriging' Library |
|---|---|
| Description: | Interface to 'libKriging' 'C++' library <https://github.com/libKriging> that should provide most standard Kriging / Gaussian process regression features (like in 'DiceKriging', 'kergp' or 'RobustGaSP' packages). 'libKriging' relies on Armadillo linear algebra library (Apache 2 license) by Conrad Sanderson, 'lbfgsb_cpp' is a 'C++' port around by Pascal Have of 'lbfgsb' library (BSD-3 license) by Ciyou Zhu, Richard Byrd, Jorge Nocedal and Jose Luis Morales used for hyperparameters optimization. |
| Authors: | Yann Richet [aut, cre] (ORCID: <https://orcid.org/0000-0002-5677-8458>), Pascal Havé [aut], Yves Deville [aut], Conrad Sanderson [ctb], Ciyou Zhu [ctb], Richard Byrd [ctb], Jorge Nocedal [ctb], Jose Luis Morales [ctb], Mike Smith [ctb] |
| Maintainer: | Yann Richet <[email protected]> |
| License: | Apache License (>= 2) |
| Version: | 1.0-0 |
| Built: | 2026-06-12 10:54:24 UTC |
| Source: | https://github.com/cran/rlibkriging |
Get activation function name
activation(object, ...)activation(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get activation function for an MLPKriging model
## S3 method for class 'MLPKriging' activation(object, ...)## S3 method for class 'MLPKriging' activation(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
km ObjectCoerce an object into an object with S4 class "km" from the
DiceKriging package.
as.km(x, ...)as.km(x, ...)
x |
Object to be coerced. |
... |
Further arguments for methods. |
Such a coercion is typically used to compare the performance of the methods implemented in the current rlibkriging package to those which are available in the DiceKriging package.
An object with S4 class "km".
Kriging object into the "km" class of the
DiceKriging package.Coerce a Kriging object into the "km" class of the
DiceKriging package.
## S3 method for class 'Kriging' as.km(x, .call = NULL, ...)## S3 method for class 'Kriging' as.km(x, .call = NULL, ...)
x |
An object with S3 class |
.call |
Force the |
... |
Not used. |
An object of having the S4 class "KM" which extends
the "km" class of the DiceKriging package and
contains an extra Kriging slot.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2") print(k) k_km <- as.km(k) print(k_km)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2") print(k) k_km <- as.km(k) print(k_km)
Kriging Object into a ListCoerce a Kriging Object into a List
## S3 method for class 'Kriging' as.list(x, ...)## S3 method for class 'Kriging' as.list(x, ...)
x |
An object with class |
... |
Ignored |
A list with its elements copying the content of the
Kriging object fields: kernel, optim,
objective, theta (vector of ranges),
sigma2 (variance), X, centerX,
scaleX, y, centerY, scaleY,
regmodel, F, T, M, z,
beta.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x ) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") l <- as.list(k) cat(paste0(names(l), " =" , l, collapse = "\n"))f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x ) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") l <- as.list(k) cat(paste0(names(l), " =" , l, collapse = "\n"))
Get trend coefficients beta
beta(object, ...)beta(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get trend coefficients beta for an MLPKriging model
## S3 method for class 'MLPKriging' beta(object, ...)## S3 method for class 'MLPKriging' beta(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get trend coefficients beta for a WarpKriging model
## S3 method for class 'WarpKriging' beta(object, ...)## S3 method for class 'WarpKriging' beta(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get input centering vector
centerX(object, ...)centerX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get input centering vector for an MLPKriging model
## S3 method for class 'MLPKriging' centerX(object, ...)## S3 method for class 'MLPKriging' centerX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get input centering vector for a WarpKriging model
## S3 method for class 'WarpKriging' centerX(object, ...)## S3 method for class 'WarpKriging' centerX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output centering value
centerY(object, ...)centerY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output centering value for an MLPKriging model
## S3 method for class 'MLPKriging' centerY(object, ...)## S3 method for class 'MLPKriging' centerY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output centering value for a WarpKriging model
## S3 method for class 'WarpKriging' centerY(object, ...)## S3 method for class 'WarpKriging' centerY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Shortcut to provide functions to the S3 class "Kriging"
classKriging(nk)classKriging(nk)
nk |
A pointer to a C++ object of class "Kriging" |
An object of class "Kriging" with methods to access and manipulate the data
Shortcut to provide functions to the S3 class "MLPKriging"
classMLPKriging(obj)classMLPKriging(obj)
obj |
A list with a |
An object of class "MLPKriging" with methods accessible via $
Shortcut to provide functions to the S3 class "WarpKriging"
classWarpKriging(obj)classWarpKriging(obj)
obj |
A list with a |
An object of class "WarpKriging" with methods accessible via $
Duplicate a model given in
object.
copy(object, ...)copy(object, ...)
object |
An object representing a fitted model. |
... |
Ignored. |
The copied object.
Duplicate a Kriging Model
## S3 method for class 'Kriging' copy(object, ...)## S3 method for class 'Kriging' copy(object, ...)
object |
An S3 Kriging object. |
... |
Not used. |
The copy of object.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) print(copy(k))f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) print(copy(k))
Deep copy of MLPKriging model
## S3 method for class 'MLPKriging' copy(object, ...)## S3 method for class 'MLPKriging' copy(object, ...)
object |
MLPKriging object |
... |
ignored |
a new independent MLPKriging object
Deep copy of WarpKriging model
## S3 method for class 'WarpKriging' copy(object, ...)## S3 method for class 'WarpKriging' copy(object, ...)
object |
WarpKriging object |
... |
ignored |
a new independent WarpKriging object
Compute the covariance matrix of a model given in object,
between given set of points.
covMat(object, ...)covMat(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function (eg. points, range). |
The covariance matrix.
Compute Covariance Matrix of Kriging Model
## S3 method for class 'Kriging' covMat(object, x1, x2, ...)## S3 method for class 'Kriging' covMat(object, x1, x2, ...)
object |
An S3 Kriging object. |
x1 |
Numeric matrix of input points. |
x2 |
Numeric matrix of input points. |
... |
Not used. |
A matrix of the covariance matrix of the Kriging model.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "gauss") x1 = runif(10) x2 = runif(10) covMat(k, x1, x2)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "gauss") x1 = runif(10) x2 = runif(10) covMat(k, x1, x2)
Get trend matrix F
F_(object, ...)F_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get trend matrix F for an MLPKriging model
## S3 method for class 'MLPKriging' F_(object, ...)## S3 method for class 'MLPKriging' F_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get trend matrix F for a WarpKriging model
## S3 method for class 'WarpKriging' F_(object, ...)## S3 method for class 'WarpKriging' F_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get feature dimensionality (d_out)
feature_dim(object, ...)feature_dim(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get feature dimensionality for an MLPKriging model
## S3 method for class 'MLPKriging' feature_dim(object, ...)## S3 method for class 'MLPKriging' feature_dim(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get feature dimensionality of warped space
## S3 method for class 'WarpKriging' feature_dim(object, ...)## S3 method for class 'WarpKriging' feature_dim(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Fit a model given in
object.
fit(object, ...)fit(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function |
No return value. Kriging object argument is modified.
Kriging object on given data.The hyper-parameters (variance and vector of correlation ranges)
are estimated thanks to the optimization of a criterion given by
objective, using the method given in optim.
## S3 method for class 'Kriging' fit( object, y, X, noise = NULL, regmodel = c("constant", "linear", "interactive", "none"), normalize = FALSE, optim = c("BFGS", "none"), objective = c("LL", "LOO", "LMP"), parameters = NULL, ... )## S3 method for class 'Kriging' fit( object, y, X, noise = NULL, regmodel = c("constant", "linear", "interactive", "none"), normalize = FALSE, optim = c("BFGS", "none"), objective = c("LL", "LOO", "LMP"), parameters = NULL, ... )
object |
S3 Kriging object. |
y |
Numeric vector of response values. |
X |
Numeric matrix of input design. |
noise |
Either a numeric vector of per-observation noise variances,
or |
regmodel |
Universal Kriging linear trend: |
normalize |
Logical. If |
optim |
Character giving the Optimization method used to fit
hyper-parameters. Possible values are: |
objective |
Character giving the objective function to
optimize. Possible values are: |
parameters |
Initial values for the hyper-parameters. When
provided this must be named list with elements |
... |
Ignored. |
No return value. Kriging object argument is modified.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue", pch = 16) k <- Kriging("matern3_2") print(k) fit(k,y,X) print(k)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue", pch = 16) k <- Kriging("matern3_2") print(k) fit(k,y,X) print(k)
(Re-)fit an already-constructed MLPKriging object on new data. The MLP architecture and kernel are kept from construction.
## S3 method for class 'MLPKriging' fit( object, y, X, regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, ... )## S3 method for class 'MLPKriging' fit( object, y, X, regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, ... )
object |
MLPKriging object |
y |
numeric vector of observations (n) |
X |
numeric matrix of inputs (n x d) |
regmodel |
trend: "constant", "linear", "quadratic" |
normalize |
logical; normalise inputs? |
optim |
optimiser |
objective |
"LL" (log-likelihood) |
parameters |
optional named list of tuning parameters |
... |
ignored |
No return value. MLPKriging object argument is modified.
(Re-)fit an already-constructed WarpKriging object on new data. The warping specification and kernel are kept from construction.
## S3 method for class 'WarpKriging' fit( object, y, X, regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, noise = NULL, ... )## S3 method for class 'WarpKriging' fit( object, y, X, regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, noise = NULL, ... )
object |
WarpKriging object (created with the constructor or an empty-kernel call) |
y |
numeric vector of observations (n) |
X |
numeric matrix of inputs (n x d) |
regmodel |
trend: "constant", "linear", "quadratic" |
normalize |
logical; normalise continuous inputs? |
optim |
optimiser |
objective |
"LL" (log-likelihood) |
parameters |
optional named list of tuning parameters |
noise |
Either a numeric vector of per-observation noise variances,
or |
... |
ignored |
No return value. WarpKriging object argument is modified.
Check if the model has been fitted
is_fitted(object, ...)is_fitted(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Check whether an MLPKriging model is fitted
## S3 method for class 'MLPKriging' is_fitted(object, ...)## S3 method for class 'MLPKriging' is_fitted(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Check whether a WarpKriging model is fitted
## S3 method for class 'WarpKriging' is_fitted(object, ...)## S3 method for class 'WarpKriging' is_fitted(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get kernel name
kernel(object, ...)kernel(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Character scalar naming the covariance kernel.
Get kernel name
## S3 method for class 'WarpKriging' kernel(object, ...)## S3 method for class 'WarpKriging' kernel(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
KM ObjectCreate an object of S4 class "KM" similar to a
km object in the DiceKriging package.
KM( formula = ~1, design, response, covtype = c("matern5_2", "gauss", "matern3_2", "exp"), coef.trend = NULL, coef.cov = NULL, coef.var = NULL, nugget = NULL, nugget.estim = FALSE, noise.var = NULL, estim.method = c("MLE", "LOO"), penalty = NULL, optim.method = "BFGS", lower = NULL, upper = NULL, parinit = NULL, multistart = 1, control = NULL, gr = TRUE, iso = FALSE, scaling = FALSE, knots = NULL, kernel = NULL, ... )KM( formula = ~1, design, response, covtype = c("matern5_2", "gauss", "matern3_2", "exp"), coef.trend = NULL, coef.cov = NULL, coef.var = NULL, nugget = NULL, nugget.estim = FALSE, noise.var = NULL, estim.method = c("MLE", "LOO"), penalty = NULL, optim.method = "BFGS", lower = NULL, upper = NULL, parinit = NULL, multistart = 1, control = NULL, gr = TRUE, iso = FALSE, scaling = FALSE, knots = NULL, kernel = NULL, ... )
formula |
R formula object to setup the linear trend in
Universal Kriging. Supports |
design |
Data frame. The design of experiments. |
response |
Vector of output values. |
covtype |
Covariance structure. For now all the kernels are tensor product kernels. |
coef.trend |
Optional value for a fixed vector of trend coefficients. If given, no optimization is done. |
coef.cov |
Optional value for a fixed correlation range value. If given, no optimization is done. |
coef.var |
Optional value for a fixed variance. If given, no optimization is done. |
nugget |
Ignored (kept for DiceKriging API compatibility). |
nugget.estim |
Logical. If |
noise.var |
Numeric vector of known per-point noise variances. |
estim.method |
Estimation criterion. |
penalty |
Not implemented yet. |
optim.method |
Optimization algorithm used in the
optimization of the objective given in
|
lower, upper
|
Not implemented yet. |
parinit |
Initial values for the correlation ranges which
will be optimized using |
multistart, control, gr, iso
|
Not implemented yet. |
scaling, knots, kernel
|
Not implemented yet. |
... |
Ignored. |
The class "KM" extends the "km" class of the
DiceKriging package, hence has all slots of "km". It
also has an extra slot "Kriging" slot which contains a copy
of the original object.
A KM object. See Details.
Yann Richet [email protected]
km in the DiceKriging
package for more details on the slots.
# a 16-points factorial design, and the corresponding response d <- 2; n <- 16 design.fact <- as.matrix(expand.grid(x1 = seq(0, 1, length = 4), x2 = seq(0, 1, length = 4))) y <- apply(design.fact, 1, DiceKriging::branin) # Using `km` from DiceKriging and a similar `KM` object # kriging model 1 : matern5_2 covariance structure, no trend, no nugget effect km1 <- DiceKriging::km(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1), control = list(trace = FALSE)) KM1 <- KM(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1))# a 16-points factorial design, and the corresponding response d <- 2; n <- 16 design.fact <- as.matrix(expand.grid(x1 = seq(0, 1, length = 4), x2 = seq(0, 1, length = 4))) y <- apply(design.fact, 1, DiceKriging::branin) # Using `km` from DiceKriging and a similar `KM` object # kriging model 1 : matern5_2 covariance structure, no trend, no nugget effect km1 <- DiceKriging::km(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1), control = list(trace = FALSE)) KM1 <- KM(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1))
"km" ClassThis class is intended to be used either by using its
own dedicated S4 methods or by using the S4 methods inherited
from the "km" class of the libKriging package.
d,n,X,y,p,FNumber of (numeric) inputs, number of observations, design matrix, response vector, number of trend variables, trend matrix.
trend.formula,trend.coefFormula used for the trend, vector
of estimated (or fixed)
trend coefficients with length .
covarianceA S4 object with class "covTensorProduct"
representing a covariance kernel.
noise.flag,noise.varLogical flag and numeric value for an optional noise term.
known.paramA character code indicating what parameters are known.
lower,upperBounds on the correlation range parameters.
method,penalty,optim.method,control,gr,parinitObjects defining the estimation criterion, the optimization.
T,M,zAuxiliary variables (matrices and vectors) that can be used in several computations.
caseThe possible concentration (a.k.a. profiling) of the likelihood.
param.estimLogical. Is an estimation used?
KrigingA copy of the Kriging object used to create
the current KM object.
Yann Richet [email protected]
km-class in the
DiceKriging package. The creator KM.
"Kriging" using
the libKriging library.The hyper-parameters (variance and vector of correlation ranges)
are estimated thanks to the optimization of a criterion given by
objective, using the method given in optim.
Kriging( y = NULL, X = NULL, kernel = NULL, noise = NULL, regmodel = c("constant", "linear", "interactive", "none"), normalize = FALSE, optim = c("BFGS", "none"), objective = c("LL", "LOO", "LMP"), parameters = NULL )Kriging( y = NULL, X = NULL, kernel = NULL, noise = NULL, regmodel = c("constant", "linear", "interactive", "none"), normalize = FALSE, optim = c("BFGS", "none"), objective = c("LL", "LOO", "LMP"), parameters = NULL )
y |
Numeric vector of response values. |
X |
Numeric matrix of input design. |
kernel |
Character defining the covariance model:
|
noise |
Either a numeric vector of per-observation noise variances,
or |
regmodel |
Universal Kriging linear trend: |
normalize |
Logical. If |
optim |
Character giving the Optimization method used to fit
hyper-parameters. Possible values are: |
objective |
Character giving the objective function to
optimize. Possible values are: |
parameters |
Initial values for the hyper-parameters. When
provided this must be named list with elements |
An object with S3 class "Kriging". Should be used
with its predict, simulate, update
methods.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) ## fit and print k <- Kriging(y, X, kernel = "matern3_2") print(k) x <- as.matrix(seq(from = 0, to = 1, length.out = 101)) p <- predict(k, x = x, return_stdev = TRUE, return_cov = FALSE) plot(f) points(X, y) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) s <- simulate(k, nsim = 10, seed = 123, x = x) matlines(x, s, col = rgb(0, 0, 1, 0.2), type = "l", lty = 1)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) ## fit and print k <- Kriging(y, X, kernel = "matern3_2") print(k) x <- as.matrix(seq(from = 0, to = 1, length.out = 101)) p <- predict(k, x = x, return_stdev = TRUE, return_cov = FALSE) plot(f) points(X, y) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) s <- simulate(k, nsim = 10, seed = 123, x = x) matlines(x, s, col = rgb(0, 0, 1, 0.2), type = "l", lty = 1)
Compute the leave-One-Out error of a model given in object.
leaveOneOut(object, ...)leaveOneOut(object, ...)
object |
An object representing a fitted model. |
... |
Ignored. |
The Leave-One-Out sum of squares.
Get leaveOneOut of Kriging Model
## S3 method for class 'Kriging' leaveOneOut(object, ...)## S3 method for class 'Kriging' leaveOneOut(object, ...)
object |
An S3 Kriging object. |
... |
Not used. |
The leaveOneOut computed for fitted
.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LOO") print(k) leaveOneOut(k)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LOO") print(k) leaveOneOut(k)
Compute the leave-One-Out error of a model given in object,
at a different value of the parameters.
leaveOneOutFun(object, ...)leaveOneOutFun(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function (eg. range). |
The Leave-One-Out sum of squares.
"Kriging" representing a kriging model.The returned value is the sum of squares where is the
prediction of based on the the observations
with .
## S3 method for class 'Kriging' leaveOneOutFun(object, theta, return_grad = FALSE, bench = FALSE, ...)## S3 method for class 'Kriging' leaveOneOutFun(object, theta, return_grad = FALSE, bench = FALSE, ...)
object |
A |
theta |
A numeric vector of range parameters at which the LOO will be evaluated. |
return_grad |
Logical. Should the gradient (w.r.t. |
bench |
Logical. Should the function display benchmarking output |
... |
Not used. |
The leave-One-Out value computed for the given vector
of correlation ranges.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective = "LOO", optim="BFGS") print(k) loo <- function(theta) leaveOneOutFun(k, theta)$leaveOneOut t <- seq(from = 0.001, to = 2, length.out = 101) plot(t, loo(t), type = "l") abline(v = k$theta(), col = "blue")f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective = "LOO", optim="BFGS") print(k) loo <- function(theta) leaveOneOutFun(k, theta)$leaveOneOut t <- seq(from = 0.001, to = 2, length.out = 101) plot(t, loo(t), type = "l") abline(v = k$theta(), col = "blue")
Compute the leave-One-Out vector error of a model given in object,
at a different value of the parameters.
leaveOneOutVec(object, ...)leaveOneOutVec(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function (eg. range). |
The Leave-One-Out errors (mean and stdev) for each conditional point.
"Kriging" representing a kriging model.The returned value is the mean and stdev of , the
prediction of based on the the observations
with .
## S3 method for class 'Kriging' leaveOneOutVec(object, theta, ...)## S3 method for class 'Kriging' leaveOneOutVec(object, theta, ...)
object |
A |
theta |
A numeric vector of range parameters at which the LOO will be evaluated. |
... |
Not used. |
The leave-One-Out vector computed for the given vector
of correlation ranges.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(c(0.0, 0.25, 0.5, 0.75, 1.0)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") print(k) x <- as.matrix(seq(0, 1, , 101)) p <- predict(k, x, TRUE, FALSE) plot(f) points(X, y) lines(x, p$mean, col = 'blue') polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) # Compute leave-one-out (no range re-estimate) on 2nd point X_no2 = X[-2,,drop=FALSE] y_no2 = f(X_no2) k_no2 = Kriging(y_no2, X_no2, "matern3_2", optim = "none", parameters = list(theta = k$theta())) print(k_no2) p_no2 <- predict(k_no2, x, TRUE, FALSE) lines(x, p_no2$mean, col = 'red') polygon(c(x, rev(x)), c(p_no2$mean - 2 * p_no2$stdev, rev(p_no2$mean + 2 * p_no2$stdev)), border = NA, col = rgb(1, 0, 0, 0.2)) # Use leaveOneOutVec to get the same loov = k$leaveOneOutVec(matrix(k$theta())) points(X[2],loov$mean[2],col='red') lines(rep(X[2],2),loov$mean[2]+2*c(-loov$stdev[2],loov$stdev[2]),col='red')f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(c(0.0, 0.25, 0.5, 0.75, 1.0)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") print(k) x <- as.matrix(seq(0, 1, , 101)) p <- predict(k, x, TRUE, FALSE) plot(f) points(X, y) lines(x, p$mean, col = 'blue') polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) # Compute leave-one-out (no range re-estimate) on 2nd point X_no2 = X[-2,,drop=FALSE] y_no2 = f(X_no2) k_no2 = Kriging(y_no2, X_no2, "matern3_2", optim = "none", parameters = list(theta = k$theta())) print(k_no2) p_no2 <- predict(k_no2, x, TRUE, FALSE) lines(x, p_no2$mean, col = 'red') polygon(c(x, rev(x)), c(p_no2$mean - 2 * p_no2$stdev, rev(p_no2$mean + 2 * p_no2$stdev)), border = NA, col = rgb(1, 0, 0, 0.2)) # Use leaveOneOutVec to get the same loov = k$leaveOneOutVec(matrix(k$theta())) points(X[2],loov$mean[2],col='red') lines(rep(X[2],2),loov$mean[2]+2*c(-loov$stdev[2],loov$stdev[2]),col='red')
Load any Kriging Model from a file storage. Back to base::load if not a Kriging object.
load(filename, ...)load(filename, ...)
filename |
A file holding any Kriging object. |
... |
Arguments used by base::load. |
The loaded "*"Kriging object, or nothing if base::load is used (update parent environment).
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) print(load(outfile))f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) print(load(outfile))
Load a Kriging Model from a file storage
load.Kriging(filename, ...)load.Kriging(filename, ...)
filename |
File name to load from. |
... |
Not used. |
The loaded Kriging object.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) print(load.Kriging(outfile)) unlink(outfile)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) print(load.Kriging(outfile)) unlink(outfile)
Load an MLPKriging model from file
## S3 method for class 'MLPKriging' load(filename, ...)## S3 method for class 'MLPKriging' load(filename, ...)
filename |
path to saved file |
... |
ignored |
MLPKriging object
Load a WarpKriging model from file
## S3 method for class 'WarpKriging' load(filename, ...)## S3 method for class 'WarpKriging' load(filename, ...)
filename |
path to saved file |
... |
ignored |
WarpKriging object
Compute the log-Likelihood of a model given in object.
logLikelihood(object, ...)logLikelihood(object, ...)
object |
An object representing a fitted model. |
... |
Ignored. |
The log-likelihood.
Get Log-Likelihood of Kriging Model
## S3 method for class 'Kriging' logLikelihood(object, ...)## S3 method for class 'Kriging' logLikelihood(object, ...)
object |
An S3 Kriging object. |
... |
Not used. |
The log-Likelihood computed for fitted
.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LL") print(k) logLikelihood(k)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LL") print(k) logLikelihood(k)
Log-likelihood of the fitted model
## S3 method for class 'WarpKriging' logLikelihood(object, ...)## S3 method for class 'WarpKriging' logLikelihood(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Compute the log-Likelihood of a model given in object,
at a different value of the parameters.
logLikelihoodFun(object, ...)logLikelihoodFun(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function (eg. range). |
The log-likelihood.
Compute Log-Likelihood of Kriging Model
## S3 method for class 'Kriging' logLikelihoodFun( object, theta, return_grad = FALSE, return_hess = FALSE, bench = FALSE, ... )## S3 method for class 'Kriging' logLikelihoodFun( object, theta, return_grad = FALSE, return_hess = FALSE, bench = FALSE, ... )
object |
An S3 Kriging object. |
theta |
A numeric vector of (positive) range parameters at which the log-likelihood will be evaluated. |
return_grad |
Logical. Should the function return the gradient? |
return_hess |
Logical. Should the function return Hessian? |
bench |
Logical. Should the function display benchmarking output? |
... |
Not used. |
The log-Likelihood computed for given
.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") print(k) ll <- function(theta) logLikelihoodFun(k, theta)$logLikelihood t <- seq(from = 0.001, to = 2, length.out = 101) plot(t, ll(t), type = 'l') abline(v = k$theta(), col = "blue")f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2") print(k) ll <- function(theta) logLikelihoodFun(k, theta)$logLikelihood t <- seq(from = 0.001, to = 2, length.out = 101) plot(t, ll(t), type = 'l') abline(v = k$theta(), col = "blue")
Evaluate log-likelihood at given GP theta
## S3 method for class 'MLPKriging' logLikelihoodFun(object, theta, return_grad = FALSE, return_hess = FALSE, ...)## S3 method for class 'MLPKriging' logLikelihoodFun(object, theta, return_grad = FALSE, return_hess = FALSE, ...)
object |
MLPKriging object |
theta |
range parameter vector |
return_grad |
return gradient? |
return_hess |
return hessian? |
... |
ignored |
Evaluate log-likelihood at given theta
## S3 method for class 'WarpKriging' logLikelihoodFun(object, theta, return_grad = FALSE, return_hess = FALSE, ...)## S3 method for class 'WarpKriging' logLikelihoodFun(object, theta, return_grad = FALSE, return_hess = FALSE, ...)
object |
WarpKriging object |
theta |
range parameter vector |
return_grad |
return gradient? |
return_hess |
return hessian? |
... |
ignored |
list with logLikelihood and optionally gradient, hessian
Compute the log-Marginal Posterior of a model given in
object.
logMargPost(object, ...)logMargPost(object, ...)
object |
An object representing a fitted model. |
... |
Ignored. |
The log-marginal posterior.
Get logMargPost of Kriging Model
## S3 method for class 'Kriging' logMargPost(object, ...)## S3 method for class 'Kriging' logMargPost(object, ...)
object |
An S3 Kriging object. |
... |
Not used. |
The logMargPost computed for fitted
.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) logMargPost(k)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) logMargPost(k)
Compute the log-Marginal Posterior of a model given in
object, at a different value of the parameters.
logMargPostFun(object, ...)logMargPostFun(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function (eg. range). |
The log-marginal posterior.
Compute the log-marginal posterior of a kriging model, using the prior XXXY.
## S3 method for class 'Kriging' logMargPostFun(object, theta, return_grad = FALSE, bench = FALSE, ...)## S3 method for class 'Kriging' logMargPostFun(object, theta, return_grad = FALSE, bench = FALSE, ...)
object |
S3 Kriging object. |
theta |
Numeric vector of correlation range parameters at which the function is to be evaluated. |
return_grad |
Logical. Should the function return the gradient (w.r.t theta)? |
bench |
Logical. Should the function display benchmarking output? |
... |
Not used. |
The value of the log-marginal posterior computed for the given vector theta.
Yann Richet [email protected]
XXXY A reference describing the model (prior, ...)
rgasp in the RobustGaSP package.
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2", objective="LMP") print(k) lmp <- function(theta) logMargPostFun(k, theta)$logMargPost t <- seq(from = 0.01, to = 2, length.out = 101) plot(t, lmp(t), type = "l") abline(v = k$theta(), col = "blue")f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2", objective="LMP") print(k) lmp <- function(theta) logMargPostFun(k, theta)$logMargPost t <- seq(from = 0.01, to = 2, length.out = 101) plot(t, lmp(t), type = "l") abline(v = k$theta(), col = "blue")
Get whitened trend matrix M
M(object, ...)M(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get whitened trend matrix M for an MLPKriging model
## S3 method for class 'MLPKriging' M(object, ...)## S3 method for class 'MLPKriging' M(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get whitened trend matrix M for a WarpKriging model
## S3 method for class 'WarpKriging' M(object, ...)## S3 method for class 'WarpKriging' M(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Kriging with a joint multi-layer perceptron applied to all inputs before the GP kernel is evaluated. The MLP weights, GP range parameters, variance and trend are jointly fitted by maximising the concentrated log-likelihood.
MLPKriging( y, X, hidden_dims, d_out = 2, activation = "selu", kernel = "gauss", regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL )MLPKriging( y, X, hidden_dims, d_out = 2, activation = "selu", kernel = "gauss", regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL )
y |
numeric vector of observations (n) |
X |
numeric matrix of inputs (n x d) |
|
integer vector of hidden layer sizes, e.g. |
|
d_out |
output feature dimensionality (default 2) |
activation |
activation function: "relu", "selu", "tanh", "sigmoid", "elu" |
kernel |
covariance kernel: "gauss", "matern3_2", "matern5_2", "exp" |
regmodel |
trend: "constant", "linear", "quadratic" |
normalize |
logical; normalise inputs? |
optim |
optimiser (default "BFGS+Adam") |
objective |
"LL" (log-likelihood) |
parameters |
optional named list of tuning parameters, e.g.
|
An S3 object of class "MLPKriging".
X <- as.matrix(seq(0.01, 0.99, length.out = 10)) f <- function(x) 1 - 1/2 * (sin(12*x)/(1+x) + 2*cos(7*x)*x^5 + 0.7) y <- f(X) k <- MLPKriging(y, X, hidden_dims = c(16, 8), d_out = 2, activation = "selu", kernel = "gauss") print(k)X <- as.matrix(seq(0.01, 0.99, length.out = 10)) f <- function(x) 1 - 1/2 * (sin(12*x)/(1+x) + 2*cos(7*x)*x^5 + 0.7) y <- f(X) k <- MLPKriging(y, X, hidden_dims = c(16, 8), d_out = 2, activation = "selu", kernel = "gauss") print(k)
This is a thin compatibility wrapper. Use
KM(noise.var = ...) directly instead.
NoiseKM(..., noise.var)NoiseKM(..., noise.var)
... |
Arguments passed to |
noise.var |
Numeric vector of known per-point noise variances. |
A KM object. See KM.
Get normalize flag
normalize(object, ...)normalize(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get normalize flag for an MLPKriging model
## S3 method for class 'MLPKriging' normalize(object, ...)## S3 method for class 'MLPKriging' normalize(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get normalize flag for a WarpKriging model
## S3 method for class 'WarpKriging' normalize(object, ...)## S3 method for class 'WarpKriging' normalize(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
This is a thin compatibility wrapper. Use
KM(nugget.estim = TRUE) directly instead.
NuggetKM(...)NuggetKM(...)
... |
Arguments passed to |
A KM object. See KM.
KM ObjectCompute predictions for the response at new given input points. These conditional mean, the conditional standard deviation and confidence limits at the 95% level. Optionnally the conditional covariance can be returned as well.
## S4 method for signature 'KM' predict( object, newdata, type = "UK", se.compute = TRUE, cov.compute = FALSE, light.return = TRUE, bias.correct = FALSE, checkNames = FALSE, ... )## S4 method for signature 'KM' predict( object, newdata, type = "UK", se.compute = TRUE, cov.compute = FALSE, light.return = TRUE, bias.correct = FALSE, checkNames = FALSE, ... )
object |
|
newdata |
Matrix of "new" input points where to perform prediction. |
type |
character giving the kriging type. For now only
|
se.compute |
Logical. Should the standard error be computed? |
cov.compute |
Logical. Should the covariance matrix between newdata points be computed? |
light.return |
Logical. If |
bias.correct |
Logical. If |
checkNames |
Logical to check the consistency of the column
names between the design stored in |
... |
Ignored. |
Without a dedicated predict method for the class
"KM", this method would have been inherited from the
"km" class. The dedicated method is expected to run faster.
A comparison can be made by coercing a KM object to a
km object with as.km before calling
predict.
A named list. The elements are the conditional mean and
standard deviation (mean and sd), the predicted
trend (trend) and the confidence limits (lower95
and upper95). Optionnally, the conditional covariance matrix
is returned in cov.
Yann Richet [email protected]
## a 16-points factorial design, and the corresponding response d <- 2; n <- 16 design.fact <- expand.grid(x1 = seq(0, 1, length = 4), x2 = seq(0, 1, length = 4)) y <- apply(design.fact, 1, DiceKriging::branin) ## library(DiceKriging) ## kriging model 1 : matern5_2 covariance structure, no trend, no nugget ## m1 <- km(design = design.fact, response = y, covtype = "gauss", ## parinit = c(.5, 1), control = list(trace = FALSE)) KM1 <- KM(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1)) Pred <- predict(KM1, newdata = matrix(.5,ncol = 2), type = "UK", checkNames = FALSE, light.return = TRUE)## a 16-points factorial design, and the corresponding response d <- 2; n <- 16 design.fact <- expand.grid(x1 = seq(0, 1, length = 4), x2 = seq(0, 1, length = 4)) y <- apply(design.fact, 1, DiceKriging::branin) ## library(DiceKriging) ## kriging model 1 : matern5_2 covariance structure, no trend, no nugget ## m1 <- km(design = design.fact, response = y, covtype = "gauss", ## parinit = c(.5, 1), control = list(trace = FALSE)) KM1 <- KM(design = design.fact, response = y, covtype = "gauss", parinit = c(.5, 1)) Pred <- predict(KM1, newdata = matrix(.5,ncol = 2), type = "UK", checkNames = FALSE, light.return = TRUE)
Kriging object.Given "new" input points, the method compute the expectation, variance and (optionnally) the covariance of the corresponding stochastic process, conditional on the values at the input points used when fitting the model.
## S3 method for class 'Kriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )## S3 method for class 'Kriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )
object |
S3 Kriging object. |
x |
Input points where the prediction must be computed. |
return_stdev |
|
return_cov |
|
return_deriv |
|
... |
Ignored. |
A list containing the element mean and possibly
stdev and cov.
The names of the formal arguments differ from those of the
predict methods for the S4 classes "km" and
"KM". The formal x corresponds to
newdata, stdev corresponds to se.compute
and cov to cov.compute. These names are chosen
Python and Octave interfaces to libKriging.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue", pch = 16) k <- Kriging(y, X, "matern3_2") x <-seq(from = 0, to = 1, length.out = 101) p <- predict(k, x) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2))f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue", pch = 16) k <- Kriging(y, X, "matern3_2") x <-seq(from = 0, to = 1, length.out = 101) p <- predict(k, x) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2))
Predict with an MLPKriging model
## S3 method for class 'MLPKriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )## S3 method for class 'MLPKriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )
object |
MLPKriging object |
x |
prediction matrix (m x d) |
return_stdev |
return standard deviations? |
return_cov |
return full covariance? |
return_deriv |
return derivatives of mean and stdev wrt x? |
... |
ignored |
list with mean, optionally stdev, cov,
mean_deriv, stdev_deriv
Predict with a WarpKriging model
## S3 method for class 'WarpKriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )## S3 method for class 'WarpKriging' predict( object, x, return_stdev = TRUE, return_cov = FALSE, return_deriv = FALSE, ... )
object |
WarpKriging object |
x |
prediction matrix (m x d) |
return_stdev |
return standard deviations? |
return_cov |
return full covariance? |
return_deriv |
return derivatives of mean and stdev wrt x? |
... |
ignored |
list with mean, optionally stdev, cov,
mean_deriv, stdev_deriv
Kriging object.Print the content of a Kriging object.
## S3 method for class 'Kriging' print(x, ...)## S3 method for class 'Kriging' print(x, ...)
x |
A (S3) |
... |
Ignored. |
String of printed object.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2") print(k) ## same thing kf <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, "matern3_2") print(k) ## same thing k
Get regression model type
regmodel(object, ...)regmodel(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get regression model type for an MLPKriging model
## S3 method for class 'MLPKriging' regmodel(object, ...)## S3 method for class 'MLPKriging' regmodel(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get regression model type for a WarpKriging model
## S3 method for class 'WarpKriging' regmodel(object, ...)## S3 method for class 'WarpKriging' regmodel(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Save a Kriging Model inside a file. Back to base::save if argument is not a Kriging object.
save(object = NULL, filename = NULL, ...)save(object = NULL, filename = NULL, ...)
object |
An object representing a model. |
filename |
A file to save the object. |
... |
Arguments used by base::save. |
Yann Richet [email protected]
Save a Kriging Model to a file storage
## S3 method for class 'Kriging' save(object, filename, ...)## S3 method for class 'Kriging' save(object, filename, ...)
object |
An S3 Kriging object. |
filename |
File name to save in. |
... |
Not used. |
The loaded Kriging object.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) unlink(outfile)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) k <- Kriging(y, X, kernel = "matern3_2", objective="LMP") print(k) outfile = tempfile("k.json") save(k,outfile) unlink(outfile)
Save an MLPKriging model to file
## S3 method for class 'MLPKriging' save(object, filename, ...)## S3 method for class 'MLPKriging' save(object, filename, ...)
object |
MLPKriging object |
filename |
path to save file |
... |
ignored |
Save a WarpKriging model to file
## S3 method for class 'WarpKriging' save(object, filename, ...)## S3 method for class 'WarpKriging' save(object, filename, ...)
object |
WarpKriging object |
filename |
path to save file |
... |
ignored |
Get input scaling vector
scaleX(object, ...)scaleX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get input scaling vector for an MLPKriging model
## S3 method for class 'MLPKriging' scaleX(object, ...)## S3 method for class 'MLPKriging' scaleX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get input scaling vector for a WarpKriging model
## S3 method for class 'WarpKriging' scaleX(object, ...)## S3 method for class 'WarpKriging' scaleX(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output scaling value
scaleY(object, ...)scaleY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output scaling value for an MLPKriging model
## S3 method for class 'MLPKriging' scaleY(object, ...)## S3 method for class 'MLPKriging' scaleY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get output scaling value for a WarpKriging model
## S3 method for class 'WarpKriging' scaleY(object, ...)## S3 method for class 'WarpKriging' scaleY(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get process variance
sigma2(object, ...)sigma2(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Numeric process variance estimate.
Get process variance (concentrated MLE)
## S3 method for class 'WarpKriging' sigma2(object, ...)## S3 method for class 'WarpKriging' sigma2(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
KM ObjectThe simulate method is used to simulate paths from the
kriging model described in object.
## S4 method for signature 'KM' simulate( object, nsim = 1, seed = NULL, newdata, cond = TRUE, nugget.sim = 0, checkNames = FALSE, ... )## S4 method for signature 'KM' simulate( object, nsim = 1, seed = NULL, newdata, cond = TRUE, nugget.sim = 0, checkNames = FALSE, ... )
object |
A |
nsim |
Integer: number of response vectors to simulate. |
seed |
Random seed. |
newdata |
Numeric matrix with it rows giving the points where the simulation is to be performed. |
cond |
Logical telling wether the simulation is conditional
or not. Only |
nugget.sim |
Numeric. A postive nugget effect used to avoid numerical instability. |
checkNames |
Check consistency between the design data
|
... |
Ignored. |
Without a dedicated simulate method for the class
"KM", this method would have been inherited from the
"km" class. The dedicated method is expected to run faster.
A comparison can be made by coercing a KM object to a
km object with as.km before calling
simulate.
A numeric matrix with nrow(newdata) rows and
nsim columns containing as its columns the simulated
paths at the input points given in newdata.
XXX method simulate KM
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(5)) y <- f(X) points(X, y, col = 'blue') k <- KM(design = X, response = y, covtype = "gauss") x <- seq(from = 0, to = 1, length.out = 101) s_x <- simulate(k, nsim = 3, newdata = x) lines(x, s_x[ , 1], col = 'blue') lines(x, s_x[ , 2], col = 'blue') lines(x, s_x[ , 3], col = 'blue')f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(5)) y <- f(X) points(X, y, col = 'blue') k <- KM(design = X, response = y, covtype = "gauss") x <- seq(from = 0, to = 1, length.out = 101) s_x <- simulate(k, nsim = 3, newdata = x) lines(x, s_x[ , 1], col = 'blue') lines(x, s_x[ , 2], col = 'blue') lines(x, s_x[ , 3], col = 'blue')
Kriging model object.This method draws paths of the stochastic process at new input points conditional on the values at the input points used in the fit.
## S3 method for class 'Kriging' simulate( object, nsim = 1, seed = 123, x, with_noise = NULL, will_update = FALSE, ... )## S3 method for class 'Kriging' simulate( object, nsim = 1, seed = 123, x, with_noise = NULL, will_update = FALSE, ... )
object |
S3 Kriging object. |
nsim |
Number of simulations to perform. |
seed |
Random seed used. |
x |
Points in model input space where to simulate. |
with_noise |
Logical or numeric specification controlling whether
observation noise is included in the simulated paths. Use
|
will_update |
Set to |
... |
Ignored. |
a matrix with nrow(x) rows and nsim
columns containing the simulated paths at the inputs points
given in x.
The names of the formal arguments differ from those of the
simulate methods for the S4 classes "km" and
"KM". The formal x corresponds to
newdata. These names are chosen Python and
Octave interfaces to libKriging.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, kernel = "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) s <- simulate(k, nsim = 3, x = x) lines(x, s[ , 1], col = "blue") lines(x, s[ , 2], col = "blue") lines(x, s[ , 3], col = "blue")f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, kernel = "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) s <- simulate(k, nsim = 3, x = x) lines(x, s[ , 1], col = "blue") lines(x, s[ , 2], col = "blue") lines(x, s[ , 3], col = "blue")
Simulate from an MLPKriging model
## S3 method for class 'MLPKriging' simulate(object, nsim = 1, seed = 123, x, will_update = FALSE, ...)## S3 method for class 'MLPKriging' simulate(object, nsim = 1, seed = 123, x, will_update = FALSE, ...)
object |
MLPKriging object |
nsim |
number of simulations |
seed |
random seed |
x |
simulation matrix (m x d) |
will_update |
logical; if |
... |
ignored |
matrix (m x nsim)
Simulate from a WarpKriging model
## S3 method for class 'WarpKriging' simulate(object, nsim = 1, seed = 123, x, will_update = FALSE, ...)## S3 method for class 'WarpKriging' simulate(object, nsim = 1, seed = 123, x, will_update = FALSE, ...)
object |
WarpKriging object |
nsim |
number of simulations |
seed |
random seed |
x |
simulation matrix (m x d) |
will_update |
logical; if TRUE, cache data for update_simulate |
... |
ignored |
matrix (m x nsim)
Get Cholesky factor T
T_(object, ...)T_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get Cholesky factor T for an MLPKriging model
## S3 method for class 'MLPKriging' T_(object, ...)## S3 method for class 'MLPKriging' T_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get Cholesky factor T for a WarpKriging model
## S3 method for class 'WarpKriging' T_(object, ...)## S3 method for class 'WarpKriging' T_(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get GP range parameters
theta(object, ...)theta(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Numeric vector of range parameters.
Get GP range parameters
## S3 method for class 'WarpKriging' theta(object, ...)## S3 method for class 'WarpKriging' theta(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Update previous simulate of a model given in
object.
update_simulate(object, ...)update_simulate(object, ...)
object |
An object representing a fitted model. |
... |
Further arguments of function |
Updated simulation of model output.
Kriging model object.This method draws paths of the stochastic process conditional on the values at the input points used in the fit, plus the new input points and their values given as argument (knonw as 'update' points).
## S3 method for class 'Kriging' update_simulate(object, y_u, ..., X_u = NULL, noise_u = NULL)## S3 method for class 'Kriging' update_simulate(object, y_u, ..., X_u = NULL, noise_u = NULL)
object |
S3 Kriging object. |
y_u |
Numeric vector of new responses (output). |
... |
Ignored. |
X_u |
Numeric matrix of new input points. |
noise_u |
Optional numeric vector of observation noise variances
attached to |
a matrix with nrow(x) rows and nsim
columns containing the simulated paths at the inputs points
given in x.
Yann Richet [email protected]
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, kernel = "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) s <- k$simulate(nsim = 3, x = x, will_update = TRUE) lines(x, s[ , 1], col = "blue") lines(x, s[ , 2], col = "blue") lines(x, s[ , 3], col = "blue") X_u <- as.matrix(runif(3)) y_u <- f(X_u) points(X_u, y_u, col = "red") su <- k$update_simulate(y_u, X_u) lines(x, su[ , 1], col = "blue", lty=2) lines(x, su[ , 2], col = "blue", lty=2) lines(x, su[ , 3], col = "blue", lty=2)f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, kernel = "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) s <- k$simulate(nsim = 3, x = x, will_update = TRUE) lines(x, s[ , 1], col = "blue") lines(x, s[ , 2], col = "blue") lines(x, s[ , 3], col = "blue") X_u <- as.matrix(runif(3)) y_u <- f(X_u) points(X_u, y_u, col = "red") su <- k$update_simulate(y_u, X_u) lines(x, su[ , 1], col = "blue", lty=2) lines(x, su[ , 2], col = "blue", lty=2) lines(x, su[ , 3], col = "blue", lty=2)
Update simulated paths with new observations (FOXY algorithm)
## S3 method for class 'MLPKriging' update_simulate(object, y_u, X_u, ...)## S3 method for class 'MLPKriging' update_simulate(object, y_u, X_u, ...)
object |
MLPKriging object (must have called simulate with will_update=TRUE) |
y_u |
new observations |
X_u |
new input matrix |
... |
ignored |
matrix (m x nsim) of updated simulated paths
Update simulated paths with new observations (FOXY algorithm)
## S3 method for class 'WarpKriging' update_simulate(object, y_u, X_u, ...)## S3 method for class 'WarpKriging' update_simulate(object, y_u, X_u, ...)
object |
WarpKriging object (must have called simulate with will_update=TRUE) |
y_u |
new observations |
X_u |
new input matrix |
... |
ignored |
matrix (m x nsim) of updated simulated paths
KM Object with New PointsThe update method is used when new observations are added
to a fitted kriging model. Rather than fitting the model from
scratch with the updated observations added, the results of the
fit as stored in object are used to achieve some savings.
## S4 method for signature 'KM' update( object, newX, newy, newX.alreadyExist = FALSE, cov.reestim = TRUE, trend.reestim = cov.reestim, nugget.reestim = FALSE, newnoise.var = NULL, kmcontrol = NULL, newF = NULL, ... )## S4 method for signature 'KM' update( object, newX, newy, newX.alreadyExist = FALSE, cov.reestim = TRUE, trend.reestim = cov.reestim, nugget.reestim = FALSE, newnoise.var = NULL, kmcontrol = NULL, newF = NULL, ... )
object |
A KM object. |
newX |
A numeric matrix containing the new design points. It
must have |
newy |
A numeric vector of new response values, in
correspondence with the rows of |
newX.alreadyExist |
Logical. If TRUE, |
cov.reestim |
Logical. If |
trend.reestim |
Logical. If |
nugget.reestim |
Logical. If |
newnoise.var |
Optional variance of an additional noise on the new response. |
kmcontrol |
A list of options to tune the fit. Not available yet. |
newF |
New trend matrix. XXXY? |
... |
Ignored. |
Without a dedicated update method for the class
"KM", this would have been inherited from the class
"km". The dedicated method is expected to run faster. A
comparison can be made by coercing a KM object to a
km object with as.km before calling
update.
The updated KM object.
Yann Richet [email protected]
as.km to coerce a KM object to the
class "km".
f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(5)) y <- f(X) points(X, y, col = "blue") KMobj <- KM(design = X, response = y,covtype = "gauss") x <- seq(from = 0, to = 1, length.out = 101) p_x <- predict(KMobj, x) lines(x, p_x$mean, col = "blue") lines(x, p_x$lower95, col = "blue") lines(x, p_x$upper95, col = "blue") newX <- as.matrix(runif(3)) newy <- f(newX) points(newX, newy, col = "red") ## replace the object by its udated version KMobj <- update(KMobj, newX=newX, newy=newy) x <- seq(from = 0, to = 1, length.out = 101) p2_x <- predict(KMobj, x) lines(x, p2_x$mean, col = "red") lines(x, p2_x$lower95, col = "red") lines(x, p2_x$upper95, col = "red")f <- function(x) 1 - 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x) * x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(5)) y <- f(X) points(X, y, col = "blue") KMobj <- KM(design = X, response = y,covtype = "gauss") x <- seq(from = 0, to = 1, length.out = 101) p_x <- predict(KMobj, x) lines(x, p_x$mean, col = "blue") lines(x, p_x$lower95, col = "blue") lines(x, p_x$upper95, col = "blue") newX <- as.matrix(runif(3)) newy <- f(newX) points(newX, newy, col = "red") ## replace the object by its udated version KMobj <- update(KMobj, newX=newX, newy=newy) x <- seq(from = 0, to = 1, length.out = 101) p2_x <- predict(KMobj, x) lines(x, p2_x$mean, col = "red") lines(x, p2_x$lower95, col = "red") lines(x, p2_x$upper95, col = "red")
Kriging model object with new pointsUpdate a Kriging model object with new points
## S3 method for class 'Kriging' update(object, y_u, ..., X_u = NULL, noise_u = NULL, refit = TRUE)## S3 method for class 'Kriging' update(object, y_u, ..., X_u = NULL, noise_u = NULL, refit = TRUE)
object |
S3 Kriging object. |
y_u |
Numeric vector of new responses (output). |
... |
Ignored. |
X_u |
Numeric matrix of new input points. |
noise_u |
Optional numeric vector of observation noise variances
attached to |
refit |
Logical. If |
No return value. Kriging object argument is modified.
The method does not return the updated
object, but instead changes the content of
object. This behaviour is quite unusual in R and
differs from the behaviour of update.km
in DiceKriging and the update method for class
KM.
Yann Richet [email protected]
f <- function(x) 1- 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x)*x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) p <- predict(k, x) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) X_u <- as.matrix(runif(3)) y_u <- f(X_u) points(X_u, y_u, col = "red") ## change the content of the object 'k' update(k, y_u, X_u) ## include design points to see interpolation x <- sort(c(X,X_u,seq(from = 0, to = 1, length.out = 101))) p2 <- predict(k, x) lines(x, p2$mean, col = "red") polygon(c(x, rev(x)), c(p2$mean - 2 * p2$stdev, rev(p2$mean + 2 * p2$stdev)), border = NA, col = rgb(1, 0, 0, 0.2))f <- function(x) 1- 1 / 2 * (sin(12 * x) / (1 + x) + 2 * cos(7 * x)*x^5 + 0.7) plot(f) set.seed(123) X <- as.matrix(runif(10)) y <- f(X) points(X, y, col = "blue") k <- Kriging(y, X, "matern3_2") x <- seq(from = 0, to = 1, length.out = 101) p <- predict(k, x) lines(x, p$mean, col = "blue") polygon(c(x, rev(x)), c(p$mean - 2 * p$stdev, rev(p$mean + 2 * p$stdev)), border = NA, col = rgb(0, 0, 1, 0.2)) X_u <- as.matrix(runif(3)) y_u <- f(X_u) points(X_u, y_u, col = "red") ## change the content of the object 'k' update(k, y_u, X_u) ## include design points to see interpolation x <- sort(c(X,X_u,seq(from = 0, to = 1, length.out = 101))) p2 <- predict(k, x) lines(x, p2$mean, col = "red") polygon(c(x, rev(x)), c(p2$mean - 2 * p2$stdev, rev(p2$mean + 2 * p2$stdev)), border = NA, col = rgb(1, 0, 0, 0.2))
Update an MLPKriging model with new observations
## S3 method for class 'MLPKriging' update(object, y_u, X_u, refit = TRUE, ...)## S3 method for class 'MLPKriging' update(object, y_u, X_u, refit = TRUE, ...)
object |
MLPKriging object |
y_u |
new observations |
X_u |
new input matrix |
refit |
Logical. If |
... |
ignored |
Update a WarpKriging model with new observations
## S3 method for class 'WarpKriging' update(object, y_u, X_u, refit = TRUE, ...)## S3 method for class 'WarpKriging' update(object, y_u, X_u, refit = TRUE, ...)
object |
WarpKriging object |
y_u |
new observations |
X_u |
new input matrix |
refit |
logical; if TRUE (default), re-optimise hyperparameters |
... |
ignored |
Affine warping: w(x) = a*x + b
warp_affine()warp_affine()
string "affine"
Categorical embedding
warp_categorical(n_levels, embed_dim = 2)warp_categorical(n_levels, embed_dim = 2)
n_levels |
number of levels (integer-coded 0..n_levels-1) |
embed_dim |
embedding dimensionality (default 2) |
string e.g. "categorical(5,2)"
Implements the non-stationary input warping of Xiong, Chen, Apley & Ding (2007,
Int. J. Numer. Meth. Engng). The domain is split into
intervals by interior knots, each carrying a learnable
positive slope . The warp is:
for , giving a continuous, monotone
piecewise-linear function with free parameters.
This is the same construction as the knots argument in
DiceKriging.
warp_knots(n_knots = 3, knot_positions = NULL)warp_knots(n_knots = 3, knot_positions = NULL)
n_knots |
number of interior knots |
knot_positions |
optional numeric vector of |
warp specification string, e.g. "knots(3)" or
"knots(0.25:0.5:0.75)"
Xiong, Y., Chen, W., Apley, D. & Ding, X. (2007). A non-stationary covariance-based Kriging method for metamodelling in engineering design. International Journal for Numerical Methods in Engineering, 71(6), 733–756.
Kumaraswamy CDF warping on [0,1]
warp_kumaraswamy()warp_kumaraswamy()
string "kumaraswamy"
Per-variable MLP warping (unconstrained, multi-dim output)
warp_mlp(hidden_dims, d_out = 2, activation = "selu")warp_mlp(hidden_dims, d_out = 2, activation = "selu")
|
integer vector of hidden layer sizes |
|
d_out |
output dimensionality (default 2) |
activation |
activation: "relu","selu","tanh","sigmoid","elu" |
string e.g. "mlp(16:8,3,selu)"
Monotone neural network warping
warp_neural_mono(n_hidden = 8)warp_neural_mono(n_hidden = 8)
|
number of hidden units (default 8) |
string e.g. "neural_mono(16)"
No warping (identity)
warp_none()warp_none()
string "none"
Ordinal warping (learned ordered positions)
warp_ordinal(n_levels)warp_ordinal(n_levels)
n_levels |
number of ordered levels (0..n_levels-1) |
string e.g. "ordinal(4)"
Get warping specifications as strings
warping(object, ...)warping(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get warping specification for a WarpKriging model
## S3 method for class 'WarpKriging' warping(object, ...)## S3 method for class 'WarpKriging' warping(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Kriging with per-variable input warping. Each input dimension is independently transformed before the GP kernel is evaluated. Supports continuous, categorical, ordinal variables and joint deep kernel learning.
WarpKriging( y, X, warping, kernel = "gauss", regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, noise = NULL )WarpKriging( y, X, warping, kernel = "gauss", regmodel = "constant", normalize = FALSE, optim = "BFGS+Adam", objective = "LL", parameters = NULL, noise = NULL )
y |
numeric vector of observations (n) |
X |
numeric matrix of inputs (n x d) |
warping |
character vector of warp specifications, one per column
of X. Use |
kernel |
covariance kernel: "gauss", "matern3_2", "matern5_2", "exp" |
regmodel |
trend: "constant", "linear", "quadratic" |
normalize |
logical; normalise continuous inputs? |
optim |
optimiser (currently only one bi-level strategy) |
objective |
"LL" (log-likelihood) |
parameters |
optional named list of tuning parameters,
e.g. |
noise |
Either a numeric vector of per-observation noise variances,
or |
An S3 object of class "WarpKriging".
# Continuous with Kumaraswamy warping X <- as.matrix(c(0.0, 0.25, 0.5, 0.75, 1.0)) f <- function(x) 1 - 1/2 * (sin(12*x)/(1+x) + 2*cos(7*x)*x^5 + 0.7) y <- f(X) k <- WarpKriging(y, X, warping = "kumaraswamy", kernel = "gauss") print(k) # Mixed: 1 continuous + 1 categorical (3 levels) n <- 15 X_mix <- cbind(runif(n), rep(0:2, length.out = n)) y_mix <- sin(2 * pi * X_mix[, 1]) * c(1, 2, 0.5)[X_mix[, 2] + 1] k2 <- WarpKriging(y_mix, X_mix, warping = c("mlp(16:8,2,selu)", "categorical(3,2)"), kernel = "matern5_2")# Continuous with Kumaraswamy warping X <- as.matrix(c(0.0, 0.25, 0.5, 0.75, 1.0)) f <- function(x) 1 - 1/2 * (sin(12*x)/(1+x) + 2*cos(7*x)*x^5 + 0.7) y <- f(X) k <- WarpKriging(y, X, warping = "kumaraswamy", kernel = "gauss") print(k) # Mixed: 1 continuous + 1 categorical (3 levels) n <- 15 X_mix <- cbind(runif(n), rep(0:2, length.out = n)) y_mix <- sin(2 * pi * X_mix[, 1]) * c(1, 2, 0.5)[X_mix[, 2] + 1] k2 <- WarpKriging(y_mix, X_mix, warping = c("mlp(16:8,2,selu)", "categorical(3,2)"), kernel = "matern5_2")
Get training input matrix
X(object, ...)X(object, ...)
object |
A fitted model object |
... |
ignored |
matrix of training inputs
Get training input matrix
## S3 method for class 'MLPKriging' X(object, ...)## S3 method for class 'MLPKriging' X(object, ...)
object |
MLPKriging object |
... |
ignored |
matrix of training inputs
Get training output vector
y(object, ...)y(object, ...)
object |
A fitted model object |
... |
ignored |
vector of training outputs
Get training output vector
## S3 method for class 'MLPKriging' y(object, ...)## S3 method for class 'MLPKriging' y(object, ...)
object |
MLPKriging object |
... |
ignored |
vector of training outputs
Get whitened residuals z
z(object, ...)z(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get whitened residuals z for an MLPKriging model
## S3 method for class 'MLPKriging' z(object, ...)## S3 method for class 'MLPKriging' z(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |
Get whitened residuals z for a WarpKriging model
## S3 method for class 'WarpKriging' z(object, ...)## S3 method for class 'WarpKriging' z(object, ...)
object |
A Kriging/MLPKriging/WarpKriging model object. |
... |
Unused. |