Name matchData()'s arguments for what they are

NEWS.md

@@ -71,6 +71,27 @@
 
 ## Breaking changes
 
+* **`matchData()`'s arguments are now `dat` and `source`**, replacing
+  `speciesDat` and `envDat`.
+
+  The function was never specific to species observations or environmental data.
+  It is a spatiotemporal nearest-feature join between two `sf` point objects, and
+  works as well for tag positions against a model field, moorings against
+  satellite retrievals, or one gridded product against another. The old names
+  described one use of it as though it were the only one.
+
+  **The old names still work**, with a warning, so existing scripts and
+  `taupatch` keep running. They will be removed in a later version.
+
+  ```r
+  matchData(observations, env)                 # positional, unchanged
+  matchData(dat = observations, source = env)  # new names
+  matchData(speciesDat = obs, envDat = env)    # still works, warns
+  ```
+
+  One related change: a `source` column whose name collides with one already in
+  `dat` is now suffixed **`.matched`** rather than `.env`, for the same reason.
+
 * **The `BigelowLab/copernicus` dependency is gone.** It was used in two places,
   both in `accessEnvDat()`, and both are now internal. This removes the
   `Remotes:` field and the hand-created `~/.copernicusdata` file that a new

R/accessEnvDat.R

@@ -347,7 +347,7 @@ read_day <- function(item, vars) {
 #' @param n_workers <integer> how many days to download at once. See the
 #'   Downloading in parallel section. Use `n_workers = 1` to download one day at
 #'   a time.
-#' @return envDat <sf object> sf object containing requested environmental data from Copernicus Marine Service
+#' @return <sf object> sf object containing requested environmental data from Copernicus Marine Service
 #' @export
 accessEnvDat <- function(product_id = NULL, dataset_id = NULL, vars,
                          years = NULL, months = NULL,

R/matchData.R

@@ -1,108 +1,158 @@
 
-#' Match environmental data to species occurrence data
+#' Match one set of points to another in space and time
 #'
-#' Joins each species observation to the nearest environmental grid point within
-#' the same time period. The time period is the environmental data's own temporal
-#' resolution: daily products match on year/month/day, monthly products (e.g.
-#' Copernicus `...P1M-m` means) on year/month, and annual products on year alone.
+#' Joins each row of `dat` to the nearest feature of `source` within the same
+#' time period, and returns `dat` with `source`'s columns added.
 #'
-#' Matching at the environmental data's native resolution matters because a
-#' day-exact join against monthly data matches nothing - a monthly product carries
-#' one time step per month, while observations fall on arbitrary days.
+#' Neither side has to be species observations or environmental data. It is a
+#' spatiotemporal nearest-feature join between two `sf` point objects that carry
+#' `YEAR`/`MONTH`/`DAY` columns, so it works equally for stations against a
+#' covariate grid, tag positions against a model field, moorings against
+#' satellite retrievals, or one gridded product against another.
 #'
-#' @param speciesDat <sf object> species observation data (e.g., presence, density, count);
-#'                               must have spatial and temporal components. Needs
-#'                               year and month columns, plus a day column when
-#'                               matching at daily resolution.
-#' @param envDat <sf object> environmental data accessed using the datamatch::accessEnvDat
-#'                           function to be matched to the species observation data
-#' @param temporal_resolution <char> one of "auto" (default), "day", "month", or
-#'                                   "year". "auto" infers the resolution from
-#'                                   `envDat`'s own time steps.
-#' @return <sf object> `speciesDat` with the matched environmental variables joined on,
-#'                     one row per input observation, plus LON/LAT coordinate columns.
-#'                     Observations in a period with no environmental data get NA
-#'                     for the environmental variables, and a warning is issued.
-#'                     An environmental variable whose name collides with a column
-#'                     already in `speciesDat` is suffixed `.env`.
+#' @section Matching in time:
+#' The time period is `source`'s own resolution: daily data matches on
+#' year/month/day, monthly data (Copernicus `...P1M-m` means, say) on
+#' year/month, and annual data on year alone.
+#'
+#' That matters because a day-exact join against monthly data matches nothing. A
+#' monthly product carries one time step per month, while observations fall on
+#' arbitrary days. `temporal_resolution` overrides the inference when the data
+#' cannot speak for itself.
+#'
+#' @section What is preserved:
+#' One row out per row of `dat`, in the same order, whatever happens. A period
+#' `source` does not cover gives `NA` for its columns and a warning naming the
+#' periods, rather than dropping those rows — a silent change in row count is a
+#' worse outcome than a visible gap.
+#'
+#' `dat` keeps its own columns. One of `source`'s that collides with a name
+#' already in `dat` is suffixed `.matched`, so nothing of `dat`'s is overwritten
+#' or renamed.
+#'
+#' @param dat <sf object> the points to add columns to: observations, stations,
+#'   tag positions, anything with coordinates and time. Needs year and month
+#'   columns, plus a day column when matching at daily resolution. Columns whose
+#'   names begin with those words are recognised, so `Year` or `obs_month` work.
+#' @param source <sf object> the points to take values from, typically a grid
+#'   from [accessEnvDat()]. Must carry `YEAR`/`MONTH`/`DAY`.
+#' @param temporal_resolution <char> one of `"auto"` (default), `"day"`,
+#'   `"month"`, or `"year"`. `"auto"` uses the step `accessEnvDat()` recorded on
+#'   `source`, or infers it from `source`'s time steps.
+#' @param speciesDat,envDat deprecated names for `dat` and `source`. Still
+#'   accepted, with a warning.
+#' @return <sf object> `dat` with `source`'s columns joined on, one row per input
+#'   row, plus `LON`/`LAT` coordinate columns.
+#' @examples
+#' \dontrun{
+#' env <- accessEnvDat(vars = "SST", years = 2010, months = 1:12, bounding_box = bb)
+#'
+#' matched <- matchData(observations, env)
+#'
+#' # Chains, so several sources land on one table
+#' matched <- matchData(matched, chlorophyll)
+#' }
+#' @seealso [accessEnvDat()] for the usual `source`, [attach_bathymetry()] and
+#'   [attach_climate_index()] for covariates that are not matched this way
 #' @export
-matchData <- function(speciesDat, envDat,
-                      temporal_resolution = c("auto", "day", "month", "year")) {
+matchData <- function(dat, source,
+                      temporal_resolution = c("auto", "day", "month", "year"),
+                      speciesDat = NULL, envDat = NULL) {
+
+  # The old names were specific to one use of a function that was never specific
+  # to it. Accepted for now because taupatch and any script written against the
+  # old signature call them by name, and breaking those silently would be worse
+  # than carrying two lines.
+  if (!is.null(speciesDat)) {
+    warning("`speciesDat` is now `dat`. The old name still works but will be ",
+            "removed.", call. = FALSE)
+    if (missing(dat)) dat <- speciesDat
+  }
+  if (!is.null(envDat)) {
+    warning("`envDat` is now `source`. The old name still works but will be ",
+            "removed.", call. = FALSE)
+    if (missing(source)) source <- envDat
+  }
+  if (missing(dat) || missing(source)) {
+    if (is.null(speciesDat) || is.null(envDat)) {
+      stop("Both `dat` and `source` are required.", call. = FALSE)
+    }
+  }
 
   temporal_resolution <- match.arg(temporal_resolution)
   if (temporal_resolution == "auto") {
-    temporal_resolution <- detect_temporal_resolution(envDat)
+    temporal_resolution <- detect_temporal_resolution(source)
   }
   match_keys <- switch(temporal_resolution,
                        day = c("YEAR", "MONTH", "DAY"),
                        month = c("YEAR", "MONTH"),
                        year = "YEAR")
 
-  speciesDat <- standardize_time_columns(speciesDat, match_keys)
+  dat <- standardize_time_columns(dat, match_keys)
 
-  env_geom <- attr(envDat, "sf_column")
-  env_vars <- setdiff(names(envDat), c("YEAR", "MONTH", "DAY", env_geom))
+  source_geom <- attr(source, "sf_column")
+  source_vars <- setdiff(names(source), c("YEAR", "MONTH", "DAY", source_geom))
 
   # Both sides need a CRS before they can be reconciled. Without one,
   # st_transform() fails with "crs not found: is it missing?", which is true but
   # does not say which object or what to do. Silently assuming a CRS would be
   # worse: coordinates would be matched as though they were degrees, and every
-  # observation would join to whichever cell happened to be nearest in a
-  # meaningless space.
-  for (side in list(list(x = speciesDat, name = "speciesDat"),
-                    list(x = envDat, name = "envDat"))) {
+  # row would join to whichever feature happened to be nearest in a meaningless
+  # space.
+  for (side in list(list(x = dat, name = "dat"),
+                    list(x = source, name = "source"))) {
     if (is.na(sf::st_crs(side$x))) {
-      stop(side$name, " has no coordinate reference system, so it cannot be ",
-           "matched.\nSet one with sf::st_crs(", side$name,
+      stop("`", side$name, "` has no coordinate reference system, so it cannot ",
+           "be matched.\nSet one with sf::st_crs(", side$name,
            ") <- 4326 for longitude/latitude,\nor the EPSG code the ",
            "coordinates are actually in.", call. = FALSE)
     }
   }
 
-  envDat <- sf::st_transform(envDat, sf::st_crs(speciesDat))
+  source <- sf::st_transform(source, sf::st_crs(dat))
 
-  # Give environmental variables that share a name with a species column an
-  # explicit ".env" suffix. Otherwise st_join() disambiguates them as ".x"/".y",
-  # which both renames the species column and makes the result's column names
+  # Give a source column that shares a name with one in `dat` an explicit
+  # ".matched" suffix. Otherwise st_join() disambiguates them as ".x"/".y",
+  # which both renames the caller's column and makes the result's column names
   # depend on whether a given period actually matched anything.
-  collisions <- intersect(env_vars, names(speciesDat))
+  collisions <- intersect(source_vars, names(dat))
   if (length(collisions) > 0) {
-    renamed <- paste0(collisions, ".env")
-    names(envDat)[match(collisions, names(envDat))] <- renamed
-    env_vars[match(collisions, env_vars)] <- renamed
+    renamed <- paste0(collisions, ".matched")
+    names(source)[match(collisions, names(source))] <- renamed
+    source_vars[match(collisions, source_vars)] <- renamed
   }
 
-  periods <- unique(sf::st_drop_geometry(speciesDat)[match_keys])
+  periods <- unique(sf::st_drop_geometry(dat)[match_keys])
   matched <- vector("list", nrow(periods))
   unmatched_periods <- character()
 
   for (i in seq_len(nrow(periods))) {
-    in_period <- rep(TRUE, nrow(speciesDat))
-    env_in_period <- rep(TRUE, nrow(envDat))
+    in_period <- rep(TRUE, nrow(dat))
+    source_in_period <- rep(TRUE, nrow(source))
     for (key in match_keys) {
-      in_period <- in_period & speciesDat[[key]] == periods[[key]][i]
-      env_in_period <- env_in_period & envDat[[key]] == periods[[key]][i]
+      in_period <- in_period & dat[[key]] == periods[[key]][i]
+      source_in_period <- source_in_period & source[[key]] == periods[[key]][i]
     }
 
-    obs <- speciesDat[in_period, ]
-    env_slice <- envDat[env_in_period, c(env_vars)]
+    rows <- dat[in_period, ]
+    source_slice <- source[source_in_period, c(source_vars)]
 
-    if (nrow(env_slice) == 0) {
-      # st_nearest_feature cannot join against an empty set, so fill the
-      # environmental columns with NA rather than dropping the observations.
-      # Dropping them would silently change the row count of the result.
-      for (v in env_vars) obs[[v]] <- NA
+    if (nrow(source_slice) == 0) {
+      # st_nearest_feature cannot join against an empty set, so fill the matched
+      # columns with NA rather than dropping the rows. Dropping them would
+      # silently change the row count of the result.
+      for (v in source_vars) rows[[v]] <- NA
       unmatched_periods <- c(unmatched_periods,
                              paste(unlist(periods[i, , drop = TRUE]), collapse = "-"))
-      matched[[i]] <- obs
+      matched[[i]] <- rows
     } else {
-      matched[[i]] <- sf::st_join(obs, env_slice, join = sf::st_nearest_feature)
+      matched[[i]] <- sf::st_join(rows, source_slice, join = sf::st_nearest_feature)
     }
   }
 
   if (length(unmatched_periods) > 0) {
-    warning("No environmental data for ", length(unmatched_periods),
-            " period(s); environmental variables set to NA for: ",
+    warning("No data in `source` for ", length(unmatched_periods),
+            " period(s); matched columns set to NA for: ",
             paste(utils::head(unmatched_periods, 5), collapse = ", "),
             if (length(unmatched_periods) > 5) ", ..." else "", call. = FALSE)
   }
@@ -120,7 +170,7 @@ matchData <- function(speciesDat, envDat,
   matched_data
 }
 
-#' Infer the temporal resolution of environmental data
+#' Infer the temporal resolution of a set of time steps
 #'
 #' Reads the resolution off the time steps actually present: more than one day
 #' within any month means daily data, and more than one month within any year
@@ -135,21 +185,22 @@ matchData <- function(speciesDat, envDat,
 #' fine leaves them unmatched and warns. Pass `temporal_resolution` explicitly to
 #' override.
 #'
-#' @param envDat <sf object> environmental data with YEAR/MONTH/DAY columns
+#' @param x <sf object> an object with YEAR/MONTH/DAY columns, typically the
+#'   `source` side of a match
 #' @return one of "day", "month", or "year"
 #' @keywords internal
-detect_temporal_resolution <- function(envDat) {
+detect_temporal_resolution <- function(x) {
   # accessEnvDat() knows which dataset it fetched, so it records the step rather
   # than leaving it to be inferred. Worth trusting over the heuristics below: a
   # `dates` request of one date per month is genuinely indistinguishable from
   # monthly data by inspection, and guessing monthly would drop the day from the
   # match.
-  recorded <- attr(envDat, "datamatch_step")
+  recorded <- attr(x, "datamatch_step")
   if (!is.null(recorded) && recorded %in% c("day", "month", "year")) {
     return(recorded)
   }
 
-  times <- unique(sf::st_drop_geometry(envDat)[c("YEAR", "MONTH", "DAY")])
+  times <- unique(sf::st_drop_geometry(x)[c("YEAR", "MONTH", "DAY")])
 
   days_per_month <- tapply(times$DAY, paste(times$YEAR, times$MONTH), function(d) length(unique(d)))
   if (any(days_per_month > 1)) return("day")
@@ -163,28 +214,28 @@ detect_temporal_resolution <- function(envDat) {
   "month"
 }
 
-#' Rename a species dataset's time columns to YEAR/MONTH/DAY
+#' Rename a table's time columns to YEAR/MONTH/DAY
 #'
 #' Only the columns needed for the requested match keys are required, so monthly
 #' matching works on data that has no day column at all.
 #'
-#' @param speciesDat <sf object> species observation data
+#' @param dat <sf object> the table being matched
 #' @param match_keys <char> the standardized time columns needed, e.g. c("YEAR", "MONTH")
-#' @return `speciesDat` with its time columns renamed to YEAR/MONTH/DAY
+#' @return `dat` with its time columns renamed to YEAR/MONTH/DAY
 #' @keywords internal
-standardize_time_columns <- function(speciesDat, match_keys) {
+standardize_time_columns <- function(dat, match_keys) {
   # sf's select() method keeps the geometry column "sticky" regardless of the
   # select criteria, so it must be excluded here - otherwise the rename() below
   # would rename the geometry column itself and corrupt the sf object's tracked
   # geometry-column name.
-  geom_col <- attr(speciesDat, "sf_column")
+  geom_col <- attr(dat, "sf_column")
 
   for (key in match_keys) {
-    if (key %in% names(speciesDat)) next
+    if (key %in% names(dat)) next
 
     prefix <- tolower(key)
     candidates <- setdiff(
-      names(speciesDat |> dplyr::select(dplyr::starts_with(prefix, ignore.case = TRUE))),
+      names(dat |> dplyr::select(dplyr::starts_with(prefix, ignore.case = TRUE))),
       geom_col
     )
     # An exact match wins over a mere prefix match, so a dataset carrying both
@@ -195,17 +246,17 @@ standardize_time_columns <- function(speciesDat, match_keys) {
     }
 
     if (length(candidates) == 0) {
-      stop("speciesDat has no column for '", key, "' (looked for names starting with '",
+      stop("`dat` has no column for '", key, "' (looked for names starting with '",
            prefix, "'). It is required to match at this temporal resolution.",
            call. = FALSE)
     }
     if (length(candidates) > 1) {
-      stop("speciesDat has multiple candidate '", key, "' columns: ",
+      stop("`dat` has multiple candidate '", key, "' columns: ",
            paste(candidates, collapse = ", "),
            ". Rename the intended one to '", key, "'.", call. = FALSE)
     }
-    names(speciesDat)[names(speciesDat) == candidates] <- key
+    names(dat)[names(dat) == candidates] <- key
   }
 
-  speciesDat
+  dat
 }

R/plot.R

@@ -223,7 +223,7 @@ plot_series <- function(env_dat, vars = NULL, fun = mean, spread = TRUE, ...) {
 #' @return the plotted values, invisibly
 #' @examples
 #' \dontrun{
-#' matched <- matchData(speciesDat = observations, envDat = env)
+#' matched <- matchData(observations, env)
 #'
 #' plot_matched(matched, "SST")
 #' # Open circles are observations that matched nothing.