#creates a vector of the list of packageslist.of.packages <-c("activity","overlap","tidyverse","dplyr","lubridate") # check if you have the list of packages made above in your library if not this will automatically install any missing packagesnew.packages <- list.of.packages[!(list.of.packages %in%installed.packages()[,"Package"])]if(length(new.packages)) install.packages(new.packages,repos ="http://cran.us.r-project.org")lapply(list.of.packages, require, character.only =TRUE)
Loading required package: activity
Warning: package 'activity' was built under R version 4.4.2
Loading required package: overlap
Warning: package 'overlap' was built under R version 4.4.2
Loading required package: suntools
Warning: package 'suntools' was built under R version 4.4.2
Loading required package: tidyverse
Warning: package 'tidyverse' was built under R version 4.4.2
Warning: package 'ggplot2' was built under R version 4.4.2
Warning: package 'tidyr' was built under R version 4.4.2
Warning: package 'readr' was built under R version 4.4.2
Warning: package 'purrr' was built under R version 4.4.2
Warning: package 'dplyr' was built under R version 4.4.2
Warning: package 'forcats' was built under R version 4.4.2
Warning: package 'lubridate' was built under R version 4.4.2
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
2 Step 2 assign Latitude and Longitude to the img data set
# Add lat/long to the camera trapping raw data pulled from the deployment dataimg <- img %>%left_join(deployment %>%select(deployment_id, latitude, longitude), by ="deployment_id") # This will add latitude and longitude columns to img based on matching deployment_id
Warning in left_join(., deployment %>% select(deployment_id, latitude, longitude), : Detected an unexpected many-to-many relationship between `x` and `y`.
ℹ Row 54464 of `x` matches multiple rows in `y`.
ℹ Row 162 of `y` matches multiple rows in `x`.
ℹ If a many-to-many relationship is expected, set `relationship =
"many-to-many"` to silence this warning.
# Convert from the existing format ("%d-%m-%Y %H:%M") to POSIXctimg$timestamp <-as.POSIXct(img$timestamp, format ="%d/%m/%Y %H:%M", tz ="America/Los_Angeles")# Confirm the structurestr(img$timestamp)
# calculate solar time tmp <-solartime ( img$timestamp, # the date time column img$latitude, # Latitude img$longitude, # Longitudetz=-8, # an offset in numeric hours to UTC (California time is UTC-8)format="%Y:%m:%d %H:%M:%S") # format now in y/m/d instead of d/m/y as the timestamp has been converted to POSIXct
# add solar and clock times to the datasetimg$solar <- tmp$solarimg$clock <- tmp$clock# plot the relationship between solar and clock timesplot(img$solar, img$clock)
4 Step 4 Single species activity for Squirrel, Bear and Deer
We will do the Squirrels first as I have done them in that order throughout the all of the predictive modelling.
# Count the number of detections of the Douglas's Squirrelsquirrel_count <-sum(img$common_name =="Douglas's Squirrel")squirrel_count # this resulted in 35703 detections
[1] 35703
# Fit an activity model based on the data from imgsquirrel_act <-fitact(img$solar[img$common_name =="Douglas's Squirrel"], sample ="model", reps =10) #typically uses 1000 bootstraps but due to rendering issues with Quarto reduced to 10 plot(squirrel_act)
now for the Black Bear
# Count the number of detections of the American Black Bearbear_count <-sum(img$common_name =="American Black Bear")bear_count # this resulted in 13502 detections
[1] 13502
# Fit an activity model based on the data from imgbear_act <-fitact(img$solar[img$common_name =="American Black Bear"], sample ="model", reps =10) plot(bear_act)
now for the Mule Deer
# Count the number of detections of the Mule Deerdeer_count <-sum(img$common_name =="Mule Deer")deer_count # this resulted in 31781 detections
[1] 31781
# Fit an activity model based on the data from imgdeer_act <-fitact(img$solar[img$common_name =="Mule Deer"], sample ="model", reps =10) plot(deer_act)
It is valuable also to understand the impact that Human activity has on the three species too. I will be using genus for detection humans as the when the primary researcher went through the data they named homo sapiens as human research or homo sapian or just homo therefore this mixture of different names could impact the validity of the model esspecially as it might not be accounting for a large majority of the time humans where detected.
# Count the number of detections of Humanshuman_count <-sum(img$genus =="Homo")human_count # this resulted in detections
[1] 50821
# Fit an activity model based on the data from imghuman_act <-fitact(img$solar[img$genus =="Homo"], sample ="model", reps =10)
Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
collapsing to unique 'x' values
plot(human_act)
5 Step 5 examine the raw data from the KDE (kernal density estimation)
The squirrels results from the code below shows that the Douglas’s Squirrel was active for
# Extract the activity information from the squirrel_act KDE to see how much of the day they are active for around 31% of the day.act_info <-slot(squirrel_act, "act")act_info
act se lcl.2.5% ucl.97.5%
0.281629284 0.003564769 0.293632766 0.304222289
The deers results from the code below shows that the Mule Deer was active for
# Extract the activity information from the bear_act KDE to see how much of the day they are active for around 51% of the day.act_info <-slot(bear_act, "act")act_info
act se lcl.2.5% ucl.97.5%
0.496017843 0.009491015 0.494572316 0.522795040
The bear results from the code below shows that the American Black Bear was active for
# Extract the activity information from the deer_act KDE to see how much of the day they are active for around 50% of the day.act_info <-slot(deer_act, "act")act_info
act se lcl.2.5% ucl.97.5%
0.526242726 0.009401863 0.543062092 0.570208716
# Extract the activity information from the human_act KDE to see how much of the day they are active for around 28% of the day.act_info <-slot(human_act, "act")act_info
act se lcl.2.5% ucl.97.5%
0.280668510 0.002273329 0.283636636 0.290044322
6 Step 6 compare species activity
# plot multiple species together# let's make sure that the figure accommodate the maximum y-values for both activity curves. Set y-axis limits between 0 and 0.20. You might need to plot this a few times until you get the optimum y max valuey_lim <-c(0, 0.25)# Plot the douglas's squirrel activityplot(squirrel_act, yunit ="density", data ="none", las =1, lwd =2,ylim = y_lim, # Set y-axis limitstline =list(lwd =2), # Thick line cline =list(lty =0)) # Suppress confidence intervals# Plot the black bear activityplot(bear_act, yunit ="density", data ="none", add =TRUE, tline =list(col ="red", lwd =2),cline =list(lty =0))# Plot the mule deer activityplot(deer_act, yunit ="density", data ="none", add =TRUE, tline =list(col ="blue", lwd =2),cline =list(lty =0))# Plot the human activityplot(human_act, yunit ="density", data ="none", add =TRUE, tline =list(col ="green", lwd =2),cline =list(lty =0))# Add a legend to the plot for the four specieslegend("topleft", c("Douglas's Squirrel", "American Black Bear", "Mule Deer", "Human"), col =c("black", "red", "blue", "green"), lty =1, lwd =2)
7 Step 7 compare the activity levels the four chosen species.
although the previous visual comparison was already produced a numerical comparison will be valuable to understanding the graph created previously.
# Douglas's Squirrel vs Black BearcompareCkern(squirrel_act, bear_act, reps =100) #typically uses 1000 bootstraps but due to rendering issues with Quarto reduced to 100
# the obs column represents the observered overlap within the data which means there is a 65% overlap in the activity of bears and squirrels#seNull shows us a low standard error of 0.0025 which means there is a lowe expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasn't by chance that this high overlap exists
# Douglas's Squirrel vs Mule DeercompareCkern(squirrel_act, deer_act, reps =10) #typically 1000 repetitions however due to rendering issues reduced to 10
# the obs column represents the observered overlap within the data which means there is a 68% overlap in the activity of deer and squirrels#seNull shows us a low standard error of 0.0020 which means there is a low expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasnt by chance that this high overlap exists
# Douglas's Squirrel vs HumancompareCkern(squirrel_act, human_act, reps =10)
## the obs column represents the observered overlap within the data which means there is a 70% overlap in the activity of humans and squirrels#seNull shows us a low standard error of 0.0015 which means there is a low expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasn't by chance that this high overlap exists
# Douglas's Squirrel vs Mule DeercompareCkern(bear_act, deer_act, reps =10)
# the obs column represents the observered overlap within the data which means there is a 85% overlap in the activity of bears and deer#seNull shows us a low standard error of 0.0026 which means there is a low expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasn't by chance that this high overlap exists
# Douglas's Squirrel vs Mule DeercompareCkern(bear_act, human_act, reps =10)
# the obs column represents the observered overlap within the data which means there is a 49% overlap in the activity of bears and humans#seNull shows us a low standard error of 0.0024 which means there is a lowe expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasn't by chance that this high overlap exists
# Douglas's Squirrel vs Mule DeercompareCkern(deer_act, human_act, reps =10)
# the obs column represents the observered overlap within the data which means there is a 54% overlap in the activity of bears and squirrels#seNull shows us a low standard error of 0.0016 which means there is a low expected variation if conducted again. #the pNull shows is significantly low the output being 0.000000000 which means that it wasn't by chance that this high overlap exists
the more statically sound test is called a wald test which can test for the difference between each of the animals temporal activity.
compareAct(list(squirrel_act,deer_act))
Difference SE W p
1v2 -0.2446134 0.01005498 591.8318 0
# The test shows a significant difference in the alignment of activity patterns between douglas's squirrel and mule deer, with a difference of -0.25 indicating that their activity times are not significantly aligned (W = 5.61096, SE = 0.01275986, p = 0.02).
compareAct(list(squirrel_act,bear_act))
Difference SE W p
1v2 -0.2143886 0.01013839 447.1624 0
#The test shows a significant difference in the alignment of activity patterns between douglas's squirrel and black bear, with a difference of -0.21 indicating that their activity times are not significantly aligned (W = 372.6577, SE = 0.01110572, p < 0.001).
compareAct(list(squirrel_act,human_act))
Difference SE W p
1v2 0.0009607743 0.004227955 0.05163949 0.8202346
# The test shows a significant difference in the alignment of activity patterns between humans and douglas's squirrel, with a difference of 0.001 indicating that their activity times are nearly perfectly aligned show no significant difference in their activity patterns. (W = 0.03447079, SE = 0.005, p = 0.9).
compareAct(list(bear_act,deer_act))
Difference SE W p
1v2 -0.03022488 0.01335943 5.118625 0.02367036
# The test shows a significant difference in the alignment of activity patterns between black bears and mule deer, with a difference of -0.03 indicating that their activity times are not aligned with only slight significance (W = 5.61096, SE = 0.01275986, p = 0.02).
compareAct(list(bear_act,human_act))
Difference SE W p
1v2 0.2153493 0.009759477 486.8935 0
#The test shows a significant difference in the alignment of activity patterns between black bears and mule deer, with a difference of 0.22 indicating that their activity times are significantly alligned (W = 428.1654, SE = 0.01040729, p < 0.001).
compareAct(list(human_act,deer_act))
Difference SE W p
1v2 -0.2455742 0.0096728 644.5567 0
# The test shows a significant difference in the alignment of activity patterns between humans and mule deer, with a difference of -0.25 indicating that their activity times are significantly not closely aligned (W = 910.207, SE = 0.00813978, p < 0.001).
8 Step 8 single species activity
img$Area <-1# the package needs an 'Area' column, so we will add this
9 Step 8.1 Douglas’s Squirrel
# plot a kernel density curvesquirrel.overlap<-img$solar[img$Area==1& img$common_name =="Douglas's Squirrel"]densityPlot(squirrel.overlap, rug=T, # rug=T tells r to make a bar along the bottom fo the graph indicating when the captures where during the dayadjust=1, # degree of smoothing of 1 to avoid jagged peaks ylim =c(0,0.11)) legend("topleft", c("Degree of smoothing = 1"), col =c("black"), lty =1, lwd =2)
10 Step 8.2 American Black Bear
# plot a kernel density curvebear.overlap<-img$solar[img$Area==1& img$common_name =="American Black Bear"]densityPlot(bear.overlap, rug=T,adjust=1,ylim =c(0,0.11)) legend("topleft", c("Degree of smoothing = 1"), col =c("black"), lty =1, lwd =2)
11 Step 8.3 Mule Deer
# plot a kernel density curvedeer.overlap<-img$solar[img$Area==1& img$common_name =="Mule Deer"]densityPlot(deer.overlap, rug=T,adjust=1,ylim =c(0,0.11)) legend("topleft", c("Degree of smoothing = 1"), col =c("black"), lty =1, lwd =2)
12 Step 8.4 Humans
# plot a kernel density curvehuman.overlap<-img$solar[img$Area==1& img$genus =="Homo"]densityPlot(human.overlap, rug=T,adjust=1,ylim =c(0,0.18)) legend("topleft", c("Degree of smoothing = 1"), col =c("black"), lty =1, lwd =2)
13 Step 9 Comparing species
14 Step 9.1 Comparing Douglas’s Squirrel and American Black Bear
# Best estimator depends on the size of the smaller of the four samples which is 13502 meaning we need to use DHAT4min(length(bear.overlap), length(squirrel.overlap), length(deer.overlap), length(human.overlap))
[1] 13502
# Calculate the overlap of both species Dhat represents the degrees of similaritysquirrel_bear_est<-overlapEst(squirrel.overlap,bear.overlap)squirrel_bear_est
Dhat1 Dhat4 Dhat5
0.6474488 0.6482735 0.6653143
# Dhat1 Dhat4 Dhat5# 0.65 0.65 0.67#generate random samples 1000 times based of the original distributionssquirrel_boot<-resample(squirrel.overlap,1000) # ran 10,000 repetitions which is industry standard however due to not being able to render quarto document, reduced for rendering.bear_boot<-resample(bear.overlap, 1000)#check the dimensions of the bootstrapping to ensure that it is correct, the dimensions should consist of the number of bootstraps and the number of captures.dim(squirrel_boot)
[1] 35703 1000
# [1] 35703 10000dim(bear_boot)
[1] 13502 1000
# [1] 13502 10000# The code below generates estimates of overlap from each pair of bootstrapped datasquirrel_bear_boot<-bootEst(squirrel_boot, bear_boot,adjust=c(NA,1,NA)) #the 1 in between the NA's indicate the Dhat that we went to use in this case Dhat 4dim(squirrel_bear_boot)
[1] 1000 3
# [1] 10000 3# bootstrap mean of coefficient of overlap BSmean<-colMeans(squirrel_bear_boot)BSmean
Dhat1 Dhat4 Dhat5
NA 0.6545558 NA
# Dhat4 0.66# confidence intervals with corrections for bootstrap biastmp<-squirrel_bear_boot[,2]# So our coefficient of overlap is:squirrel_bear_est
Dhat1 Dhat4 Dhat5
0.6474488 0.6482735 0.6653143
# Dhat4 0.65# The 95% confidence interval is taken from basic0:bootCI(squirrel_bear_est[2],tmp)
overlapPlot(squirrel.overlap,bear.overlap, rug=T, main=NULL)legend('topright', c("Douglas's Squirrel (n = 35703)", "American Black Bear (n = 13502)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.648 (95% CI 0.64 - 0.657)"))))
Squirrel vs Deer
overlapPlot(squirrel.overlap,deer.overlap, rug=T, main=NULL)legend('topright', c("Douglas's Squirrel (n = 35703)", "Mule Deer (n = 31781)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.68 (95% CI 0.671 - 0.686)"))))
Squirrel vs Human
overlapPlot(squirrel.overlap,human.overlap, rug=T, main=NULL)legend('topright', c("Douglas's Squirrel (n = 35703)", "Human (n = 50821)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.704 (95% CI 0.699 - 0.71)"))))
Bear vs Deer
overlapPlot(deer.overlap,bear.overlap, rug=T, main=NULL)legend('topright', c("Mule Deer (n = 31781)", "American Black Bear (n = 13502)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.848 (95% CI 0.839 - 0.857)"))))
Bear vs Humans
overlapPlot(human.overlap,bear.overlap, rug=T, main=NULL)legend('topright', c("American Black Bear (n = 13502), Human (n = 50821)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.492 (95% CI 0.485 - 0.499)"))))
Deer vs human
overlapPlot(deer.overlap,human.overlap, rug=T, main=NULL)legend('topright', c("Mule Deer (n = 31781)", "Human (n = 50821)"), lty=c(1,2), col=c(1,4), bty='n', cex=0.85)# add text for overlap estimate and test for significance (done later on)text(x=6,y=0.12, # coordinates of this textcex=0.9,label=c(expression(paste(hat(Delta)," = 0.535 (95% CI 0.53 - 0.541)"))))
