mleHomTP | R Documentation |
Student-t process regression under homoskedastic noise based on maximum likelihood estimation of the hyperparameters. This function is enhanced to deal with replicated observations.
mleHomTP(
X,
Z,
lower = NULL,
upper = NULL,
known = list(beta0 = 0),
noiseControl = list(g_bounds = c(sqrt(.Machine$double.eps), 10000), nu_bounds = c(2 +
0.001, 30), sigma2_bounds = c(sqrt(.Machine$double.eps), 10000)),
init = list(nu = 3),
covtype = c("Gaussian", "Matern5_2", "Matern3_2"),
maxit = 100,
eps = sqrt(.Machine$double.eps),
settings = list(return.Ki = TRUE, factr = 1e+09)
)
X |
matrix of all designs, one per row, or list with elements:
|
Z |
vector of all observations. If using a list with |
lower , upper |
bounds for the |
known |
optional list of known parameters, e.g., |
noiseControl |
list with element,
|
init |
list specifying starting values for MLE optimization, with elements:
|
covtype |
covariance kernel type, either 'Gaussian', 'Matern5_2' or 'Matern3_2', see |
maxit |
maximum number of iteration for L-BFGS-B of |
eps |
jitter used in the inversion of the covariance matrix for numerical stability |
settings |
list with argument |
The global covariance matrix of the model is parameterized as K = sigma2 * C + g * diag(1/mult)
,
with C
the correlation matrix between unique designs, depending on the family of kernel used (see cov_gen
for available choices).
It is generally recommended to use find_reps
to pre-process the data, to rescale the inputs to the unit cube and to normalize the outputs.
a list which is given the S3 class "homGP
", with elements:
theta
: maximum likelihood estimate of the lengthscale parameter(s),
g
: maximum likelihood estimate of the nugget variance,
trendtype
: either "SK
" if beta0
is given, else "OK
"
beta0
: estimated trend unless given in input,
sigma2
: maximum likelihood estimate of the scale variance,
nu2
: maximum likelihood estimate of the degrees of freedom parameter,
ll
: log-likelihood value,
X0
, Z0
, Z
, mult
, eps
, covtype
: values given in input,
call
: user call of the function
used_args
: list with arguments provided in the call
nit_opt
, msg
: counts
and msg
returned by optim
Ki
, inverse covariance matrix (if return.Ki
is TRUE
in settings
)
time
: time to train the model, in seconds.
M. Binois, Robert B. Gramacy, M. Ludkovski (2018), Practical heteroskedastic Gaussian process modeling for large simulation experiments,
Journal of Computational and Graphical Statistics, 27(4), 808–821.
Preprint available on arXiv:1611.05902.
A. Shah, A. Wilson, Z. Ghahramani (2014), Student-t processes as alternatives to Gaussian processes, Artificial Intelligence and Statistics, 877–885.
M. Chung, M. Binois, RB Gramacy, DJ Moquin, AP Smith, AM Smith (2019).
Parameter and Uncertainty Estimation for Dynamical Systems Using Surrogate Stochastic Processes.
SIAM Journal on Scientific Computing, 41(4), 2212-2238.
Preprint available on arXiv:1802.00852.
predict.homTP
for predictions.
summary
and plot
functions are available as well.
##------------------------------------------------------------
## Example 1: Homoskedastic Student-t modeling on the motorcycle data
##------------------------------------------------------------
set.seed(32)
## motorcycle data
library(MASS)
X <- matrix(mcycle$times, ncol = 1)
Z <- mcycle$accel
plot(X, Z, ylim = c(-160, 90), ylab = 'acceleration', xlab = "time")
noiseControl = list(g_bounds = c(1e-3, 1e4))
model <- mleHomTP(X = X, Z = Z, lower = 0.01, upper = 100, noiseControl = noiseControl)
summary(model)
## Display averaged observations
points(model$X0, model$Z0, pch = 20)
xgrid <- matrix(seq(0, 60, length.out = 301), ncol = 1)
preds <- predict(x = xgrid, object = model)
## Display mean prediction
lines(xgrid, preds$mean, col = 'red', lwd = 2)
## Display 95% confidence intervals
lines(xgrid, preds$mean + sqrt(preds$sd2) * qt(0.05, df = model$nu + nrow(X)), col = 2, lty = 2)
lines(xgrid, preds$mean + sqrt(preds$sd2) * qt(0.95, df = model$nu + nrow(X)), col = 2, lty = 2)
## Display 95% prediction intervals
lines(xgrid, preds$mean + sqrt(preds$sd2 + preds$nugs) * qt(0.05, df = model$nu + nrow(X)),
col = 3, lty = 2)
lines(xgrid, preds$mean + sqrt(preds$sd2 + preds$nugs) * qt(0.95, df = model$nu + nrow(X)),
col = 3, lty = 2)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.