library(devtools)
library(knitr)
library(ezknitr)
library(dplyr)
library(gdata)
library(sp)
library(lubridate)
library(rgeos)
library(rgdal)
# Clear environment
remove(list=ls())
Import processed fall15 dataset
fall.df<-read.csv("Processed Data/fall15.df.csv")
Split dataset by each population
ep.fall<-fall.df[fall.df$pop=="EP",]
mcp.fall<-fall.df[fall.df$pop=="MCP",]
ep.sp<-ep.fall
mcp.sp<-mcp.fall
Create SpatialPointsDataFrames
coordinates(ep.sp)<-c("x","y")
coordinates(mcp.sp)<-c("x","y")
Define projection
proj4string(ep.sp)<-CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5
+lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ")
proj4string(mcp.sp)<-CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5
+lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ")
Buffer each point by 3 kilometers
ep.buff<-gBuffer(ep.sp, width=3000)
mcp.buff<-gBuffer(mcp.sp, width=3000)
plot(ep.buff, col="red" )
plot(mcp.buff, col="blue", add=T)
both.pop<-gIntersection(ep.buff, mcp.buff)
plot(both.pop, col="brown", add=T)
Convert from SpatialPolygons to SpatialPolygonsDataFrame so it can be exported as a shapefile. Create a dataframe and display default rownames
buff.df<-data.frame(ID=1:length(both.pop))
both.pop<- SpatialPolygonsDataFrame(both.pop, buff.df)
Can easily be exported as shapefile for use in GIS
# writeOGR(obj=both.pop, dsn="2015 overlap/overlap_polys15/pops",
#layer="fall_pop_overlap15_no_2flyway_cranes", driver="ESRI Shapefile")
Turn points layer into spatial first
fall.sp<-fall.df
coordinates(fall.sp)<-c("x", "y")
proj4string(fall.sp)<-CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5
+lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ")
Extract points within the overlap polygons
over.pts<-over(fall.sp, both.pop)
ind<-which(!is.na(over.pts))
tot.overlap<-fall.df[ind,] #this is a dataframe of all the points that fall into any of the 3 overlap regions
Split overlap SpatialPolygonsDataFrame into each individual area of overlap
test<-SpatialPolygons(both.pop@polygons, proj4string = both.pop@proj4string) #convert to spatialPolygon
ptest<-disaggregate(test) #split apart from list of one into list of many
Create SpatialPolygonDataframe for each area of overlap, which can be exported as a shapefile
for (i in 1:3){
poly<-ptest@polygons[[i]] #extract each polygon object
poly<-SpatialPolygons(list(poly),proj4string=CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5
+lon_0=-96 +x_0=0 +y_0=0 +ellps=GRS80 +datum=NAD83 +units=m +no_defs ")) #convert to SpatialPolygon object
poly.id <- sapply(slot(poly, "polygons"), function(x) slot(x, "ID")) # pull info on ID slot
poly.df<-data.frame(ID=1:length(poly), row.names=poly.id) # use that info to create dataframe
poly<-SpatialPolygonsDataFrame(poly, poly.df) # create SpatialPolygonsDataFrame object
assign(paste("stage",i,sep=""),poly) # give unique name for each area of overlap
#export here as a shapefile if desired
}
Extract crane observations from periods of overlap
st1<-over(fall.sp, stage1) # finds matching geometry between polygons and points
ind1<-which(!is.na(st1$ID)) # creates index for points that were in overlap polygons
st1.pts<-fall.df[ind1,] # subsets points layer for locations that were in area of overlap
st1.pts<-droplevels(st1.pts) # removes id labels from cranes that aren't present
# This gives overlap, but doesn't explicity say when cranes from diff pop's were there at the same time...
# I'll run a loop to determine which days have multiple cranes, and
# also which days have cranes from different populations present simultaneously
Create dataframe of overlap information
#Empty objects to fill in loop
dat.summary.1<-NULL
tempdat<-NULL
pops<-c("EP", "MCP")
st1alt<-st1.pts[order(st1.pts$day, st1.pts$id),] # order points by ascending day and id
st1days<-as.list(unique(st1alt$day)) # pull off days that had cranes present
days.vec<-sort(unlist(st1days)) # convert from list to vector
firstday<-days.vec[1] # pull off value for first day
lastday<-last(days.vec) # pull off value for last day
for(i in firstday:lastday){
tempdat<-st1alt[which(st1alt$day==i),] # pull off info for each day
if(nrow(tempdat)>0){
tempdat$mult<-ifelse(length(unique(tempdat$id))>1,"T","F") # assign if there were multiple cranes present on that day
ifelse(pops%in%tempdat$pop,tempdat$pop.over<-"T",tempdat$pop.over<-"F") # assign if there were cranes from different populations present on that day
dat.summary.1<-rbind(dat.summary.1, tempdat) # store results
}}
#This process can be automated to do all areas at once, but focusing on each area one at a
#time makes it easier to focus on the results of a specific area of overlap
Calculate the number of days visited by each crane
aggregate(dat.summary.1$day, by=list(id=dat.summary.1$id), FUN=function(x){length(unique(x))})
## id x
## 1 0C (Callaway colt 2014) 4
## 2 0J (Waubun adult #2) 50
## 3 4K (Waubun #1) 62
## 4 4M (Ogema adult) 3
## 5 9J (Waubun colt) 13
Create indicator variable of days when multiple cranes were present on a staging area
st1mult<-dat.summary.1[dat.summary.1$mult=="T",]
length(unique(st1mult$day)) #52 days
## [1] 52
Create indicator variable of days in which multiple cranes from different population were present
st1pov<-dat.summary.1[dat.summary.1$pop.over=="T",]
length(unique(st1pov$day)) #50 days
## [1] 50
Recreate previous steps for the second staging area
st2<-over(fall.sp, stage2) # finds matching geometry between polygons and points
ind2<-which(!is.na(st2$ID)) # creates index for points that were in overlap polygons
st2.pts<-fall.df[ind2,] # subsets points layer for locations that were in area of overlap
st2.pts<-droplevels(st2.pts) # removes id labels from cranes that aren't present
#Empty objects to fill in loop
dat.summary.2<-NULL
tempdat<-NULL
st2alt<-st2.pts[order(st2.pts$day, st2.pts$id),] # order points by ascending day and id
st2days<-as.list(unique(st2alt$day)) # pull off days that had cranes present
days.vec<-sort(unlist(st2days)) # convert from list to vector
firstday<-days.vec[1] # pull off value for first day
lastday<-last(days.vec) # pull off value for last day
for(i in firstday:lastday){
tempdat<-st2alt[which(st2alt$day==i),] # pull off info for each day
if(nrow(tempdat)>0){
tempdat$mult<-ifelse(length(unique(tempdat$id))>1,"T","F") # assign if there were multiple cranes present on that day
ifelse(pops%in%tempdat$pop,tempdat$pop.over<-"T",tempdat$pop.over<-"F") # assign if there were cranes from different populations present on that day
dat.summary.2<-rbind(dat.summary.2, tempdat) # store results
}}
Calculate the number of days visited by each crane
aggregate(dat.summary.2$day, by=list(id=dat.summary.2$id), FUN=function(x){length(unique(x))})
## id x
## 1 2M (Bagley) 38
## 2 3C (Larson adult #2) 19
## 3 3E (Wambach adult) 18
## 4 3K (Rasmussen adult) 4
## 5 4E Larson adult #1 7
## 6 5C (Stockyard colt) 2
## 7 6M (Larson colt #2) 18
## 8 7E (Nora Township colt) 12
## 9 7K (Larson colt #1) 18
Create indicator variable of days when multiple cranes were present on a staging area
st2mult<-dat.summary.2[dat.summary.2$mult=="T",]
length(unique(st2mult$day)) #20 days
## [1] 20
Create indicator variable of days in which multiple cranes from different population were present
st2pov<-dat.summary.2[dat.summary.2$pop.over=="T",]
length(unique(st2pov$day)) #18 days
## [1] 18
The third area of overlap
st3<-over(fall.sp, stage3) # finds matching geometry between polygons and points
ind3<-which(!is.na(st3$ID)) # creates index for points that were in overlap polygons
st3.pts<-fall.df[ind3,] # subsets points layer for locations that were in area of overlap
st3.pts<-droplevels(st3.pts) # removes id labels from cranes that aren't present
#Empty objects to fill in loop
dat.summary.3<-NULL
tempdat<-NULL
st3alt<-st3.pts[order(st3.pts$day, st3.pts$id),] # order points by ascending day and id
st3days<-as.list(unique(st3alt$day)) # pull off days that had cranes present
days.vec<-sort(unlist(st3days)) # convert from list to vector
firstday<-days.vec[1] # pull off value for first day
lastday<-last(days.vec) # pull off value for last day
for(i in firstday:lastday){
tempdat<-st3alt[which(st3alt$day==i),] # pull off info for each day
if(nrow(tempdat)>0){
tempdat$mult<-ifelse(length(unique(tempdat$id))>1,"T","F") # assign if there were multiple cranes present on that day
ifelse(pops%in%tempdat$pop,tempdat$pop.over<-"T",tempdat$pop.over<-"F") # assign if there were cranes from different populations present on that day
dat.summary.3<-rbind(dat.summary.3, tempdat) # store results
}
}
Calculate the number of days visited by each crane
aggregate(dat.summary.3$day, by=list(id=dat.summary.3$id), FUN=function(x){length(unique(x))})
## id x
## 1 0C (Callaway colt 2014) 34
## 2 2A (Helliksen July 2015 colt) 6
## 3 5K (Helliksen July adult 2015) 9
Create indicator variable of days when multiple cranes were present on a staging area
st3mult<-dat.summary.3[dat.summary.3$mult=="T",]
length(unique(st3mult$day)) #9 days
## [1] 9
Create indicator variable of days in which multiple cranes from different population were present
st3pov<-dat.summary.3[dat.summary.3$pop.over=="T",]
length(unique(st3pov$day)) #9 days
## [1] 9
Combine all dataframes together for export
dat.summary.1$area<-"1"
dat.summary.2$area<-"2"
dat.summary.3$area<-"3"
overlap.df<-rbind(dat.summary.1, dat.summary.2, dat.summary.3)
stage.summary<-aggregate(overlap15.df$day, by=list(id=overlap15.df$id, pop=overlap15.df$pop,
age=overlap15.df$age, gender=overlap15.df$gender),
FUN=function(x){length(unique(x))})
## Error in aggregate(overlap15.df$day, by = list(id = overlap15.df$id, pop = overlap15.df$pop, : object 'overlap15.df' not found
stage.summary[order(stage.summary$pop, stage.summary$x),]
## id pop age gender x
## 13 7K (Larson colt #1) EP colt male 2
## 3 8M (Stockyard adult) EP adult female 13
## 10 3K (Rasmussen adult) EP adult male 16
## 12 7E (Nora Township colt) EP colt male 17
## 8 7J (Melby colt #1) EP colt female 18
## 7 6M (Larson colt #2) EP colt female 19
## 1 2M (Bagley) EP adult female 33
## 2 4E Larson adult #1 EP adult female 35
## 6 5C (Stockyard colt) EP colt female 44
## 9 0J (Waubun adult #2) EP adult male 52
## 11 3E (Wambach adult) MCP adult male 10
## 5 4M (Ogema adult) MCP adult female 17
## 14 2A (Helliksen July 2015 colt) MCP colt male 17
## 15 5K (Helliksen July adult 2015) MCP colt male 22
## 4 4K (Waubun #1) MCP adult female 57
# So 16 cranes were in the overlapped staging areas, and here are the number of days each spend in any one of them.
Export for use in Event Chart script
# write.csv(x=overlap.df, file="overlap15.df.csv", row.names=FALSE)
Session Info
devtools::session_info()
## Session info -------------------------------------------------------------
## setting value
## version R version 3.4.0 (2017-04-21)
## system x86_64, mingw32
## ui RStudio (1.0.143)
## language (EN)
## collate English_United States.1252
## tz America/Denver
## date 2017-06-01
## Packages -----------------------------------------------------------------
## package * version date source
## assertthat 0.2.0 2017-04-11 CRAN (R 3.4.0)
## base * 3.4.0 2017-04-21 local
## compiler 3.4.0 2017-04-21 local
## datasets * 3.4.0 2017-04-21 local
## DBI 0.6-1 2017-04-01 CRAN (R 3.4.0)
## devtools * 1.13.1 2017-05-13 CRAN (R 3.4.0)
## digest 0.6.12 2017-01-27 CRAN (R 3.4.0)
## dplyr * 0.5.0 2016-06-24 CRAN (R 3.4.0)
## evaluate 0.10 2016-10-11 CRAN (R 3.4.0)
## ezknitr * 0.6 2016-09-16 CRAN (R 3.4.0)
## gdata * 2.17.0 2015-07-04 CRAN (R 3.4.0)
## graphics * 3.4.0 2017-04-21 local
## grDevices * 3.4.0 2017-04-21 local
## grid 3.4.0 2017-04-21 local
## gtools 3.5.0 2015-05-29 CRAN (R 3.4.0)
## highr 0.6 2016-05-09 CRAN (R 3.4.0)
## hms 0.3 2016-11-22 CRAN (R 3.4.0)
## knitr * 1.15.1 2016-11-22 CRAN (R 3.4.0)
## lattice 0.20-35 2017-03-25 CRAN (R 3.4.0)
## lazyeval 0.2.0 2016-06-12 CRAN (R 3.4.0)
## lubridate * 1.6.0 2016-09-13 CRAN (R 3.4.0)
## magrittr 1.5 2014-11-22 CRAN (R 3.4.0)
## markdown 0.8 2017-04-20 CRAN (R 3.4.0)
## memoise 1.1.0 2017-04-21 CRAN (R 3.4.0)
## methods * 3.4.0 2017-04-21 local
## mime 0.5 2016-07-07 CRAN (R 3.4.0)
## R.methodsS3 1.7.1 2016-02-16 CRAN (R 3.4.0)
## R.oo 1.21.0 2016-11-01 CRAN (R 3.4.0)
## R.utils 2.5.0 2016-11-07 CRAN (R 3.4.0)
## R6 2.2.1 2017-05-10 CRAN (R 3.4.0)
## Rcpp 0.12.10 2017-03-19 CRAN (R 3.4.0)
## readr * 1.1.0 2017-03-22 CRAN (R 3.4.0)
## rgdal * 1.2-7 2017-04-25 CRAN (R 3.4.0)
## rgeos * 0.3-23 2017-04-06 CRAN (R 3.4.0)
## rstudioapi 0.6 2016-06-27 CRAN (R 3.4.0)
## sp * 1.2-4 2016-12-22 CRAN (R 3.4.0)
## stats * 3.4.0 2017-04-21 local
## stringi 1.1.5 2017-04-07 CRAN (R 3.4.0)
## stringr 1.2.0 2017-02-18 CRAN (R 3.4.0)
## tibble 1.3.0 2017-04-01 CRAN (R 3.4.0)
## tools 3.4.0 2017-04-21 local
## utils * 3.4.0 2017-04-21 local
## withr 1.0.2 2016-06-20 CRAN (R 3.4.0)
spun with: ezknitr::ezspin(“Overlap Scripts/Fall 2015 Overlap.R”, out_dir=“output”)