R/edge_delay.R

In ropensci/spatsoc: Group Animal Relocation Data by Spatial and Temporal Relationship

#' Directional correlation delay based edge lists
#'
#' \code{edge_delay} returns edge lists defined by the directional correlation
#' delay between individuals. The function expects a distance based edge list
#' generated by \code{edge_dist} or \code{edge_nn}, a \code{data.table} with
#' relocation data, individual identifiers and a window argument. The window
#' argument is used to specify the temporal window within which to measure the
#' directional correlation delay. Relocation data should be in two columns
#' representing the X and Y coordinates.
#'
#' The \code{edges} and \code{DT} must be \code{data.table}s. If your data is a
#' \code{data.frame}, you can convert it by reference using
#' \code{\link[data.table:setDT]{data.table::setDT}}.
#'
#' The \code{edges} argument expects a distance based edge list generated with
#' \code{edge_nn} or \code{edge_dist}. The \code{DT} argument expects relocation
#' data with a timegroup column generated with \code{group_times}.
#'
#' The rows in \code{edges} and \code{DT} are internally matched in
#' \code{edge_delay} using the columns \code{timegroup} (from
#' \code{group_times}) and \code{ID1} and \code{ID2} (in \code{edges}, from
#' \code{dyad_id}) with \code{id} (in \code{DT}). This function expects a
#' \code{fusionID} present, generated with the \code{fusion_id} function, and a
#' \code{dyadID} present, generated with the \code{dyad_id} function. The
#' \code{id}, and \code{direction} arguments expect the names of a column in
#' \code{DT} which correspond to the id, and direction columns.
#'
#' @inheritParams centroid_fusion
#' @inheritParams direction_group
#' @param window temporal window in unit of timegroup column generated with
#'   \code{group_times}, eg. \code{window = 4} corresponds to the 4 timegroups
#'   before and after the focal observation
#'
#' @return \code{edge_delay} returns the input \code{edges} appended with
#'   a 'dir_corr_delay' column indicating the temporal delay (in units of
#'   timegroups) at which ID1's direction of movement is most similar to
#'   ID2's direction of movement, within the temporal window defined. For
#'   example, if focal individual 'A' moves in a 45 degree direction at time 2
#'   and individual 'B' moves in a most similar direction within the window
#'   at time 5, the directional correlation delay between A and B is 3. Positive
#'   values of directional correlation delay indicate a directed leadership
#'   edge from ID1 to ID2.
#'
#' @export
#'
#' @references
#'
#' The directional correlation delay is defined in Nagy et al. 2010
#' (<https://doi.org/10.1038/nature08891>).
#'
#' See examples of measuring the directional correlation delay:
#'  * <https://doi.org/10.1016/j.anbehav.2013.07.005>
#'  * <https://doi.org/10.1073/pnas.1305552110>
#'  * <https://doi.org/10.1111/jfb.15315>
#'  * <https://doi.org/10.1371/journal.pcbi.1003446>
#'
#' @family Edge-list generation
#' @family Direction functions
#'
#' @examples
#' # Load data.table
#' library(data.table)
#' \dontshow{data.table::setDTthreads(1)}
#'
#' # Read example data
#' DT <- fread(system.file("extdata", "DT.csv", package = "spatsoc"))
#'
#' # Select only individuals A, B, C for this example
#' DT <- DT[ID %in% c('A', 'B', 'C')]
#'
#' # Cast the character column to POSIXct
#' DT[, datetime := as.POSIXct(datetime, tz = 'UTC')]
#'
#' # Temporal grouping
#' group_times(DT, datetime = 'datetime', threshold = '20 minutes')
#'
#' # Calculate direction
#' direction_step(
#'   DT = DT,
#'   id = 'ID',
#'   coords = c('X', 'Y'),
#'   projection = 32736
#' )
#'
#' # Distance based edge list generation
#' edges <- edge_dist(
#'   DT,
#'   threshold = 100,
#'   id = 'ID',
#'   coords = c('X', 'Y'),
#'   timegroup = 'timegroup',
#'   returnDist = TRUE,
#'   fillNA = FALSE
#' )
#'
#' # Generate dyad id
#' dyad_id(edges, id1 = 'ID1', id2 = 'ID2')
#'
#' # Generate fusion id
#' fusion_id(edges, threshold = 100)
#'
#' # Directional correlation delay
#' delay <- edge_delay(
#'   edges = edges,
#'   DT = DT,
#'   window = 3,
#'   id = 'ID'
#' )
#'
#' delay[, mean(dir_corr_delay, na.rm = TRUE), by = .(ID1, ID2)][V1 > 0]
edge_delay <- function(
    edges,
    DT,
    window = NULL,
    id = NULL,
    direction = 'direction') {
  # due to NSE notes in R CMD check
  . <- timegroup <- fusionID <- min_timegroup <- max_timegroup <-
    delay_timegroup <- ID1  <- ID2 <- dir_corr_delay <- NULL

  if (is.null(DT)) {
    stop('input DT required')
  }

  if (is.null(edges)) {
    stop('input edges required')
  }

  if (is.null(id)) {
    stop('id column name required')
  }

  check_cols_edges <- c('ID1', 'ID2', 'timegroup')
  if (any(!(check_cols_edges %in% colnames(edges)))) {
    stop(paste0(
      as.character(paste(setdiff(
        check_cols_edges,
        colnames(edges)
      ), collapse = ', ')),
      ' field(s) provided are not present in input DT'
    ))
  }

  check_cols_DT <- c(id, 'timegroup', direction)
  if (any(!(check_cols_DT %in% colnames(DT)
  ))) {
    stop(paste0(
      as.character(paste(setdiff(
        check_cols_DT,
        colnames(DT)
      ), collapse = ', ')),
      ' field(s) provided are not present in input DT'
    ))
  }

  if (is.null(window)) {
    stop('window is required')
  }

  if (!is.numeric(window)) {
    stop('window should be a numeric, in the units of timegroup')
  }

  if (!is.integer(DT$timegroup) || !is.integer(edges$timegroup)) {
    stop('timegroup should be an integer, did you use group_times?')
  }

  if (!'fusionID' %in% colnames(edges)) {
    stop('fusionID field not present in edges, did you run fusion_id?')
  }

  if (!'dyadID' %in% colnames(edges)) {
    stop('dyadID field not present in edges, did you run dyad_id?')
  }

  # "Forward": all edges ID1 -> ID2
  forward <- edges[!is.na(fusionID),
                   data.table::first(.SD),
                   by = .(fusionID, timegroup)]

  forward[, min_timegroup :=
            data.table::fifelse(timegroup - window < min(timegroup),
                                min(timegroup),
                                timegroup - window),
          by = fusionID,
          env = list(window = window)]

  forward[, max_timegroup :=
            data.table::fifelse(timegroup + window > max(timegroup),
                                max(timegroup),
                                timegroup + window),
          by = fusionID,
          env = list(window = window)]

  forward[, delay_timegroup := {
    focal_direction <- DT[timegroup == .BY$timegroup &
                            id == ID1, direction]
    DT[between(timegroup, min_timegroup, max_timegroup) & id == ID2,
       timegroup[which.min(diff_rad(focal_direction, direction))]]
  },
  by = c('timegroup',  'dyadID'),
  env = list(id = id, direction = direction)]

  forward[, dir_corr_delay := delay_timegroup - timegroup]

  data.table::set(forward,
                  j = c('min_timegroup', 'max_timegroup','delay_timegroup'),
                  value = NULL)

  # "Reverse": replicate forward but reverse direction ID1 <- ID2
  reverse <- data.table::copy(forward)
  setnames(reverse, c('ID1', 'ID2'), c('ID2', 'ID1'))
  reverse[, dir_corr_delay := - dir_corr_delay]

  out <- data.table::rbindlist(list(
    forward,
    reverse
  ), use.names = TRUE)

  data.table::setorder(out, timegroup)
  data.table::setcolorder(
    out,
    c('timegroup', 'ID1', 'ID2', 'dyadID', 'fusionID')
  )

  return(out)
}