-
Notifications
You must be signed in to change notification settings - Fork 0
/
fit_basecase.R
51 lines (43 loc) · 1.82 KB
/
fit_basecase.R
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
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
## File to run the fits to the real data
chains <- 6
options(mc.cores = chains)
td <- 17
ad <- .85
iter <- 800
warmup <- 300
## Base case for paper: combined
control <- list(model='combined', n_x=400,
n_eps1="IID", n_eps2="IID", n_omega2="IID", n_omega1="IID",
make_plots=TRUE)
savedir <- paste0(getwd(), '/mcmcfit_400_combined')
source("prepare_inputs.R")
fit <- tmbstan(Obj, lower=TmbList$Lower, upper=TmbList$Upper, chains=chains,
iter=iter, open_progress=FALSE, warmup=warmup,
init=prior.fn, seed=12512,
control=list(max_treedepth=td, adapt_delta=ad))
saveRDS(object = fit, file=paste0(savedir,'/mcmcfit.RDS'))
plot.mcmc(Obj, savedir, fit)
## Base case for paper: ATS
control <- list(model='ats', n_x=400,
n_eps1=1, n_eps2=1, n_omega2=1, n_omega1=1,
make_plots=TRUE)
savedir <- paste0(getwd(), '/mcmcfit_400_ats')
source("prepare_inputs.R")
fit <- tmbstan(Obj, lower=TmbList$Lower, upper=TmbList$Upper, chains=chains,
iter=iter, open_progress=FALSE, warmup=warmup,
init=prior.fn, seed=1245123,
control=list(max_treedepth=td, adapt_delta=ad))
saveRDS(object = fit, file=paste0(savedir,'/mcmcfit.RDS'))
plot.mcmc(Obj, savedir, fit)
## Base case for paper: BTS
control <- list(model='bts', n_x=400,
n_eps1=1, n_eps2=1, n_omega2=1, n_omega1=1,
make_plots=TRUE)
savedir <- paste0(getwd(), '/mcmcfit_400_bts')
source("prepare_inputs.R")
fit <- tmbstan(Obj, lower=TmbList$Lower, upper=TmbList$Upper, chains=chains,
iter=iter, open_progress=FALSE, warmup=warmup,
init=prior.fn, seed=653,
control=list(max_treedepth=td, adapt_delta=ad))
saveRDS(object = fit, file=paste0(savedir,'/mcmcfit.RDS'))
plot.mcmc(Obj, savedir, fit)