## ----setup, include=FALSE----------------------------------------------------- knitr::opts_chunk$set(echo = TRUE, error = FALSE, fig.retina = 1, dpi = 80) ## ----packages, message=FALSE, warning=FALSE----------------------------------- library(voluModel) # Because of course library(ggplot2) # For fancy plotting library(rangeBuilder) # To compare marineBackground to getDynamicAlphaHull library(dplyr) # To filter data library(terra) # Now being transitioned in library(sf) # Now being transitioned in ## ----occurrence dataset, message=FALSE, warning = FALSE----------------------- # Get points occs <- read.csv(system.file("extdata/Steindachneria_argentea.csv", package='voluModel')) occurrences <- occs %>% dplyr::select(decimalLongitude, decimalLatitude, depth) %>% dplyr::distinct() %>% dplyr::filter(depth %in% 1:2000) ## ----occurrence dataset plotting code, message=FALSE, warning = FALSE, eval=TRUE, echo = FALSE---- # Map occurrences land <- rnaturalearth::ne_countries(scale = "small", returnclass = "sf")[1] pointMap(occs = occurrences, ptCol = "orange", landCol = "black", spName = "Steindachneria argentea", ptSize = 3, land = land) ## ----environmental data loading, eval=T, asis=T, message = FALSE, warning=FALSE---- # Temperature td <- tempdir() unzip(system.file("extdata/woa18_decav_t00mn01_cropped.zip", package = "voluModel"), exdir = paste0(td, "/temperature"), junkpaths = T) temperature <- vect(paste0(td, "/temperature/woa18_decav_t00mn01_cropped.shp")) # Looking at the dataset head(temperature) ## ----temperature to compatible RasterBrick------------------------------------ # Creating a SpatRaster vector template <- centerPointRasterTemplate(temperature) tempTerVal <- rasterize(x = temperature, y = template, field = names(temperature)) # Get names of depths envtNames <- gsub("[d,M]", "", names(temperature)) envtNames[[1]] <- "0" names(tempTerVal) <- envtNames temperature <- tempTerVal # Here's a sampling of depth plots from the 102 depth layers available plot(temperature[[c(1, 50)]]) rm(tempTerVal) ## ----column interpretations, message=TRUE, warning=TRUE----------------------- occsTest <- occurrences[19:24,] xyzSample(occs = occsTest, envBrick = temperature) colnames(occsTest) <- c("x", "y", "z") xyzSample(occs = occsTest, envBrick = temperature) rm(occsTest) ## ----downsample to voxel, eval=TRUE, warning=FALSE, message=FALSE------------- # Gets the layer index for each occurrence by matching to depth layerNames <- as.numeric(names(temperature)) occurrences$index <- unlist(lapply(occurrences$depth, FUN = function(x) which.min(abs(layerNames - x)))) indices <- unique(occurrences$index) downsampledOccs <- data.frame() for(i in indices){ tempPoints <- occurrences[occurrences$index==i,] tempPoints <- downsample(tempPoints, temperature[[1]]) tempPoints$depth <- rep(layerNames[[i]], times = nrow(tempPoints)) downsampledOccs <- rbind(downsampledOccs, tempPoints) } occurrences <- downsampledOccs head(occurrences) print(paste0("Original number of points: ", nrow(occs), "; number of downsampled occs: ", nrow(occurrences))) ## ----plot downsample, warning=FALSE------------------------------------------- pointCompMap(occs1 = occs, occs2 = occurrences, occs1Name = "original", occs2Name = "cleaned", spName = "Steindachneria argentea", land = land, verbose = FALSE) ## ----temperature extraction--------------------------------------------------- # Extraction occurrences$temperature <- xyzSample(occs = occurrences, envBrick = temperature) # Add "response" column for modeling occurrences$response <- rep(1, times = nrow(occurrences)) occurrences <- occurrences[complete.cases(occurrences),] head(occurrences) ## ----alpha hull demonstration, message=FALSE, warning=FALSE, eval=FALSE------- # trainingRegion <- marineBackground(occurrences, # fraction = 1, partCount = 1, buff = 1000000, # clipToOcean = F) # plot(trainingRegion, border = F, col = "gray", # main = "100% Points, Max 1 Polygon Permitted, 100 km Buffer", # axes = T) # plot(land, col = "black", add = T) # points(occurrences[,c("decimalLongitude", "decimalLatitude")], # pch = 20, col = "red", cex = 1.5) ## ----plot alpha hull demonstration, echo=FALSE, fig.width=7------------------- knitr::include_graphics("alphaHullDemonstration-1.png", ) ## ----clipToOcean demo, message=FALSE, warning=FALSE, eval=F------------------- # trainingRegion <- marineBackground(occurrences, # buff = 1000000, # clipToOcean = T) # plot(trainingRegion, border = F, col = "gray", # main = "100 km Buffer, Clipped to Occupied Polygon", # axes = T) # plot(land, col = "black", add = T) # points(occurrences[,c("decimalLongitude", "decimalLatitude")], # pch = 20, col = "red", cex = 1.5) ## ----plot clipToOcean demo, echo=FALSE---------------------------------------- trainingRegion <- vect(system.file("extdata/backgroundSamplingRegions.shp", package='voluModel')) knitr::include_graphics("clipToOceanDemo-1.png") ## ----meridian wrap demo, warning=FALSE, message=FALSE, eval=F----------------- # # Fictional example occurrences # pacificOccs <- occurrences # pacificOccs$decimalLongitude <- pacificOccs$decimalLongitude - 100 # for (i in 1:length(pacificOccs$decimalLongitude)){ # if (pacificOccs$decimalLongitude[[i]] < -180){ # pacificOccs$decimalLongitude[[i]] <- pacificOccs$decimalLongitude[[i]] + 360 # } # } # # # marine Background # pacificTrainingRegion <- marineBackground(pacificOccs, # fraction = 0.95, partCount = 3, # clipToOcean = T) # plot(pacificTrainingRegion, border = F, col = "gray", # main = "marineBackground Antimeridian Wrap", # axes = T) # plot(land, col = "black", add = T) # points(pacificOccs[,c("decimalLongitude", "decimalLatitude")], # pch = 20, col = "red", cex = 1.5) ## ----plot meridian wrap demo, echo=FALSE-------------------------------------- knitr::include_graphics("meridianWrapDemo-1.png") ## ----training points---------------------------------------------------------- # Background backgroundVals <- mSampling3D(occs = occurrences, envBrick = temperature, mShp = trainingRegion, depthLimit = c(50, 1500)) backgroundVals$temperature <- xyzSample(occs = backgroundVals, temperature) #Remove incomplete cases backgroundVals <- backgroundVals[complete.cases(backgroundVals),] ## ----------------------------------------------------------------------------- # Add "response" column for modeling backgroundVals$response <- rep(0, times = nrow(backgroundVals)) head(backgroundVals) ## ----glm data prep part 2----------------------------------------------------- # Sample background points weighted by distance from mean of occurrence temperatures meanTemp <- mean(occurrences[,c("temperature")]) backgroundVals$distance <- abs(backgroundVals[,"temperature"] - meanTemp) backgroundVals$sampleWeight <- (backgroundVals$distance - min(backgroundVals$distance))/(max(backgroundVals$distance) - min(backgroundVals$distance)) sampleForAbsence <- sample(x = rownames(backgroundVals), size = nrow(occurrences) * 100, prob = backgroundVals$backgroundVals) backgroundVals <- backgroundVals[match(sampleForAbsence, rownames(backgroundVals)),] backgroundVals$response <- as.factor(backgroundVals$response) # Unite datasets and see how things look dataForModeling <- rbind(occurrences, backgroundVals[,colnames(occurrences)]) ggplot(dataForModeling, aes(x = temperature, fill = response, color = response)) + geom_density(alpha = .6) + scale_color_manual(values=c("#999999", "#E69F00"), labels = c("Pseudoabsence", "Presence")) + scale_fill_manual(values=c("#999999", "#E69F00"), labels = c("Pseudoabsence", "Presence")) + labs(title="Temperature Sampling Density,\nOccurrences vs. Pseudoabsences", x="Temperature (C)", y = "Density")+ theme_classic() ## ----cleanup temporary directory---------------------------------------------- unlink(td, recursive = T)