randPointsRespectingSelfOther1: Randomizes the location of two sets of geographic points...

In adamlilith/enmSdm: Tools for Modeling Niches and Distributions of Species

Randomizes the location of two sets of geographic points while respecting spatial autocorrelation

Description

THIS FUNCTION PROBABLY WONT WORK AS-IS–IT NEEDS A NEW RANDOMIZATION PROCEDURE TO ACHIEVE CONVERGENCE. This function randomizes the location of two sets of geographic points with respect to one another retaining (more or less) the same distribution of pairwise distances between points with and between sets (plus or minus a user-defined tolerance).

Usage

randPointsRespectingSelfOther1(
  x1,
  x2,
  rast,
  tol1 = NULL,
  tol12 = NULL,
  distFunct = NULL,
  keepData = FALSE,
  verbose = FALSE,
  ...
)

Arguments

x1	Matrix, data frame, SpatialPoints, or SpatialPointsDataFrame object representing geographic points to be shuffled. If this is a matrix or data frame, the first two columns must represent longitude and latitude (in that order). If x is a matrix or data frame, the coordinates are assumed to be unprojected (WGS84) (a coordinate reference system proj4 string or CRS object can be passed into the function using ...). If x is a SpatialPoints or SpatialPointsDataFrame and not in WGS84 or NAD83, then coordinates are projected to WGS84 (with a warning).
x2	As x1. Note that only points in x1. Points in x2 are kept fixed on the landscape.
rast	Raster, RasterStack, or RasterBrick used to locate presences randomly. If this is a RasterStack or a RasterBrick then the first layer will be used (i.e., so cells with NA will not have points located within them).
tol1	Numeric >0, maximum root-mean-square distance allowed between the set of observed pairwise distances between points in x1 and the set of randomized pairwise distances between points simulating x1. The algorithm will shuffle points until the calculated difference is <= this number. Units are the same as units used by the coordinate reference system of x (usually meters).
tol12	As tol1 but for the root-mean-square deviation between points in x1 and x2.
distFunct	Either a function or NULL. If NULL then distGeo is used to calculate distances. Other "dist" functions (e.g., distGeo) can be used. Alternatively, a custom function can be used so long as its first argument is a 2-column numeric matrix with one row for the x- and y-coordinates of a single point and its second argument is a two-column numeric matrix with one or more rows of other points.
keepData	Logical, if TRUE then the original data in x (i.e., columns that do not represent coordinates) will be retained in the output but the coordinates will be shuffled. If FALSE (default) then the returned value will have just shuffled coordinates.
verbose	Logical, if FALSE (default) show no progress indicator. If TRUE then display updates and graph.
...	Arguments to pass to distGeo or sampleRast. Note that if x is a matrix or data frame a coordinate reference system may be passed using crs = <proj4 string code> or crs = <object of class CRS> (see sp package). Otherwise the coordinates are assumed to be unprojected (WGS84).

Value

Object of the same class as x1 but with coordinates randomized.

See Also

randPointsRespectingSelf, randPointsRespectingSelfOther2, randPointsBatch, randPointsBatchSampled, randPointsBatchExtract, randPointsBatchNicheOverlap

Examples

library(dismo)
library(raster)

data(lemurs, package='enmSdm')
longLat <- c('decimalLongitude', 'decimalLatitude')

mad <- raster::getData('GADM', country='MDG', level=0)
elev <- raster::getData('alt', country='MDG', mask=TRUE, res=2.5)

# plot data as-is
plot(mad)
species <- sort(unique(lemurs$species))

for (i in seq_along(species)) {
	
	thisLemur <- lemurs[lemurs$species == species[i], longLat]
	points(thisLemur, pch=i, col=i)
	
}

legend('bottomleft', legend=species, pch=seq_along(species), col=seq_along(species))

# geographically thin presences of each species
thinLemurs <- data.frame()

for (i in seq_along(species)) {
	
	thisLemur <- lemurs[lemurs$species == species[i], ]
	thinned <- geoThin(thisLemur, minDist=10000, longLat=longLat)
	thinLemurs <- rbind(thinLemurs, thinned)
	
}

# plot geographically thinned data
plot(mad)

for (i in seq_along(species)) {
	
	thisLemur <- thinLemurs[thinLemurs$species == species[i], longLat]
	points(thisLemur, pch=i, col=i)
	
}

