GPdiffMesh {diffMeshGP} | R Documentation |
This R function implements the nonstationary Kriging model proposed by Tuo, Wu and Yu (2014) for analyzing multi-fidelity computer outputs. This function computes the maximum likelihood estimates for the model parameters as well as the predictive means and variances of the exact solution (i.e., the conceptually highest fidelity).
GPdiffMesh(x,X,meshT,Y,regFunX = function(x){ return(0*matrix(x[,1]))},
regFunT = function(x){ return(1*matrix(x[,1]))},
phi1 = 1,
sigma12 = 1,
sigma22 = 1,
phi2 = 1,
mybeta = FALSE,
l = 4)
x |
An |
X |
An |
meshT |
An |
Y |
An |
regFunX |
A scalar or a vector-valued regression function |
regFunT |
A scalar or a vector-valued regression function |
phi1,sigma12,sigma22,phi2,mybeta |
Initial values of the parameters for the maximum likelihood estimation. The default values are |
l |
The value of parameter |
This R function implements the nonstationary Kriging model proposed by Tuo, Wu and Yu (2014) for the modeling and analysis of multi-fidelity computer experiments. Denote an input-output pair from the computer simulation by (x, t, y)
, where x
is the vector of input variables, t
is the mesh density, and y
is the corresponding computer output. Tuo, Wu and Yu (2014) use the following Gaussian process model:
y = \eta(x, t) = \eta(x, 0) + \delta(x, t),
where \eta(x, 0)
and \delta(x, t)
are realizations of two mutually independent Gaussian processes V(x)
and Z(x, t)
. Assume
E(V(x))=f_{1}^{T}(x)\beta_1, E(Z(x, t))=f^{T}_{2}(x,t)\beta_2,
where f_{1}(x)
can be set in regFunX, and f_{2}(x,t)
can be set in regFunT. Both f_{1}(x)
and f_{2}(x,t)
can be vector-valued function. f_{2}(x,t)
should satisfy for any x
, \lim_{t\rightarrow 0}f_{2}(x,t) = 0
. The default functions are f_{1}(x) = 1
and f_{2}(x,t) = t
.
Assume the covariance of V(x)
is
Cov(V(x_{1}),V(x_{2}))=\sigma_{1}^{2}e^{-\phi^2_1\|x_{1}-x_{2}\|^{2}_{2}},
and the covariance of Z(x, t)
is
Cov(Z(x_{1}, t_{1}),Z(x_{2}, t_{2}))=\sigma_{2}^{2}\min(t_{1},t_{2})^{l}e^{-\phi_{2}^{2}\|x_{1}-x_{2}\|^{2}_{2}},
where l
is a fixed parameter and is not estimated in this function.
The goal of the Tuo-Wu-Yu model is to predict for the exact solution
\varphi(.):=\eta(.,0)
.
For a set of control points x=(x_{1},\ldots,x_{n_0})
, this R function predicts \eta(x_{i},0)=\varphi(x_{i})
for i=1,\ldots,n_{0}
.
This R function uses maximum likelihood method to estimate the model parameters. Nelder-Mead method is used to maximize the likelihood function. The solution may depend on the choice of initial values. Users may specify the initial values or use the default values.
outy |
An |
sigy |
An |
estipar |
A list of parameter estimates. |
hig02 <- function(s)
{
#The test function is from [2].
y <- s*sin(s) / 10
return(y)
}
myX <- matrix(c(seq(from = 0,to = 10, by = 1),
seq(from = 0,to = 10, by = 1)),ncol = 2)
myy <- hig02(matrix(myX[,1]))
myT <- matrix(c(0.01,0.5,0.01,0.02,0.02,0.01,0.01,0.02,0.002,
0.003,0.03))
myregf <- function(x){
return(x)
}
myregfn <- function(s){
return(cbind((matrix(s[,1])^2*matrix(s[,2])),
(matrix(s[,1])*matrix(s[,2]))))
}
#Here s=cbind(t,x), where x is a matrix of input variables and
# t is the corresponding mesh density matrix.
x <- matrix(c(seq(from = 0,to = 10, by = 0.1),
seq(from = 0,to = 10, by = 0.1)),ncol = 2)
myploty <- hig02(matrix(x[,1]))
y <- GPdiffMesh(x, myX, myT, myy, regFunX = myregf, regFunT = myregfn)
# The regression function is beta_0 + beta_1 x + (t^2x, tx)' beta_2
y$outy
y$sigy
y$estipar
plot(x[,1], myploty,"l")
lines(x[,1],y$outy, type="o", pch=22, lty=2, col="red")