R/tfd-class.R

In tf: S3 Classes and Methods for Tidy Functional Data

warn_na_entries_created <- function(na_indices) {
  n_entries <- length(na_indices)
  entry_label <- if (n_entries == 1) "entry" else "entries"
  function_label <- if (n_entries == 1) "function" else "functions"
  index_label <- if (n_entries == 1) "index" else "indices"
  shown_indices <- head(na_indices, 10)
  shown_string <- paste(shown_indices, collapse = ", ")
  if (n_entries > 10) {
    shown_string <- paste0(shown_string, ", ...")
  }
  cli::cli_warn(c(
    "{n_entries} {.code NA} {entry_label} (empty {function_label}) created.",
    i = "Affected {index_label}: {shown_string}"
  ))
}

new_tfd <- function(
  arg = NULL,
  datalist = NULL,
  regular = TRUE,
  domain = NULL,
  evaluator
) {
  # evaluator argument parsing needs to deal with indirection/lazy evals:
  # 1) evaluator given as bare/quoted function name:
  evaluator_f <- try(
    get(evaluator, mode = "function", envir = parent.frame()),
    silent = TRUE
  )
  if (!inherits(evaluator_f, "try-error")) {
    # turn bare into quoted name if necessary
    if (!is.character(evaluator)) evaluator <- deparse(evaluator)
  } else {
    # 2) given as string (x <- "tf_approx_bla") that's name of a function
    evaluator <- get(evaluator, parent.frame())
    evaluator_f <- get(evaluator, mode = "function", envir = parent.frame())
  }

  if (
    vec_size(datalist) == 0 || allMissing(unlist(datalist, use.names = FALSE))
  ) {
    arg <- arg %||% list(numeric())
    domain <- domain %||% numeric(2)
    subclass <- if (regular) "tfd_reg" else "tfd_irreg"
    datalist <- list()
    # message("empty or missing input `data`; returning prototype of length 0")
    ret <- new_vctr(
      datalist,
      arg = arg,
      domain = domain,
      evaluator = evaluator_f,
      evaluator_name = evaluator,
      class = c(subclass, "tfd", "tf")
    )
    return(ret)
  }

  assert_string(evaluator)
  assert_function(evaluator_f, args = c("x", "arg", "evaluations"), nargs = 3)

  # sort args and values by arg:
  arg_o <- map(arg, order)
  arg <- map2(arg, arg_o, \(x, y) x[y])
  datalist <- map2(datalist, arg_o, \(x, y) unname(x[y]))

  domain <- domain %||% range(arg, na.rm = TRUE)
  assert_numeric(
    domain,
    finite = TRUE,
    any.missing = FALSE,
    sorted = TRUE,
    len = 2,
    unique = TRUE
  )
  u_args <- unlist(arg, use.names = FALSE)
  if (domain[1] > min(u_args) || max(u_args) > domain[2]) {
    cli::cli_abort("Evaluations must be inside the domain.")
  }

  if (!regular) {
    datalist <- map2(
      datalist,
      arg,
      function(x, y) {
        if (is.null(x) || allMissing(x)) return(NULL)
        this_arg <- unname(y[!is.na(x)])
        list(arg = this_arg, value = unname(x[!is.na(x)]))
      }
    )
    nas <- map_lgl(datalist, is.null)
    n_null <- sum(nas)
    if (n_null > 0) {
      warn_na_entries_created(which(nas))
    }
    arg <- numeric(0)
    class <- "tfd_irreg"
  } else {
    nas <- map_lgl(datalist, \(x) is.null(x) || allMissing(x))
    if (any(nas)) {
      warn_na_entries_created(which(nas))
    }
    datalist <- map_if(datalist, nas, \(x) NULL)
    arg <- list(arg[[1]])
    class <- "tfd_reg"
  }

  if (!is.null(names(datalist))) {
    # ensure "unique" names
    names(datalist) <- vec_as_names(names(datalist), repair = "unique")
  }

  ret <- new_vctr(
    datalist,
    arg = arg,
    domain = domain,
    evaluator = evaluator_f,
    evaluator_name = evaluator,
    class = c(class, "tfd", "tf")
  )
  assert_arg(tf_arg(ret), ret)
  ret
}

