S3-eu-west-1.amazonaws.com
Supplementary InformationTable S1: Sample size and efforts in four study sites, namely Bhadra, Nagarahole, Bandipur and Bilgiri Ranganaswamy Temple, in Western Ghats, India.BhadraBiligiriNagaraholeBandipurSurvey duration29 Jan - 28 Feb 201318 May - 17 Jun 201308 Mar - 07 Apr 201313 Apr - 13 May 2013Sampling occasions30303030Camera trap stations12299161180Trap array polygon area (km2)571.88351.7503.77647.77Average Euclidean inter-trap distance in km 1.59 1.55 1.53 1.70 Sampling effort (number of trap-nights)3629296647605381Photo-encounters of carnivoresDhole: 56Leopard: 316Tiger: 82Dhole: 112Leopard: 53Tiger: 172Dhole: 145Leopard: 391Tiger: 326Dhole: 190Leopard: 258Tiger: 350Table S2. Model comparisons to test for independence in site-use between carnivore species pairs in the four parks. For each species-pair we compared (1) a full model where parameter ? is estimated and (2) a reduced model where ?=1. K= number of parameters. ParkSpeciesModelAICΔAICKDeviance? (SE)BhadraDhole-LeopardSIF(.)1523.1051513.11.03 (0.10)SIF=11523.20.151513.2-Dhole-TigerSIF(.)1523.1051513.11.14 (0.25)SIF=11523.20.151513.2-Leopard-TigerSIF(.)1700.41051690.411.14 (0.07)*SIF=11704.454.0451694.45-NagaraholeDhole-LeopardSIF(.)2124.73052114.731.07 (0.09)SIF=12125.410.6852115.41-Dhole-TigerSIF(.)1999.44051989.440.99 (0.07)*SIF=12125.41125.9752115.41-Leopard-TigerSIF(.)2850.72052840.721.08 (0.04)*SIF=12856.475.7552846.47-Table S2: contd.ParkSpeciesModelAICΔAICKDeviance? (SE)BandipurDhole-LeopardSIF(.)2166.63052156.631.13 (0.08)*SIF=12169.212.5852159.21-Dhole-TigerSIF(.)2456.15052446.151.03 (0.05)SIF=12456.460.3152446.46-Leopard-TigerSIF(.)2719.43052709.431.10 (0.05)*SIF=12725.145.7152715.14-BiligiriDhole-LeopardSIF(.)845.0605835.061.06 (0.17)SIF=1845.20.145835.2-Dhole-TigerSIF(.)1262.85051252.851.13 (0.08)*SIF=11265.642.7951255.64-Leopard-TigerSIF(.)1129.21051119.210.97 (0.10)SIF=11129.30.0951119.3-* indicate models where SIF ≠ 1 are statistically significant. SIF(.) implies full model with ψA(.), ψB(.), ?AB(.), pA=rA, pB=rB, δAB=1. SIF=1 implies a reduced model with ψA(.), ψB(.), ?AB=1, pA=rA, pB=rB, δAB=1. Table S3: Parameter estimates for species-wise site-use probability ψ and detectability p, and, ? for species-pairs based on multi-species occupancy models. Estimates reflect values from the best model, chosen based on AIC values. Values in parentheses are estimated Standard Errors for the parameters. The subscripts ‘D’, ‘L’ and ‘T’ for the parameters indicate dhole, leopard and tiger respectively??BhadraBiligiriNagaraholeBandipurDhole-LeopardpD0.11 (0.02)0.11 (0.02)0.14 (0.02)0.13 (0.01)?pL0.25 (0.02)0.08 (0.02)0.26 (0.01)0.18 (0.01)?ψD0.32 (0.07)0.61 (0.09)0.38 (0.05)0.60 (0.06)?ψL0.79 (0.04)0.58 (0.13)0.72 (0.04)0.62 (0.05)??DL1.03 (0.10)1.06 (0.17)1.07 (0.09)1.13 (0.08)Dhole-TigerpD0.11 (0.02)0.11 (0.02)0.14 (0.02)0.13 (0.01)?pT0.13 (0.02)0.15 (0.02)0.18 (0.01)0.18 (0.01)?ψD0.32 (0.09)0.61 (0.09)0.38 (0.05)0.60 (0.06)?ψT0.44 (0.06)0.86 (0.07)0.85 (0.05)0.83 (0.04)??DT1.14 (0.25)1.13 (0.08)0.99 (0.07)1.03 (0.05)??Leopard-TigerpL0.25 (0.02)0.08 (0.02)0.26 (0.01)0.18 (0.01)?pT0.13 (0.02)0.15 (0.02)0.18 (0.01)0.18 (0.01)?ψL0.79 (0.04)0.58 (0.13)0.72 (0.04)0.62 (0.05)?ψT0.44 (0.06)0.86 (0.07)0.85 (0.05)0.83 (0.04)??LT1.14 (0.07)0.97 (0.10)1.08 (0.04)1.10 (0.05)Figure S1. Locations of camera trap stations in (a) Bhadra, (b) Nagarahole, (c) Biligiri and (d) Bandipur for photo-capture of large carnivores.center23558500Figure S2. Probability distributions for temporal activities of the three carnivores in (a) Bhadra, (b) Nagarahole, (c) Bandipur and (d) Biligiri. The distributions are plotted around a 24-hour circular diel. The arrows within the plot indicate the mean value (for the most likely distribution for mixture distributions) for the three species in the corresponding parks.20574003175000 R Code ADDIN EN.REFLIST used for assessing temporal and spatio-temporal overlap# Set the working directorysetwd("/RWorking/Sympatry")# Data Input and formatting########################### Read in file; for one example site# Note: I first formatted the txt file to be somewhat similar; please see files# Note: all spaces in location names were removed in the txt files. # Note: All column names are kept as the same. dhole=read.table("Dhole.txt", header=T, colClasses=rep("character", 3))tiger=read.table("Tiger.txt", header=T, colClasses=rep("character", 3))leopard=read.table("Leopard.txt", header=T, colClasses=rep("character", 3))# Adding a species columndhole$Spp=rep("Dhole", nrow(dhole))tiger$Spp=rep("Tiger", nrow(tiger))leopard$Spp=rep("Leopard", nrow(leopard))################ Creating a time axis############################################ Install library chron#install.packages("chron") # for the first timelibrary(chron)# Forming a time column in time format# Note: format in code below should be as per format of the file. dhole$timeaxis=chron(times=dhole$Time, format=c(times="h:m:s"))# Converting it into a number, and scaling it to a 24 hr clockdhole$timenum=as.numeric(dhole$timeaxis) * 24# Same for tigertiger$timeaxis=chron(times=tiger$Time, format=c(times="h:m:s"))tiger$timenum=as.numeric(tiger$timeaxis) * 24# And for leopardleopard$timeaxis=chron(times=leopard$Time, format=c(times="h:m:s"))leopard$timenum=as.numeric(leopard$timeaxis) * 24########################### Creating a circular object for circular statistics############################## Installing package circular#install.packages("circular") # required the first timelibrary(circular)dhole$Circular=circular(dhole$timenum, units="hours", template="clock24")# Checking#summary(dhole$Circular)plot(dhole$Circular, pch=16, col="blue", stack=T, shrink=1, bins=24*4, ticks=F)rose.diag(dhole$Circular, bins=12, col="darkgrey", cex=1.5, prop=1.3, add=T, axes=F) # Oriana# for tigertiger$Circular=circular(tiger$timenum, units="hours", template="clock24")# Checkingsummary(tiger$Circular)plot(tiger$Circular, pch=16, col="blue", stack=T, shrink=1, bins=24*4, ticks=F)rose.diag(tiger$Circular, bins=12, col="darkgrey", cex=1.5, prop=1.3, add=TRUE, axes=F) # Oriana# for leopardleopard$Circular=circular(leopard$timenum, units="hours", template="clock24")# Checkingsummary(leopard$Circular)plot(leopard$Circular, pch=16, col="blue", stack=T, shrink=1, bins=24*4, ticks=F)rose.diag(leopard$Circular, bins=16, col="darkgrey", cex=1.5, prop=1.3, add=TRUE, axes=F) # Oriana########################## Objective 1: Bimodal mean. subobjective: obtaining distributions for all################################### dhole# Watsons GOF test for von mises distributionwatson.test(dhole$Circular, dist="vonmises")# Test Statistic: U2 1.1752 , P-value < 0.01# Not a von mises distributionwatson.test(dhole$Circular, dist="uniform")#Test Statistic: 2.7404 P-value < 0.01 # Not a uniform distributionwatson.test(tiger$Circular, dist="vonmises")# Test Statistic: 0.521, P-value < 0.01 # Not a von mises distribution at 0.05 alphawatson.test(leopard$Circular, dist="vonmises")# Test Statistic: 0.4262 P-value < 0.01 # Not a von mises distbrition at 0.05 alpha#install.packages("movMF")library(movMF)# requires points on a unit circle# See: Unit circle is circle of radius 1, so x = cos, y = sin# Conversion of a 24-scale to a 2pi-scale = x * 2pi / 24dhole$unitcirclex=cos(dhole$timenum * pi / 12) # dhole$unitcircley=sin(dhole$timenum * pi / 12)AIC(movMF(cbind(dhole$unitcirclex, dhole$unitcircley), 1)) # -62.10AIC(movMF(cbind(dhole$unitcirclex, dhole$unitcircley), 2)) # -243.43AIC(movMF(cbind(dhole$unitcirclex, dhole$unitcircley), 3)) # -238.71# Choose model 2dhole.model=movMF(cbind(dhole$unitcirclex, dhole$unitcircley), 2)dhole.mu = atan2(dhole.model$theta[,2], dhole.model$theta[,1]) # conversion back# back transformationdhole.mu.hr=dhole.mu / pi * 12dhole.mu.hr[dhole.mu.hr<0] = 24+dhole.mu.hr[dhole.mu.hr<0]# 1. 7.378733 2. 18.034436 dhole.kappa = sqrt(rowSums(dhole.model$theta^2)) # is this correct?dhole.kappa.hr= dhole.kappa / pi * 12dhole.alphas = dhole.model$alphadhole.fn = function (x) { (dhole.alphas[1] * dvonmises(x, mu=circular(dhole.mu[1], units="radians"), kappa=dhole.kappa[1]) + dhole.alphas[2] * dvonmises(x, mu=circular(dhole.mu[2], units="radians"), kappa=dhole.kappa[2]))}curve.circular(dvonmises(x, mu=circular(dhole.mu[1], units="radians"), kappa=dhole.kappa[1]), to=2*pi, shrink=1.8)curve.circular(dhole.fn, from=0, to=2*pi, shrink=1.5, axes=F)rose.diag(12-dhole$Circular-6, add=T, bins=12, axes=F)# tigertiger$unitcirclex=cos(tiger$timenum * pi / 12) # tiger$unitcircley=sin(tiger$timenum * pi / 12)AIC(movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 1)) # -195.72AIC(movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 2)) # -245.72AIC(movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 3)) # -238.32AIC(movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 4))# Choose model 2tiger.model=movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 3)tiger.mu = atan2(tiger.model$theta[,2], tiger.model$theta[,1]) # conversion backtiger.mu=tiger.mu # back transformationtiger.mu.hr=tiger.mu / pi * 12tiger.mu.hr[tiger.mu.hr<0] = 24+tiger.mu.hr[tiger.mu.hr<0]#20.156616 2.282643 tiger.kappa = sqrt(rowSums(tiger.model$theta^2)) # is this correct?tiger.kappa.hr= tiger.kappa / pi * 12tiger.alphas = tiger.model$alphatiger.model2=movMF(cbind(tiger$unitcirclex, tiger$unitcircley), 2)tiger.mu2 = atan2(tiger.model2$theta[,2], tiger.model2$theta[,1]) # conversion backtiger.mu2=tiger.mu2 # back transformationtiger.mu2.hr=tiger.mu2 / pi * 12tiger.mu2.hr[tiger.mu2.hr<0] = 24+tiger.mu2.hr[tiger.mu2.hr<0]#20.156616 2.282643 tiger.kappa2 = sqrt(rowSums(tiger.model2$theta^2)) # is this correct?tiger.kappa2.hr= tiger.kappa2 / pi * 12tiger.alphas2 = tiger.model2$alphatiger.alphas2 = c(tiger.alphas2, 0.000000001)tiger.fn = function (x) { (tiger.alphas[1] * dvonmises(x, mu=circular(tiger.mu[1], units="radians"), kappa=tiger.kappa[1]) + tiger.alphas[2] * dvonmises(x, mu=circular(tiger.mu[2], units="radians"), kappa=tiger.kappa[2])) + tiger.alphas[3] * dvonmises(x, mu=circular(tiger.mu[3], units="radians"), kappa=tiger.kappa[3])}tiger.fn2 = function (x) { (tiger.alphas2[1] * dvonmises(x, mu=circular(tiger.mu2[1], units="radians"), kappa=tiger.kappa2[1]) + tiger.alphas2[2] * dvonmises(x, mu=circular(tiger.mu2[2], units="radians"), kappa=tiger.kappa2[2]))}curve.circular(tiger.fn, from=0, to=2*pi, shrink=1.5, axes=T)curve.circular(tiger.fn2, from=0, to=2*pi, shrink=1.5, axes=T, col=2, add=T)rose.diag(12-tiger$Circular-6, add=T, bins=12, axes=F)lines(density.circular(6-tiger$Circular, bw=25), shrink=1.2, lty=2, lwd=2)# leopardleopard$unitcirclex=cos(leopard$timenum * pi / 12) # leopard$unitcircley=sin(leopard$timenum * pi / 12)AIC(movMF(cbind(leopard$unitcirclex, leopard$unitcircley), 1)) # -220.88AIC(movMF(cbind(leopard$unitcirclex, leopard$unitcircley), 2)) # -267.15AIC(movMF(cbind(leopard$unitcirclex, leopard$unitcircley), 3)) # -266.56# Choose model 2leopard.model=movMF(cbind(leopard$unitcirclex, leopard$unitcircley), 2)leopard.mu = atan2(leopard.model$theta[,2], leopard.model$theta[,1]) # conversion backleopard.mu=leopard.mu # back transformationleopard.mu.hr=leopard.mu / pi * 12leopard.mu.hr[leopard.mu.hr<0] = 24+leopard.mu.hr[leopard.mu.hr<0]#4.642433 22.720060 leopard.kappa = sqrt(rowSums(leopard.model$theta^2)) # is this correct?leopard.kappa.hr= leopard.kappa / pi * 12leopard.alphas = leopard.model$alphaleopard.fn = function (x) { (leopard.alphas[1] * dvonmises(x, mu=circular(leopard.mu[1], units="radians"), kappa=leopard.kappa[1]) + leopard.alphas[2] * dvonmises(x, mu=circular(leopard.mu[2], units="radians"), kappa=leopard.kappa[2])) }##################### Objective 2: Density diagrams# plotted in window########################################## Option 1: Kernel density estimator from the data itselfplot(density.circular(dhole$Circular, bw=25), shrink=1.2, lty=1, lwd=2, main="", xlab="Kernel density (hours)", ylab="", axes=F)lines(density.circular(tiger$Circular, bw=25), shrink=1.2, lty=2, lwd=2)lines(density.circular(leopard$Circular, bw=25), shrink=1.2, lty=3, lwd=2)legend("topright", c("Dhole", "Tiger", "Leopard"), lty=1:3, lwd=2, cex=0.75)# Option 2: Use the distributionscurve.circular(dhole.fn, from=0, to=2*pi, shrink=1.5, axes=F, rotation="clock", zero=pi/2, lwd=2)circ.ticks=c(0, pi/2, pi, 3*pi/2)circ.labels=seq(0,18,by=6)axis.circular(at=circ.ticks, labels=sim.labels, cex=0.9, zero=pi/2, rotation="clock", units="hours", lwd=2)curve.circular(tiger.fn, from=0, to=2*pi, shrink=1.5, axes=T, add=T, lty=2, rotation="clock", zero=pi/2, lwd=2)curve.circular(leopard.fn, from=0, to=2*pi, shrink=1.5, axes=T, add=T, lty=3, rotation="clock", zero=pi/2, lwd=2)legend("topleft", c("Dhole", "Tiger", "Leopard"), lty=1:3, lwd=2, cex=0.7)arrows.1 <- function(length.ar, angle.ar, ...){ angle.ar=(2*pi) - angle.ar + (pi/2) ab <- cos(angle.ar) * length.ar bc <- sign(sin(angle.ar)) * sqrt(length.ar^2 - ab^2) x1 <- ab y1 <- bc arrows(0, 0, x1, y1, ...)}# dhole.alphas # Choose the largestarrows.1(length.ar = dhole.alphas[1], angle.ar = dhole.mu[1], lwd=2)#arrows.1(length.ar = dhole.alphas[2], angle.ar = dhole.mu[2])# tiger.alphasarrows.1(length.ar = tiger.alphas[2], angle.ar = tiger.mu[2], lty=2, lwd=2)# leopard.alphasarrows.1(length.ar = leopard.alphas[2], angle.ar = leopard.mu[2], lty=3, lwd=2)#################### Objective 3: Watsons U across spp, within site. # Saved as dataframe Watson.df#################################wats.dt=watson.two.test(x=dhole$Circular, y=tiger$Circular)wats.dt$prange="p < 0.001"wats.dl=watson.two.test(x=dhole$Circular, y=leopard$Circular)wats.dl$prange="p < 0.001"wats.lt=watson.two.test(x=leopard$Circular, y=tiger$Circular)wats.lt$prange="0.05 < p < 0.10"# Creating a matrixWatson.df=data.frame(Species = c("Dhole","", "Tiger", "", "Leopard", ""), Watson.Test = rep(c("U2 Statistic", "p val"), 3), Tiger = c(round(wats.dt$statistic,2), wats.dt$prange, "-", "-", round(wats.lt$statistic,2), wats.lt$prange), Leopard = c(round(wats.dl$statistic,2), wats.dl$prange, rep("-", 4)))##################################### Watson's U Across sites within spp is from Oriana....but you can modify above code# Overlap with density kernel and bootstrapped CI# # Point 3 in DJ notes###################################################library(overlap)getBandWidth(dhole$timerad) # 41.57987getBandWidth(leopard$timerad) # 20.44067getBandWidth(tiger$timerad) # 21.78904# Modification of function overlapEst for 3 spp. overlapEst1=function (A, B, C, kmax = 3, adjust = c(1, NA, NA), n.grid = 128) { if (length(adjust) == 1) adjust <- rep(adjust, 3) grid <- seq(0, 2 * pi, length = n.grid) dA.current <- NA out1 <- rep(NA, 7) names(out1) <- c("All", "AB", "BC", "AC", "OnlyA", "OnlyB", "OnlyC") bw0.A <- getBandWidth(A, kmax = kmax) bw0.B <- getBandWidth(B, kmax = kmax) bw0.C <- getBandWidth(C, kmax = kmax) if (!is.na(bw0.A) && !is.na(bw0.B) && !is.na(bw0.C)) { if (!is.na(adjust[1])) { dA.current <- 1 dA <- densityFit(A, grid, bw0.A/adjust[1]) dB <- densityFit(B, grid, bw0.B/adjust[1]) dC <- densityFit(C, grid, bw0.C/adjust[1]) dA1 <- dA[-1]/sum(dA[-1]) # why -1 - to eliminate the 0, and just take the bands dB1 <- dB[-1]/sum(dB[-1]) dC1 <- dC[-1]/sum(dC[-1]) out1[1] <- sum(pmin(dA1, dB1, dC1)) # All three species overlap out1[2] = sum(pmin(dA1, dB1)) # D and L overlap (ignoring T) out1[3] = sum(pmin(dB1, dC1)) # L and T overlap (ignoring D) out1[4] = sum(pmin(dA1, dC1)) # T and D overlap (ignoring L) out1[5]=sum(apply(cbind(dA1, dB1, dC1), 1, function(x) {max(x[1]-max(x[2],x[3]),0)})) # Only D out1[6]=sum(apply(cbind(dA1, dB1, dC1), 1, function(x) {max(x[2]-max(x[1],x[3]),0)})) # Only L out1[7]=sum(apply(cbind(dA1, dB1, dC1), 1, function(x) {max(x[3]-max(x[2],x[1]),0)})) # Only T } } return(out1)}dhl_nh_boot<-resample(dhole$timerad, 100, smooth=TRUE) # Please increase no. of sims after testing for no errorslpd_nh_boot<-resample(leopard$timerad, 100, smooth=TRUE)tgr_nh_boot<-resample(tiger$timerad, 100, smooth=TRUE)nh.est=overlapEst1(dhole$timerad, leopard$timerad, tiger$timerad)nh.est# modification for 3 spp. bootEst1 = function (Amat, Bmat, Cmat, kmax = 3, adjust = c(1, NA, NA), n.grid = 128) { nboot <- min(ncol(Amat), ncol(Bmat), ncol(Cmat)) bsamp <- matrix(NA, nboot, 7) colnames(bsamp) <- c("All", "AB", "BC", "AC", "OnlyA", "OnlyB", "OnlyC") for (i in 1:nboot) bsamp[i, ] <- overlapEst1(Amat[,i], Bmat[,i], Cmat[,i], kmax = kmax, adjust = adjust, n.grid = n.grid) return(bsamp)}nhboot<-bootEst1(dhl_nh_boot, lpd_nh_boot, tgr_nh_boot)nhMean<-colMeans(nhboot)# modification for 3 spp. bootCIlogit1 = function(est, boot) { out = list() for (l in 1:length(est)) { out[[l]] = bootCIlogit(est[l],boot[,l]) } return(out)}bootCIlogit1(nh.est,nhboot)####################################################### MRPP############################# summary(factor(dhole$Location))# summary(factor(tiger$Location))# summary(factor(leopard$Location))dhole$datetime=chron(dates=dhole$Date, times=dhole$Time, format=c(dates="d-month-Y", times="h:m:s"))dhole$datetime.num=as.numeric(dhole$datetime)tiger$datetime=chron(dates=tiger$Date, times=tiger$Time, format=c(dates="d-month-Y", times="h:m:s"))tiger$datetime.num=as.numeric(tiger$datetime)leopard$datetime=chron(dates=leopard$Date, times=leopard$Time, format=c(dates="d-month-Y", times="h:m:s"))leopard$datetime.num=as.numeric(leopard$datetime)dt.times.to=rep(NA, nrow(dhole))for (i in 1:nrow(dhole)) { # subsetting those caps from the same location tiger.temp=tiger[tiger$Location==dhole[i,"Location"],] if (nrow(tiger.temp)>0) { #min time between capture of dhole and tiger dt.times.to[i]=min(abs(as.numeric(tiger.temp$datetime - dhole[i,"datetime"]))) } else {dt.times.to[i] = "Nodata"}}# proportion of sites that dholes were there, where tigers were notprop.dh.exclusive=length(which(dt.times.to=="Nodata"))/length(dt.times.to)# 26.21%# histogram of the minimum times to capture (between dhole and tiger)hist(as.numeric(dt.times.to[dt.times.to!="Nodata"]))med.dt.time=median(as.numeric(dt.times.to[dt.times.to!="Nodata"]), na.rm=T)# median of these times = 3 days and 15 and a half hourschron(med.dt.time) # (01/04/70 17:52:12) - from origin 01/01/70# Randomizationsims=1000meds.dt.sims=rep(NA, sims)prop.excl.sims=rep(NA, sims)for (sim in 1:sims) { # shuffling the locations, keeping the datetime the same. sim.tiger=sample(tiger$Location, length(tiger$Location), replace=F) times.sims=rep(NA, nrow(dhole)) for (i in 1:nrow(dhole)) { # subsetting based on the shuffled locations tiger.temp=tiger[sim.tiger==dhole[i,"Location"],] if (nrow(tiger.temp)>0) { # minimum time to capture times.sims[i]=min(abs(as.numeric(tiger.temp$datetime - dhole[i,"datetime"]))) } else {times.sims[i] = "Nodata"} } prop.excl.sims[sim]=length(which(times.sims=="Nodata"))/length(times.sims) meds.dt.sims[sim]=median(as.numeric(times.sims[times.sims!="Nodata"]), na.rm=T)}plot(density(meds.dt.sims, na.rm=T))abline(v=med.dt.time, col="red")# Probability of getting a greater valuelength(meds.dt.sims[meds.dt.sims>med.dt.time])/1000 # 0.183# leopard and tigerlt.times.to=rep(NA, nrow(leopard))for (i in 1:nrow(leopard)) { # subsetting those caps from the same location tiger.temp=tiger[tiger$Location==leopard[i,"Location"],] if (nrow(tiger.temp)>0) { #min time between capture of dhole and tiger lt.times.to[i]=min(abs(as.numeric(tiger.temp$datetime - leopard[i,"datetime"]))) } else {lt.times.to[i] = "Nodata"}}# proportion of sites that dholes were there, where tigers were notprop.lt.exclusive=length(which(lt.times.to=="Nodata"))/length(lt.times.to)# 20.46%# histogram of the minimum times to capture (between dhole and tiger)hist(as.numeric(lt.times.to[lt.times.to!="Nodata"]))med.lt.time=median(as.numeric(lt.times.to[lt.times.to!="Nodata"]), na.rm=T)chron(med.lt.time) # (01/05/70 16:16:09) - from origin 01/01/70# Randomizationsims=1000meds.lt.sims=rep(NA, sims)prop.excl.simslt=rep(NA, sims)for (sim in 1:sims) { # shuffling the locations, keeping the datetime the same. sim.tiger=sample(tiger$Location, length(tiger$Location), replace=F) times.sims=rep(NA, nrow(leopard)) for (i in 1:nrow(leopard)) { # subsetting based on the shuffled locations tiger.temp=tiger[sim.tiger==leopard[i,"Location"],] if (nrow(tiger.temp)>0) { # minimum time to capture times.sims[i]=min(abs(as.numeric(tiger.temp$datetime - leopard[i,"datetime"]))) } else {times.sims[i] = "Nodata"} } prop.excl.simslt[sim]=length(which(times.sims=="Nodata"))/length(times.sims) meds.lt.sims[sim]=median(as.numeric(times.sims[times.sims!="Nodata"]), na.rm=T)}plot(density(meds.lt.sims, na.rm=T))abline(v=med.lt.time, col="red")# Probability of getting a greater valuelength(meds.lt.sims[meds.lt.sims>med.lt.time])/1000 # 0.024# leopard and dholeld.times.to=rep(NA, nrow(leopard))for (i in 1:nrow(leopard)) { # subsetting those caps from the same location dhole.temp=dhole[dhole$Location==leopard[i,"Location"],] if (nrow(dhole.temp)>0) { #min time between capture of dhole and tiger ld.times.to[i]=min(abs(as.numeric(dhole.temp$datetime - leopard[i,"datetime"]))) } else {ld.times.to[i] = "Nodata"}}# proportion of sites that dholes were there, where tigers were notprop.ld.exclusive=length(which(ld.times.to=="Nodata"))/length(ld.times.to)# 67.78%# histogram of the minimum times to capture (between dhole and tiger)hist(as.numeric(ld.times.to[ld.times.to!="Nodata"]))med.ld.time=median(as.numeric(ld.times.to[ld.times.to!="Nodata"]), na.rm=T)chron(med.ld.time) # (01/06/70 13:04:45) - from origin 01/01/70# Randomizationsims=1000meds.ld.sims=rep(NA, sims)prop.excl.simsld=rep(NA, sims)for (sim in 1:sims) { # shuffling the locations, keeping the datetime the same. sim.dhole=sample(dhole$Location, length(dhole$Location), replace=F) times.sims=rep(NA, nrow(leopard)) for (i in 1:nrow(leopard)) { # subsetting based on the shuffled locations dhole.temp=dhole[sim.dhole==leopard[i,"Location"],] if (nrow(dhole.temp)>0) { # minimum time to capture times.sims[i]=min(abs(as.numeric(dhole.temp$datetime - leopard[i,"datetime"]))) } else {times.sims[i] = "Nodata"} } prop.excl.simsld[sim]=length(which(times.sims=="Nodata"))/length(times.sims) meds.ld.sims[sim]=median(as.numeric(times.sims[times.sims!="Nodata"]), na.rm=T)}plot(density(meds.ld.sims, na.rm=T))abline(v=med.ld.time, col="red")# Probability of getting a greater valuelength(meds.ld.sims[meds.ld.sims>med.ld.time])/1000 # 0.065# For plot of all locationstimes.to.all=data.frame(Spp = rep(rep(c("LT", "DT", "LD"), each=sims), 4), Locn = rep(c("Nagarahole", "Bandipur", "Bhadra", "BRT"), each=3*sims), Times = c(matrix(meds.lt.sims.nh, ncol=1), matrix(meds.dt.sims.nh, ncol=1), matrix(meds.ld.sims.nh, ncol=1), matrix(meds.lt.sims.bp, ncol=1), matrix(meds.dt.sims.bp, ncol=1), matrix(meds.ld.sims.bp, ncol=1), matrix(meds.lt.sims.bd, ncol=1), matrix(meds.dt.sims.bd, ncol=1), matrix(meds.ld.sims.bd, ncol=1), matrix(meds.lt.sims.br, ncol=1), matrix(meds.dt.sims.br, ncol=1), matrix(meds.ld.sims.br, ncol=1)))write.table(x=times.to.all, file="SimTimes.txt", row.names=F)med.times.all = data.frame(Spp = rep(c("LT", "DT", "LD"), 4), Locn = rep(c("Nagarahole", "Bandipur", "Bhadra", "BRT"), each=3), Med.time = c(med.lt.time.nh, med.dt.time.nh, med.ld.time.nh, med.lt.time.bp, med.dt.time.bp, med.ld.time.bp, med.lt.time.bd, med.dt.time.bd, med.ld.time.bd, med.lt.time.br, med.dt.time.br, med.ld.time.br))write.table(x=med.times.all, file="MedTimes.txt", row.names=F)times.to.all = read.table("SimTimes.txt", header=T)med.times = read.table("MedTimes.txt", header=T)times.to.all$Spp = factor(times.to.all$Spp, levels=c("LD", "DT", "LT"), labels = c("Dhole - Leopard", "Dhole - Tiger", "Leopard - Tiger"))med.times$Spp = factor(med.times$Spp, levels=c("LD", "DT", "LT"), labels = c("Dhole - Leopard", "Dhole - Tiger", "Leopard - Tiger"))times.to.all$Locn = factor(times.to.all$Locn, levels=c("Bhadra", "Nagarahole", "Bandipur", "BRT"))times.to.all$Locn = factor(times.to.all$Locn, levels=c("Bhadra", "Nagarahole", "Bandipur", "BRT"))probs.all = rep(NA, nrow(med.times))sims=10000for (i in 1:nrow(med.times)) { probs.all[i] = nrow(times.to.all[times.to.all$Locn==med.times$Locn[i] & times.to.all$Spp==med.times$Spp[i] & times.to.all$Times>med.times$Med.time[i],])/sims}probs.all=round(probs.all, 3)med.times$prob=probs.alllibrary(ggplot2)ggplot(times.to.all, aes(x = Times)) + geom_density(adjust=4, fill="grey90") + geom_vline(data=med.times, aes(xintercept=Med.time), linetype="dashed") + xlab("Days") + ylab("Probability") + ylim(c(0,1)) + geom_text(data = med.times, x = 14, y = 0.8, aes(label=paste("p =", prob)), size = 4, fontface="italic") + facet_grid(Locn~Spp) + theme(text=element_text(size=15)) + theme_bw()plot(density(times.to.all[times.to.all$Locn=="Nagarahole" & times.to.all$Spp=="Leopard - Tiger", "Times"]))location.latlong=read.csv("LocationCoord.csv", header=T)location.latlong=location.latlong[order(location.latlong$Location),]dist.locations=dist(location.latlong)dist.locn=as.matrix(dist.locations)dist.locn[lower.tri(dist.locn, diag=F)]=NAmean(c(dist.locn), na.rm=T)median(c(dist.locn), na.rm=T)100*(sd(dist.locn,na.rm=T)/mean(dist.locn,na.rm=T)) # CV# less than some numberdist.locn[dist.locn>10000]=NAmean(c(dist.locn), na.rm=T)median(c(dist.locn), na.rm=T)100*(sd(dist.locn,na.rm=T)/mean(dist.locn,na.rm=T)) # CVcolnames(dist.locn) = location.latlong$Locationrownames(dist.locn) = location.latlong$Location ................
................
In order to avoid copyright disputes, this page is only a partial summary.
To fulfill the demand for quickly locating and searching documents.
It is intelligent file search solution for home and business.