SCRi.fn SCRi.fn <- function(scrobj, ni=1100,burn=100,skip=2,nz=200,theta=NA, Msigma=1,Mb=0,Msex=0,Msexsigma = 0,Xeff=NULL,Xsex=NULL, ss.prob=NULL, coord.scale=5000,area.per.pixel=1,thinstatespace=1,maxNN=20,dumprate=1000){ # Added input vector ss.prob, which gives a proportional (RSF-like) weight to each node in the statespace call <- match.call() traps<-scrobj$traps captures<-scrobj$captures statespace<-scrobj$statespace alive<-scrobj$alive Xd<- scrobj$Xd ## current SCRh and prior versions: trapid, individual, occasion ## ## EDF: session, individual, occasion, trapID ## Now we assume EDF as input and need to reorganize to the old format so the rest of the code works captures<- cbind(captures[,4],captures[,2],captures[,3]) if( length(unique(captures[,2])) != length(min(captures[,2]):max(captures[,2])) ) { cat("Error: individuals not numbered sequentially, renumbering them now",fill=TRUE) captures[,2]<- as.numeric(factor(captures[,2])) } # ni = number of iterations, total # burn = number to discard # skip = thin rate (i.e., keep one in skip iterations) # Msigma = 0 fits non-spatial model # statespace (formerly GRID) is the discrete grid of the study area # nz = number of "all zero" encounter histories to add ##traps=datafile containing the XY locations of the grid center points ###captures=data file containing the ID of the animal, the Sampling occassion the animals was captured on, and the number of the trap the animal was captured in. ###statespace=data file containing a fine grid of points ## SCRf.fn expands code to include 3 types of "fixed" covariates that are easy to deal with ## (1) trap-specific only ## (2) time-specific only ## (3) trap and time-specific (e.g., effort) ## all of these can be put into a trap x time matrix "Xeff" ## Other types: vary by individual (sex) ## Spatial: habitat or things that vary over the state-space ## ## "effort" = Xeff = ntraps x ndays matrix or NULL. ## Xeff1 = test data, X traps x Y periods (random numbers) assume log(eff) Y<-captures traplocs<-traps[,2:3] MASK<-as.matrix(traps[,4:ncol(traps)]) nind<-max(Y[,2]) T<-dim(MASK)[2] M<-nind+nz # total size of data set ntraps<-nrow(traplocs) ################## ################## ################## ################## totalarea<- nrow(statespace)*area.per.pixel thinned<- seq(1,nrow(statespace),thinstatespace) statespace<-statespace[thinned,] Xd<-Xd[thinned] goodbad<-statespace[,3] G<-statespace[,1:2] G<-G[goodbad==1,] Gunscaled<-G Xd<- Xd[goodbad==1] nG<-nrow(G) new.area.per.pixel<- totalarea/nG # The following lines thin the statespace probability weights to match the thinned statespace ######################################################### if(!is.null(ss.prob)){ ss.prob = ss.prob[thinned] ss.prob = ss.prob[goodbad==1] } ########################################################## ### ### ## following lines scale coordinate system. This is kind of arbitrary. ### ### mgx<-min(traplocs[,1]) mgy<-min(traplocs[,2]) traplocs[,1]<-(traplocs[,1]-min(traplocs[,1]))/coord.scale traplocs[,2]<-(traplocs[,2]-min(traplocs[,2]))/coord.scale G[,1]<-(G[,1]-mgx)/coord.scale G[,2]<-(G[,2]-mgy)/coord.scale ### # create "Data" vector but with trap mask information ### msk2<-array(NA,c(nind+nz,T,ntraps)) for(i in 1:(nind+nz)){ msk2[i,1:T,1:ntraps]<-t(MASK[1:ntraps,1:T]) } msk2<-as.vector(msk2) ### ### # expand the data to a 3-d array Ynew<-array(0,dim=c(nind,T,ntraps)) Ynew[cbind(Y[,2],Y[,3],Y[,1])]<-1 Y<-Ynew ### ### data augmentation ### Yaug<-array(0,dim=c(nind+nz,T,ntraps)) for(j in 1:nind){ Yaug[j,1:T,1:ntraps]<-Y[j,1:T,1:ntraps] } if(!is.null(Xeff)){ Xeffnew<-array(0,dim=c(nind+nz,T,ntraps)) for(j in 1:M){ Xeffnew[j,1:T,1:ntraps]<-t(Xeff) } Xeff.tf<-TRUE } if(is.null(Xeff)){ Xeffnew<-array(0,dim=c(nind+nz,T,ntraps)) Xeff.tf<-FALSE } Xeff<-Xeffnew if(!is.null(Xsex)){ Xsexnew<-c(Xsex,rep(NA,nz)) } if(is.null(Xsex)){ Xsexnew<-rep(0,nind+nz) } Xsex<-Xsexnew sex.naflag<-is.na(Xsex) ### # create covariate of previous capture ### prevcap<-array(0,c(nind+nz,T,ntraps)) for(i in 1:(nind)){ for(j in 1:ntraps){ tmp<-Yaug[i,1:T,j] if(any(tmp==1)){ fst<- min( (1:T)[tmp==1] ) if(fst0) Xsex[sex.naflag]<-rbinom(sum(sex.naflag),1,psi.sex) # modified 10/15/13 JFG - changed mean binom prob to mean of known pop. instead of .65 if(is.null(Xd)){ Xd<-rep(1,nG) } beta.den<-0 ### ### ## initializing things and making utility funcs ### ### lik.fn<-function(lpc,y1){ llvector.new<- -1*exp(lpc) part2<- exp(exp(lpc[y1])) - 1 part2[part2==0]<-.Machine$double.eps llvector.new[y1]<- llvector.new[y1] + log(part2) llvector.new } trapgridbig<-traplocs[trapid,] # streteches out the trap coord matrix y1<-y==1 c1<- (S[indid,1]-trapgridbig[,1])^2 c2<- (S[indid,2]-trapgridbig[,2])^2 gof.new<-gof.data<-rep(NA,(ni-burn)/skip) out<-matrix(NA,nrow=(ni-burn)/skip,ncol=14) dimnames(out)<-list(NULL,c("bsigma","sigma","bsigma2","sigma2","lam0","beta.behave","beta1(effort)","beta.sex","psi","psi.sex","Nsuper","theta","beta.density","D")) zout<-matrix(NA,nrow=(ni-burn)/skip,ncol=M) Sout<-matrix(NA,nrow=(ni-burn)/skip,ncol=M) LLout<-rep(NA, (ni-burn)/skip) m<-1 # edits 8/14/2013 #LM1<-LM2<-matrix(NA,nrow=M,ncol=length(indid)/M) LM1<-LM2<-matrix(0,nrow=M,ncol=length(indid.LM)/M) ones<-rep(1,ncol(LM1)) if(Msexsigma==0) lp.sigma<-Msigma*bsigma*(c1+c2)^theta if(Msexsigma==1) lp.sigma<-bsigma[Xsex[indid]+1]*(c1+c2)^theta acc.count<-0 delta<-.05 for(i in 1:ni){ cat("iter: ",i,fill=TRUE) ######################## ######################## ######################## ######################## # PART 1 OF THE MCMC ALGORITHM UPDATES THE REGRESSION PARAMETERS. FOR THIS MODEL # THE REGRESSION PARAMETERS ARE (1) INTERCEPT (2) EFFECT OF PREVIOUS CAPTURE # (BEHAVIORAL RESPONSE) (3) THE SPATIAL PARAMETER "sigma" ### Updating parameters here should only involve guys with z = 1 (i.e., members of the population) ### update loglam0 lp<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] loglam0c<-rnorm(1,loglam0,.1) lpc<- loglam0c + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector<-lik.fn(lp,y1) llvector.new<-lik.fn(lpc,y1) #LM1[1:length(LM1)]<-llvector.new #LM2[1:length(LM2)]<- llvector ### July 2013 LM1[aliveid==1]<-llvector.new LM2[aliveid==1]<- llvector if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ loglam0<-loglam0c lam0<-exp(loglam0) llvector<-llvector.new lp<-lpc LM2<-LM1 } ## update theta if(update.theta){ lp<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] thetac<-rnorm(1,theta,.02) # theta must be between exponential (.5) and gaussian (1) if(thetac>=0.5 & thetac<=1){ if(Msexsigma==0) lp.sigmac<-Msigma*bsigma*(c1+c2)^thetac if(Msexsigma==1) lp.sigmac<-bsigma[Xsex[indid]+1]*(c1 + c2)^thetac lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector<-lik.fn(lp,y1) llvector.new<-lik.fn(lpc,y1) #LM1[1:length(LM1)]<-llvector.new #LM2[1:length(LM2)]<- llvector # Robin Russell # July 2013 LM1[aliveid==1]<-llvector.new LM2[aliveid==1]<- llvector if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ theta<-thetac lp.sigma<-lp.sigmac llvector<-llvector.new lp<-lpc LM2<-LM1 } } } ### Update bsigma ## do this differently depending on if sigma depends on sex if(Msexsigma==0){ bsigmac<-exp(rnorm(1,log(bsigma),delta)) lp.sigmac<- Msigma*bsigmac*(c1+c2)^theta lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ##LM1[1:length(LM1)]<- llvector.new LM1[aliveid==1]<- llvector.new ## July 2013 if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ lp.sigma<-lp.sigmac bsigma<-bsigmac llvector<-llvector.new lp<-lpc LM2<-LM1 acc.count<- acc.count - 1 } else{ acc.count<-acc.count+1 } } if(Msexsigma==1){ bsigmac<-c(exp(rnorm(1,log(bsigma[1]),2*delta)),bsigma[2]) lp.sigmac<- bsigmac[Xsex[indid]+1]*(c1+c2)^theta lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ###LM1[1:length(LM1)]<- llvector.new # July 2013 LM1[aliveid==1]<- llvector.new if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ lp.sigma<-lp.sigmac bsigma<-bsigmac llvector<-llvector.new lp<-lpc LM2<-LM1 # delta<-delta*1.2 } else{ # delta<- delta*.8 } bsigmac<-c(bsigma[1],exp(rnorm(1,log(bsigma[2]),2*delta))) lp.sigmac<- bsigmac[Xsex[indid]+1]*(c1+c2)^theta lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ####LM1[1:length(LM1)]<- llvector.new # July 2013 LM1[aliveid==1]<- llvector.new if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ lp.sigma<-lp.sigmac bsigma[2]<-bsigmac[2] llvector<-llvector.new lp<-lpc LM2<-LM1 } } cat("DELTA: ",delta,fill=TRUE) if(any(bsigma<0) ){ cat("negative bsigma....",fill=TRUE) return(0) } if(Msex==1){ beta.sexc<- rnorm(1,beta.sex,.1) lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sexc*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ##LM1[1:length(LM1)]<- llvector.new # July 2013 LM1[aliveid==1]<- llvector.new if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ beta.sex<- beta.sexc llvector<-llvector.new lp<-lpc LM2<-LM1 # LM2 is current value } } ### This code allows a "local behavioral effect" as opposed to a conventional ### "global effect" which would not be trap specific. if(Mb==1){ beta.behave.c<- rnorm(1,beta.behave,.1) lpc<- loglam0 + Mb*beta.behave.c*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ##LM1[1:length(LM1)]<- llvector.new ## July 2013 LM1[aliveid==1]<- llvector.new if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ beta.behave<- beta.behave.c llvector<-llvector.new lp<-lpc LM2<-LM1 # LM2 is current value } } if(Xeff.tf){ ### This block of code below deals with effort. beta1c<- rnorm(1,beta1,.1) lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1c*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ##LM1[1:length(LM1)]<- llvector.new ## July 2013 LM1[aliveid==1]<- llvector.new if(runif(1)< exp(sum(( (LM1[z==1,]-LM2[z==1,])%*%ones)))){ beta1<- beta1c llvector<-llvector.new lp<-lpc LM2<-LM1 # LM2 is current value } } ############################## ############################## ######################## ######################## ######################## # PART 2 OF THE MCMC ALGORITHM UPDATES THE DATA AUGMENTATION PARAMETERS # THIS INCLUDES THE LATENT "z" VARIABLES AS WELL AS THE CRITICAL # PARAMETER "psi" ######################## ######################## # This is the data augmentation part. This code updates each z[i] # for all i=1,2,...M individuals. z[i]=1 is a "real" animal that lives # in S, whereas z[i]=0 are excess zeros in the data set # this is an application of Bayes rule to get Pr(z=1| y[i,,]=0) probz<- exp( rowsum(llvector[indid>nind],indid[indid>nind]) ) # only for nind+1...M probz<- (probz*psi )/(probz*psi + (1-psi)) z[(nind+1):M]<-rbinom(M-nind,1,probz) psi<-rbeta(1,1+sum(z),1+M-sum(z)) ### ### This updates SEX identifier and sex ratio ### if(!is.null(Xsex)){ # here we only switch sex state if sex is missing for an individual tmp.sex<-Xsex tmp.sex[sex.naflag]<- 1-Xsex[sex.naflag] if(Msexsigma==0) lp.sigmac<-Msigma*bsigma*(c1+c2)^theta if(Msexsigma==1) lp.sigmac<-bsigma[tmp.sex[indid]+1]*(c1+c2)^theta lpc<- loglam0 + Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*tmp.sex[indid] ## + sex contribution if sex is in the model! llvector.new<-lik.fn(lpc,y1) lik.othersex<- exp(rowsum(llvector.new,indid)) lik.sex<- exp( rowsum(llvector,indid) ) # need to do this only for MISSING SEX is done here for all animals prior.curr<-(psi.sex^Xsex)*((1-psi.sex)^(1-Xsex)) prior.cand<-(psi.sex^tmp.sex)*((1-psi.sex)^(1-tmp.sex)) swtch<- sex.naflag & (runif(M,0,1)< ( (lik.othersex*prior.cand)/(lik.sex*prior.curr))) Xsex[swtch]<- 1-Xsex[swtch] psi.sex<-rbeta(1,.1+sum(Xsex[z==1]),.1+sum(z)-sum(Xsex[z==1])) if(Msexsigma==0) lp.sigma<-Msigma*bsigma*(c1+c2)^theta if(Msexsigma==1) lp.sigma<-bsigma[Xsex[indid]+1]*(c1+c2)^theta lp<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lp,y1) ####LM1[1:length(LM1)]<- llvector.new ## July 2013 LM1[aliveid==1]<- llvector.new llvector<-llvector.new LM2<-LM1 } ######################## ######################## ## PART III OF THE ALGORITHM -- THIS BLOCK OF CODE UPDATES THE ## ACTIVITY CENTERS. ######################## ######################## ################################################################## ### Gets new centers with probability weightings - 10/9/13 JFG if (!is.null(ss.prob)){ newcenters = rep(NA, M) for (gc in 1:M){ newcenters[gc] <- sample(1:numnn[centers[gc]],size=1,replace=TRUE, prob=ss.prob[na.omit(unique(NN[centers[gc],]))]) } } if (is.null(ss.prob)){ newcenters <- trunc(runif(M, 0, numnn[centers])) + 1 } ################################################################### newcenters<- NN[cbind(centers,newcenters)] # these probabilities are needed in the Metrpolis "acceptance probability" calculation # since the proposal is not symmetric (i.e., near the boundary of the state-space) qnew<- 1/numnn[centers] qold<- 1/numnn[newcenters] Sc<-G[newcenters,] c1c<- (Sc[indid,1]-trapgridbig[,1])^2 c2c<- (Sc[indid,2]-trapgridbig[,2])^2 if(Msexsigma==0) lp.sigmac<-Msigma*bsigma*(c1c+c2c)^theta if(Msexsigma==1) lp.sigmac<-bsigma[Xsex[indid]+1]*(c1c+c2c)^theta lpc<- loglam0+ Mb*beta.behave*prevcap - lp.sigmac + beta1*Xeff + Msex*beta.sex*Xsex[indid] llvector.new<- lik.fn(lpc,y1) ####LM1[1:length(LM1)]<- llvector.new LM1[aliveid==1]<- llvector.new likdiff<- (LM1-LM2)%*%ones likdiff[z==0]<-0 # all other things equal, this lines sets acceptance prob to 1 for z=0 guys ### edits 7/23/2013 logprior.new<- Xd[newcenters]*beta.den logprior.old<- Xd[centers]*beta.den #likdiff<-likdiff + log(qold/qnew) ### I think dimension of these is wrong...... likdiff<-likdiff + log(qold/qnew) + (logprior.new-logprior.old) accept<- runif(M)burn) & (i%%skip == 0) ){ sigma<- sqrt(1/(2*bsigma)) if(Msexsigma == 0){ sigmatmp<-c(sigma,sigma) bsigmatmp<-c(bsigma,bsigma) } else{ sigmatmp<-sigma bsigmatmp<-bsigma } ##### ## Jan 25 GoF stuff right here for JWM paper ## 09/19/2011 ##### # mean of "real individuals" (not fixed zeros) # realguys<- z[indid]==1 # cat(sum(realguys)," real observations for GoF",fill=TRUE) logmu<- loglam0 + Mb*beta.behave*prevcap - lp.sigma + beta1*Xeff + Msex*beta.sex*Xsex[indid] mu<- ( 1-exp(-exp(logmu)))*z[indid] # zeros out the z=0 guys so they contribute nothing # JFG Added Oct. '13 these lines calculate the probability of observing the data given the parameter values logmu.array <- array(log(exp(-exp(logmu)))*z[indid], dim=c(max(indid),max(repid), max(trapid))) #logmu.array[as.matrix(captures[,c(2,3,1)])] <- array(log(1-exp(-exp(logmu)))*z[indid], dim=c(max(indid),max(repid), max(trapid)))[as.matrix(captures[,c(2,3,1)])] #logmu.array[is.infinite(logmu.array)]=.Machine$double.min.exp mu.array <- array(1-mu, dim=c(max(indid),max(repid), max(trapid))) mu.array[as.matrix(captures[,c(2,3,1)])] <- 1 - mu.array[as.matrix(captures[,c(2,3,1)])] logmu.array = log(mu.array) logmu.array[is.infinite(logmu.array)]=.Machine$double.min.exp which(is.infinite(logmu.array),arr.ind=TRUE) newy<-rbinom(length(mu),1,mu) gof.stats<-cbind(y,newy,mu) gof.stats<-aggregate(gof.stats,list(indid),sum) gof.data[m]<- sum( (sqrt(gof.stats[,2])-sqrt(gof.stats[,4]))[z==1]^2) gof.new[m]<- sum( (sqrt(gof.stats[,3])-sqrt(gof.stats[,4]))[z==1]^2) ##### ##### ##### ##### density<- sum(z)/totalarea zout[m,]<-z Sout[m,]<- centers # Compute the likelihood of the data given the model/parameter values LLout[m]<- sum(logmu.array)#prod(c(mu.array)) out[m,]<-c(bsigmatmp[1],sigmatmp[1],bsigmatmp[2],sigmatmp[2], lam0, beta.behave, beta1,beta.sex,psi,psi.sex,sum(z),theta,beta.den,density) print(out[m,]) if(m%%dumprate==0){ #write a file here not implemented yet } m<-m+1 } } #### vector of model effects parms.2.report<- c(TRUE,TRUE,Msexsigma==1,Msexsigma==1, TRUE, Mb==1, Xeff.tf, Msex==1, TRUE, (Msex==1|Msexsigma==1), TRUE, update.theta, sum(Xd)>0, TRUE ) # JFG 10/15/13 edited so that psi.sex is included if Msexsigma = 1. #model # ni=1100,burn=100,skip=2,nz=200,theta=NA, # Msigma=1,Mb=0,Msex=0,Msexsigma = 0,Xeff=NULL,Xsex=NULL, out<- list(mcmchist=out,G=G,Gunscaled=Gunscaled,traplocs=traplocs,Sout=Sout,zout=zout, likelihood=LLout,statespace=statespace,gof.data=gof.data,gof.new=gof.new,call=call,parms2report=parms.2.report) # JFG edited to return the likelihood class(out) <- c("scrfit","list") return(out) } Model fitting: traps.x = c(217416.10928451, 222416.10928451, 227416.10928451, 217416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 237416.10928451, 217416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 212416.10928451, 217416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 217416.10928451, 222416.10928451, 227416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 227416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 272416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 272416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 272416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451, 232416.10928451, 237416.10928451, 242416.10928451, 247416.10928451, 252416.10928451, 257416.10928451, 262416.10928451, 267416.10928451) traps.y = c(153560.81238841, 153560.81238841, 153560.81238841, 158560.81238841, 158560.81238841, 158560.81238841, 158560.81238842, 163560.81238841, 163560.81238841, 163560.81238841, 168560.81238841, 168560.81238841, 168560.81238841, 168560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 173560.81238841, 178560.81238841, 178560.81238841, 178560.81238842, 178560.81238842, 178560.81238841, 178560.81238841, 178560.81238841, 178560.81238841, 178560.81238841, 178560.81238842, 178560.81238842, 183560.81238842, 183560.81238841, 183560.81238841, 183560.81238841, 183560.81238841, 183560.81238841, 183560.81238841, 183560.81238842, 183560.81238842, 183560.81238841, 183560.81238841, 188560.81238841, 188560.81238841, 188560.81238841, 188560.81238842, 188560.81238842, 188560.81238841, 188560.81238841, 188560.81238841, 188560.81238842, 188560.81238842, 193560.81238841, 193560.81238841, 193560.81238842, 193560.81238842, 193560.81238841, 193560.81238841, 193560.81238841, 193560.81238842, 193560.81238842, 198560.81238841, 198560.81238841, 198560.81238841, 198560.81238841, 198560.81238842, 198560.81238842, 198560.81238841, 198560.81238841, 198560.81238841, 203560.81238841, 203560.81238841, 203560.81238841, 203560.81238842, 203560.81238842, 203560.81238842, 203560.81238841, 203560.81238841, 208560.81238841, 208560.81238841, 208560.81238842, 208560.81238842, 208560.81238841, 208560.81238841, 208560.81238841, 208560.81238841, 213560.81238841, 213560.81238841, 213560.81238842, 213560.81238842, 213560.81238841, 213560.81238841, 213560.81238841, 213560.81238841, 218560.81238841, 218560.81238841, 218560.81238842, 218560.81238842, 218560.81238841, 218560.81238841, 218560.81238841, 218560.81238841) #Create trapping grid with all traps open in all time periods. traps.data<-data.frame(grid_id=1:105, x=traps.x, y=traps.y, V4=matrix(1,nrow=105,ncol=4)) #Create state space, buffer= buffer + 1/2 grid in meters, entire study area considered habitat ss.animal<-make.ss(traps.data, 12500,2000)$ss.animal ss.animal$habitat = 1 ss.animal = ss.animal[c(147, 148, 149, 150, 151, 152, 159, 160, 194, 195, 196, 197, 198, 199, 203, 204, 205, 206, 207, 239, 240, 241, 242, 243, 244, 245, 246, 247, 249, 250, 251, 252, 253, 254, 255, 285, 286, 287, 288, 289, 290, 291, 292, 293, 294, 295, 296, 297, 298, 299, 300, 301, 302, 312, 317, 318, 330, 331, 332, 333, 334, 335, 336, 337, 338, 339, 340, 341, 342, 343, 344, 345, 346, 347, 348, 356, 357, 358, 359, 362, 363, 364, 365, 374, 375, 376, 377, 378, 379, 380, 381, 382, 383, 384, 385, 386, 387, 388, 389, 390, 391, 392, 393, 394, 401, 402, 403, 404, 405, 406, 407, 408, 409, 410, 411, 420, 421, 422, 423, 424, 425, 426, 427, 428, 429, 430, 431, 432, 433, 434, 435, 436, 437, 438, 439, 440, 444, 445, 447, 448, 449, 450, 451, 452, 453, 454, 455, 456, 457, 458, 459, 460, 464, 465, 466, 467, 468, 469, 470, 471, 472, 473, 474, 475, 476, 477, 478, 479, 480, 481, 482, 483, 484, 485, 486, 487, 488, 489, 490, 491, 492, 493, 494, 495, 496, 497, 498, 499, 500, 501, 502, 503, 504, 505, 506, 510, 511, 512, 513, 514, 515, 516, 517, 518, 519, 520, 521, 522, 523, 524, 525, 526, 527, 528, 529, 530, 531, 532, 533, 534, 535, 536, 537, 538, 539, 540, 541, 542, 543, 544, 545, 546, 547, 548, 549, 550, 551, 552, 556, 557, 558, 559, 560, 561, 562, 563, 564, 565, 566, 567, 568, 569, 570, 571, 572, 573, 574, 575, 576, 577, 578, 579, 580, 581, 582, 583, 584, 585, 586, 587, 588, 589, 590, 591, 592, 593, 594, 595, 596, 597, 598, 603, 604, 605, 606, 607, 608, 609, 610, 611, 612, 613, 614, 615, 616, 617, 618, 619, 620, 621, 622, 623, 624, 625, 626, 627, 628, 629, 630, 631, 632, 633, 634, 635, 636, 637, 638, 639, 640, 641, 642, 643, 644, 649, 650, 651, 652, 653, 654, 655, 656, 657, 658, 659, 660, 661, 662, 663, 664, 665, 666, 667, 668, 669, 670, 671, 672, 673, 674, 675, 676, 677, 678, 679, 680, 681, 682, 683, 684, 685, 686, 687, 688, 689, 690, 697, 698, 699, 700, 701, 702, 703, 704, 705, 706, 707, 708, 709, 710, 711, 712, 713, 714, 715, 716, 717, 718, 719, 720, 721, 722, 723, 724, 725, 726, 727, 728, 729, 730, 731, 732, 733, 734, 735, 736, 744, 745, 746, 747, 748, 749, 750, 751, 752, 753, 754, 755, 756, 757, 758, 759, 760, 761, 762, 763, 764, 765, 766, 767, 768, 769, 770, 771, 772, 773, 774, 775, 776, 777, 778, 779, 780, 781, 782, 789, 790, 791, 792, 793, 794, 795, 796, 797, 798, 799, 800, 801, 802, 803, 804, 805, 806, 807, 808, 809, 810, 811, 812, 813, 814, 815, 816, 817, 818, 819, 820, 821, 822, 823, 824, 825, 826, 827, 828, 836, 837, 838, 839, 840, 841, 842, 843, 844, 845, 846, 847, 848, 849, 850, 851, 852, 853, 854, 855, 856, 857, 858, 859, 860, 861, 862, 863, 864, 865, 866, 867, 868, 869, 870, 871, 872, 873, 874, 882, 883, 884, 885, 886, 887, 888, 889, 890, 891, 892, 893, 894, 895, 896, 897, 898, 899, 900, 901, 902, 903, 904, 905, 906, 907, 908, 909, 910, 911, 912, 913, 914, 915, 916, 917, 918, 919, 920, 929, 930, 931, 932, 933, 934, 935, 936, 937, 938, 939, 940, 941, 942, 943, 944, 945, 946, 947, 948, 949, 950, 951, 952, 953, 954, 955, 956, 957, 958, 959, 960, 961, 962, 963, 964, 965, 966, 976, 977, 978, 979, 980, 981, 982, 983, 984, 985, 986, 987, 988, 989, 990, 991, 992, 993, 994, 995, 996, 997, 998, 999, 1000, 1001, 1002, 1003, 1004, 1005, 1006, 1007, 1008, 1009, 1010, 1011, 1012, 1024, 1025, 1026, 1027, 1028, 1029, 1030, 1031, 1032, 1033, 1034, 1035, 1036, 1037, 1038, 1039, 1040, 1041, 1042, 1043, 1044, 1045, 1046, 1047, 1048, 1049, 1050, 1051, 1052, 1053, 1054, 1055, 1056, 1057, 1058, 1071, 1072, 1073, 1074, 1075, 1076, 1077, 1078, 1079, 1080, 1081, 1082, 1083, 1084, 1085, 1086, 1087, 1088, 1089, 1090, 1091, 1092, 1093, 1094, 1095, 1096, 1097, 1098, 1099, 1100, 1101, 1102, 1103, 1104, 1117, 1118, 1119, 1120, 1121, 1122, 1123, 1124, 1125, 1126, 1127, 1128, 1129, 1130, 1131, 1132, 1133, 1134, 1135, 1136, 1137, 1138, 1139, 1140, 1141, 1142, 1143, 1144, 1145, 1146, 1147, 1148, 1149, 1150, 1166, 1167, 1168, 1169, 1170, 1171, 1172, 1173, 1174, 1175, 1176, 1177, 1178, 1179, 1180, 1181, 1182, 1183, 1184, 1185, 1186, 1187, 1188, 1189, 1190, 1191, 1192, 1193, 1194, 1195, 1196, 1212, 1213, 1214, 1215, 1216, 1217, 1218, 1219, 1220, 1221, 1222, 1223, 1224, 1225, 1226, 1227, 1228, 1229, 1230, 1231, 1232, 1233, 1234, 1235, 1236, 1237, 1238, 1239, 1240, 1241, 1242, 1257, 1258, 1259, 1260, 1261, 1262, 1263, 1264, 1265, 1266, 1267, 1268, 1269, 1270, 1271, 1272, 1273, 1274, 1275, 1276, 1277, 1278, 1279, 1280, 1281, 1282, 1283, 1284, 1285, 1286, 1287, 1288, 1301, 1302, 1303, 1304, 1305, 1306, 1307, 1308, 1309, 1310, 1311, 1312, 1313, 1314, 1315, 1316, 1317, 1318, 1319, 1320, 1321, 1322, 1323, 1324, 1325, 1326, 1327, 1328, 1329, 1330, 1331, 1332, 1333, 1334, 1345, 1346, 1347, 1348, 1349, 1350, 1351, 1352, 1353, 1354, 1355, 1356, 1357, 1358, 1359, 1360, 1361, 1362, 1363, 1364, 1365, 1366, 1367, 1368, 1369, 1370, 1371, 1372, 1373, 1374, 1375, 1376, 1377, 1378, 1379, 1380, 1391, 1392, 1393, 1394, 1395, 1396, 1397, 1398, 1399, 1400, 1401, 1402, 1403, 1404, 1405, 1406, 1407, 1408, 1409, 1410, 1411, 1412, 1413, 1414, 1415, 1416, 1417, 1418, 1419, 1420, 1421, 1422, 1423, 1424, 1425, 1426, 1432, 1433, 1437, 1438, 1439, 1440, 1441, 1442, 1443, 1444, 1445, 1446, 1447, 1448, 1449, 1450, 1451, 1452, 1453, 1454, 1455, 1456, 1457, 1458, 1459, 1460, 1461, 1462, 1463, 1464, 1465, 1466, 1467, 1468, 1469, 1470, 1471, 1472, 1478, 1479, 1480, 1481, 1482, 1483, 1484, 1485, 1486, 1487, 1488, 1489, 1490, 1491, 1492, 1493, 1494, 1495, 1496, 1497, 1498, 1499, 1500, 1501, 1502, 1503, 1504, 1505, 1506, 1507, 1508, 1509, 1510, 1511, 1512, 1513, 1514, 1515, 1516, 1517, 1518, 1519, 1520, 1521, 1523, 1524, 1525, 1526, 1527, 1528, 1529, 1530, 1531, 1532, 1533, 1534, 1535, 1536, 1537, 1538, 1539, 1540, 1541, 1542, 1543, 1544, 1545, 1546, 1547, 1548, 1549, 1550, 1551, 1552, 1553, 1554, 1555, 1556, 1557, 1558, 1559, 1560, 1561, 1562, 1563, 1564, 1565, 1566, 1567, 1568, 1569, 1570, 1571, 1572, 1573, 1574, 1575, 1576, 1577, 1578, 1579, 1580, 1581, 1582, 1583, 1584, 1585, 1586, 1587, 1588, 1589, 1590, 1591, 1592, 1593, 1594, 1595, 1596, 1597, 1598, 1599, 1600, 1601, 1602, 1603, 1604, 1605, 1606, 1607, 1608, 1609, 1610, 1611, 1612, 1613, 1614, 1615, 1616, 1617, 1618, 1619, 1620, 1621, 1622, 1623, 1624, 1625, 1626, 1627, 1628, 1629, 1630, 1631, 1632, 1633, 1634, 1635, 1636, 1637, 1638, 1639, 1640, 1641, 1642, 1643, 1644, 1645, 1646, 1647, 1648, 1649, 1650, 1651, 1652, 1653, 1654, 1655, 1656, 1657, 1658, 1659, 1660, 1661, 1662, 1663, 1664, 1665, 1666, 1667, 1668, 1669, 1670, 1671, 1672, 1673, 1674, 1675, 1676, 1677, 1678, 1679, 1680, 1681, 1682, 1683, 1684, 1685, 1686, 1687, 1688, 1689, 1690, 1691, 1692, 1693, 1694, 1695, 1696, 1697, 1698, 1699, 1700, 1701, 1702, 1703, 1704, 1705, 1706, 1707, 1708, 1709, 1710, 1711, 1712, 1713, 1714, 1715, 1716, 1717, 1718, 1719, 1720, 1721, 1722, 1723, 1724, 1725, 1726, 1727, 1728, 1729, 1730, 1731, 1732, 1733, 1734, 1735, 1736, 1737, 1738, 1739, 1740, 1741, 1742, 1743, 1744, 1745, 1746, 1747, 1748, 1749, 1750, 1751, 1752, 1753, 1754, 1755, 1756, 1757, 1758, 1759, 1760, 1761, 1762, 1763, 1764, 1765, 1766, 1767, 1768, 1769, 1770, 1771, 1772, 1773, 1774, 1775, 1776, 1777, 1778, 1779, 1780, 1781, 1782, 1783, 1784, 1785, 1786, 1787, 1788, 1789, 1790, 1791, 1792, 1793, 1794, 1795, 1796, 1797, 1798, 1799, 1800, 1801, 1802, 1803, 1804, 1805, 1806, 1807, 1808, 1809, 1810, 1811, 1812, 1813, 1814, 1815, 1816, 1817, 1818, 1819, 1820, 1821, 1822, 1823, 1824, 1825, 1826, 1827, 1828, 1829, 1830, 1831, 1832, 1833, 1834, 1835, 1836, 1837, 1838, 1839, 1840, 1841, 1842, 1843, 1844, 1845, 1846, 1847, 1848, 1849, 1850, 1851, 1852, 1853, 1854, 1855, 1856, 1857, 1858, 1859, 1860, 1861, 1862, 1863, 1864, 1865, 1866, 1867, 1868, 1869, 1870, 1871, 1872, 1873, 1874, 1875, 1876, 1877, 1878, 1879, 1880, 1881, 1882, 1883, 1884, 1885, 1886, 1887, 1888, 1889, 1890, 1891, 1892, 1893, 1894, 1895, 1896, 1897, 1898, 1899, 1900, 1901, 1902, 1903, 1904, 1905, 1906, 1907, 1908, 1909, 1910, 1911, 1912, 1913, 1914, 1915, 1916, 1917, 1918, 1919, 1920, 1921, 1922, 1923, 1924, 1925, 1926, 1927, 1928, 1929, 1930, 1931, 1932, 1933, 1934, 1935, 1936, 1937, 1938, 1939, 1940, 1941, 1942, 1943, 1944, 1945, 1946, 1947, 1948, 1949, 1950, 1951, 1952, 1953, 1954, 1955, 1956, 1957, 1958, 1959, 1960, 1961, 1962, 1963, 1964, 1965, 1966, 1967, 1968, 1969, 1970, 1971, 1972, 1973, 1974, 1975, 1976, 1977, 1978),] ss.rsf = c(0.00051837, 0.002439, 0.0007725, 0.00982532, 0.00188936, 0.00722496, 0.02093748, 0.0683585, 0.00494246, 0.00477991, 0.01257677, 0.01377301, 0.01297616, 0.05618855, 0.08489683, 0.03425201, 0.10684451, 0.00933296, 0.00192994, 0.0015077, 0.00342773, 0.00946239, 0.00422736, 0.01067345, 0.10886528, 0.00129805, 0.17532538, 0.08282711, 0.03713504, 0.16266094, 0.03150069, 0.01258559, 0.00771276, 0.0003596, 0.03336592, 0.00560588, 0.01464138, 0.00231042, 0.00106127, 0.00638394, 0.11116429, 0.03202997, 0.09493393, 0.05373818, 0.01500099, 0.19816311, 0.17920354, 0.05953353, 0.03831936, 0.0005269, 0.00319136, 0.00774046, 0.00381688, 0.03234502, 0.00809918, 0.00231309, 0.01385049, 0.00024425, 0.01260602, 0.01684313, 0.00284751, 0.00012306, 0.00055104, 0.05475359, 0.1450506, 0.06539202, 0.00930316, 0.25532186, 0.26189253, 0.03646847, 0.05966365, 0.29076564, 0.00407024, 0.00027964, 0.01829401, 0.00126835, 0.03178028, 0.01785823, 0.01653123, 0.12497701, 0.00407442, 0.07600817, 0.0722881, 0.01693118, 0.02522674, 0.00488651, 0.00484634, 7.749e-05, 0.03389356, 0.0019124, 0.00064859, 6.866e-05, 0.01792942, 0.02645094, 0.15779705, 0.05989114, 0.1875805, 0.2381698, 0.28758672, 0.08864473, 0.00029737, 0.00776061, 0.00616372, 5.399e-05, 0.00638243, 0.03805692, 0.03840059, 0.01251155, 0.00458311, 0.00030468, 0.02522291, 0.00047852, 0.03264622, 0.18983965, 0.00249177, 0.00410083, 0.0040837, 0.0549717, 0.05374086, 0.0191049, 0.08682147, 0.00364726, 0.01603135, 0.01105204, 0.00691787, 0.04642892, 0.40217087, 0.00219798, 0.04606065, 0.5164578, 0.23537523, 0.05260616, 0.07961183, 0.00038173, 0.01755701, 0.00085494, 0.00013072, 0.00254843, 0.01194531, 0.09634438, 0.08116564, 0.02531559, 0.0043448, 0.00399964, 0.03316554, 0.00295903, 0.07010128, 0.24671568, 0.00716672, 0.01862798, 0.09389009, 0.01005894, 0.00654108, 0.00017046, 0.00178282, 0.0919202, 0.00316829, 0.13595727, 0.4028123, 0.00241184, 0.01165937, 0.00568923, 0.09335816, 0.00352651, 0.00423098, 0.5520342, 0.04756167, 0.23794669, 0.15871897, 0.28164956, 0.36248216, 0.01686679, 0.0695897, 1.416e-05, 0.00075993, 0.02902362, 5.753e-05, 0.02790798, 0.06400511, 0.00013674, 0.00024427, 0.02264175, 0.0003521, 0.00966394, 0.19543175, 0.00099746, 0.00082647, 0.02493562, 0.00105458, 0.17944403, 0.00767188, 0.04630907, 0.00329996, 0.00025006, 0.01882491, 0.00363853, 0.01276501, 0.01860438, 0.26611674, 0.1011148, 0.21013588, 0.04843203, 0.1899333, 0.10726625, 0.00249138, 0.26701576, 0.14935936, 0.06581633, 0.58310795, 0.24206303, 0.53714049, 0.490051, 0.29856661, 0.53439826, 0.13214129, 0.00015424, 0.00109202, 0.01027097, 0.00538442, 0.01343166, 0.0016132, 0.00268828, 0.15551706, 0.00182742, 0.00592641, 0.00014441, 0.00079113, 0.59094417, 0.00148027, 0.00014474, 0.22892582, 0.00318985, 0.52397716, 0.05252479, 0.00064694, 0.02794711, 0.16612564, 0.01639798, 0.00027456, 0.00135967, 0.02428773, 0.38749257, 0.38248172, 0.4828375, 0.15927719, 0.22676654, 0.02977854, 0.24564289, 0.23457864, 0.27192256, 0.3656795, 0.48368037, 0.57878041, 0.09579825, 0.41062152, 0.57671881, 0.68443853, 0.02785188, 0.00705036, 1.8e-07, 0.01237652, 0.00183097, 0.00858618, 0.01250069, 1.69e-06, 0.00121177, 0.22172022, 0.0055714, 0.0064593, 0.00043537, 0.67528898, 0.01440276, 0.00115747, 0.16125694, 0.01629554, 0.57574683, 0.00653235, 0.00191184, 0.00391056, 0.39901429, 0.36006299, 0.00326621, 0.01227972, 0.03923592, 0.16614343, 0.08688933, 0.18015087, 0.14338176, 0.28088108, 0.32490629, 0.07439461, 0.0502125, 0.4659532, 0.26903161, 0.32523328, 0.28303674, 0.21102376, 0.50895625, 0.24755454, 0.14785637, 0.01052557, 0.00010894, 0.1158487, 1.928e-05, 0.00322833, 0.00196946, 3.731e-05, 1.77e-06, 0.34006488, 0.02045398, 0.0003265, 0.02136431, 0.59888113, 0.00580173, 0.00436943, 0.39525324, 0.41343412, 0.006214, 0.04651783, 5.652e-05, 0.03301104, 0.07249092, 0.37128943, 0.00104292, 0.00314388, 0.01431255, 0.0064971, 0.00659814, 0.00344528, 0.19511709, 0.1110094, 0.0229234, 0.10987342, 0.51015854, 0.24228434, 0.66311193, 0.48046312, 0.38107541, 0.36943355, 0.64794701, 0.51364887, 0.62199616, 0.06436446, 0.17523985, 0.02332958, 1.97e-06, 0.00343012, 0.01420863, 0.00011917, 0.00734235, 0.00399406, 0.57653517, 0.00050643, 0.07376879, 0.2343642, 0.00256185, 0.03338887, 0.17188449, 0.62224102, 0.00803295, 0.0105026, 0.0019362, 0.09211858, 0.01046318, 0.22485372, 0.00014439, 0.00068581, 0.01366691, 0.05675692, 0.30725637, 0.03258862, 0.00294969, 0.16155565, 0.04673055, 0.35110945, 0.25485125, 0.15936878, 0.07474511, 0.38601995, 0.62750095, 0.58751917, 0.52297515, 0.42609754, 0.14128946, 0.00294597, 1.975e-05, 0.00392142, 0.04482024, 0.00189545, 0.02876142, 5.807e-05, 0.67283529, 0.53284907, 0.1293893, 0.07329269, 0.03482775, 0.20789868, 0.78339499, 0.55983526, 0.00097728, 0.00023596, 1.868e-05, 0.40392569, 0.14959808, 0.00717241, 0.01129587, 0.03066332, 0.08215461, 0.08707474, 0.35048547, 0.01574548, 0.2962622, 0.00373275, 0.01253806, 0.31449401, 0.13341719, 0.02904445, 0.26592991, 0.21247531, 0.39327294, 0.7207489, 0.37535065, 0.41229379, 0.19448458, 0.03527429, 0.04165455, 0.34010229, 0.00471307, 0.14530738, 0.05430017, 0.01157587, 0.33601511, 0.25365183, 0.3961384, 0.4862054, 0.56816888, 0.32443428, 0.15185493, 0.05089233, 0.00206592, 0.00011649, 0.59137082, 0.35337877, 0.00072477, 0.11145443, 0.0083195, 0.06423467, 0.02296656, 0.07767079, 0.07032613, 0.00109319, 0.32095051, 0.1038546, 0.17219327, 0.00178753, 0.01035311, 0.00114371, 0.033267, 0.03813501, 0.23436227, 0.1814748, 0.59349579, 0.61916721, 0.63244861, 0.43746561, 0.34004769, 0.43549329, 0.27582544, 0.43457878, 0.30629578, 0.6604501, 0.28021693, 0.38090974, 0.4503141, 0.29839152, 0.61921489, 0.28659898, 0.46081221, 0.46498084, 0.28857189, 0.05356539, 0.36444858, 0.46832874, 0.0598658, 0.06354035, 0.13067289, 0.34127256, 0.35773325, 0.10397971, 0.02108067, 0.10280392, 0.02953533, 0.00151822, 0.00436264, 0.00061064, 0.0127612, 0.00844276, 0.01099218, 0.00648713, 0.0342974, 0.35773492, 0.23578493, 0.15864788, 0.55545259, 0.40297186, 0.40217051, 0.28260505, 0.45816326, 0.41518942, 0.04364638, 0.02497335, 0.24536288, 0.23197113, 0.37829942, 0.31933871, 0.33813781, 0.25107384, 0.23524222, 0.44671303, 0.33354655, 0.46987355, 0.62141073, 0.42798072, 0.01529346, 0.01992041, 0.0446452, 0.15001401, 0.03179132, 0.00718612, 0.18581195, 0.00694876, 0.00511268, 0.0034352, 0.00064527, 6.768e-05, 0.00439249, 0.25277922, 0.4237099, 0.24289045, 0.16137549, 0.35306129, 0.50996554, 0.3943058, 0.22622588, 0.15972453, 0.42286244, 0.22521809, 0.32828611, 0.19010344, 0.30263218, 0.09388892, 0.00899358, 0.02214231, 0.11264347, 0.11378901, 0.1150786, 0.01519787, 0.00726789, 0.20509537, 0.18536842, 0.03016328, 0.07955509, 0.00595266, 0.00329946, 0.04950253, 0.10969663, 0.00099748, 0.00035859, 0.00927942, 0.01591198, 0.00500895, 0.00102485, 0.03714934, 0.00227313, 0.032241, 0.08424335, 0.10002245, 0.0624181, 0.1772895, 0.19421247, 0.13955307, 0.60709202, 0.35642022, 0.07302589, 0.14297314, 0.1753957, 0.19219594, 0.13697855, 0.23859563, 0.62250006, 0.6292811, 0.35074759, 0.58700746, 0.2279363, 0.14730711, 0.1291972, 0.12356707, 0.08063962, 0.12034382, 0.11541807, 0.00078819, 0.00035737, 0.00243089, 0.00620542, 0.07597866, 0.24804522, 0.07942418, 0.00265157, 0.00884698, 0.10301261, 0.10800795, 9.797e-05, 0.00063017, 0.00900309, 0.02905369, 0.07256779, 0.07573431, 0.38148594, 0.32644993, 0.26826647, 0.00565511, 0.1829246, 0.12078883, 0.40281668, 0.4822228, 0.09349564, 0.54326892, 0.44623953, 0.42206639, 0.50731087, 0.24674952, 0.08811926, 0.04469164, 0.09852679, 0.0009022, 0.00394577, 0.0009279, 7.28e-06, 0.0089742, 0.00607559, 0.0502856, 0.00425242, 0.00733682, 0.02034996, 0.00947734, 0.11497666, 0.11611257, 0.04800113, 0.09803711, 0.27950794, 0.31305736, 0.32771447, 0.50716585, 0.12116278, 0.31450856, 0.2420862, 0.60398859, 0.13283984, 0.48676684, 0.43833059, 0.48210463, 0.35220146, 0.16399701, 0.47402796, 0.67290246, 0.53980893, 0.2140926, 0.0081316, 0.00420095, 0.14545411, 0.00230563, 0.04827971, 6.042e-05, 0.00801454, 0.00736888, 0.00501986, 0.00346068, 0.00429833, 0.18615369, 0.06689611, 0.06677036, 0.10989636, 0.57041132, 0.10010898, 0.4095211, 0.25550875, 0.33856335, 0.39039302, 0.43389046, 0.31485695, 0.41323107, 0.01527414, 0.19019346, 0.19715965, 0.23901303, 0.07048763, 0.29048944, 0.40209037, 0.63278204, 0.22348914, 0.19758056, 0.03522246, 0.22754815, 0.05192252, 0.01939473, 2.93e-06, 0.00063103, 0.00281222, 0.0012697, 0.01476145, 0.01235264, 0.04998837, 0.10692106, 0.00440892, 0.00317288, 0.0126139, 0.13730232, 0.58722335, 0.33227172, 0.23131269, 0.27748075, 0.4531132, 0.41082793, 0.21615732, 0.39503399, 0.05459786, 0.49504495, 0.01713016, 0.05289591, 0.02040648, 0.45581287, 0.40204743, 0.3457875, 0.7192179, 0.25393364, 0.60863698, 0.32967046, 0.19541363, 0.4142653, 0.00169988, 0.00049405, 0.00139842, 0.0103935, 0.00805256, 0.00204634, 0.01525829, 0.00417099, 0.49372581, 0.52872318, 0.54453498, 0.16443953, 0.10605877, 0.09167064, 0.17716533, 0.38009685, 0.58460838, 0.22193104, 0.6761384, 0.44576329, 0.2735064, 0.19633043, 0.08383237, 0.22815916, 0.28506428, 0.44346419, 0.28842995, 0.59366739, 0.03281879, 0.12229277, 0.08444994, 0.00144954, 0.00048261, 0.00951155, 0.01583645, 0.09401568, 0.1014349, 0.24966459, 0.29307476, 0.38185677, 0.41026607, 0.2854169, 0.35012543, 0.36791363, 0.27671242, 0.31473938, 0.10586109, 0.13662858, 0.03954377, 0.16564262, 0.08437948, 0.19116783, 0.11732311, 0.05682075, 0.23908637, 0.96113318, 0.25277236, 0.49891892, 0.48015833, 0.20695388, 0.28121397, 0.25083604, 0.23269804, 0.06274977, 0.01596376, 0.018668, 0.11613376, 0.28223711, 0.55329704, 0.38496995, 0.18685478, 0.16296792, 0.01048083, 0.01707985, 0.12625405, 0.04697806, 0.1862691, 0.48331112, 0.22472344, 0.35174614, 0.07578761, 0.0632691, 0.04685754, 0.06514346, 0.16209877, 0.46965012, 0.22511195, 0.50376368, 0.05060122, 0.65663773, 0.42001373, 0.21406189, 0.48952129, 0.40664849, 0.32934472, 0.13707264, 0.00335513, 0.01391848, 0.00477599, 0.09564293, 0.0620766, 0.07828538, 0.03886888, 0.41603601, 0.16095221, 0.16547997, 0.01794108, 0.01938532, 0.1464514, 0.10331356, 0.24680895, 0.5316602, 0.23270842, 0.11317386, 0.22024608, 0.05675948, 0.01191306, 0.00196164, 0.32270548, 0.03719775, 0.01479611, 0.28264952, 0.00642169, 0.46389952, 0.36107296, 0.35725966, 0.37923548, 0.5166204, 0.7502073, 0.00714276, 0.0125056, 0.00343323, 0.01522753, 0.01548862, 0.01069178, 0.00243785, 0.00518907, 0.01074473, 0.02766695, 0.15718792, 0.14910427, 0.50312978, 0.17427346, 0.1839951, 0.07203588, 0.22170104, 0.47532332, 0.18990181, 0.1765305, 0.10783837, 0.20126756, 0.03069184, 0.00192554, 0.29908124, 0.08973527, 0.04155784, 0.06076574, 0.0528445, 0.13736549, 0.28589422, 0.55859822, 0.28644699, 0.36812788, 0.0240329, 0.28439063, 0.10114127, 0.01287542, 0.02468057, 0.03918234, 0.0176722, 0.04489787, 0.02375951, 0.01568362, 0.00297496, 0.00188131, 0.00533253, 0.09653932, 0.24406262, 0.34918219, 0.45539209, 0.12566592, 0.21652393, 0.11863341, 0.09463897, 0.12838502, 0.06856698, 0.1522377, 0.0441996, 0.09588458, 0.01467962, 0.04805065, 0.04524963, 0.02085913, 0.00205178, 0.04974868, 0.0960291, 0.30771416, 0.07032807, 0.00221404, 0.02622068, 0.03690849, 0.00318223, 0.07953928, 0.00232473, 0.01351681, 0.00915872, 0.07736257, 0.04381923, 0.0202934, 0.00534093, 0.00267488, 0.00125153, 0.01264697, 0.08270916, 0.22022697, 0.01149347, 0.03562329, 0.52898151, 0.43436095, 0.32633719, 0.06633555, 0.009545, 0.06644872, 0.0474646, 0.04198, 0.01621062, 0.03947908, 0.01444426, 0.03440441, 0.02508589, 0.09421494, 0.0347643, 0.10676348, 0.65539265, 0.3001093, 0.61690319, 0.00764088, 0.01641656, 0.0007436, 0.00101684, 0.01600028, 0.01350147, 0.04845988, 0.01594205, 0.02587838, 0.00530368, 0.03554783, 0.07465969, 0.05802345, 0.01755941, 0.00342191, 0.00190626, 0.0009388, 0.15634128, 0.02309357, 0.02021657, 0.04897255, 0.11360943, 0.1283358, 0.48394084, 0.26510182, 0.10259545, 0.17778072, 0.01811686, 0.01627871, 0.03741498, 0.00831906, 0.01805909, 0.01033184, 0.00896127, 0.02560114, 0.05440479, 0.31046405, 0.28837034, 0.22853026, 0.19794178, 0.28250638, 0.40134937, 0.02476933, 0.01096141, 1.74e-06, 3.32e-06, 8.7e-06, 0.04363589, 0.00889725, 0.01170753, 0.02108573, 0.0206585, 0.00390028, 0.001687, 0.01871311, 0.05793261, 0.03495154, 0.02112458, 0.00630004, 0.00052214, 0.00651057, 0.00310923, 0.00121098, 0.00126056, 0.00266441, 0.06556749, 0.0720272, 0.1373495, 0.0425417, 0.40061009, 0.30952734, 0.1592018, 0.00901301, 0.04312877, 0.01309928, 0.01328351, 0.00072813, 0.00229221, 0.00776999, 0.00122279, 0.0690735, 0.06790443, 0.26487923, 0.07896442, 0.0527431, 0.0189077, 0.04986068, 0.03985619, 0.26498541, 0.00025943, 3.39e-06, 5.391e-05, 3.35e-06, 0.00444011, 0.02218192, 0.01672925, 0.0192097, 0.01114734, 0.00195905, 0.03152582, 0.00143268, 0.05305094, 0.00482321, 0.00043609, 0.00086901, 0.00519172, 0.00610492, 0.00034096, 0.00187562, 0.00017255, 0.00353173, 0.00275077, 0.02267857, 0.04713988, 0.1871983, 0.56789047, 0.23902343, 0.1132102, 0.33175316, 0.07402591, 0.01456622, 0.00292405, 0.00012068, 0.0002391, 6.074e-05, 0.00704391, 0.0672012, 0.05703699, 0.0330767, 0.00457497, 0.01247387, 0.00688501, 0.00699177, 0.02219804, 0.02684083, 0.00030009, 7.021e-05, 0.00223949, 8.806e-05, 0.00623299, 0.01634803, 0.03372898, 0.05304797, 0.02647677, 0.00415373, 0.05686059, 0.00912895, 0.05990744, 0.02725546, 0.00301108, 0.00082126, 0.02528892, 0.00573057, 0.00523498, 0.00707193, 0.00090024, 2.063e-05, 0.00715144, 0.01199026, 0.11883367, 0.21559137, 0.14890514, 0.39989871, 0.08511005, 0.27120021, 0.01932714, 0.02860415, 0.00011919, 0.00017926, 0.01577118, 0.00073013, 8.96e-06, 0.00019243, 0.00407561, 0.0007868, 0.00692403, 0.01357324, 0.03501772, 0.01893584, 0.00342706, 0.05595878, 0.00023469, 0.0051245, 0.00312469, 0.0023634, 0.00569196, 0.01797792, 0.02051551, 0.05380342, 0.02134339, 0.04325941, 0.07303385, 0.0505375, 0.01731098, 0.0180342, 0.02669911, 0.01484313, 0.00800216, 0.00432449, 0.03385837, 0.00043631, 0.00019977, 2.306e-05, 0.00045186, 0.00466469, 0.01639639, 0.08150879, 0.4270899, 0.13607731, 0.39306679, 0.09888619, 0.07944791, 0.00102573, 0.00199084, 0.00427397, 0.00016176, 0.0293911, 0.00604749, 0.00074651, 0.00552509, 0.01646228, 0.01192578, 0.03421257, 0.01652478, 0.01343166, 0.01295767, 0.00263111, 0.001321, 0.00747661, 0.00186066, 0.00263931, 0.01002125, 0.0161701, 0.02685163, 0.03122711, 0.03363749, 0.0357973, 0.01171081, 0.03175202, 0.02382034, 0.02513663, 0.07247423, 0.01195222, 0.00329627, 0.06035299, 0.00655462, 0.00659015, 1.366e-05, 0.00118536, 0.00018012, 0.00240617, 0.01802137, 0.17884699, 0.31048214, 0.20118399, 0.19688839, 0.17499346, 0.30365005, 0.10214488, 0.06057713, 0.00354225, 0.00383303, 0.00056849, 0.05329061, 0.00670427, 0.01067279, 0.00915228, 0.00156638, 0.02107044, 0.01546381, 0.01957427, 0.00683843, 0.00031586, 0.00747788, 0.00759157, 0.00313472, 0.00160422, 0.01829003, 0.02346847, 0.01289519, 0.01507582, 0.03613751, 0.0006288, 0.01224184, 0.0119888, 0.00952176, 0.06833741, 0.04213351, 0.02481868, 0.04523927, 0.01628467, 0.03412702, 0.0029772, 0.00070974, 0.00748885, 8.893e-05, 0.00534198, 0.05379217, 0.0738226, 0.33232194, 0.05205332, 0.1031408, 0.09105534, 0.07866833, 0.00750966, 0.01153013, 0.00568382, 0.01780549, 0.012556, 0.00588575, 0.09201591, 0.02818937, 0.01774773, 0.00644869, 0.02109465, 0.03981993, 0.04385727, 0.06002031, 0.01243629, 0.01289333, 0.01338403, 0.01911651, 0.01920028, 0.01411846, 0.01983492, 0.00514464, 0.01329644, 0.00078475, 0.01738276, 0.02030691, 0.03194029, 0.06958617, 0.08209866, 0.04209353, 0.05314489, 0.02331794, 0.0356955, 0.01041537, 0.00312601, 0.04598991, 0.00518782, 0.00254647, 0.0036626, 0.03323975, 0.04137051, 0.1732351, 0.04340437, 0.0049354, 0.00093796, 0.00071838, 0.0012601, 0.003029, 0.00328953, 0.00445964, 0.02353463, 0.03200113, 0.07532023, 0.01966049, 0.02131355, 0.02299212, 0.02318804, 0.02747566, 0.00659459, 0.15180896, 0.00951052, 0.02074206, 0.00918815, 0.03205645, 0.0303946, 0.01004658, 0.02900779, 0.00198184, 0.01108289, 0.03381174, 0.01006765, 0.01362867, 0.01327477, 0.02120636, 0.02378243, 0.0350005, 0.01256214, 0.00632067, 0.00807559, 0.01427796, 0.02010233, 0.04933847, 0.00201583, 0.00053546, 0.00284748, 0.02019239, 0.07620189, 0.10143526, 0.09152082, 0.04376992, 0.00183931, 0.00869186, 0.02256092, 0.01740519, 0.00553328, 0.0032313, 0.00036876, 0.00506874, 0.06036168, 0.14467335, 0.06299875, 0.01706719, 0.00516626, 0.00069836, 0.00030053, 0.05139191, 0.02682956, 0.0105174, 0.02611497, 0.0193756, 0.02970197, 0.00442629, 0.00306888, 4.602e-05, 0.00829043, 0.03600598, 0.00474385, 0.01241115, 0.01247483, 0.05543686, 0.04855109, 0.03709855, 0.04095275, 0.03015478, 0.00587251, 0.05602265, 0.04157021, 0.00231775, 0.004744, 9.228e-05, 0.00029994, 0.00061877, 0.00149206, 0.11955388, 0.01028467, 9.69e-05, 0.0026549, 0.03791618, 0.06233439, 0.06160427, 0.06397662, 0.03065981, 0.00405603, 0.02579539, 0.01008702, 0.13764885, 0.12259026, 0.0205679, 0.003442, 0.00157734, 0.01305576, 0.05414071, 0.16115214, 0.00203389, 0.03887318, 0.04657627, 0.05059889, 0.03320333, 0.00585435, 0.00389705, 0.04900789, 0.05665816, 0.05675859, 0.02077645, 0.01771012, 0.02411255, 0.04686045, 0.05992252, 0.02970779, 0.06963541, 0.01531503, 0.0731689, 0.04976667, 0.00903265, 0.00064113, 0.00165458, 0.00033849, 0.00039663, 0.00837803, 0.08101855, 0.00897081, 0.00508225, 0.02137074, 0.00579307, 0.00226102, 0.0106742, 0.07919265, 0.05517551, 0.01707546, 0.03207049, 0.07848138, 0.11634813, 0.15153126, 0.07922855, 0.14834045, 0.0544363, 0.05369232, 0.1146599, 0.1407354) #Read in capture data, set date format, assign captures to sampling period #Sampling periods, Dec (Dec 1-Dec 31), Jan, Feb, and March until end. April captures assigned to March period. captures = data.frame(“Status”= c(Dead, Dead, Dead, Dead, Dead, Live, Live, Live, Dead, Live, Live, Live, Live, Dead, Live, Live, Live, Dead, Live, Live, Live, Live, Live, Live, Live, Live, Live, Live, Live, Dead, Live, Dead, Live, Live, Dead, Dead, Dead, Live, Dead, Dead, Dead, Dead, Dead, Dead, Live, Dead, Live, Live, Live, Live, Live, Live, Dead, Live, Live, Live, Live, Live, Live, Live, Live, Live, Dead, Live, Live, Dead, Dead, Dead, Dead, Dead, Dead, Live, Live, Live, Live, Live, Live, Live, Dead, Dead), “Sex”= c(Male, Male, Male, Male, Male, Male, Female, Female, Male, Male, Male, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Female, Male, Male, Male, Male, Male, Male, Male, Female, Male, Female, Female, Male, Female, Female, Female, Female, Female, Male, Male, Female, Female, Female, Female, Male, Female, Female, Female, Male, Female, Male, Female, Female, Male, Female, Female, Female, Female, Female, Female, Female, Male, Female, Female, Female, Female, Female, Female, Male, Female, Male, Male, Female, Male)) scrCaptures = data.frame(session=rep(1,80), indid = c(1, 2, 3, 4, 5, 6, 7, 7, 8, 8, 9, 10, 11, 11, 12, 13, 14, 14, 15, 15, 15, 15, 16, 16, 17, 17, 18, 19, 20, 20, 21, 21, 22, 23, 24, 25, 26, 26, 27, 28, 29, 30, 31, 32, 33, 33, 34, 35, 35, 36, 37, 38, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 47, 47, 48, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62), SO= c(2, 1, 1, 1, 1, 2, 2, 2, 3, 2, 1, 2, 2, 3, 2, 2, 1, 3, 1, 1, 2, 3, 1, 3, 1, 3, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 2, 2, 3, 2, 2, 2, 2, 4, 2, 2, 2, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 2, 2, 3, 3, 3, 3, 3, 3, 3, 4, 3, 4, 3, 3, 4, 3, 3), trapid= c(103, 17, 9, 98, 98, 52, 12, 17, 57, 74, 32, 98, 66, 66, 74, 66, 56, 91, 37, 17, 28, 37, 77, 77, 76, 67, 9, 63, 48, 17, 66, 75, 63, 66, 98, 43, 17, 12, 102, 5, 12, 57, 98, 102, 100, 98, 76, 92, 92, 48, 67, 66, 58, 27, 13, 100, 67, 50, 63, 92, 50, 43, 32, 43, 43, 51, 9, 63, 93, 11, 48, 55, 92, 50, 62, 78, 75, 57, 78, 93)) #Create the alive matrix. Assign individuals as live or dead at the start (ie. Alive Jan 1? Alive Feb 1?) of each sampling period alive = matrix(1,nrow=max(scrCaptures$indid),ncol=max(scrCaptures$SO)) dead.col = scrCaptures$SO[captures$Status=="Dead"]+1 dead.row = scrCaptures$indid[(captures$Status=="Dead")] rep.deadr = rep(dead.row[dead.col<4],times = (5-dead.col[dead.col<4])) rep.deadc = unlist(lapply(dead.col[dead.col<4],function(x){seq(x,4,1)})) alive[cbind(rep.deadr,rep.deadc)]=0 Xsex<-(as.numeric(captures$Sex)-1)[-which(duplicated(scrCaptures$indid))] eff = data.frame(“trap”=1:105, “dec”= c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 12.176037, 4.123253, 0, 0, 0, 8.581865, 0, 0, 0, 0, 9.827748, 0, 0, 0, 0, 82.500455, 10.380116, 0, 0, 0, 15.914408, 6.81028, 2.712217, 9.563802, 6.50441, 35.541576, 9.611951, 0, 0, 0.167156, 23.10318, 27.380993, 4.759285, 22.966487, 12.413805, 6.915492, 9.519142, 25.028365, 0, 8.638221, 27.152321, 0.250346, 0, 0, 0, 20.261722, 25.036613, 34.739985, 0.2417, 33.474248, 18.497867, 21.052492, 14.231936, 2.660814, 2.64201, 8.49353, 24.667454, 3.629889, 1.299645, 24.981959, 33.435883, 28.144104, 10.602757, 0, 0, 0, 0, 0, 0, 17.686969, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), “jan” = c(2.479411, 5.346401, 0, 12.094471, 41.121339, 6.420272, 2.754368, 9.16501, 89.768888, 31.939385, 60.890072, 133.466275, 20.949944, 0, 2.05324, 69.519967, 144.442112, 5.33025, 0, 0, 0, 35.235273, 16.373733, 0, 5.185546, 38.664427, 172.528213, 22.985612, 0, 0.449955, 0, 159.023606, 32.474534, 0, 0, 27.545649, 123.458867, 59.097595, 2.77506, 19.024039, 33.078396, 271.949888, 55.61538, 0, 0, 0, 89.501959, 149.760569, 28.428335, 128.162348, 63.05252, 65.439603, 49.788398, 8.979387, 0, 11.151276, 122.116898, 22.251089, 15.865045, 21.195684, 42.798451, 197.457101, 97.281799, 59.652021, 9.568614, 208.404454, 83.277058, 82.001942, 35.532041, 27.962273, 2.043634, 0.410362, 11.253626, 32.689559, 45.423118, 60.202177, 60.086692, 47.987505, 0, 0, 0, 19.336273, 88.412034, 12.163104, 1.658244, 2.173222, 0, 0, 0, 159.939888, 89.408594, 47.575146, 8.762239, 2.357185, 0, 0, 0, 50.403225, 61.651879, 31.557358, 35.68458, 45.314474, 0, 7.306049, 0), “feb” = c(0, 4.44917, 0, 0.779841, 52.527146, 22.74438, 0, 5.807116, 123.559905, 26.649062, 33.2356, 58.868443, 11.028419, 3.128674, 0, 53.564053, 82.695568, 4.784824, 0, 0, 0, 27.954774, 0, 0, 0, 43.982829, 126.08563, 0, 0, 0, 0, 104.032325, 8.073885, 0, 0, 21.057896, 97.682692, 61.282897, 3.273549, 6.943018, 11.569427, 67.700385, 0, 0.468639, 0, 0.870467, 54.104924, 47.907733, 13.484995, 36.514056, 48.0149, 105.697379, 77.179344, 55.428426, 5.607452, 3.964262, 28.993429, 0.289605, 0, 0, 6.094131, 141.408556, 200.523656, 135.392175, 0, 67.348065, 18.661182, 14.435283, 75.500426, 28.750326, 0, 0.410362, 45.903894, 0, 22.620082, 24.745581, 70.753574, 95.00799, 46.271961, 0, 0, 0, 15.27311, 0, 5.791374, 12.098667, 0, 0, 0, 34.860559, 18.921576, 60.474506, 0, 0, 0, 0, 0, 45.156849, 56.510393, 6.066692, 14.362176, 12.645611, 0, 0, 0), “mar.apr”= c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 16.804316, 10.043035, 0, 0, 0, 9.76843, 9.93064, 0, 0, 0.912497, 6.176141, 0, 0, 11.282083, 15.968469, 21.657361, 13.914059, 0, 0, 7.545615, 19.381007, 5.515278, 17.009144, 0, 6.270923, 124.813088, 14.338875, 1.451798, 0, 38.482777, 12.549118, 14.04417, 0, 0, 0, 0, 0, 0, 0, 43.15425, 7.973873, 3.846206, 0, 0, 0, 0, 10.811539, 57.664686, 0, 0, 0, 0, 0, 18.32501, 18.910477, 42.577476, 0, 0, 0, 0, 0, 10.025186, 8.196708, 0, 0, 0, 0, 0, 0)) avg.eff = apply(eff[,2:5],2,function(x){mean(x[x!=0])}) for (c in 2:5){ eff[eff[,c]==0,c]=avg.eff[c-1] } Xeff =log(as.matrix(eff[,2:5])) scrobj<-list(traps=traps.data, captures=scrCaptures, statespace=ss.animal, alive = alive, Xd= log(ss.rsf/(1-ss.rsf))) ######################################### ###Define and run models of interest. This code is for the top model, Effort + RSF + Sex + ssex ######################################### model_example <-SCRi.fn(ni=30000, burn=10000, skip=1,nz=1000,scrobj=scrobj,Msigma=1, Mb=0, Msex=1, Msexsigma=1, Xeff=Xeff, Xsex=Xsex, area.per.pixel=4)