ARCokrig {ARCokrig} R Documentation

## Fit the AR-Cokriging model and make predictions

### Description

This is a simple and high-level funciton to fit autoregressive cokriging models to multifidelity computer model outputs.

### Usage

ARCokrig(
formula = list(~1, ~1),
output,
input,
cov.model = "matern_5_2",
nugget.est = FALSE,
input.new,
prior = list(),
opt = list(),
NestDesign = TRUE,
tuning = list(),
info = list()
)


### Arguments

 formula a list of s elements, each of which contains the formula to specify fixed basis functions or regressors. output a list of s elements, each of which contains a matrix of computer model outputs. input a list of s elements, each of which contains a matrix of inputs. cov.model a string indicating the type of covariance function in AR-cokriging models. Current covariance functions include expproduct form of exponential covariance functions. matern_3_2product form of Matern covariance functions with smoothness parameter 3/2. matern_5_2product form of Matern covariance functions with smoothness parameter 5/2. Gaussianproduct form of Gaussian covariance functions. powexpproduct form of power-exponential covariance functions with roughness parameter fixed at 1.9. nugget.est a logical value indicating whether nugget parameter is included or not. Default value is FALSE. input.new a matrix including new inputs for making prediction prior a list of arguments to setup the prior distributions namethe name of the prior. Current implementation includes JR, Reference, Jeffreys, Ind_Jeffreys hyperparamhyperparameters in the priors. For jointly robust (JR) prior, three parameters are included: a refers to the polynomial penalty to avoid singular correlation matrix with a default value 0.2; b refers to the exponenetial penalty to avoid diagonal correlation matrix with a default value 1; nugget.UB is the upper bound of the nugget variance with default value 1, which indicates that the nugget variance has support (0, 1). opt a list of arguments to setup the optim routine. NestDesign a logical value indicating whether the experimental design is hierarchically nested within each level of the code. tuning a list of arguments to control the MCEM algorithm for non-nested design. It includes the arguments maxitthe maximum number of MCEM iterations. tola tolerance to stop the MCEM algorithm. If the parameter difference between any two consecutive MCEM algorithm is less than this tolerance, the MCEM algorithm is stopped. n.samplethe number of Monte Carlo samples in the MCEM algorithm. verbosea logical value to show the MCEM iterations if it is true. info a list that contains iternumber of iterations used in the MCEM algorithm epsparameter difference after the MCEM algorithm stops

### Value

The main call inside ARCokrig consists of cokm, cokm.fit, and cokm.predict. Thus, the function returns the cokm object and predictions over new inputs.

### Author(s)

Pulong Ma <mpulong@gmail.com>

### References

• Ma, P. (2019). “Objective Bayesian Analysis of a Cokriging Model for Hierarchical Multifidelity Codes." arXiv:1910.10225. https://arxiv.org/abs/1910.10225.

• Ma, P., Karagiannis, G., Konomi, B., Asher, T., Toro, G., and Cox, A. (2019) “Multifidelity Computer Model Emulation with High-Dimensional Output: An Application to Storm Surge." arXiv:1909.01836. https://arxiv.org/abs/1909.01836.

cokm, cokm.param, cokm.fit, cokm.predict

### Examples


##############################################################
##############################################################
############ Example
Funcc = function(x){
return(0.5*(6*x-2)^2*sin(12*x-4)+10*(x-0.5)-5)
}

Funcf = function(x){
z1 = Funcc(x)
z2 = 2*z1-20*x+20 + sin(10*cos(5*x))
return(z2)
}

#####################################################################
###### Nested design
#####################################################################
Dc <- seq(-1,1,0.1)
indDf <- c(1, 3, 6, 8, 10, 13, 17, 21)
zc <- Funcc(Dc)
Df <- Dc[indDf]
zf <- Funcf(Df)

input.new = as.matrix(seq(-1,1,length.out=200))

## fit and predict with the AR-Cokriging model

out = ARCokrig(formula=list(~1,~1+x1), output=list(c(zc), c(zf)),
input=list(as.matrix(Dc), as.matrix(Df)),
cov.model="matern_5_2",
input.new=input.new)

## plot results

library(ggplot2)
cokrig = out$cokrig df.l1 = data.frame(x=c(Dc), y=c(zc)) df.l2 = data.frame(x=c(Df), y=c(zf)) CI.lower = cokrig$lower95[[2]]
CI.upper = cokrig$upper95[[2]] df.CI = data.frame(x=c(input.new),lower=CI.lower, upper=CI.upper) df.pred = data.frame(x=c(input.new), y=cokrig$mu[[2]])

g = ggplot(data.frame(x=c(-1,1)), aes(x)) +
stat_function(fun=Funcc, geom="line", aes(colour="level 1"), n=500) +
stat_function(fun=Funcf, geom="line", aes(colour="level 2"), n=500)

g = g + geom_point(data=df.l1, mapping=aes(x=x, y=y), shape=16, size=2, color="black") +
geom_point(data=df.l2, mapping=aes(x=x, y=y), shape=17, size=2, color="black")

g = g + geom_line(data=df.pred, aes(x=x, y=y, colour="cokriging"), inherit.aes=FALSE) +
geom_ribbon(data=df.CI, mapping=aes(x=x,ymin=lower, ymax=upper), fill="gray40",
alpha=0.3, inherit.aes=FALSE)
g = g + scale_colour_manual(name=NULL, values=c("red","blue", "turquoise3"),
breaks=c("cokriging","level 1", "level 2"))

g = g + ggtitle("A Two-Level Example") +
theme(plot.title=element_text(size=14),
axis.title.x=element_text(size=14),
axis.text.x=element_text(size=14),
axis.title.y=element_text(size=14),
axis.text.y=element_text(size=14),
legend.text = element_text(size=12),
legend.direction = "horizontal",
legend.position = c(0.6, 0.1)) + xlab("") + ylab("")
print(g)



[Package ARCokrig version 0.1.2 Index]