-
Notifications
You must be signed in to change notification settings - Fork 5
/
Copy pathJARA_7examples.R
73 lines (60 loc) · 2.42 KB
/
JARA_7examples.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
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
##############################################################################
# JARA: Just Another Redlisting Assessment
# Bayesian state-space redlisting tool
# Prime File (Settings)
# Developed by: Henning Winker* Nathan Pacoureau & Richard Sherley**
# * JRC - European Commission, Ispra, Italy
# Email: [email protected]
# **University of Exeter, UK ([email protected])
# Email: [email protected]
#
##############################################################################
rm(list=ls())
cat("\014") # clear console# set working directory
library(JARA)
data("jaradata")
dat = jrdat
# Set Working directory file, where assessments are stored
File = "C:/Work/Research/Github/JARAruns"
# Read assessment species information file
sp.assess = jara.examples
# Check species
print(data.frame(Assessment.list=sp.assess[,1]))
# Select assessment species by row numbers
specs = 1:nrow(sp.assess) # select all species
specs=1
# Loop over selected species
for(sp in specs){
spsel= sp.assess[sp,]
#-----------------------------------------------------------------
# # Extract info for build_jara() for each example
#-----------------------------------------------------------------
ts = dat[[paste0(spsel$assessment)]] # get time series
assessment = spsel$assessment # assessment name
run = spsel$run
abundance = spsel$abundance # c("census","relative")
GL = spsel$generation.length # numeric generation length (years)
pk.i = NULL
if(spsel[,1] =="Mountain_zebra") pk.i = 1/c(1.55,1.01,0.85,1.75,1.18,1.02,1.27,0.88,0.95)
fixed.obsE = spsel$sigma.obs.add # numeric, fixed component of observation error
# Organise folders
output = file.path(File,assessment)
dir.create(output,showWarnings = F)
#--------------------------------------------------
# Build JARA model
#--------------------------------------------------
jara.input = build_jara(I=ts$I,se=ts$SE,model.type = abundance,
assessment = assessment,scenario = "noK",
GL=GL,fixed.obsE=fixed.obsE,
pk.prior = c(0.0001,0.1),
pk.i = NULL # Carrying capacitiy example for Mountain Zebra
)
# Check input
jrplot_indices(jara.input,as.png = T,output.dir = output)
# Run JARA with csv file output
fit = fit_jara(jarainput =jara.input,save.csvs = T,output.dir = output)
# plot all routine output
jara_plots(fit,output.dir = output)
}
#-----------------
# The End