#------------------------------------------------------------------------------

#' Constructors for vectors of "raw" functional data
#'
#' Various constructor methods for `tfd`-objects.\cr
#' `tfd` objects contain vectors of function evaluations at observed `arg`-values,
#' either all at the same `arg`-values (`tfd_reg`) or at different `arg`-values (`tfd_irreg`).
#' `NA`-functions are represented by `NULL`-entries in that list.
#'
#' @details
#' `tfd`-objects are list-`vctrs` of numeric vectors containing function
#' evaluations.
#'
#' **`evaluator`**: must be the (quoted or bare) name of a
#' function with signature `function(x, arg, evaluations)` that returns
#' the functions' (approximated/interpolated) values at locations `x` based on
#' the function `evaluations` available at locations `arg`.\cr
#' Available `evaluator`-functions:
#' - `tf_approx_linear` for linear interpolation without extrapolation (i.e.,
#' [zoo::na.approx()] with `na.rm = FALSE`)  -- this is the default,
#' - `tf_approx_spline` for cubic spline interpolation, (i.e., [zoo::na.spline()]
#' with `na.rm = FALSE`),
#' - `tf_approx_fill_extend` for linear interpolation and constant extrapolation
#' (i.e., [zoo::na.fill()] with `fill = "extend"`)
#' - `tf_approx_locf` for "last observation carried forward" (i.e.,
#' [zoo::na.locf()] with `na.rm = FALSE`)
#' - `tf_approx_nocb` for "next observation carried backward" (i.e.,
#' [zoo::na.locf()] with `na.rm = FALSE, fromLast = TRUE`).
#' See `tf:::zoo_wrapper` and `tf:::tf_approx_linear`, which is simply
#' `zoo_wrapper(zoo::na.approx, na.rm = FALSE)`, for examples of
#' implementations of this.
#'
#'
#' @param data a `matrix`, `data.frame` or `list` of suitable shape, or another
#'   `tf`-object. when this argument is `NULL` (i.e. when calling `tfd()`) this
#'   returns a prototype of class `tfd`.
#' @param ... not used in `tfd`, except for `tfd.tf` -- specify `arg` and
#'   `interpolate = TRUE` to turn an irregular `tfd` into a regular one, see
#'   examples.
#' @returns a `tfd`-object (or a `data.frame`/`matrix` for the conversion
#'   functions, obviously).
#' @family tfd-class
#' @export
tfd <- function(data, ...) UseMethod("tfd")

#' @export
#' @rdname tfd
#' @description `tfd.matrix` accepts a numeric matrix with one function per
#'   *row* (!). If `arg` is not provided, it tries to guess `arg` from the
#'   column names and falls back on `1:ncol(data)` if that fails.
#' @param arg For the `list`- and `matrix`-methods:
#'   `numeric`, or list of `numeric`s. The evaluation grid.
#'   For the `data.frame`-method: the
#'   name/number of the column defining the evaluation grid.
#'   The `matrix` method
#'   will try to guess suitable `arg`-values from the column names of `data` if
#'   `arg` is not supplied. Other methods fall back on integer sequences
#'   (`1:<length of data>`) as the default if not provided.
#' @param domain range of the `arg`.
#' @param evaluator a function accepting arguments `x, arg, evaluations`. See
#'   details for [tfd()].
tfd.matrix <- function(
  data,
  arg = NULL,
  domain = NULL,
  evaluator = tf_approx_linear,
  ...
) {
  assert_numeric(data)
  evaluator <- as_name(enexpr(evaluator))
  arg <- find_arg(data, arg) # either arg or numeric colnames or 1:ncol
  id <- unique_id(rownames(data) %||% seq_len(nrow(data)))
  # make factor conversion explicit to avoid reordering
  datalist <- split(data, factor(id, unique(as.character(id))))
  names(datalist) <- rownames(data)
  # don't count as irregular if entire rows are NA and nowhere else:
  irregular <- {
    na_rows <- rowSums(is.na(data)) == ncol(data)
    data_ <- data[!na_rows, ]
    anyNA(data_)
  }
  new_tfd(arg, datalist, !irregular, domain, evaluator)
}

