Mark-recapture with occasion and individual effects: Abundance estimation through Bayesian model selection in a fixed dimensional parameter space J.W. DURBAN AND D.A. ELSTON ********************************************************************************** WINBUGS v1.4 code to perform mark-recapture analyses described in manuscript in press with Journal of Agricultural, Biologial and Encironmental Statistics (JABES; http://www.amstat.org/publications/jabes/). Created by John Durban (John.Durban@noaa.gov) and David Elston (d.elston@macaulay.ac.uk) WinBUGS v1.4 software and instruction for running this code can be downloaded from the web at http://www/mrc-bsu.cam.ac.uk/bugs/ ********************************************************************************** DATA list(n=68,J=6,M=200) # change M for different numbers of candidate unobserved individuals #Also load binary matrix (Y) for individual by occasion capture matrix #Also load binary matrix (X) for M-rows and J-columns of zeros for unobserved #individuals INITIAL VALUES # set initial values for each MCMC chain list(mulogit = -2, n0=10) list(mulogit = -2, n0=25) list(mulogit = -2, n0=50) MODEL{ ###observed individuals### for (i in 1 : n) {for (j in 1 : J) { logit(p[i, j]) <- mulogit + theta[i] + alpha[j] Y[i, j] ~ dbern(p[i, j])}} ###Priors### mulogit~dnorm(0,0.1) meanprob <- exp(mulogit)/(1+exp(mulogit)) for(i in 1:n){ theta[i] ~ dnorm(0, tau.th)} tau.th ~ dgamma(a.th,b.th) a.th <- v0.th/2 b.th <- v0.th*priorvar.th/2 priorvar.th <- 1 v0.th<-1 #priorvar.th is our prior belief about the variance of theta #degree of strength of the prior beliefs is contained in the degrees of freedom #parameter v0.th sigma2.th <- 1/tau.th sigma.th <- sqrt(sigma2.th) for (j in 1:J) { alpha[j] ~ dnorm(0, tau.a)} tau.a ~ dgamma(a.a,b.a) a.a <- v0.a/2 b.a <- v0.a*priorvar.a/2 priorvar.a <- 1 v0.a<-1 sigma2.a <- 1/tau.a sigma.a <- sqrt(sigma2.a) ###allow for the fraction of the population seen### d <- 0.000001 n ~ dbin(d,N) ###prior for N### N <- n + n0 n0 ~ dcat(prob[]) for(j in 1: M){ prob[j] <- weight[j] / sum(weight[]) weight[j] <- pow(n+j, -delta)} delta <- 0 #for uniform prior # delta <- 1 for Jeffreys prior ###unobserved individuals### for (i in 1:M){ for (j in 1:J){ logit(p0[i,j])<- mulogit + theta0[i] + alpha[j] -100000*step(i-n0-0.5) X[i,j] ~dbern(p0[i,j]) #X to be read in as matrix of zeros } inmod[i] <- step(n0+0.5-i) mu.th0[i]<-equals(inmod[i],1)*mu.th0.in + equals(inmod[i],0)*mu.th0.out #mu.th0[i] <- 0 #for pseudoprior calculation tau.th0[i]<-equals(inmod[i],1)*tau.th + equals(inmod[i],0)*tau.out theta0[i]~ dnorm(mu.th0[i],tau.th0[i]) } ###calculating pseudoprior by monitoring an included individual### th0.1.in <- equals(inmod[1],1)*theta0[1] + equals(inmod[1],0)*mu.th0.out n.1.in <- equals(inmod[1],1) ### posterior means (sd) from pilot run of th0_1_in are -0.84 (0.90)### mu.th0.out <- -0.84 sd.th0.out <- 0.90 sdnew <- G*sd.th0.out G <- 1 tau.out <- 1/(sdnew*sdnew) #tau.out <- 0.001 # for psuedoprior calculation mu.th0.in <- 0 } -- Biomathematics and Statistics Scotland (BioSS) is formally part of The James Hutton Institute (JHI), a registered Scottish charity No. SC041796 and a company limited by guarantee No. SC374831