RobustCalibration-package {RobustCalibration} | R Documentation |
Robust Calibration of Imperfect Mathematical Models
Description
Implements full Bayesian analysis for calibrating mathematical models with new methodology for modeling the discrepancy function. It allows for emulation, calibration and prediction using complex mathematical model outputs and experimental data. See the reference: Mengyang Gu and Long Wang, 2018, Journal of Uncertainty Quantification; Mengyang Gu, Fangzheng Xie and Long Wang, 2022, Journal of Uncertainty Quantification; Mengyang Gu, Kyle Anderson and Erika McPhillips, 2023, Technometrics.
Details
The DESCRIPTION file:
Package: | RobustCalibration |
Type: | Package |
Title: | Robust Calibration of Imperfect Mathematical Models |
Version: | 0.5.5 |
Date: | 2024-05-29 |
Authors@R: | c(person(given="Mengyang",family="Gu",role=c("aut","cre"),email="mengyang@pstat.ucsb.edu")) |
Maintainer: | Mengyang Gu <mengyang@pstat.ucsb.edu> |
Author: | Mengyang Gu [aut, cre] |
Description: | Implements full Bayesian analysis for calibrating mathematical models with new methodology for modeling the discrepancy function. It allows for emulation, calibration and prediction using complex mathematical model outputs and experimental data. See the reference: Mengyang Gu and Long Wang, 2018, Journal of Uncertainty Quantification; Mengyang Gu, Fangzheng Xie and Long Wang, 2022, Journal of Uncertainty Quantification; Mengyang Gu, Kyle Anderson and Erika McPhillips, 2023, Technometrics. |
License: | GPL (>= 2) |
Depends: | methods |
Imports: | Rcpp (>= 0.12.3), RobustGaSP (>= 0.6.4), nloptr (>= 1.0.4) |
LinkingTo: | Rcpp, RcppEigen |
NeedsCompilation: | yes |
Repository: | CRAN |
Packaged: | 2018-10-07 20:06:35 UTC; gumengyang |
RoxygenNote: | 5.0.1 |
Date/Publication: | 2018-05-14 04:26:37 UTC |
Index of help topics:
RobustCalibration-package Robust Calibration of Imperfect Mathematical Models predict Prediction for the robust calibration model predict_MS Prediction for the robust calibration model for multiple sources predictobj.rcalibration-class Predictive results for the Robust Calibration class predictobj.rcalibration_MS-class Predictive results for the Robust Calibration class rcalibration Setting up the robust Calibration model rcalibration-class Robust Calibration class rcalibration_MS Setting up the robust Calibration model for multiple sources data rcalibration_MS-class Robust Calibration for multiple sources class show Show an Robust Calibration object.
Robust calibration of imperfect mathematical models and prediction using experimental data
Author(s)
Mengyang Gu [aut, cre]
Maintainer: Mengyang Gu <mengyang@pstat.ucsb.edu>
References
A. O'Hagan and M. C. Kennedy (2001), Bayesian calibration of computer models, Journal of the Royal Statistical Society: Series B (Statistical Methodology, 63, 425-464.
Bayarri, Maria J and Berger, James O and Paulo, Rui and Sacks, Jerry and Cafeo, John A and Cavendish, James and Lin, Chin-Hsu and Tu, Jian (2007) A framework for validation of computer models. Technometrics. 49, 138–154.
M. Gu (2016), Robust Uncertainty Quantification and Scalable Computation for Computer Models with Massive Output, Ph.D. thesis., Duke University.
M. Gu and L. Wang (2017) Scaled Gaussian Stochastic Process for Computer Model Calibration and Prediction. arXiv preprint arXiv:1707.08215.
M. Gu (2018) Jointly Robust Prior for Gaussian Stochastic Process in Emulation, Calibration and Variable Selection . arXiv preprint arXiv:1804.09329.
See Also
Examples
##---------------------------------------------------
##A simple example where the math model is not biased
##---------------------------------------------------
## the reality
test_funct_eg1<-function(x){
sin(pi/2*x)
}
## obtain 25 data from the reality plus a noise
set.seed(1)
## 10 data points are very small, one may want to add more data
n=15
input=seq(0,4,4/(n-1))
input=as.matrix(input)
output=test_funct_eg1(input)+rnorm(length(input),mean=0,sd=0.2)
## plot input and output
#plot(input,output)
#num_obs=n=length(output)
## the math model
math_model_eg1<-function(x,theta){
sin(theta*x)
}
##fit the S-GaSP model for the discrepancy
##one can choose the discrepancy_type to GaSP, S-GaSP or no discrepancy
##p_theta is the number of parameters to calibrate and user needs to specifiy
##one may also want to change the number of posterior samples by change S and S_0
##one may change sd_proposal for the standard derivation of the proposal distribution
## one may also add a mean by setting X=... and have_trend=TRUE
model_sgasp=rcalibration(design=input, observations=output, p_theta=1,simul_type=1,
math_model=math_model_eg1,theta_range=matrix(c(0,3),1,2)
,S=10000,S_0=2000,discrepancy_type='S-GaSP')
##posterior samples of calibration parameter and value
## the value is
plot(model_sgasp@post_sample[,1],type='l',xlab='num',ylab=expression(theta))
plot(model_sgasp@post_value,type='l',xlab='num',ylab='posterior value')
show(model_sgasp)
#------------------------------------------------------------------------------
# Example: an example used in Susie Bayarri et. al. 2007 Technometrics paper
#------------------------------------------------------------------------------
##reality
test_funct_eg1<-function(x){
3.5*exp(-1.7*x)+1.5
}
##math model
math_model_eg1<-function(x,theta){
5*exp(-x*theta)
}
## noise observations (sampled from reality + independent Gaussian noises)
## each has 3 replicates
input=c(rep(.110,3),rep(.432,3),rep(.754,3),rep(1.077,3),rep(1.399,3),rep(1.721,3),
rep(2.043,3),rep(2.366,3),rep(2.688,3),rep(3.010,3))
output=c(4.730,4.720,4.234,3.177,2.966,3.653,1.970,2.267,2.084,2.079,2.409,2.371,1.908,1.665,1.685,
1.773,1.603,1.922,1.370,1.661,1.757,1.868,1.505,1.638,1.390,1.275,1.679,1.461,1.157,1.530)
n_stack=length(output)/3
output_stack=rep(0,n_stack)
input_stack=rep(0,n_stack)
for(j in 1:n_stack){
output_stack[j]=mean(output[ ((j-1)*3+1):(3*j)])
input_stack[j]=mean(input[ ((j-1)*3+1):(3*j)])
}
output_stack=as.matrix(output_stack)
input_stack=as.matrix(input_stack)
## plot the output and stack
#plot(input,output,pch=16,col='red')
#lines(input_stack,output_stack,pch=16,col='blue',type='p')
## fit the model with S-GaSP for the discrepancy
model_sgasp=rcalibration(design=input_stack, observations=output_stack, p_theta=1,simul_type=1,
math_model=math_model_eg1,theta_range=matrix(c(0,10),1,2),S=10000,
S_0=2000,discrepancy_type='S-GaSP')
#posterior
plot(model_sgasp@post_sample[,1],type='l',xlab='num',ylab=expression(theta))
show(model_sgasp)