#' @rdname tfd
#' @export
tfd.numeric <- function(
  data,
  arg = NULL,
  domain = NULL,
  evaluator = tf_approx_linear,
  ...
) {
  evaluator <- as_name(enexpr(evaluator))
  data <- t(as.matrix(data))
  # dispatch to matrix method
  args <- list(data, arg = arg, domain = domain, evaluator = evaluator)
  do.call(tfd, args)
}

#' @description `tfd.data.frame` uses the first 3 columns of `data` for
#'   `id` (function ID), `arg` (argument value) and `value` (function value)
#'   by default.
#' @export
#' @rdname tfd
#' @param id The name or number of the column defining which data belong to
#'   which function.
#' @param value The name or number of the column containing the function
#'   evaluations.
tfd.data.frame <- function(
  data,
  id = 1,
  arg = 2,
  value = 3,
  domain = NULL,
  evaluator = tf_approx_linear,
  ...
) {
  assert_numeric(data[[arg]])
  assert_numeric(data[[value]])

  evaluator <- as_name(enexpr(evaluator))
  # keep observations with NA values -- otherwise this never
  # creates NA-functions and risks dropping entire id-levels!
  keep <- which(!(is.na(data[, id]) | is.na(data[, arg])))
  data <- data[keep, c(id, arg, value)]

  # make factor conversion explicit to avoid reordering
  id <- factor(data[[1]], levels = as.factor(unique(data[[1]])))
  datalist <- split(data[[3]], id)
  arg <- split(data[[2]], id)

  # regular data always has non-NA values at all the same args:
  regular <- {
    data_ <- na.omit(data)
    arg_ <- split(data_[[2]], factor(data_[[1]])) # drop missing id levels
    length(arg_) == 1 || all(duplicated(arg_)[-1])
  }

  new_tfd(arg, datalist, regular, domain, evaluator)
}

# TODO this will break for multivariate data!
#' @description `tfd.list` accepts a list of vectors of identical lengths
#' containing evaluations or a list of 2-column matrices/data.frames with
#' `arg` in the first and evaluations in the second column
#' @export
#' @rdname tfd
tfd.list <- function(
  data,
  arg = NULL,
  domain = NULL,
  evaluator = tf_approx_linear,
  ...
) {
  evaluator <- as_name(enexpr(evaluator))
  vectors <- map_lgl(data, \(x) is.null(x) || (is.numeric(x) & !is.array(x)))
  if (all(vectors)) {
    where_na <- map(data, is.na)
    data <- map2(data, where_na, \(x, y) x[!y])
    lens <- lengths(data)
    empty <- lens != 0
    regular <- all(lens[!empty] == lens[!empty][1]) &
      (is.numeric(arg) || all(duplicated(arg)[-1]))
    # duplicated(NULL) == TRUE!
    if (!regular) {
      if (is.null(arg)) {
        cli::cli_abort("{.arg arg} cannot be NULL")
      }
      if (length(arg) != length(data)) {
        cli::cli_abort(
          "Length of {.arg arg} list does not match {.arg data} list."
        )
      }
      if (any(lengths(arg)[!empty] != lengths(where_na)[!empty])) {
        cli::cli_abort(
          "Lengths of {.arg arg} vectors do not match lengths of {.arg data} list entries."
        )
      }
      arg <- map2(arg, where_na, \(x, y) x[!y])
    } else {
      if (is.null(arg)) {
        cli::cli_warn("No {.arg arg} values supplied, using index positions.")
        arg <- map(data, seq_along)
      }
      arg <- ensure_list(arg)
      assert_numeric(
        arg[[1]],
        finite = TRUE,
        any.missing = FALSE,
        sorted = TRUE
      )
    }
  }
  if (!any(vectors)) {
    dims <- map(data, dim)
    if (any(lengths(dims) != 2) || any(map_int(dims, 2) != 2)) {
      cli::cli_abort("{.arg data} cannot be formatted into dimension 2.")
    }
    if (!all(rapply(data, is.numeric))) {
      cli::cli_abort("{.arg data} must be numeric.")
    }
    arg <- map(data, \(x) unlist(x[, 1], use.names = FALSE))
    data <- map(data, \(x) unlist(x[, 2], use.names = FALSE))
    regular <- length(data) == 1 || all(duplicated(arg)[-1])
  }
  new_tfd(arg, data, regular = regular, domain = domain, evaluator = evaluator)
}