legend('bottomleft', legend=species, pch=seq_along(species), col=seq_along(species))

# randomize one species with respect to itself 
x <- thinLemurs[thinLemurs$species == 'Eulemur fulvus', longLat]

set.seed(123)
x1rand <- randPointsRespectingSelf(x=x, rast=elev, tol=24000, verbose=TRUE)

# plot observed and randomized occurrences
plot(mad)
points(x, pch=16)
points(x1rand, col='red')

# randomize two species with respect to selves and others
species1 <- species[1]
species2 <- species[3]

x1 <- thinLemurs[thinLemurs$species == species1, longLat]
x2 <- thinLemurs[thinLemurs$species == species2, longLat]

set.seed(123)
tol1 <- tol2 <- tol12 <- 16000
x12rand <- randPointsRespectingSelfOther2(x1=x1, x2=x2, rast=elev,
	tol1=tol1, tol2=tol2, tol12=tol12, verbose=TRUE)

# plot geographically thinned data
plot(mad)
points(x1, pch=21, bg='cornflowerblue')
points(x2, pch=24, bg='cornflowerblue')
points(x12rand$x1rand, pch=1, col='red')
points(x12rand$x2rand, pch=2, col='red')

legend('bottomleft', legend=c(species1, species2,
	legend=paste('rand', species1), paste('rand', species2)),
	pch=c(21, 24, 1, 2), col=c('black', 'black', 'red', 'red'),
	pt.bg=c('cornflowerblue', 'cornflowerblue', NA, NA))

### batch mode
## Not run: 

# download climate data
clim <- raster::getData('worldclim', var='bio', res=2.5)

# lemur data
data(lemurs, package='enmSdm')
longLat <- c('decimalLongitude', 'decimalLatitude')

# geographically thin presences of each species
thinLemurs <- data.frame()

for (i in seq_along(species)) {
	
	thisLemur <- lemurs[lemurs$species == species[i], ]
	thinned <- geoThin(thisLemur, minDist=10000, longLat=longLat)
	thinLemurs <- rbind(thinLemurs, thinned)
	
}

# randomize two species with respect to selves and others
species1 <- species[1]
species2 <- species[3]

x1 <- thinLemurs[thinLemurs$species == species1, longLat]
x2 <- thinLemurs[thinLemurs$species == species2, longLat]

# create null distributions
set.seed(123)
tol1 <- tol2 <- tol12 <- 24000
iterations <- 100 # for analysis set this to 100 or more
# for testing use a small number!

x12rand <- randPointsBatch('randPointsRespectingSelfOther2', x1=x1, x2=x2,
	rast=clim[[1]], tol1=tol1, tol2=tol2, tol12=tol12, iterations=iterations,
	verbose=TRUE)

# get environment that was sampled to use as background
bg <- randPointsBatchSampled(x12rand)
bgEnv <- raster::extract(clim, bg)

# create PCA of environmental space
vars <- paste0('bio', 1:19)
bgPca <- princomp(bgEnv[ , vars], cor=TRUE)

x1env <- raster::extract(clim, x1)
x2env <- raster::extract(clim, x2)

nas1 <- omnibus::naRows(x1env)
nas2 <- omnibus::naRows(x2env)

if (length(nas1) > 0) x1env <- x1env[-nas1, ]
if (length(nas2) > 0) x2env <- x2env[-nas2, ]

# observed niche overlap
obsOverlap <- enmSdm::nicheOverlap(
	x1=x1env,
	x2=x2env,
	env=bgPca,
	vars=vars,
	bins=100,
	cor=TRUE
)

# extract climate at randomized sites
x12rand <- randPointsBatchExtract(x12rand, clim, verbose=TRUE)

# null niche overlap
nullOverlap <- randPointsBatchNicheOverlap(
	rands=x12rand,
	env=bgPca,
	vars=vars,
	bins=100,
	cor=TRUE
)

hist(nullOverlap$d, 20, main='Niche Overlap',
	xlab='Schoener\'s D', xlim=c(0, 1))
abline(v=obsOverlap[['d']], col='blue', lwd=3)
legend('topright', legend='Observed', lwd=3, col='blue')


## End(Not run)

adamlilith/enmSdm documentation built on Jan. 6, 2023, 11 a.m.