Metrics of environmental and community dissimilarity between regions

knitr::opts_chunk$set(cache = TRUE, dev = c("png", "pdf"))
library(here)
library(eucs)
library(gaussDiff)
library(raster)
library(rgeos)
library(vegan)

Calculate the geographic distance between the centroid of the Greater Grampians IBRA subregion and the centroids of the other South-eastern subregions.

rgn_names <- sort(unique(southeast$ibra_subregion))
ibra_subs <- aggregate(map_ibra_southeast, by = "IBRA_SUB_N")
centroids <- gCentroid(
  ibra_subs, byid = TRUE, id = ibra_subs@data$IBRA_SUB_N
)
dists     <- spDistsN1(centroids, centroids["Greater Grampians"], TRUE)
dists     <- dists[rownames(centroids@coords) %in% rgn_names]
names(dists) <- rgn_names

Calculate the Kullback-Leibler distance between the Grampians and the other regions based on the environment defined by relief, topographic wetness, moisture index and sand.

vars <- c(
  "mean_moisture_index_of_low_qtr", "topographic_wetness_index_3sec",
  "relief_1000m_radius", "total_nitrogen"
)

kldists <- vapply(
  rgn_names,
  function(x) {
    r <- covariates_southeast$ibra_subregion == x
    f <- function(x, y) {
      x <- na.omit(as.matrix(x))
      for (i in y) {
        x[, i] <- log(x[, i])
      }
      x
    }
    s <- f(covariates_southeast[r, vars], c(1, 3))
    g <- f(covariates_grampians[vars], c(1, 3))
    normdiff(colMeans(g), cov(g), colMeans(s), cov(s), method = "KL")
  },
  numeric(1)
)

cat(
  round(cor(dists, kldists), 1),
  file = here(
    "notebooks/environmental_distance_metrics_files/figure-html/kldistcor.tex"
  )
)

kldists_univariate <- lapply(
  rgn_names,
  function(x) {
    r <- covariates_southeast$ibra_subregion == x
    f <- function(x, y) {
      x <- na.omit(as.matrix(x))
      for (i in y) {
        x[, i] <- log(x[, i])
      }
      x
    }
    s <- f(covariates_southeast[r, vars], c(1, 3))
    g <- f(covariates_grampians[vars], c(1, 3))
    ans <- sapply(
      1:4,
      function(i) {
        normdiff(
          mean(g[, i]), matrix(var(g[, i])), mean(s[, i]), matrix(var(s[, i])),
          method = "KL"
        )
      }
    )
    ans <- rbind(ans)
    rownames(ans) <- x
    colnames(ans) <- vars
    ans
  }
)

kldists_univariate <- do.call(rbind, kldists_univariate)
colnames(kldists_univariate) <- c("mlq", "twi", "r1k", "tn")

kldists_univariate <- reshape(
  data.frame(region = rownames(kldists_univariate), kldists_univariate),
  idvar = "region",
  varying = colnames(kldists_univariate),
  v.names = "kldist_var",
  timevar = "var",
  times = colnames(kldists_univariate),
  direction = "long",
  new.row.names = 1:length(kldists_univariate)
)

Calculate the Jaccard dissimilarity between the Grampians and the Southeast IBRA subregions.

subreg_prev <- do.call(
  rbind,
  lapply(
    split(wide_southeast[-1:-6], wide_southeast$ibra_subregion), colMeans
  )
)
grampians_prev <- colMeans(wide_grampians[-1:-5])
grampians_only <- setdiff(names(grampians_prev), colnames(subreg_prev))
southeast_only <- setdiff(colnames(subreg_prev), names(grampians_prev))
grampians_extra <- rep(0, length(southeast_only))
southeast_extra <- matrix(0, nrow = NROW(subreg_prev), length(grampians_only))
names(grampians_extra) <- southeast_only
colnames(southeast_extra) <- grampians_only
rownames(southeast_extra) <- rownames(subreg_prev)
subreg_prev <- cbind(subreg_prev, southeast_extra)
grampians_prev <- c(grampians_prev, grampians_extra)
subreg_prev <- subreg_prev[, sort(colnames(subreg_prev))]
grampians_prev <- grampians_prev[sort(names(grampians_prev))]
subreg_prev <- rbind(grampians_prev, subreg_prev)
community_dist <- as.matrix(
  vegdist(subreg_prev, method = "jaccard")
)[rgn_names, "grampians_prev"]

Plot of geographic and Jaccard distances of IBRA subregions from the Grampians.

par(mfrow = c(1, 2), oma = c(5, 0, 1, 1), mar = c(0, 5, 0, 0))

plot(
  dists, community_dist, las = 1, xlab = "",
  ylab = "Jaccard dissimilarity", pch = 19, xlim = c(0, 1000),
  ylim = c(.8, 1.01)
)
text(
  dists, community_dist, rgn_names,cex = .75,
  pos =
    ifelse(
      rgn_names %in% c(
        "Strzelecki Ranges", "East Gippsland Lowlands", "Monaro",
        "Highlands-Northern Fall", "Snowy Mountains",
        "Lower Slopes", "Bateman"
      ),
      2,
      ifelse(
        rgn_names == "Jervis",
        3,
        ifelse(rgn_names %in% c("Victorian Alps", "Ettrema"), 1, 4)
      )
    )
)

plot(
  dists, kldists, las = 1, xlab = "Distance from Grampians (km)",
  ylab = "Kullback Leibler Divergence", pch = 19, xlim = c(0, 1000),
  ylim = c(5, 35)
)
text(
  dists, kldists, rgn_names,
  pos =
    ifelse(
      rgn_names %in% c(
        "Strzelecki Ranges", "Highlands-Southern Fall",
        "Bondo", "Monaro", "Jervis", "East Gippsland Lowlands"
      ),
      2, 4
    ),
  cex = .75
)
mtext("Distance from Grampians (km)", 1, 3, outer = TRUE)