## MOVEMENT SUMMARIES

# loading packages
pcks <- list("dplyr", "lubridate", "adehabitatLT", "ggplot2", "move") #, "sf")
sapply(pcks, require, char = TRUE)

# loading caribou data (a data frame and a 'move' object)
load("./_data/caribou.df.rda")
load("./_data/caribou.utm.rda")

# calculating step lengths for an individual
FA1403 <- subset(caribou.df, id == "FA1403")
Z <- FA1403$x + (0+1i) * FA1403$y
StepLength <- Mod(diff(Z))

# calculating time lags between locations
dT <- as.numeric(diff(FA1403$time))
table(dT)

# calculating speeds and looking at a histogram of speeds in km/hr
Speed <- StepLength/dT
hist(Speed/1000, breaks = 50, col = "grey", xlab="km/hour", main="Speeds of caribou FA1403")
summary(Speed/1000)

# creating a function to calculate these variables and add them to a data frame for multiple animals
computeDD <- function(data) {
    Z <- data$y + (0+1i) * data$y
    StepLength <- c(NA, Mod(diff(Z)))
    dT <- c(NA, as.numeric(diff(data$time)))
    Speed <- (StepLength/dT)/1000
    return(data.frame(data, StepLength, dT, Speed))
}

# actually adding the variables to the data frame
caribou.df <- caribou.df %>%
  group_by(id) %>% 
  computeDD

# a boxplot of speed versus hour of day (looks weird because this is in UTC!)
boxplot((Speed) ~ hour(time), data = caribou.df, 
    ylim = c(0, 4), ylab = "km/day", xlab="hour of day", pch = 19, cex = 0.5, col = rgb(0, 0, 0, 0.5))

# ggplot of speed versus hour of day (but this is in UTC!))
ggplot(subset(caribou.df, Speed>0), aes(x = hour(time), y = Speed)) + 
  facet_grid(~id) + geom_point(alpha = 0.3) + 
  stat_smooth() + scale_y_log10() +
  ylab("Speed (km/hr)") +  xlab("Hour of day")

# fixing the problem so we plot local times
caribou.df$time.ak <- caribou.df$time - dhours(9)
tz(caribou.df$time.ak) <- "US/Alaska"

# Corrected plot with local times
ggplot(subset(caribou.df, Speed>0), aes(x = hour(time.ak), y = Speed)) + 
  facet_grid(~id) + geom_point(alpha = 0.3) + stat_smooth() +
  scale_y_log10() +  ylab("Speed (km/hr)") + xlab("Hour of day")


# Calculating distances, time lags and speeds in the move package
distances.move <- distance(caribou.utm)
timeLag.move <- timeLag(caribou.utm, units="hours")
speeds.move <- move::speed(caribou.utm)

# Transforming the data from UTMs back to lat longs (required for 'move' to calculate turn angles)
# Then calculating turn angles and looking at a histogram
caribou.move <- spTransform(caribou.utm, "+init=epsg:4326")
turn.angles.move <- turnAngleGc(caribou.move)
hist(unlist(turn.angles.move), main="Histogram of relative turn angles of four caribou", xlab="Relative turn angles", breaks=40)


# Creating an ltraj object in adehabitatLT, which calculates all these variables automatically
caribou.ltraj <- as.ltraj(xy = caribou.df[,c("x","y")], date = caribou.df$time, id = caribou.df$id)
caribou.ltraj
head(caribou.ltraj[[1]])


# loading the RasterStack from the SpatialData lab
load("./_data/fmch.habitat.rda")
elevation <- fmch.habitat[[1]]

# extracting elevation data for each caribou location and adding it to the caribou data frame
caribou.df$elevation.dem <- raster::extract(elevation, caribou.df[, c("x", "y")]) 

# making a boxplot of elevation versus week of year, colored by animal id
ggplot(caribou.df, aes(as.factor(week(time)), elevation, col=id)) + 
  geom_boxplot() +
  xlab("week of year")

# a simple scatterplot comparing elevation collected by the collar versus elevation from the DEM
with(caribou.df, plot(elevation, elevation.dem))



# creating an 'sf' points object from the caribou data frame
require(sf)
caribou.sf <- st_as_sf(caribou.df, 
                        coords = c("x", "y"),
                        crs = "+proj=utm +zone=6 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")

# loading in the fire data shapefile and projecting it to UTM
fires <- st_read("./_data/fmch_fire_tws.shp")
fires.utm <- st_transform(fires, "+proj=utm +zone=6 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")

# overlaying the caribou locations with the fire polygons
# Then assigning locations not in a burn as 0 and those in a burn as 1
caribou.sf <- st_join(caribou.sf, dplyr::select(fires.utm, FireYear))
caribou.sf$fire <- as.factor(ifelse(is.na(caribou.sf$FireYear), 0, 1))
summary(caribou.sf$fire)

# Plotting caribou speed within burns and outside of burns
ggplot(caribou.sf, aes(id, Speed, col=fire)) + 
  geom_boxplot()

# plotting GPS location elevation in burns and outside of burns
ggplot(caribou.sf, aes(id, elevation.dem, col=fire)) + 
  geom_boxplot()


# loading the roads polyline shapefile and transforming to UTM
roads <- st_read("./_data/fmch_roads_tws.shp")
roads.utm <- st_transform(roads, "+proj=utm +zone=6 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")

# Calculating the distance from each GPS location to the nearest road
dist.roads <- st_distance(caribou.sf, roads.utm)
caribou.sf$dist.roads <- as.numeric(apply(dist.roads,1,min))

# Plot of distance to road by week of year, colored by animal ID
ggplot(caribou.sf, aes(as.factor(week(time)), dist.roads, col=id)) + 
  geom_boxplot() +
  xlab("Week of year") +
  ylab("distance to roads")