#' @export
#' @examples
#' # turn irregular to regular tfd by evaluating on a common grid:
#'
#' f <- c(
#'   tf_rgp(1, arg = seq(0, 1, length.out = 11)),
#'   tf_rgp(1, arg = seq(0, 1, length.out = 21))
#' )
#' tfd(f, arg = seq(0, 1, length.out = 21))
#'
#' set.seed(1213)
#' f <- tf_rgp(3, arg = seq(0, 1, length.out = 51)) |> tf_sparsify(0.9)
#' # does not yield regular data because linear extrapolation yields NAs
#' #   outside observed range:
#' tfd(f, arg = seq(0, 1, length.out = 101))
#' # this "works" (but may not yield sensible values..!!) for
#' #   e.g. constant extrapolation:
#' tfd(f, evaluator = tf_approx_fill_extend, arg = seq(0, 1, length.out = 101))
#' plot(f, col = 2)
#' tfd(f,
#'   arg = seq(0, 1, length.out = 151), evaluator = tf_approx_fill_extend
#' ) |> lines()
#' @rdname tfd
tfd.tf <- function(data, arg = NULL, domain = NULL, evaluator = NULL, ...) {
  evaluator_name <- enexpr(evaluator)
  evaluator <- if (is_tfd(data) && is.null(evaluator)) {
    attr(data, "evaluator_name")
  } else if (is.null(evaluator)) {
    "tf_approx_linear"
  } else {
    as_name(evaluator_name)
  }
  domain <- (domain %||% unlist(arg, use.names = FALSE) %||% tf_domain(data)) |>
    range()
  re_eval <- !is.null(arg)
  na_mask <- is.na(data)
  arg <- ensure_list(arg %||% tf_arg(data))

  if (any(na_mask) && re_eval) {
    # process only non-NA entries, keep NULLs for NA entries
    evaluator_f <- get(evaluator, mode = "function", envir = parent.frame())
    if (all(na_mask)) {
      evaluations <- vector("list", length(data))
      evaluations[] <- list(NULL)
    } else {
      # subset arg for per-function arg lists (length > 1) to match non-NA data
      arg_for_eval <- if (length(arg) > 1) arg[!na_mask] else arg
      non_na_evals <- tf_evaluate(
        data[!na_mask],
        arg = arg_for_eval,
        evaluator = evaluator_f
      )
      evaluations <- vector("list", length(data))
      evaluations[!na_mask] <- non_na_evals
      evaluations[na_mask] <- list(NULL)
    }
  } else if (re_eval) {
    evaluator_f <- get(evaluator, mode = "function", envir = parent.frame())
    evaluations <- tf_evaluate(data, arg = arg, evaluator = evaluator_f)
  } else {
    evaluations <- tf_evaluations(data)
  }
  # handle NAs within non-NULL evaluations (e.g. from interpolation outside domain)
  non_null <- !map_lgl(evaluations, is.null)
  nas <- map(evaluations[non_null], \(x) which(is.na(x)))
  if (re_eval && any(lengths(nas))) {
    n <- length(evaluations)
    was_shared <- length(arg) == 1

    # normalize arg to per-function (length n) for uniform processing
    if (was_shared) {
      full_arg <- vector("list", n)
      full_arg[non_null] <- list(arg[[1]])
      full_arg[!non_null] <- list(numeric(0))
    } else {
      full_arg <- arg
    }

    # extract NA arg positions for warning (before pruning)
    na_arg_vals <- map2(full_arg[non_null], nas, \(x, y) x[y])
    same_nas <- all(duplicated(na_arg_vals)[-1])
    n_na <- length(unlist(nas, use.names = FALSE))

    if (!same_nas) {
      cli::cli_warn(c(
        i = "{n_na} evaluations were {.code NA}",
        x = "Returning irregular {.cls tfd}."
      ))
    } else {
      na_arg_string <- prettyNum(na_arg_vals[[1]]) |> paste(collapse = ", ")
      if (nchar(na_arg_string) > options()$width) {
        na_arg_string <- substr(na_arg_string, 1, options()$width - 15) |>
          paste0("[... truncated]")
      }
      cli::cli_warn(c(
        i = "All {n_na} evaluations on {.code arg = ({na_arg_string})} were {.code NA}",
        x = "Returning regular data {.cls tfd_reg} on the reduced grid."
      ))
    }

    # prune NAs from arg and evaluations (non-null entries only)
    full_arg[non_null] <- map2(
      full_arg[non_null],
      nas,
      \(x, y) if (length(y)) x[-y] else x
    )
    evaluations[non_null] <- map2(
      evaluations[non_null],
      nas,
      \(x, y) if (length(y)) x[-y] else x
    )

    # collapse back to shared arg only if originally shared and NAs were uniform
    arg <- if (same_nas && was_shared) full_arg[non_null][1] else full_arg
  }
  names(evaluations) <- names(data)
  new_tfd(
    arg,
    evaluations,
    regular = (length(arg) == 1),
    domain = domain,
    evaluator = evaluator
  )
}

