Description Usage Arguments Details Author(s) References See Also Examples
Perform O'Hagan's three stage optimization on the toy dataset. Function
stage1()
and stage2()
find the optimal values for
the hyperparameters and stage3()
finds the optimal values for
the three parameters.
1 2 3 4 5 6 7 8 9 10 11 12  stage1(D1, y, H1, maxit, trace=100, method="NelderMead",
directory = ".", do.filewrite=FALSE, do.print=TRUE,
phi.fun, lognormally.distributed=FALSE, include.prior=TRUE, phi)
stage2(D1, D2, H1, H2, y, z, maxit, trace=100, method = "NelderMead",
directory = ".", do.filewrite=FALSE, do.print=TRUE, extractor,
phi.fun, E.theta, Edash.theta, isotropic=FALSE,
lognormally.distributed = FALSE, include.prior = TRUE,
use.standin = FALSE, rho.eq.1 = TRUE, phi)
stage3(D1, D2, H1, H2, d, maxit, trace=100, method="NelderMead",
directory = ".", do.filewrite=FALSE, do.print=TRUE,
include.prior = TRUE, lognormally.distributed=FALSE,
theta.start=NULL, phi)

maxit 
Maximum number of iterations as passed to 
trace 
Amount of information displayed, as passed to 
D1 
Matrix whose rows are points at which code output is known 
D2 
Matrix whose rows are points at which observations were made 
H1,H2 
Regressor basis functions for 
y 
Code outputs. Toy example is 
z 
Observations. Toy example is 
d 
Data vector consisting of the code runs and observations 
extractor 
extractor function for 
E.theta,Edash.theta 
Expectation WRT theta, and dashed theta.
Toy examples are 
phi.fun 
Function to create hyperparameters; passed (in

method 
Method argument passed to 
include.prior 
Boolean variable with default 
lognormally.distributed 
Boolean with 
do.filewrite 
Boolean, with 
directory 
The directory to write files to; only matters if

isotropic 
In function 
do.print 
Boolean, with default 
use.standin 
In 
rho.eq.1 
Boolean, with default 
theta.start 
In 
phi 
Hyperparameters. Used as initial values for the hyperparameters in the optimization routines 
The three functions documented here carry out the multistage
optimization detailed in KOH2001 (actually, KOH2001 only defined stage
1 and stage 2, which estimated the hyperparameters. What is here
called “stage3()
” estimates the true value of
theta given the hyperparameters).
stage1()
carries out stage 1 of KOH2001 which is used to
estimate psi1 using optimization.
In function stage2()
, setting argument isotropic
to
TRUE
will force phi$omegastar_x
to be a function of a
length one scalar. The value of phi$omegastar_x
used will
depend on pdm.maker.psi2()
(an internal function appearing in
hpa.fun.toy()
). In stage2()
, several kludges are made.
The initial conditions are provided by argument phi
. The
relevant part of this is phi$psi2
.
Function stage2()
estimates psi2 and
rho and lambda, using
optimization. Note that psi2 includes
sigma2squared in addition to omegastar_X
(in
the toy case, psi2 has three elements: the first two are
the diagonal of omegastar_x
and the third is
sigma2squared but this information is
encoded in phi.fun.toy()
, which changes from application to
application).
Function stage3()
attempts to find the maximum likelihood
estimate of theta, given hyperparameters and
observations, using optimization
Robin K. S. Hankin
M. C. Kennedy and A. O'Hagan 2001. Bayesian calibration of computer models. Journal of the Royal Statistical Society B, 63(3) pp425464
M. C. Kennedy and A. O'Hagan 2001. Supplementary details on Bayesian calibration of computer models, Internal report, University of Sheffield. Available at http://www.tonyohagan.co.uk/academic/ps/calsup.ps
R. K. S. Hankin 2005. Introducing BACCO, an R bundle for Bayesian analysis of computer code output, Journal of Statistical Software, 14(16)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26  data(toys)
stage1(D1=D1.toy,y=y.toy,H1=H1.toy, maxit=5, phi.fun=phi.fun.toy, phi=phi.toy)
##now try with a slightly bigger dataset:
##Examples below take a few minutes to run:
set.seed(0)
data(toys)
jj < create.new.toy.datasets(D1.toy , D2.toy)
y.toy < jj$y.toy
z.toy < jj$z.toy
d.toy < jj$d.toy
phi.toy.stage1 < stage1(D1=D1.toy, y=y.toy, H1=H1.toy, maxit=10, phi.fun=phi.fun.toy, phi=phi.toy)
phi.toy.stage2 < stage2(D1=D1.toy, D2=D2.toy, H1=H1.toy, H2=H2.toy,
y=y.toy, z=z.toy, extractor=extractor.toy,
phi.fun=phi.fun.toy, E.theta=E.theta.toy, Edash.theta=Edash.theta.toy,
maxit=3, phi=phi.toy.stage1)
stage3(D1=D1.toy, D2=D2.toy, H1=H1.toy, H2=H2.toy, d=d.toy, maxit=3, phi=phi.toy.stage2)
# Now try with the true values of the hyperparameters:
phi.true < phi.true.toy(phi=phi.toy)
stage3(D1=D1.toy, D2=D2.toy, H1=H1.toy, H2=H2.toy, d=d.toy, maxit=3, phi=phi.true)

Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.