R/raster_resample.R

In cszang/raresa: Spatial resampling of clustered data.

##' Resample spatial configurations within a regular grid
##'
##' Divide region into regular grid cells and sample sites from each grid cell.
##'
##' A regular grid of given size is superimposed on sampling design, and the
##' individual plots are attributed to specific grid cells. Then, from each grid
##' cell a poisson distributed number of sites is selected randomly, and a list
##' of new resampled configurations of the original data is returned.
##' @param x a \code{data.frame} with as least lon and lat as decimal degrees as
##'   variables
##' @param size how big (in degrees) should the grid cells for resampling be?
##' @param n how many resampling configurations shall be computed?
##' @param lambda lambda parameter for poisson distribution that gives count of
##'   sites per grid cells
##'
##' @import dplyr
##' @export
raster_resample <- function(x, size = 5, n = 5, lambda = 3) {

  precision <- 1

  expand_extreme <- function(co, mode = "min") {
    ee <- round(co, precision)
    correction <- 10^(-precision)
    if ((ee - co) >= 0) {
      if (mode == "min") {
        ee - correction
      } else {
        ee
      }
    } else {
      if (mode == "min") {
        ee
      } else {
        ee + correction
      }
    }
  }

  draw_grid <- function(min_lon, max_lon, min_lat, max_lat,
                        size, offset_x = 0, offset_y = 0) {
    xs <- seq(min_lon + offset_x, max_lon + offset_x, by = size)
    xs <- c(xs, tail(xs, 1) + size)
    ys <- seq(min_lat + offset_y, max_lat + offset_y, by = size)
    ys <- c(ys, tail(ys, 1) + size)
    fgrid <- expand.grid(xs, ys)
    names(fgrid) <- c("lon_min", "lat_min")
    fgrid %>% mutate(lat_max = lat_min + size,
                     lon_max = lon_min + size) -> fgrid
    fgrid$gridid <- 1:nrow(fgrid)
    fgrid
  }

  in_which_cell <- function(lon, lat, fgrid) {
    fgrid %>% filter(lon_min < lon, lon_max >= lon,
                     lat_min < lat, lat_max >= lat) %>%
      .$gridid
  }

  min_lon <- expand_extreme(min(x$lon))
  max_lon <- expand_extreme(max(x$lon))
  min_lat <- expand_extreme(min(x$lat))
  max_lat <- expand_extreme(max(x$lat))

  selections <- list()

  for (i in 1:n) {
    .x <- x
    offsets <- -1 * sample(seq(0, size, by = 10^(-precision)), 2)
    the_grid <- draw_grid(min_lon, max_lon, min_lat, max_lat, size = size,
                          offset_x = offsets[1], offset_y = offsets[2])
    .x$cell <- mapply(function(x, y) in_which_cell(x, y, the_grid), x$lon, x$lat)
    ncells <- length(unique(.x$cell))

    counts <- data.frame(cell = unique(.x$cell),
                         count = rpois(ncells, lambda = lambda))

    selection <- NULL
    for (j in 1:nrow(counts)) {
      subs <- .x %>% filter(cell == counts$cell[j])
      nsubs <- nrow(subs)
      if (nsubs >= counts$count[j]) {
        selection <- rbind(selection, subs[sample(1:nsubs, counts$count[j]),])
      } else {
        selection <- rbind(selection, subs)
      }
    }
    selections[[i]] <- selection
  }
  class(selections) <- list("raster_resample", "list")
  selections
}