## ## BUGS code for the site occupancy model ## var ## Data Ns, # Number of species Nq, # Number of quadrats Y1[Ns, Nq], # Detection/non-detection of each species # in each quadrat in 1992 Y2[Ns, Nq], # Detection/non-detection of each species # in each quadrat in 2014 Gap[Nq], # Is in quadrats affected by oak dieback? Nadj, # Data for GeoBUGS Adj[Nadj], Weights[Nadj], Num[Nq], ## Parameters psi1[Ns, Nq], # Probability of presence # in 1992 for each species psi2[Ns, Nq], # Probability of presence # in 2014 for each species phi[Ns, Nq], # Survival probability gamma[Ns, Nq], # Colonization probability beta.o, # Intercept of occurrence beta.s, # Intercept of survival beta.sg, # Fixed effect of gap on survival beta.c, # Intercept of colonization beta.cg, # Fixed effect of gap on colonization r[Nq], # Spatial random effect epsilon.o # Random spcies effect on occurrence intercept epsilon.s # Random spcies effect on survival intercept epsilon.sg # Random spcies effect on survival coefficient of the gap epsilon.c # Random spcies effect on colonization intercept epsilon.cg # Random spcies effect on colonization coefficient of the gap sigma[6], tau[6]; model { ## Priors beta.o ~ dnorm(0, 1.0E-4) I(-10, 10); beta.s ~ dnorm(0, 1.0E-4) I(-10, 10); beta.sg ~ dnorm(0, 1.0E-4) I(-10, 10); beta.c ~ dnorm(0, 1.0E-4) I(-10, 10); beta.cg ~ dnorm(0, 1.0E-4) I(-10, 10); for (i in 1:Ns) { epsilon.o[i] ~ dnorm(0, tau[1]) I(-10, 10); epsilon.s[i] ~ dnorm(0, tau[2]) I(-10, 10); epsilon.sg[i] ~ dnorm(0, tau[3]) I(-10, 10); epsilon.c[i] ~ dnorm(0, tau[4]) I(-10, 10); epsilon.cg[i] ~ dnorm(0, tau[5]) I(-10, 10); } r[1:Nq] ~ car.normal(Adj[], Weights[], Num[], tau[6]); for (i in 1:6) { tau[i] <- 1 / (sigma[i] * sigma[i]); sigma[i] ~ dunif(0, 1.0E+2); } for (i in 1:Ns) { for (j in 1:Nq) { Y1[i, j] ~ dbern(psi1[i, j]); Y2[i, j] ~ dbern(psi2[i, j]); logit(psi1[i, j]) <- beta.o + epsilon.o[i] + r[j]; psi2[i, j] <- Y1[i, j] * phi[i, j] + (1 - Y1[i, j]) * gamma[i, j]; logit(phi[i, j]) <- beta.s + epsilon.s[i] + (beta.sg+ epsilon.sg[i]) * Gap[j]; logit(gamma[i, j]) <- beta.c + epsilon.c[i] + (beta.cg + epsilon.cg[i]) * Gap[j]; } } }