BclimMCMC(Bclimdata,chron.loc,nchron=10000,control.mcmc=list(iterations=100000,
burnin=20000,thinby=40,report=100),control.chains=list(v.mh.sd=2,phi1.mh.sd=1,
phi2.mh.sd=10,vstart=statmod::rinvgauss(Bclimdata$n-1,2,1),
Zstart=sample(1:Bclimdata$G,Bclimdata$n,replace=TRUE),phi1start=rep(3,Bclimdata$m),
phi2start=rep(20,Bclimdata$m)),control.priors=list(phi1dlmean=rep(1.275,Bclimdata$m),
phi1dlsd=rep(0.076,Bclimdata$m),phi2dlmean=rep(4.231,Bclimdata$m),
phi2dlsd=rep(0.271,Bclimdata$m)))
BclimMixPar
or BclimMixSer
BclimMixPar
or BclimMixSer
, and produces estimated climate and climate volatilities using an MCMC algorithm. The full details are in the Arxiv paper referenced below. The options listed above allow quite a detailed level of control over the behaviour of the algorithm, and convergence should be checked using suitable means (see e.g. the R package boa or coda).BclimRun
. See also the other 3 stages: BclimInterp
, BclimMixSer
(or BclimMixPar
), and BclimLayer
,# Set the working directory using setwd (not shown)
# Download and load in the response surfaces:
url1 <- 'http://mathsci.ucd.ie/~parnell_a/media/requireddata3D.RData'
download.file(url1,'required_data3D.RData')
# and now the pollen
url2 <- 'http://mathsci.ucd.ie/~parnell_a/media/SlugganPollen.txt'
download.file(url2,'SlugganPollen.txt')
# and finally the chronologies
url3 <- 'http://mathsci.ucd.ie/~parnell_a/media/Sluggan_2chrons.txt'
download.file(url3,'Slugganchrons.txt')
# Create variables which state the locations of the pollen and chronologies
pollen.loc <- paste(getwd(),'/SlugganPollen.txt',sep='')
chron.loc <- paste(getwd(),'/Slugganchrons.txt',sep='')
# Load in the response surfaces
load('required_data3D.RData')
## note that all of these functions have further options you can change with
step1 <- BclimLayer(pollen.loc,required.data3D=required.data3D)
step2 <- BclimMixSer(step1)
# See also the parallelised version BclimMixPar if you have doMC and foreach installed
step3 <- BclimMCMC(step2,chron.loc)
# You should probably do some convergence checking after this step
step4 <- BclimInterp(step2,step3)
results <- BclimCompile(step1,step2,step3,step4,core.name="Sluggan Moss")
# Create a plot of MTCO (dim=2)
plotBclim(results,dim=2)
# Create a volatility plot
plotBclimVol(results,dim=2)
Run the code above in your browser using DataLab