-
Notifications
You must be signed in to change notification settings - Fork 2
/
EGVU_IPM_2020_v1_chicksupplementation.jags
356 lines (255 loc) · 17.6 KB
/
EGVU_IPM_2020_v1_chicksupplementation.jags
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
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
model {
#-------------------------------------------------
# integrated population model for the balkan population
# - age structured model with 6 age classes:
# - age-specific probability to recruit at ages 4 or 5
# - age-specific survival derived from tracking data
# - adult survival based on territory occupancy; CHANGES in 2016
# - pre breeding census, female-based assuming equal sex ratio & survival
# - FUTURE PROJECTION COMBINES SCENARIOS: baseline scenario (=do nothing, no chick removal)
# - productivity supplemented by captive-bred juveniles (0-10 per year)
# - fecundity reduced by harvesting second eggs ONLY FOR SCENARIOS OF CAPTIVE RELEASES
# - rescue of wild chicks optional (currently disabled)
#-------------------------------------------------
#-------------------------------------------------
# 1. PRIORS FOR ALL DATA SETS
#-------------------------------------------------
# Priors and constraints FOR FECUNDITY
mu.fec ~ dunif(0,2) # Priors on fecundity can range from 0- 2 chicks per pair (uninformative)
# Priors and constraints FOR POPULATION COUNTS OBSERVATION
sigma.obs.count ~ dunif(0,100) #Prior for SD of observation process (variation in detectability)
tau.obs.count<-pow(sigma.obs.count,-2)
# Priors and constraints FOR JUVENILE SURVIVAL FROM TELEMETRY
#### MONTHLY SURVIVAL PROBABILITY
for (i in 1:nind.telemetry){
for (t in f.telemetry[i]:(n.occasions.telemetry)){
logit(phi.telemetry[i,t]) <- lp.mean.telemetry + ### intercept for mean survival
b.phi.capt*(capt.telemetry[i]) + ### survival dependent on captive-release (captive-raised or other)
b.phi.mig*(mig.telemetry[i,t]) + ### survival dependent on first migration across the sea
b.phi.age*(age.telemetry[i,t]) ### survival dependent on age (juvenile or other)
} #t
} #i
#### BASELINE FOR SURVIVAL PROBABILITY (wild adult stationary from east)
mean.phi.telemetry ~ dunif(0.9, 1) # uninformative prior for all MONTHLY survival probabilities
lp.mean.telemetry <- log(mean.phi.telemetry/(1 - mean.phi.telemetry)) # logit transformed survival intercept
#### SLOPE PARAMETERS FOR SURVIVAL PROBABILITY
b.phi.age ~ dnorm(0, 0.001) # Prior for COST OF MIGRATION migration on survival probability on logit scale
b.phi.capt ~ dnorm(0, 0.001) # Prior for captive release on survival probability on logit scale
b.phi.mig ~ dunif(-2,0) # Prior for captive release on survival probability on logit scale
#### TAG FAILURE AND LOSS PROBABILITY
for (i in 1:nind.telemetry){
for (t in f.telemetry[i]:(n.occasions.telemetry)){
logit(p.obs.telemetry[i,t]) <- base.obs.telemetry
logit(tag.fail.telemetry[i,t]) <- base.fail.telemetry
logit(p.found.dead.telemetry[i,t]) <- base.recover.telemetry
} #t
} #i
##### SLOPE PARAMETERS FOR OBSERVATION PROBABILITY
base.obs.telemetry ~ dnorm(0, 0.001) # Prior for intercept of observation probability on logit scale
base.fail.telemetry ~ dnorm(0, 0.001) # Prior for intercept of tag failure probability on logit scale
base.recover.telemetry ~ dnorm(0, 0.001) # Prior for intercept of tag failure probability on logit scale
# Priors and constraints FOR ADULT SURVIVAL FROM TERRITORY MONITORING
## Priors for detection probability
# Det prob varies for each site, det prob on logit scale, around 0.05 - 0.95 on normal scale
for (ipterr in 1:nsite.terrvis){
lmu.p.terrvis[ipterr] ~ dunif(-3, 2)
} ## close k loop for det prob prior
# Det prob relationship with observation effort
beta.obs.eff.terrvis ~ dnorm(0, 0.0001)
# Priors for survival
## annual survival must be greater than 80% on logit scale
for (nypterr in 1:2){ ## only 2 survival periods
lm.phi.terrvis[nypterr] <- log(mean.phi.terrvis[nypterr]/(1 - mean.phi.terrvis[nypterr])) # logit transformed survival intercept
mean.phi.terrvis[nypterr] ~ dunif(0.75, 1) # informative prior for annual survival probabilities BEFORE 2016
}
### RANDOM OBSERVATION EFFECT FOR EACH TERRITORY AND YEAR
for (nyRpterr in 1:nprim.terrvis){
rand.obs.terrvis[nyRpterr] ~ dnorm(0, tau.obs.terrvis)
}
tau.obs.terrvis <- 1 / (sd.obs.terrvis * sd.obs.terrvis)
sd.obs.terrvis ~ dunif(0, 3)
# Priors for population process and immigration
# Initial population sizes for first year of monitoring
nestlings[1] ~ dunif(50, 70) ##changed from JUV
N1[1] ~ dunif(15, 35)
N2[1] ~ dunif(10, 32)
N3[1] ~ dunif(5, 20)
N4[1] ~ dunif(5, 17)
N5[1] ~ dunif(5, 15)
N6[1] ~ dunif(100, 140)
#-------------------------------------------------
# 2. LIKELIHOODS AND ECOLOGICAL STATE MODEL
#-------------------------------------------------
# -------------------------------------------------
# 2.1. System process: female based matrix model
# -------------------------------------------------
for (tt in 2:T.count){
nestlings[tt] <- mu.fec * 0.5 * Nterr[tt] ### number of local recruits
N1[tt] ~ dbin(ann.phi.juv.telemetry, round(nestlings[tt-1])) ### number of 1-year old survivors - add CAPT.ADD in here
N2[tt] ~ dbin(ann.phi.sec.telemetry, round(N1[tt-1])) ### number of 2-year old survivors
N3[tt] ~ dbin(ann.phi.third.telemetry, round(N2[tt-1])) ### number of 3-year old survivors
N4[tt] ~ dbin(mean.phi.terrvis[phase[tt]], round(N3[tt-1])) ### number of 4-year old survivors
N5[tt] ~ dbin(mean.phi.terrvis[phase[tt]], round(N4[tt-1])) ### number of 5-year old survivors
N6[tt] ~ dbin(mean.phi.terrvis[phase[tt]], round((N5[tt-1]+N6[tt-1]))) ### number of 6-year or older (adult) birds
} # tt
# -------------------------------------------------
# 2.2. Observation process for population counts: state-space model of annual counts
# -------------------------------------------------
for (tlc in 1:T.count){
Nterr[tlc] <- N4[tlc] * 0.024 + N5[tlc] * 0.124 + N6[tlc] ### number of observable territorial birds
y.count[tlc] ~ dnorm(Nterr[tlc], tau.obs.count) # Distribution for random error in observed numbers (counts)
#}
# -------------------------------------------------
# 2.3. Likelihood for fecundity: Poisson regression from the number of surveyed broods (or binom for single broods after egg harvest)
# -------------------------------------------------
#for (t in 1:nyear.fec){
J.fec[tlc] ~ dpois(rho.fec[tlc])
rho.fec[tlc] <- R.fec[tlc]*mu.fec
} # close loop over every year in which we have count and fecundity data
# -------------------------------------------------
# 2.4. Likelihood for juvenile survival from telemetry
# -------------------------------------------------
# -------------------------------------------------
# Define state-transition and observation matrices
# -------------------------------------------------
for (i in 1:nind.telemetry){
for (t in f.telemetry[i]:(n.occasions.telemetry-1)){
# Define probabilities of state S(t+1) [last dim] given S(t) [first dim]
ps[1,i,t,1]<-1 ## dead birds stay dead
ps[1,i,t,2]<-0
ps[1,i,t,3]<-0
ps[2,i,t,1]<-(1-phi.telemetry[i,t])
ps[2,i,t,2]<-phi.telemetry[i,t] * (1-tag.fail.telemetry[i,t])
ps[2,i,t,3]<-phi.telemetry[i,t] * tag.fail.telemetry[i,t]
ps[3,i,t,1]<-(1-phi.telemetry[i,t])
ps[3,i,t,2]<-0
ps[3,i,t,3]<-phi.telemetry[i,t]
# Define probabilities of O(t) [last dim] given S(t) [first dim]
po[1,i,t,1]<-0
po[1,i,t,2]<-p.obs.telemetry[i,t] * (1-tag.fail.telemetry[i,t]) * (1-p.found.dead.telemetry[i,t])
po[1,i,t,3]<-0
po[1,i,t,4]<-p.found.dead.telemetry[i,t]
po[1,i,t,5]<-(1-p.obs.telemetry[i,t]) * tag.fail.telemetry[i,t] * (1-p.found.dead.telemetry[i,t])
po[2,i,t,1]<-p.obs.telemetry[i,t] * (1-tag.fail.telemetry[i,t])
po[2,i,t,2]<-0
po[2,i,t,3]<-0
po[2,i,t,4]<-0
po[2,i,t,5]<-(1-p.obs.telemetry[i,t]) * tag.fail.telemetry[i,t]
po[3,i,t,1]<-0
po[3,i,t,2]<-0
po[3,i,t,3]<-0
po[3,i,t,4]<-0
po[3,i,t,5]<-1
} #t
} #i
# Likelihood
for (i in 1:nind.telemetry){
# Define latent state at first capture
z.telemetry[i,f.telemetry[i]] <- 2 ## alive when first marked
for (t in (f.telemetry[i]+1):n.occasions.telemetry){
# State process: draw S(t) given S(t-1)
z.telemetry[i,t] ~ dcat(ps[z.telemetry[i,t-1], i, t-1,])
# Observation process: draw O(t) given S(t)
y.telemetry[i,t] ~ dcat(po[z.telemetry[i,t], i, t-1,])
} #t
} #i
# -------------------------------------------------
# 2.5. Likelihood for adult survival from territory monitoring
# -------------------------------------------------
### ECOLOGICAL STATE MODEL WITH ESTIMATE OF SURVIVAL
for (ilterr in 1:nsite.terrvis){
#z.terrvis[ilterr,f.obsvis[ilterr]] <- firstobs[ilterr]
for (klterr in (f.obsvis[ilterr]+1):nprim.terrvis){
z.terrvis[ilterr,klterr] ~ dbin(mean.phi.terrvis[phase[klterr]],z.terrvis[ilterr,klterr-1])
} # close klterr loop over primary period - years
} # close ilterr loop over sites
### OBSERVATION MODEL WITH RANDOM EFFECT OF YEAR
## observation prob depends on monthly observation effort
for (iobs in 1:nsite.terrvis){
for (kobs in f.obsvis[iobs]:nprim.terrvis){
for (jobs in 1:nrep.terrvis){
linpred.terrvis[iobs,jobs,kobs] <- lmu.p.terrvis[iobs]+beta.obs.eff.terrvis*eff.terrvis[iobs,jobs,kobs] + rand.obs.terrvis[kobs]
p.terrvis[iobs,jobs,kobs] <- 1 / (1 + exp(-linpred.terrvis[iobs,jobs,kobs]))
y.terrvis[iobs,jobs,kobs] ~ dbin(p.terrvis[iobs,jobs,kobs], z.terrvis[iobs,kobs])
} # close kobs loop over monthly surveys
} # close jobs loop over year
} # close iobs loop over sites
# -------------------------------------------------
# 3. DERIVED PARAMETERS
# -------------------------------------------------
### 3.1 TELEMETRY DERIVED SURVIVAL ESTIMATES
## for WILD BIRDS
for (ageprog in 1:36){
logit(phi.wild.telemetry[ageprog]) <- lp.mean.telemetry + ### intercept for mean survival
b.phi.mig*(migprog.age[ageprog]) + ### survival dependent on migration (first-time crossing of sea)
b.phi.age*ageprog ### survival dependent on age (juvenile or other)
}
ann.phi.juv.telemetry<-phi.wild.telemetry[1]*phi.wild.telemetry[2]*phi.wild.telemetry[3]*phi.wild.telemetry[4]*phi.wild.telemetry[5]*phi.wild.telemetry[6]*phi.wild.telemetry[7]*phi.wild.telemetry[8]*phi.wild.telemetry[9]*phi.wild.telemetry[10]*phi.wild.telemetry[11]*phi.wild.telemetry[12] ### multiply monthly survival from age 1-12
ann.phi.sec.telemetry<-phi.wild.telemetry[13]*phi.wild.telemetry[14]*phi.wild.telemetry[15]*phi.wild.telemetry[16]*phi.wild.telemetry[17]*phi.wild.telemetry[18]*phi.wild.telemetry[19]*phi.wild.telemetry[20]*phi.wild.telemetry[21]*phi.wild.telemetry[22]*phi.wild.telemetry[23]*phi.wild.telemetry[24] ### multiply monthly survival from age 13-24
ann.phi.third.telemetry<-phi.wild.telemetry[25]*phi.wild.telemetry[26]*phi.wild.telemetry[27]*phi.wild.telemetry[28]*phi.wild.telemetry[29]*phi.wild.telemetry[30]*phi.wild.telemetry[31]*phi.wild.telemetry[32]*phi.wild.telemetry[33]*phi.wild.telemetry[34]*phi.wild.telemetry[35]*phi.wild.telemetry[36] ### multiply monthly survival from age 25-36
## for CAPTIVE-REARED DELAYED RELEASE BIRDS
for (captageprog in 1:36){
logit(phi.capt.telemetry[captageprog]) <- lp.mean.telemetry + ### intercept for mean survival
b.phi.capt + ### survival dependent on captive-release (captive-raised or other)
b.phi.mig*(captprog.age[captageprog]) + ### survival dependent on migration (first-time crossing of sea)
b.phi.age*captageprog ### survival dependent on age (juvenile or other)
}
ann.phi.capt.rel.first.year<-phi.capt.telemetry[8]*phi.capt.telemetry[9]*phi.capt.telemetry[10]*phi.capt.telemetry[11]*phi.capt.telemetry[12]*phi.capt.telemetry[13]*phi.capt.telemetry[14]*phi.capt.telemetry[15]*phi.capt.telemetry[16]*phi.capt.telemetry[17]*phi.capt.telemetry[18]*phi.capt.telemetry[19]*phi.capt.telemetry[20]*phi.capt.telemetry[21]*phi.capt.telemetry[22]*phi.capt.telemetry[23]*phi.capt.telemetry[24] ### first year for delayed-release bird is longer, from month 8 to 24
### 3.2 POPULATION GROWTH DERIVED FROM COUNTS
# Annual population growth rate
for (ipop in 1:(T.count-1)){
lambda.t[ipop] <- Nterr[ipop+1] / max(Nterr[ipop],1) ### inserted max of 1 to avoid this becoming an invalid parent error
loglambda.t[ipop]<-log(lambda.t[ipop])## for calculating geometric mean of overall population growth rate
}
#### OVERALL POPULATION GROWTH RATE #########
mean.lambda<-exp((1/(T.count-1))*sum(loglambda.t[1:(T.count-1)])) # Geometric mean
# -------------------------------------------------
# 4. PREDICTION INTO THE FUTURE
# -------------------------------------------------
### INCLUDE SCENARIOS FOR N CAPTIVE BIRDS AND SURVIVAL IMPROVEMENT
# CAPTIVE RELEASE OF JUVENILE BIRDS
for (ncr in 1:scen.capt.release){
### FUTURE FECUNDITY IS BASED ON mu.fec - removed the chick removal scenario
for (fut in 1:PROJECTION){
fut.fec[fut,ncr] <-mu.fec ## no chicks taken after 5 years anymore
}
# SPECIFY IMPROVEMENT OF SURVIVAL
for (is in 1:scen.imp.surv){
## POPULATION PROCESS
### DOES NOT WORK WITH SAME ARRAY OF POP TRAJECTORY AS 3 dimensions required
## need to copy previous array elements
nestlings.f[ncr,is,1] <- round(nestlings[T.count]) ##JUV[T.count]
N1.f[ncr,is,1] <- N1[T.count]
N2.f[ncr,is,1] <- N2[T.count]
N3.f[ncr,is,1] <- N3[T.count]
N4.f[ncr,is,1] <- N4[T.count]
N5.f[ncr,is,1] <- N5[T.count]
N6.f[ncr,is,1] <- N6[T.count]
Nterr.f[ncr,is,1] <- Nterr[T.count]
N2wild.f[ncr,is,1] <- 0
N2released.f[ncr,is,1] <- 0
need.to.breed[ncr,is,1] <- max(0,capt.release[1,ncr]-nestlings.f[ncr,is,1])
for (fut in 2:PROJECTION){
fut.survival[ncr,is,fut] <-min(imp.surv[fut,is]*mean.phi.terrvis[2],1) ### invalid parent error if survival>1
### probabilistic formulation
nestlings.f[ncr,is,fut] <- max(0,((fut.fec[fut,ncr] * 0.5 * Nterr.f[ncr,is,fut])-capt.release[fut,ncr])) ### number of local recruits calculated as local production minus those taken into captivity
need.to.breed[ncr,is,fut] <- max(0,(capt.release[fut,ncr]-(fut.fec[fut,ncr] * 0.5 * Nterr.f[ncr,is,fut])))
N1.f[ncr,is,fut] ~ dbin(min((imp.surv[fut,is]*ann.phi.juv.telemetry),1),round(nestlings.f[ncr,is,fut-1])) ### +rescued[ncr,is,fut] number of 1-year old survivors
N2wild.f[ncr,is,fut] ~ dbin(min((imp.surv[fut,is]*ann.phi.sec.telemetry),1),round(N1.f[ncr,is,fut-1])) ### number of 2-year old wild survivors
N2released.f[ncr,is,fut] ~ dbin(min((imp.surv[fut,is]*ann.phi.capt.rel.first.year),1),round(capt.release[fut,ncr])) ### +rescued[ncr,is,fut] number of 1-year old survivors
N2.f[ncr,is,fut] <- N2wild.f[ncr,is,fut] + N2released.f[ncr,is,fut] ### sum of the N1 cohort derived from wild and captive-bred juveniles
N3.f[ncr,is,fut] ~ dbin(min((imp.surv[fut,is]*ann.phi.third.telemetry),1),round(N2.f[ncr,is,fut-1])) ### number of 3-year old survivors
N4.f[ncr,is,fut] ~ dbin(fut.survival[ncr,is,fut],round(N3.f[ncr,is,fut-1])) ### number of 4-year old survivors
N5.f[ncr,is,fut] ~ dbin(fut.survival[ncr,is,fut],round(N4.f[ncr,is,fut-1])) ### number of 5-year old survivors
N6.f[ncr,is,fut] ~ dbin(fut.survival[ncr,is,fut],round((N5.f[ncr,is,fut-1]+N6.f[ncr,is,fut-1]))) ### number of 6-year or older (adult) birds
Nterr.f[ncr,is,fut] <- round((N4.f[ncr,is,fut] * 0.024) + (N5.f[ncr,is,fut] * 0.124) + (N6.f[ncr,is,fut]))
} # fut
for (fut2 in 1:(PROJECTION-1)){
lambda.t.f[ncr,is,fut2] <- Nterr.f[ncr,is,fut2+1] / max(Nterr.f[ncr,is,fut2],1)
loglambda.t.f[ncr,is,fut2]<-log(lambda.t.f[ncr,is,fut2])## for calculating geometric mean of overall population growth rate
} # fut2
#### FUTURE POPULATION GROWTH RATE #########
fut.lambda[ncr,is]<-exp((1/(PROJECTION-1))*sum(loglambda.t.f[ncr,is,1:(PROJECTION-1)])) # Geometric mean
} # end scenario of imp.surv
} # end scenario of n capt released
} # END of the JAGS model