bremla_simulationsummarizer: Simulation-summarizer

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

Simulation-summarizer

Description

Computes posterior marginal mean and uncertainty intervals from simulations.

Usage

bremla_simulationsummarizer(object, sync = TRUE, print.progress = FALSE)

Arguments

object	List object which is the output of function bremla_chronology_simulation
sync	boolean. Set TRUE if summary statistics should be computed for syncrhonized samples (if any). FALSE for unsynchronized.
print.progress	Boolean describing whether or not progress should be printed on screen.

Value

Returns the object list from the input and appends additional summary statistics: posterior marginal mean and uncertainty interval-

Author(s)

Eirik Myrvoll-Nilsen, [e-mail]

See Also

bremla_chronology_simulation

Examples

if(inlaloader()){
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))
control.fit = list(ncores=2,noise="ar1")
control.sim=list(synchronized=2,
                 summary=list(compute=TRUE))

object = bremla_prepare(formula,data,nsims=5000,reference.label="simulated timescale",
                        events = events,
                        control.fit=control.fit,
                        control.sim=control.sim)
object = bremla_modelfitter(object)
object = bremla_chronology_simulation(object)
object = bremla_simulationsummarizer(object,sync=FALSE,print.progress=TRUE)
summary(object)
plot(object)
}

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