scbrown86

## GISFrag metric
## https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19950017676.pdf
## 1) Produce a proximity (distance) map between 'patches'
## 2) GISFrag == mean of all values on the proximity map
## 3) Large GISFrag values reflect low forest fragmentation, low values == high fragmentation

gisFrag <- function(x, ...) {
  ## x needs to be a raster where cells with suitable habitat are coded as 1
  ## unsuitable cells coded with 0
  ## extract cell numbers for suitable cells

scbrown86

## Function for returning weighted centroid location
## currently hard-coded to return mean, q10, q20, q80 & q90.
weightedCentre <- function(x, y, z) {
  require(matrixStats); require(Hmisc)
  if (anyNA(c(x, y, z))) {
    stop("There are missing values present in x, y, or z")
  }
  if (length(z) != length(x)) {
    stop("Number of weights supplied not equal to number of coordinates")
  }

scbrown86

# Convert raster to polygon using GDAL
## edited 2019-03-27 to use sf::st_read
polygonise <- function(x, outshape=NULL, gdalformat = 'ESRI Shapefile', 
                        pypath=NULL, readpoly=TRUE, quietish=TRUE) {
  ## x: an R Raster layer, or the file path to a raster file recognised by GDAL
  ## outshape: the path to the output shapefile (if NULL, a temporary file will be created)
  ## gdalformat: the desired OGR vector format. Currently only supports ESRI Shapefile!
  ## pypath: the path to gdal_polygonize.py (if NULL, an attempt will be made to determine the location
  ## readpoly: should the polygon shapefile be read back into R, and returned by this function? (logical)
  ## quietish: should (some) messages be suppressed? (logical)

scbrown86

# Use the Cumulative Relative Frequency Curve hotspot method
# to identify hotspots of high values

# Ref:
## Bartolino, V., Maiorano, L., & Colloca, F. (2011).
## A frequency distribution approach to hotspot identification.
## Population Ecology, 53(2), 351-359. doi:10.1007/s10144-010-0229-2


library(raster)

scbrown86

# Use the Cumulative Relative Frequency Curve hotspot method
# to identify hotspots of high values

# Ref:
## Bartolino, V., Maiorano, L., & Colloca, F. (2011).
## A frequency distribution approach to hotspot identification.
## Population Ecology, 53(2), 351-359. doi:10.1007/s10144-010-0229-2

# example array data
# could do with raster/terra, but wanted example in base r

scbrown86

library(raster)
library(ncdf4)

prec <- getData("worldclim", res = 10, var = "prec")
tmax <- getData("worldclim", res = 10, var = "tmax")

prec
tmax
tmax[] <- tmax[]/10

scbrown86

# The function that produces the colour matrix
colmat <- function(nbreaks = 3, breakstyle = "quantile",
                   upperleft = "#0096EB", upperright = "#820050", 
                   bottomleft = "#BEBEBE", bottomright = "#FFE60F",
                   xlab = "x label", ylab = "y label", plotLeg = TRUE,
                   saveLeg = FALSE) {
   # TODO - replace any tidyr, dplyr etc. functions with data.table #
  library(tidyverse)
  require(ggplot2)
  require(classInt)

scbrown86

# Code from https://www.uva.nl/en/profile/l/o/e.e.vanloon/e.e.vanloon.html

# Implementation of the MapCurves algorithm (DOI 10.1007/s10109-006-0025-x)

mapcurves <- function(A, B, plotGOF = TRUE) {
  ## R implementation of the mapcurves goodness of fit measure
  ## (Hargrove et al. 20006, see full reference below) for comparing
  ## two categorical maps
  ##
  ## usage:

scbrown86

# Convert a GeoJSON object to a geopackage
# Add a second geoJSON object in as a second layer

# https://gis.stackexchange.com/q/223240/24249

library(rgdal)
library(sf)

# Here I open a shapefile, but readOGR can be used to read geoJSON files
nc <- st_as_sf(readOGR(system.file("shape/nc.shp", package = "sf")))

scbrown86

# Functions for rescaling layers/stacks
## Cleaned up/edited the source code for sdvmspecies library (https://cran.r-project.org/web/packages/sdmvspecies/index.html)

rescaleLayer <- function(raster.layer) {
  new.min <- 0; new.max <- 1 ## change for different output range (e.g. 0-255)
  min.value <- cellStats(raster.layer, min, na.rm = TRUE)
  max.value <- cellStats(raster.layer, max, na.rm = TRUE)
  raster.layer <- new.min + (raster.layer - min.value) * ((new.max - new.min) / (max.value - min.value))
  return(raster.layer)
}