#' @rdname tfd
#' @description `tfd.default` returns class prototype when argument to tfd() is
#'   `NULL` or not a recognised class.
#' @export
tfd.default <- function(
  data,
  arg = NULL,
  domain = NULL,
  evaluator = tf_approx_linear,
  ...
) {
  if (!missing(data)) {
    cli::cli_warn(
      "Input {.arg data} not a recognized class; returning prototype of length 0."
    )
  }
  datalist <- list()
  evaluator <- as_name(enexpr(evaluator))
  new_tfd(
    arg = arg,
    datalist = datalist,
    domain = domain,
    regular = TRUE,
    evaluator = evaluator
  )
}

#-------------------------------------------------------------------------------

#' @rdname tfd
#' @export
as.tfd <- function(data, ...) UseMethod("as.tfd")

#' @export
as.tfd.default <- function(data, ...) {
  tfd(data, ...)
}

#' @rdname tfd
#' @description `as.tfd_irreg` converts regular `tfd` or `tfb` objects into
#' irregular ones. Mainly used internally for `tf_rebase` operations etc.
#' @export
as.tfd_irreg <- function(data, ...) UseMethod("as.tfd_irreg")

#' @export
as.tfd_irreg.tfd_reg <- function(data, ...) {
  arg <- ensure_list(tf_arg(data))
  ret <- map2(tf_evaluations(data), arg, \(x, y) {
    if (is.null(x)) return(NULL)
    list(arg = y, value = x)
  })
  attributes(ret) <- attributes(data)
  attr(ret, "arg") <- numeric(0)
  class(ret)[1] <- "tfd_irreg"
  ret
}

#' @export
as.tfd_irreg.tfd_irreg <- function(data, ...) {
  data
}

#' @export
as.tfd_irreg.tfb <- function(data, ...) {
  tfd(data) |> as.tfd_irreg()
}