numOSL (version 2.0)

mcMAM: Minimum age model optimization (using a Markov chain Monte Carlo method)

Description

Sampling from the joint-likelihood function of the minimum age model using a Markov chain Monte Carlo (MCMC) method .

Usage

mcMAM(EDdata, ncomp = -1, addsigma = 0, iflog = TRUE,  nsim = 50000, inis = list(), control.args = list())

Arguments

EDdata

matrix(required): a two-column matrix (i.e., equivalent dose values and associated standard errors)

ncomp

integer(with default): number of components, -1=MAM3, -2=MAM4

addsigma

double(with default): additional uncertainty

iflog

logical(with default): transform equivalent dose values to log-scale or not

nsim

integer(with default): deseried number of iterations

inis

list(with default): initial state of parameters. Example: inis=list(p=0.1,gamma=20,sigma=0.5) in MAM3

control.args

list(with default): arguments used by the Slice Sampling algorithm, see function mcFMM for details

Value

Return an invisible list of S3 class "mcAgeModels". See mcFMM for details.

References

Galbraith RF, Roberts RG, Laslett GM, Yoshida H, Olley JM, 1999. Optical dating of single grains of quartz from Jinmium rock shelter, northern Australia. Part I: experimental design and statistical models. Archaeometry, 41(2): 339-364.

Neal RM, 2003. "Slice sampling" (with discussion). Annals of Statistics, 31(3): 705-767. Software is freely available at http://www.cs.utoronto.ca/~radford/slice.software.html.

See Also

mcFMM; reportSAM; RadialPlotter

Examples

  # Do not run.
  # data(EDdata)
  # Construct a MCMC chain for MAM3.
  # obj<-mcMAM(EDdata$al3,ncomp=-1,addsigma=0.1,nsim=5000)
  # reportSAM(obj,burn=1e3,thin=2)
  #
  # The convergence of the simulations may be diagnosed with 
  # the Gelman and Rubin's convergence diagnostic.
  # library(coda)   # Only if package "coda" has been installed.
  # args<-list(nstart=50)
  # inis1<-list(p=0.01,gamma=26,mu=104,sigma=0.01)
  # inis2<-list(p=0.99,gamma=100,mu=104,sigma=4.99)
  # obj1<-mcMAM(EDdata$al3,ncomp=-2,nsim=3000,inis=inis1,control.args=args)
  # obj2<-mcMAM(EDdata$al3,ncomp=-2,nsim=3000,inis=inis2,control.args=args)
  # chain1<-mcmc(obj1$chains)
  # chain2<-mcmc(obj2$chains)
  # chains<-mcmc.list(chain1,chain2)
  # gelman.plot(chains)