tiepointsimmer: Tie-point simulation

In eirikmn/bremla: Bayesian Regression Modeling of Layer-Counted Archives

Tie-point simulation

Description

Wrapper function which calls specified functions to produce samples from tie-point distributions and places them in bremla object.

Usage

tiepointsimmer(object, synchronization, print.progress = FALSE, ...)

Arguments

object	Bremla list object to store samples in
synchronization	List containing specifications for generating tie-points. Must include synchronization\$locations. Includes different methods for producing tie-points specified by synchronization\$method, and can also import pre-computed samples in synchronization\$samples. See synchronization.default for more details.
print.progress	Boolean. If TRUE progress will be printed to screen.
...	Further arguments to be passed down to the selected tie-point simulation function.

Value

returns the same object from input but appends tie-point samples and information in object\$tie_points.

Author(s)

Eirik Myrvoll-Nilsen, [e-mail]

Examples

require(stats)
set.seed(1)
n <- 1000
phi <- 0.8
sigma <- 1.2
a_lintrend <- 0.3; a_proxy = 0.8
dy_noise <- as.numeric(arima.sim(model=list(ar=c(phi)),n=n,sd=sqrt(1-phi^2)))
lintrend <- seq(from=10,to=15,length.out=n)

proxy <- as.numeric(arima.sim(model=list(ar=c(0.9)),n=n,sd=sqrt(1-0.9^2)))
dy <- a_lintrend*lintrend + a_proxy*proxy + sigma*dy_noise

y0 = 11700;z0=1200
age = y0+cumsum(dy)
depth = 1200 + 1:n*0.05
depth2 = depth^2/depth[1]^2 #normalize for stability


formula = dy~-1+depth2 + proxy

data = data.frame(age=age,dy=dy,proxy=proxy,depth=depth,depth2=depth2)
data = rbind(c(y0,NA,NA,z0,NA),data) #First row is only used to extract y0 and z0.

events=list(locations=c(1210,1220,1240))
synchronization=list(locations=depth[c(100,400,700)],method="gauss",
                           params=list(mean=c(age[c(100,400,700)]+c(30,-100,50)),
                                       sd=c(50,20,100)
                                       )
                       )
object = bremla_prepare(formula,data,nsims=5000,
                        reference.label="simulated timescale",
                        events = events,
                        synchronization=synchronization)
object = tiepointsimmer(object)
summary(object)
plot(object)

eirikmn/bremla documentation built on Jan. 25, 2025, 4:41 a.m.