From 0479b9204911c61e9e30aade03c86e56c83d2730 Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Tue, 22 Sep 2026 12:47:04 -0400 Subject: [PATCH 1/8] local implementation --- DESCRIPTION | 1 + R/NSR_local.R | 474 +++++++++++++++++++++++++++ R/NSR_local_set.R | 293 +++++++++++++++++ R/local_build.R | 239 ++++++++++++++ R/local_cache.R | 133 ++++++++ R/local_import.R | 219 +++++++++++++ R/local_regions.R | 199 ++++++++++++ dev_notes/01-offline-nsr-design.md | 494 +++++++++++++++++++++++++++++ tests/testthat/test-local-nsr.R | 138 ++++++++ 9 files changed, 2190 insertions(+) create mode 100644 R/NSR_local.R create mode 100644 R/NSR_local_set.R create mode 100644 R/local_build.R create mode 100644 R/local_cache.R create mode 100644 R/local_import.R create mode 100644 R/local_regions.R create mode 100644 dev_notes/01-offline-nsr-design.md create mode 100644 tests/testthat/test-local-nsr.R diff --git a/DESCRIPTION b/DESCRIPTION index 3c00a79..019ef2b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -17,6 +17,7 @@ Imports: httr Suggests: knitr, + data.table, rmarkdown, testthat, devtools, diff --git a/R/NSR_local.R b/R/NSR_local.R new file mode 100644 index 0000000..c711254 --- /dev/null +++ b/R/NSR_local.R @@ -0,0 +1,474 @@ +#' Determine native status without an internet connection +#' +#' Offline Native Status Resolver. For each taxon and place, the status every built +#' checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s columns, +#' with four added. +#' +#' \strong{Place.} Give coordinates (\code{latitude}, \code{longitude}), political +#' division names (\code{country}, \code{state_province}, \code{county_parish}), or both. +#' Coordinates are the better input: each source is consulted in the geography it +#' publishes against (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil +#' against GADM states), with no crosswalk between them. Names are resolved to GADM +#' units through the GNRS backbone and then carried to other geographies by the spatial +#' link table. +#' +#' \strong{Precedence.} If any source says native, the answer is native: asserting +#' nativity is a positive claim, whereas "introduced" and "absent" are often artefacts of +#' a list's scope, age or purpose. Disagreement is recorded, not hidden, in +#' \code{native_status_conflict} and \code{native_status_opinions}. +#' +#' \strong{Which polygons answer.} A place is judged by the polygons it lies IN. An +#' opinion about a polygon containing the place applies to it (POWO's finest statement +#' about Guadeloupe is "native in the Leeward Islands"), and among those any native +#' opinion wins. Polygons INSIDE the place describe only parts of it, so they are not its +#' status: they answer only when they all agree, and otherwise the answer is \code{P} with +#' the reason saying the status varies and how many sub-polygons say what +#' (\code{n_subpolygons_native}, \code{n_subpolygons_introduced}). Give coordinates and +#' the question does not arise: the record is judged on the ground it sits on. +#' \code{native_status_scope} records which of these produced the answer. +#' +#' \strong{Endemism is the exception}, deliberately: \code{Ne} and \code{Ie} are claims +#' about the taxon's whole range rather than about one polygon, so they draw on evidence +#' from elsewhere. A record of a taxon confined to California, found in Michigan, is +#' introduced there however little Michigan's own checklists say. +#' @param occurrence_dataframe A data.frame with \code{species}, and either coordinates +#' or political division names (see Place). +#' @param dir Cache directory, shared with GNRS and GVS. +#' @param resolve_names Resolve submitted names against WCVP with +#' \code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none. +#' @param min_overlap Ignore region links covering less than this share of a region. +#' @param quiet Suppress progress messages? +#' @return A data.frame, one row per input row. +#' @export +NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names = TRUE, + min_overlap = 0.01, quiet = FALSE) { + if (!inherits(occurrence_dataframe, "data.frame")) { + stop("occurrence_dataframe should be a data.frame", call. = FALSE) + } + x <- occurrence_dataframe + col <- function(nm) if (nm %in% names(x)) as.character(x[[nm]]) else rep(NA_character_, nrow(x)) + num <- function(nm) if (nm %in% names(x)) suppressWarnings(as.numeric(x[[nm]])) else rep(NA_real_, nrow(x)) + species <- col("species") + country <- col("country") + state <- col("state_province") + county <- col("county_parish") + lon <- num("longitude") + lat <- num("latitude") + n <- nrow(x) + db <- nsr_local_db(dir) + + # ---- taxa ------------------------------------------------------------------------- + # WCVP records distributions for genera too, so a genus query is answerable directly. + # What must not happen is a bare genus being RESOLVED to some species of that genus, + # which would answer confidently about the wrong taxon. + above_species <- !is.na(species) & nzchar(species) & + lengths(strsplit(trimws(species), "[[:space:]]+")) < 2 + taxon_id <- db$name_index$taxon_id[match(species, db$name_index$species_name)] + need <- unique(species[is.na(taxon_id) & !is.na(species) & nzchar(species) & !above_species]) + if (resolve_names && length(need)) { + if (!quiet) message("Resolving ", format(length(need), big.mark = ","), " names ...") + r <- nsr_resolve_names(need, quiet = quiet, species_index = db$name_index) + taxon_id[is.na(taxon_id)] <- r$taxon_id[match(species[is.na(taxon_id)], r$name)] + } + + # ---- place ------------------------------------------------------------------------ + places <- nsr_query_regions(lon, lat, country, state, county, dir, db) + + # ---- opinions --------------------------------------------------------------------- + use_set <- !is.null(db$chk_dt) && nrow(occurrence_dataframe) >= getOption("NSR.set_min", 200) + res <- if (use_set) nsr_resolve_status_set(taxon_id, places, db, min_overlap) else + nsr_resolve_status(taxon_id, places, db, min_overlap) + + out <- data.frame( + family = db$taxa$family[match(taxon_id, db$taxa$taxon_id)], + genus = db$taxa$genus[match(taxon_id, db$taxa$taxon_id)], + species = species, + country = country, state_province = state, county_parish = county, + latitude = lat, longitude = lon, + poldiv_full = places$label, + poldiv_type = places$level, + native_status_country = res$country_code, + native_status_state_province = res$state_code, + native_status_county_parish = rep(NA_character_, n), + native_status = res$code, + native_status_reason = res$reason, + native_status_sources = res$sources, + isIntroduced = as.integer(res$code %in% c("I", "Ie")), + isEndemic = as.integer(res$code == "Ne"), + isCultivatedNSR = res$cultivated, + # NA, not production's 0: the flag is dropped (a taxon-level use flag with no + # polygon attached says nothing about the record in hand) and a never-populated + # column should not read as a real negative. See dev_notes/01-offline-nsr-design.md, + # open question 3. isCultivatedNSR above is the version that has a geography. + is_cultivated_taxon = NA_integer_, + native_status_conflict = res$conflict, + native_status_conflict_type = res$conflict_type, + native_status_scope = res$scope, + n_subpolygons_native = res$n_sub_native, + n_subpolygons_introduced = res$n_sub_introduced, + native_status_opinions = res$opinions, + regions_matched = places$matched, + taxon_evaluable = !is.na(taxon_id) & taxon_id %in% db$evaluable, + stringsAsFactors = FALSE + ) + rownames(out) <- NULL + out +} + +nsr_session <- new.env(parent = emptyenv()) + +#' Load the local tables once per call +#' @keywords internal +#' @noRd +nsr_local_db <- function(dir) { + stamp <- paste(normalizePath(dir), file.mtime(nsr_table_path("checklist", dir)), + file.mtime(nsr_table_path("region-links", dir))) + hit <- nsr_session$db + if (!is.null(hit) && identical(hit$stamp, stamp)) return(hit$value) + need <- c("sources", "taxa", "regions", "checklist") + miss <- vapply(need, function(t) !file.exists(nsr_table_path(t, dir)), logical(1)) + if (any(miss)) { + stop("The local NSR is not built in ", dir, ".\nRun NSR_local_build().", call. = FALSE) + } + db <- lapply(need, function(t) as.data.frame(nanoparquet::read_parquet(nsr_table_path(t, dir)))) + names(db) <- need + lf <- nsr_table_path("region-links", dir) + db$links <- if (file.exists(lf)) as.data.frame(nanoparquet::read_parquet(lf)) else + data.frame(from_region = character(0), to_region = character(0), relation = character(0), + fraction = numeric(0), stringsAsFactors = FALSE) + db$link_idx <- split(seq_len(nrow(db$links)), db$links$from_region) + db$chk_idx <- split(seq_len(nrow(db$checklist)), + paste(db$checklist$taxon_id, db$checklist$region_key)) + db$chk_env <- list2env(db$chk_idx, envir = new.env(hash = TRUE, parent = emptyenv())) + db$link_env <- list2env(db$link_idx, envir = new.env(hash = TRUE, parent = emptyenv())) + if (requireNamespace("data.table", quietly = TRUE)) { + db$chk_dt <- data.table::as.data.table(db$checklist) + data.table::setkeyv(db$chk_dt, c("taxon_id", "region_key")) + db$links_dt <- data.table::as.data.table(db$links) + data.table::setkeyv(db$links_dt, "from_region") + nat <- unique(db$chk_dt[status == "native", list(taxon_id, region_key)]) + db$native_dt <- nat + data.table::setkeyv(db$native_dt, "taxon_id") + db$confined_dt <- nsr_confined_ranges(nat, db) + } + # WCVP sometimes carries the same name at two ranks (a variety row also called + # "Pinus ponderosa"), and only one of them holds the distributions: index names to the + # species-rank id, and among those to the one with opinions + tx <- db$taxa + tx$is_species <- tolower(tx$rank) %in% c("species", "") + tx$has_rows <- tx$taxon_id %in% db$checklist$taxon_id + tx <- tx[order(!tx$has_rows, !tx$is_species), , drop = FALSE] + db$name_index <- tx[!duplicated(tx$species_name), c("species_name", "taxon_id")] + db$chk_status <- db$checklist$status + db$chk_source <- db$checklist$source_name + db$chk_cult <- db$checklist$is_cultivated + db$chk_region <- db$checklist$region_key + db$link_to <- db$links$to_region + db$link_rel <- db$links$relation + db$link_frac <- db$links$fraction + db$chk_idx_taxon <- split(seq_len(nrow(db$checklist)), db$checklist$taxon_id) + db$evaluable <- unique(db$checklist$taxon_id) + db$comprehensive <- db$sources$source_name[db$sources$is_comprehensive %in% TRUE] + db$covered <- unique(db$checklist$region_key[db$checklist$source_name %in% db$comprehensive]) + nsr_session$db <- list(stamp = stamp, value = db) + db +} + +#' The regions a query refers to, in every system +#' +#' Internal. With coordinates, each system is looked up directly. With names, the GADM +#' unit comes from the GNRS backbone and other systems are reached through the link +#' table. Returns, per row, the regions to consult and how each relates to the query. +#' @keywords internal +#' @noRd +nsr_query_regions <- function(lon, lat, country, state, county, dir, db) { + n <- length(lon) + direct <- vector("list", n) + fine <- vector("list", n) + ctry <- vector("list", n) + label <- rep(NA_character_, n) + level <- rep("country", n) + matched <- rep("none", n) + + loc <- if (any(is.finite(lon) & is.finite(lat))) nsr_locate_regions(lon, lat, dir) else NULL + bb <- try(nsr_gnrs_backbone(dir), silent = TRUE) + pd <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, state, bb) else NULL + pd0 <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, NULL, bb) else NULL + + for (i in seq_len(n)) { + keys <- character(0) + if (!is.null(loc)) keys <- c(keys, stats::na.omit(c(loc$wgsrpd3[i], loc$gadm[i], loc$gadm0[i]))) + if (!is.null(pd)) { + gid1 <- if (!is.na(pd$state_province_id[i])) + bb$state$gid_1[match(pd$state_province_id[i], bb$state$state_province_id)] else NA_character_ + gid0 <- if (!is.na(pd0$country_id[i])) + bb$country$gid_0[match(pd0$country_id[i], bb$country$country_id)] else NA_character_ + if (!is.na(gid1)) keys <- c(keys, paste0("gadm1:", gid1)) + if (!is.na(gid0)) keys <- c(keys, paste0("gadm0:", gid0)) + } + keys <- unique(keys[!is.na(keys)]) + # the finest place the query actually names, kept apart from its country: a question + # about Amazonas must not inherit Brazil's answer + fine[[i]] <- grep("^gadm0:", keys, value = TRUE, invert = TRUE) + ctry[[i]] <- grep("^gadm0:", keys, value = TRUE) + direct[[i]] <- keys + has_state <- any(grepl("^gadm1:", keys)) || (!is.null(pd) && !is.na(pd$state_province_id[i])) + level[i] <- if (has_state) "state_province" else "country" + matched[i] <- if (!length(keys)) "none" else + if (!is.null(loc) && !is.na(loc$gadm[i])) "coordinates" else "names" + label[i] <- paste(stats::na.omit(c(country[i], if (!is.na(state[i]) && nzchar(state[i])) state[i])), + collapse = ":") + if (!nzchar(label[i]) && length(keys)) label[i] <- paste(keys, collapse = " + ") + } + list(direct = direct, fine = fine, country = ctry, label = label, level = level, + matched = matched) +} + +#' Gather and reduce every opinion bearing on each row +#' @keywords internal +#' @noRd +nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { + n <- length(taxon_id) + blank <- rep(NA_character_, n) + out <- list(code = blank, reason = blank, sources = blank, opinions = blank, + country_code = blank, state_code = blank, scope = blank, + conflict = rep(FALSE, n), conflict_type = rep(NA_character_, n), + cultivated = rep(NA_integer_, n), n_sub_native = rep(NA_integer_, n), + n_sub_introduced = rep(NA_integer_, n)) + for (i in seq_len(n)) { + tid <- taxon_id[i] + fine <- places$fine[[i]] + ctry <- places$country[[i]] + keys <- if (length(fine)) fine else ctry # answer at the finest level named + if (is.na(tid) || !length(c(fine, ctry))) { + out$code[i] <- "UNK" + out$reason[i] <- if (is.na(tid)) "Taxon not matched to a species in the backbone" else + "Place not matched to any region" + next + } + op <- nsr_opinions_for(tid, keys, db, min_overlap) + # a coarser polygon the place sits in (its country) also contains it, so opinions + # recorded ON that polygon apply; its OTHER sub-polygons do not, so only direct rows + # are taken, never its links + anc <- setdiff(ctry, keys) + if (length(anc)) { + rows <- unlist(mget(paste(tid, anc), db$chk_env, ifnotfound = list(NULL)), use.names = FALSE) + if (length(rows)) { + add <- list(status = db$chk_status[rows], source_name = db$chk_source[rows], + is_cultivated = db$chk_cult[rows], relation = rep("within", length(rows))) + op <- if (is.null(op)) add else Map(c, op, add) + } + } + consulted <- nsr_consulted_regions(keys, db) + ev <- tid %in% db$evaluable + r <- nsr_reduce(op, consulted, ev, db) + if (identical(r$code, "N") && nsr_is_endemic(tid, keys, db, min_overlap)) { + r$code <- "Ne" + r$reason <- paste0(sub("^Native", "Native and endemic", r$reason)) + } + if (identical(r$code, "A")) { + ee <- nsr_endemic_elsewhere(tid, keys, db, min_overlap) + if (!is.na(ee)) { + r$code <- "Ie" + r$reason <- paste0("Absent from this region and endemic to ", ee, + ", so introduced here (inferred)") + } + } + for (f in names(r)) out[[f]][i] <- r[[f]] + # per-level codes, for the service's columns + out$country_code[i] <- if (length(ctry)) { + nsr_reduce(nsr_opinions_for(tid, ctry, db, min_overlap), nsr_consulted_regions(ctry, db), ev, db)$code + } else NA_character_ + out$state_code[i] <- if (length(fine)) { + nsr_reduce(nsr_opinions_for(tid, fine, db, min_overlap), nsr_consulted_regions(fine, db), ev, db)$code + } else NA_character_ + } + out +} + +#' Opinions about one taxon bearing on a set of regions +#' +#' Internal. Direct opinions, plus opinions about regions the query's regions are inside +#' of or contain, with the relation recorded so the reducer can apply the one-way +#' inheritance rules. +#' @keywords internal +#' @noRd +nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01) { + if (!length(keys)) return(NULL) + rows <- unlist(mget(paste(taxon_id, keys), db$chk_env, ifnotfound = list(NULL)), use.names = FALSE) + rel <- if (length(rows)) rep("same", length(rows)) else character(0) + li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) + if (length(li)) { + to <- db$link_to[li]; rl <- db$link_rel[li]; fr <- db$link_frac[li] + known <- unique(sub(":.*", "", keys)) # geographies this place is located in + good <- !(to %in% keys) & fr >= min_overlap & !(sub(":.*", "", to) %in% known) + if (any(good)) { + to <- to[good]; rl <- rl[good] + idx <- mget(paste(taxon_id, to), db$chk_env, ifnotfound = list(NULL)) + n <- lengths(idx) + if (sum(n)) { + rows <- c(rows, unlist(idx, use.names = FALSE)) + rel <- c(rel, rep(rl, n)) + } + } + } + if (!length(rows)) return(NULL) + list(status = db$chk_status[rows], source_name = db$chk_source[rows], + is_cultivated = db$chk_cult[rows], relation = rel) +} + +#' Is the taxon native only inside the queried place? +#' +#' Internal. Endemism (\code{Ne}) is read off the checklist: every region any source +#' calls it native must be the queried region, inside it, or overlapping it - never a +#' region lying outside. Only claimed when at least one source has a native opinion. +#' @keywords internal +#' @noRd +nsr_is_endemic <- function(taxon_id, keys, db, min_overlap = 0.01) { + rows <- db$chk_idx_taxon[[as.character(taxon_id)]] + if (is.null(rows)) return(FALSE) + nat <- db$chk_region[rows][db$chk_status[rows] == "native"] + if (!length(nat)) return(FALSE) + li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) + inside <- if (!length(li)) keys else + c(keys, db$link_to[li][db$link_rel[li] %in% c("same", "contains") & db$link_frac[li] >= min_overlap]) + all(nat %in% inside) +} + +#' Is the taxon endemic to somewhere else? +#' +#' Internal. The service's `Ie`: a taxon absent from the queried region but whose whole +#' native range lies elsewhere and is confined - endemic to one region, or to one country - +#' cannot be native here, so an occurrence must be introduced. A widely native species +#' merely unrecorded here stays `A`: absence is not evidence of introduction unless the +#' species could not have been native in the first place. +#' @return The name of the region it is endemic to, or NA. +#' @keywords internal +#' @noRd +nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { + rows <- db$chk_idx_taxon[[as.character(taxon_id)]] + if (is.null(rows)) return(NA_character_) + nat <- unique(db$chk_region[rows][db$chk_status[rows] == "native"]) + if (!length(nat)) return(NA_character_) + if (any(nat %in% nsr_consulted_regions(keys, db))) return(NA_character_) + nm <- function(k) { + i <- match(k, db$regions$region_key) + if (is.na(i)) k else db$regions$region_name[i] + } + if (length(nat) == 1) return(nm(nat)) + containers <- lapply(nat, function(r) { + li <- db$link_env[[r]] + if (is.null(li)) return(character(0)) + keep <- db$link_rel[li] %in% c("same", "within") & db$link_frac[li] >= min_overlap & + startsWith(db$link_to[li], "gadm0:") + db$link_to[li][keep] + }) + common <- Reduce(intersect, containers) + if (length(common)) nm(common[1]) else NA_character_ +} + +#' Every region an answer may draw on: the query's own, and those linked to them +#' +#' Internal. Used to decide whether absence is interpretable: a query about California +#' finds no source publishing against GADM California, but WGSRPD's CAL is the same place +#' and is comprehensively listed, so absence there does mean something. +#' @keywords internal +#' @noRd +nsr_consulted_regions <- function(keys, db) { + if (!length(keys)) return(character(0)) + li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) + if (!length(li)) return(keys) + unique(c(keys, db$link_to[li][db$link_rel[li] %in% c("same", "within", "contains")])) +} + +#' What kind of disagreement is this? +#' +#' Internal. Two quite different things look like conflict. Sources can genuinely +#' disagree about a polygon (VASCAN calls a plant introduced in Ontario where POWO calls it +#' native). Or one source can hold both statuses for the same polygon because WCVP records +#' one subspecies as native and another as introduced, and we roll infraspecific taxa up to +#' the species - no one disagrees with anyone there. The flag says which. +#' @return "none", "sources", "within source" or "both". +#' @keywords internal +#' @noRd +nsr_conflict_type <- function(status, source) { + if (!length(status)) return("none") + by_src <- split(status, source) + within <- any(vapply(by_src, function(z) length(unique(z)) > 1, logical(1))) + sets <- vapply(by_src, function(z) paste(sort(unique(z)), collapse = "+"), character(1)) + between <- length(unique(sets)) > 1 + if (within && between) "both" else if (between) "sources" else if (within) "within source" else "none" +} + +#' Reduce opinions to one status code +#' +#' Internal. A place is judged by the polygons it lies IN: an opinion about a polygon +#' containing it applies to it, and among those, any native opinion wins (BM's +#' precedence). Polygons INSIDE the place describe only parts of it, so they are not +#' its status: they answer only when they all agree, and otherwise the answer is that it +#' depends where the record falls - report `P` and say so. Polygons merely overlapping +#' describe neither. This keeps a record's verdict tied to the evidence for the ground +#' it actually sits on, which is what coordinates give directly. +#' +#' Internal. Inheritance: an opinion about a region the query sits INSIDE ("within" +#' from the query's side) carries introduced downward but not native; an opinion about a +#' region inside the query ("contains") carries native upward but not introduced; +#' partial overlap carries neither. Then: any native wins, else introduced, else +#' present; else absent where a comprehensive source covers the region and the taxon is +#' evaluable, else unknown. +#' @keywords internal +#' @noRd +nsr_reduce <- function(op, keys, evaluable, db) { + if (is.null(op)) { + st <- src <- rel <- character(0); cult <- integer(0) + } else { + st <- op$status; src <- op$source_name; rel <- op$relation; cult <- op$is_cultivated + } + opinions <- if (!length(st)) NA_character_ else + paste(paste0(src, ":", st, ifelse(rel == "same", "", paste0("(", rel, ")"))), collapse = "; ") + cultivated <- if (!length(cult)) NA_integer_ else as.integer(any(cult %in% c(1, "1", TRUE))) + inc <- rel %in% c("same", "within") # polygons the place lies in + sub <- rel == "contains" # polygons inside the place + ovl <- rel == "overlaps" + n_sub_n <- sum(sub & st == "native"); n_sub_i <- sum(sub & st == "introduced") + out <- function(code, reason, scope, srcs = character(0), conflict = FALSE, + conflict_type = "none") { + list(code = code, reason = reason, scope = scope, + sources = if (!length(srcs)) NA_character_ else paste(sort(unique(srcs)), collapse = ", "), + opinions = opinions, conflict = conflict, conflict_type = conflict_type, + cultivated = cultivated, n_sub_native = n_sub_n, n_sub_introduced = n_sub_i) + } + if (any(inc)) { + si <- st[inc] + ct <- nsr_conflict_type(si, src[inc]) + code <- if ("native" %in% si) "N" else if ("introduced" %in% si) "I" else "P" + want <- c(N = "native", I = "introduced", P = "present")[[code]] + stated <- any(rel == "same" & st == want) + word <- c(N = "Native", I = "Introduced", P = "Present")[[code]] + reason <- if (ct != "none" && code == "N") { + paste0("Native to this polygon as per checklist (", ct, " disagree; any native opinion is taken)") + } else if (stated) paste0(word, " in this polygon, as per checklist") + else paste0(word, " in a polygon containing this place, as per checklist") + return(out(code, reason, if (stated) "polygon" else "containing polygon", src[inc], + ct != "none", ct)) + } + if (any(sub)) { + ss <- unique(st[sub]) + if (length(ss) == 1) { + code <- c(native = "N", introduced = "I", present = "P")[[ss]] + word <- c(native = "Native", introduced = "Introduced", present = "Present")[[ss]] + return(out(code, paste0(word, " in every listed polygon within this place, as per checklist"), + "sub-polygons agree", src[sub])) + } + return(out("P", paste0("Status varies among the polygons within this place (", n_sub_n, + " native, ", n_sub_i, " introduced); give coordinates or a finer division"), + "sub-polygons differ", src[sub])) + } + if (any(ovl)) { + return(out("UNK", "Only polygons partly overlapping this place have a status; none covers it", + "overlapping only", src[ovl])) + } + if (!evaluable) return(out("UNK", "No source holds native status information for this taxon", "none")) + if (any(keys %in% db$covered)) return(out("A", "Absent from the comprehensive checklists for this polygon", "polygon")) + out("UNK", "No comprehensive checklist covers this polygon", "none") +} diff --git a/R/NSR_local_set.R b/R/NSR_local_set.R new file mode 100644 index 0000000..f160563 --- /dev/null +++ b/R/NSR_local_set.R @@ -0,0 +1,293 @@ +# data.table is used inside the package without being attached; this tells it the +# package knows what it is doing, so [.data.table keeps its own semantics here. +.datatable.aware <- TRUE + +# Set-based resolution: every query answered in a handful of joins rather than one pass +# per row. Same semantics as the row path in NSR_local.R - a place is judged by the +# polygons containing it, sub-polygons answer only when unanimous, endemism is the +# exception that may look outside - but it scales to the millions of taxon x polygon +# pairs an occurrence pipeline asks about. + +#' Which taxa have a native range confined to one region or one country +#' +#' Internal. Computed once per build, for the `Ie` rule: a taxon absent from the place +#' queried but confined elsewhere cannot be native there. A taxon native in a single +#' region is confined to it; otherwise it is confined to a country if every region it is +#' native in lies inside that one country. +#' @keywords internal +#' @noRd +nsr_confined_ranges <- function(nat, db) { + dt <- data.table::data.table + cnt <- nat[, list(n_nat = .N), by = "taxon_id"] + single <- merge(nat, cnt[n_nat == 1L], by = "taxon_id") + nm <- function(k) { + i <- match(k, db$regions$region_key) + data.table::fifelse(is.na(i), k, db$regions$region_name[i]) + } + out <- dt(taxon_id = single$taxon_id, confined_to = nm(single$region_key)) + multi <- cnt[n_nat > 1L] + if (nrow(multi)) { + lk <- db$links_dt[relation %in% c("same", "within") & + startsWith(to_region, "gadm0:"), + list(region_key = from_region, country = to_region)] + m <- merge(nat[taxon_id %in% multi$taxon_id], lk, by = "region_key", allow.cartesian = TRUE) + per <- m[, list(n_in = data.table::uniqueN(region_key)), by = c("taxon_id", "country")] + per <- merge(per, multi, by = "taxon_id") + conf <- per[n_in == n_nat] + conf <- conf[!duplicated(conf$taxon_id)] + if (nrow(conf)) out <- data.table::rbindlist(list(out, + dt(taxon_id = conf$taxon_id, confined_to = nm(conf$country)))) + } + data.table::setkeyv(out, "taxon_id") + out +} + +#' Native status for taxa and polygons you have already resolved +#' +#' A pipeline that has located its records itself - in the WGSRPD raster and the GADM +#' index, as an occurrence workflow does once for all of its coordinates - already knows +#' the polygons each record sits in. This takes those directly, skipping name resolution +#' and point-in-polygon, and answers with the same rules as \code{\link{NSR_local}}. +#' +#' @param taxon_id WCVP accepted ids (character), one per query. +#' @param region_keys The polygons each record sits in, as a list of character vectors +#' (\code{"wgsrpd3:BZL"}, \code{"gadm1:BRA.25_1"}), or a single character vector for +#' one polygon each. +#' @param country_keys Optional country polygons (\code{"gadm0:BRA"}), same shape. +#' @param dir Cache directory. +#' @param min_overlap Ignore region links covering less than this share of a region. +#' @return A data.frame with \code{native_status} and the same companion columns +#' \code{NSR_local()} returns. +#' @export +NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, + dir = nsr_cache_dir(), min_overlap = 0.01) { + if (!is.list(region_keys)) region_keys <- as.list(region_keys) + if (is.null(country_keys)) country_keys <- vector("list", length(taxon_id)) + if (!is.list(country_keys)) country_keys <- as.list(country_keys) + stopifnot(length(taxon_id) == length(region_keys), + length(taxon_id) == length(country_keys)) + db <- nsr_local_db(dir) + clean <- function(z) { z <- z[!is.na(z) & nzchar(z)]; if (!length(z)) character(0) else z } + places <- list(fine = lapply(region_keys, clean), country = lapply(country_keys, clean)) + res <- if (!is.null(db$chk_dt)) nsr_resolve_status_set(as.character(taxon_id), places, db, min_overlap) + else nsr_resolve_status(as.character(taxon_id), places, db, min_overlap) + data.frame(taxon_id = as.character(taxon_id), + native_status = res$code, native_status_reason = res$reason, + native_status_sources = res$sources, native_status_opinions = res$opinions, + native_status_scope = res$scope, native_status_conflict_type = res$conflict_type, + isIntroduced = as.integer(res$code %in% c("I", "Ie")), + isEndemic = as.integer(res$code == "Ne"), + isCultivatedNSR = res$cultivated, + n_subpolygons_native = res$n_sub_native, + n_subpolygons_introduced = res$n_sub_introduced, + stringsAsFactors = FALSE) +} + +#' Resolve every query at once, with joins +#' +#' Internal. Needs \code{data.table}; \code{\link{NSR_local}} falls back to the row path +#' without it. Returns the same list of per-row vectors as \code{nsr_resolve_status()}. +#' @keywords internal +#' @noRd +nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { + dt <- function(...) data.table::data.table(...) + n <- length(taxon_id) + blank <- rep(NA_character_, n) + out <- list(code = blank, reason = blank, sources = blank, opinions = blank, + country_code = blank, state_code = blank, scope = blank, + conflict = rep(FALSE, n), conflict_type = rep(NA_character_, n), + cultivated = rep(NA_integer_, n), n_sub_native = rep(NA_integer_, n), + n_sub_introduced = rep(NA_integer_, n)) + + fine <- places$fine + ctry <- places$country + nf <- lengths(fine) + nc <- lengths(ctry) + has_fine <- nf > 0 + # the polygons the query names: its own ("same"), and its country, which contains it + Q <- data.table::rbindlist(list( + dt(qid = rep(seq_len(n), nf), region_key = unlist(fine), relation = "same"), + dt(qid = rep(seq_len(n), nc), region_key = unlist(ctry), + relation = "same"))) # fixed below for rows that have a finer key + if (!nrow(Q)) return(nsr_set_finish(out, taxon_id, db, n)) + Q[, relation := data.table::fifelse(has_fine[qid] & region_key %in% unlist(ctry), "within", relation)] + Q[, taxon_id := taxon_id[qid]] + Q <- Q[!is.na(taxon_id) & !is.na(region_key)] + + # polygons related to the ones named: only the query's OWN polygons follow links, so a + # country's other sub-polygons never reach a state-level query + own <- Q[relation == "same"] + L <- db$links_dt[own, on = c(from_region = "region_key"), nomatch = 0L, + allow.cartesian = TRUE] + if (nrow(L)) { + known <- own[, list(sys = unique(sub(":.*", "", region_key))), by = "qid"] + L[, to_sys := sub(":.*", "", to_region)] + L <- L[!paste(qid, to_sys) %in% paste(known$qid, known$sys)] + L <- L[fraction >= min_overlap & to_region != from_region, + .(qid, region_key = to_region, relation = relation, taxon_id)] + Q <- data.table::rbindlist(list(Q, L), use.names = TRUE) + } + # keep the strongest relation per (query, polygon) + ord <- c(same = 1L, within = 2L, contains = 3L, overlaps = 4L) + Q[, r_ord := ord[relation]] + data.table::setorder(Q, qid, region_key, r_ord) + Q <- unique(Q, by = c("qid", "region_key")) + + # the opinions themselves + O <- db$chk_dt[Q, on = c("taxon_id", "region_key"), nomatch = 0L, allow.cartesian = TRUE] + if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q)) + + O[, `:=`(inc = relation %in% c("same", "within"), sub = relation == "contains", + ovl = relation == "overlaps")] + + # --- verdict from the polygons the place lies in ------------------------------------- + inc <- O[inc == TRUE] + agg <- inc[, .( + has_nat = any(status == "native"), has_int = any(status == "introduced"), + stated_nat = any(relation == "same" & status == "native"), + stated_int = any(relation == "same" & status == "introduced"), + stated_pre = any(relation == "same" & status == "present"), + n_src_sets = data.table::uniqueN(.SD[, paste(sort(unique(status)), collapse = "+"), by = source_name]$V1), + within_src = any(.SD[, data.table::uniqueN(status), by = source_name]$V1 > 1), + srcs = paste(sort(unique(source_name)), collapse = ", "), + cult = as.integer(any(is_cultivated %in% c(1, "1", TRUE))) + ), by = qid] + + sub <- O[sub == TRUE] + subagg <- sub[, .(n_sub_nat = sum(status == "native"), n_sub_int = sum(status == "introduced"), + sub_sets = data.table::uniqueN(status), sub_status = status[1], + sub_srcs = paste(sort(unique(source_name)), collapse = ", ")), by = qid] + ovlagg <- O[ovl == TRUE, .(ovl_srcs = paste(sort(unique(source_name)), collapse = ", ")), by = qid] + allop <- O[, .(opinions = paste(paste0(source_name, ":", status, + data.table::fifelse(relation == "same", "", + paste0("(", relation, ")"))), + collapse = "; "), + cult_any = as.integer(any(is_cultivated %in% c(1, "1", TRUE)))), by = qid] + + res <- dt(qid = seq_len(n)) + res <- merge(res, agg, by = "qid", all.x = TRUE) + res <- merge(res, subagg, by = "qid", all.x = TRUE) + res <- merge(res, ovlagg, by = "qid", all.x = TRUE) + res <- merge(res, allop, by = "qid", all.x = TRUE) + + res[, conflict_type := data.table::fcase( + is.na(has_nat), NA_character_, + within_src & n_src_sets > 1, "both", + n_src_sets > 1, "sources", + within_src == TRUE, "within source", + default = "none")] + res[, code := data.table::fcase( + !is.na(has_nat) & has_nat, "N", + !is.na(has_int) & has_int, "I", + !is.na(has_nat), "P", + !is.na(sub_sets) & sub_sets == 1L, + c(native = "N", introduced = "I", present = "P")[sub_status], + !is.na(sub_sets), "P", + !is.na(ovl_srcs), "UNK", + default = NA_character_)] + res[, scope := data.table::fcase( + !is.na(has_nat) & ((code == "N" & stated_nat) | (code == "I" & stated_int) | + (code == "P" & stated_pre)), "polygon", + !is.na(has_nat), "containing polygon", + !is.na(sub_sets) & sub_sets == 1L, "sub-polygons agree", + !is.na(sub_sets), "sub-polygons differ", + !is.na(ovl_srcs), "overlapping only", + default = NA_character_)] + res[, reason := data.table::fcase( + code == "N" & conflict_type != "none", + paste0("Native to this polygon as per checklist (", conflict_type, + " disagree; any native opinion is taken)"), + code == "N" & scope == "polygon", "Native in this polygon, as per checklist", + code == "N" & scope == "containing polygon", + "Native in a polygon containing this place, as per checklist", + code == "N", "Native in every listed polygon within this place, as per checklist", + code == "I" & scope == "polygon", "Introduced in this polygon, as per checklist", + code == "I" & scope == "containing polygon", + "Introduced in a polygon containing this place, as per checklist", + code == "I", "Introduced in every listed polygon within this place, as per checklist", + code == "P" & scope == "sub-polygons differ", + paste0("Status varies among the polygons within this place (", n_sub_nat, " native, ", + n_sub_int, " introduced); give coordinates or a finer division"), + code == "P" & scope == "polygon", "Present in this polygon, as per checklist", + code == "P", "Present in a polygon containing this place, as per checklist", + code == "UNK", "Only polygons partly overlapping this place have a status; none covers it", + default = NA_character_)] + res[, srcs_any := data.table::fcase(!is.na(srcs), srcs, !is.na(sub_srcs), sub_srcs, + !is.na(ovl_srcs), ovl_srcs, default = NA_character_)] + + for (f in c("code", "reason", "scope", "conflict_type")) out[[f]] <- res[[f]] + out$sources <- res$srcs_any + out$opinions <- res$opinions + out$cultivated <- res$cult_any + out$n_sub_native <- data.table::fifelse(is.na(res$n_sub_nat), 0L, as.integer(res$n_sub_nat)) + out$n_sub_introduced <- data.table::fifelse(is.na(res$n_sub_int), 0L, as.integer(res$n_sub_int)) + out$conflict <- !is.na(res$conflict_type) & res$conflict_type != "none" + nsr_set_finish(out, taxon_id, db, n, Q) +} + +#' Fill in the answers that need no opinions: absence, unknowns, and endemism +#' @keywords internal +#' @noRd +nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL) { + dt <- function(...) data.table::data.table(...) + evaluable <- !is.na(taxon_id) & taxon_id %in% db$evaluable + # the regions each query consulted, for the coverage test + covered <- rep(FALSE, n) + if (!is.null(Q) && nrow(Q)) { + cv <- Q[relation %in% c("same", "within", "contains") & region_key %in% db$covered, + .(any = TRUE), by = qid] + covered[cv$qid] <- TRUE + } + no_answer <- is.na(out$code) + out$code[no_answer] <- data.table::fcase( + !evaluable[no_answer], "UNK", + covered[no_answer], "A", + default = "UNK") + out$reason[no_answer] <- data.table::fcase( + !evaluable[no_answer], "No source holds native status information for this taxon", + covered[no_answer], "Absent from the comprehensive checklists for this polygon", + default = "No comprehensive checklist covers this polygon") + out$scope[no_answer & out$code == "A"] <- "polygon" + out$scope[no_answer & out$code == "UNK"] <- "none" + out$conflict_type[no_answer] <- "none" + out$n_sub_native[is.na(out$n_sub_native)] <- 0L + out$n_sub_introduced[is.na(out$n_sub_introduced)] <- 0L + out$cultivated[out$code == "A"] <- 0L + + # rows with no place at all, or no taxon + none <- is.na(taxon_id) + if (any(none)) { + out$code[none] <- "UNK" + out$reason[none] <- "Taxon not matched to a species in the backbone" + out$scope[none] <- "none" + } + + # --- endemism, the one rule allowed to look outside the queried polygon -------------- + if (!is.null(Q) && nrow(Q)) { + inside <- Q[relation %in% c("same", "contains"), .(qid, region_key)] + nat <- db$native_dt + tq <- dt(qid = seq_len(n), taxon_id = taxon_id)[!is.na(taxon_id)] + # how many of the taxon's native regions fall inside the queried place + m <- nat[tq, on = "taxon_id", nomatch = 0L, allow.cartesian = TRUE] + if (nrow(m)) { + m[, in_place := paste(qid, region_key) %in% paste(inside$qid, inside$region_key)] + per <- m[, .(n_nat = .N, n_in = sum(in_place)), by = qid] + endemic <- per[n_nat > 0 & n_in == n_nat, qid] + i <- endemic[out$code[endemic] == "N"] + out$code[i] <- "Ne" + out$reason[i] <- sub("^Native", "Native and endemic", out$reason[i]) + # Ie: absent here, and the whole native range lies elsewhere and is confined + consulted <- Q[relation %in% c("same", "within", "contains"), .(qid, region_key)] + m[, consulted_here := paste(qid, region_key) %in% paste(consulted$qid, consulted$region_key)] + elsewhere <- m[, .(n_nat = .N, n_here = sum(consulted_here)), by = qid][n_nat > 0 & n_here == 0] + conf <- db$confined_dt + elsewhere <- merge(elsewhere, tq, by = "qid") + e <- merge(elsewhere, conf, by = "taxon_id", all.x = FALSE, sort = FALSE) + j <- e$qid[out$code[e$qid] == "A" & !is.na(e$confined_to)] + out$code[j] <- "Ie" + out$reason[j] <- paste0("Absent from this region and endemic to ", + e$confined_to[match(j, e$qid)], ", so introduced here (inferred)") + } + } + out +} diff --git a/R/local_build.R b/R/local_build.R new file mode 100644 index 0000000..cd5e437 --- /dev/null +++ b/R/local_build.R @@ -0,0 +1,239 @@ +#' Build the local native-status reference +#' +#' Downloads (or reads) each checklist source, resolves its names against WCVP with +#' \code{TNRS::TNRS_local()} and its political divisions against the GNRS backbone, and +#' writes the shared cache tables. Sources are fetched on the user's machine; nothing +#' derived ships with the package. +#' +#' Tables written: \code{nsr-sources} (one per source, with licence and whether it is +#' comprehensive), \code{nsr-taxa} (one per taxon, keyed on the WCVP accepted id), +#' \code{nsr-regions} (one per region, in the geography its source publishes against) and +#' \code{nsr-checklist} (one per taxon x region x source). \code{nsr-region-links}, the +#' spatial relations between geographies, is built too when the GADM index from GVS is in +#' the cache. +#' +#' @param sources Which sources to build. See \code{NSR_local_status()}. +#' @param dir Cache directory, shared with GNRS and GVS. +#' @param files Named list of local archives to read instead of downloading, e.g. +#' \code{list(flbr = "flbr_dwca.zip")}. For \code{powo}, the WCVP zip; if absent, the +#' copy in the TNRS cache is used when there is one. +#' @param overwrite Rebuild sources that are already built? +#' @param quiet Suppress progress messages? +#' @return Invisibly, \code{NSR_local_status()}. +#' @export +NSR_local_build <- function(sources = c("powo", "vascan", "flbr"), + dir = nsr_cache_dir(create = TRUE), + files = list(), overwrite = FALSE, quiet = FALSE) { + for (pkg in c("nanoparquet", "TNRS")) { + if (!requireNamespace(pkg, quietly = TRUE)) { + stop("Building the local NSR needs the '", pkg, "' package.", call. = FALSE) + } + } + reg <- nsr_builtin_registry() + sources <- match.arg(sources, names(reg), several.ok = TRUE) + msg <- function(...) if (!quiet) message(...) + + existing <- lapply(nsr_tables(), function(t) { + p <- nsr_table_path(t, dir) + if (file.exists(p)) as.data.frame(nanoparquet::read_parquet(p)) else NULL + }) + names(existing) <- nsr_tables() + done <- if (is.null(existing$sources)) character(0) else existing$sources$source_name + todo <- setdiff(sources, if (overwrite) character(0) else done) + if (!length(todo)) { + msg("Nothing to build (", paste(sources, collapse = ", "), " already built).") + return(invisible(NSR_local_status(dir))) + } + + bb <- nsr_gnrs_backbone(dir) + raw <- list() + for (s in todo) { + msg("Importing ", s, " ...") + t0 <- Sys.time() + raw[[s]] <- switch(s, + powo = nsr_import_powo(files$powo, bb, quiet), + vascan = nsr_import_vascan(files$vascan, bb, quiet), + flbr = nsr_import_flbr(files$flbr, bb, quiet) + ) + msg(" ", format(nrow(raw[[s]]), big.mark = ","), " records (", + round(as.numeric(difftime(Sys.time(), t0, units = "mins")), 1), " min)") + } + + # ---- names -------------------------------------------------------------------- + # sources keyed on WCVP already carry the taxon id (POWO is the backbone); the rest + # are resolved once, together, against the same backbone + need <- unique(unlist(lapply(raw, function(x) x$taxon_name[is.na(x$taxon_id)]))) + need <- need[!is.na(need) & nzchar(need)] + res <- if (length(need)) { + msg("Resolving ", format(length(need), big.mark = ","), " distinct names against WCVP ...") + # POWO's species table, when it is being built in this run, keys the others + sp_index <- do.call(rbind, lapply(raw, function(d) { + if (!all(c("taxon_id", "species_name") %in% names(d))) return(NULL) + x <- unique(d[!is.na(d$taxon_id), c("taxon_id", "species_name")]) + x[!duplicated(x$species_name), , drop = FALSE] + })) + if (is.null(sp_index) && !is.null(existing$taxa)) { + sp_index <- existing$taxa[, c("taxon_id", "species_name")] + } + r <- nsr_resolve_names(need, quiet = quiet, species_index = sp_index) + msg(" matched ", format(sum(!is.na(r$taxon_id)), big.mark = ","), " (", + round(100 * mean(!is.na(r$taxon_id)), 1), "%)") + r + } else { + data.frame(name = character(0), taxon_id = character(0), species_name = character(0), + family = character(0), genus = character(0), rank = character(0), + stringsAsFactors = FALSE) + } + + # ---- assemble --------------------------------------------------------------------- + keep_old <- function(x, col) if (is.null(x)) NULL else x[!x[[col]] %in% todo, , drop = FALSE] + chk <- lapply(names(raw), function(s) { + d <- raw[[s]] + fill <- is.na(d$taxon_id) + d$taxon_id[fill] <- res$taxon_id[match(d$taxon_name[fill], res$name)] + d <- d[!is.na(d$taxon_id) & !is.na(d$region_key), , drop = FALSE] + data.frame(taxon_id = d$taxon_id, region_key = d$region_key, status = d$status, + is_cultivated = d$is_cultivated, source_name = s, stringsAsFactors = FALSE) + }) + checklist <- unique(do.call(rbind, c(list(keep_old(existing$checklist, "source_name")), chk))) + + taxa_res <- res[!is.na(res$taxon_id), c("taxon_id", "family", "genus", "species_name", "rank")] + taxa_src <- do.call(rbind, lapply(raw, function(d) { + if (!all(c("family", "genus", "species_name", "rank") %in% names(d))) return(NULL) + x <- d[!is.na(d$taxon_id), c("taxon_id", "family", "genus", "species_name", "rank")] + x[!duplicated(x$taxon_id), , drop = FALSE] + })) + taxa_new <- rbind(taxa_res, taxa_src) + taxa <- unique(rbind(existing$taxa, taxa_new)) + taxa <- taxa[!duplicated(taxa$taxon_id), , drop = FALSE] + + reg_new <- unique(do.call(rbind, lapply(raw, function(d) { + x <- d[!is.na(d$region_key), c("region_key", "region_name")] + x[!duplicated(x$region_key), , drop = FALSE] + }))) + reg_new$system <- sub(":.*", "", reg_new$region_key) + reg_new$code <- sub("^[^:]*:", "", reg_new$region_key) + regions <- unique(rbind(existing$regions, reg_new)) + regions <- regions[!duplicated(regions$region_key), , drop = FALSE] + + src_new <- do.call(rbind, lapply(names(raw), function(s) { + r <- reg[[s]] + n <- sum(checklist$source_name == s) + data.frame(source_name = s, source_name_full = r$source_name_full, version = r$version, + url = r$url, licence = r$licence, citation = r$citation, + is_comprehensive = r$is_comprehensive, region_system = r$region_system, + coverage = r$coverage, date_accessed = as.character(Sys.Date()), + n_records = n, n_taxa = length(unique(checklist$taxon_id[checklist$source_name == s])), + stringsAsFactors = FALSE) + })) + srcs <- rbind(keep_old(existing$sources, "source_name"), src_new) + + for (nm in nsr_tables()) { + x <- switch(nm, sources = srcs, taxa = taxa, regions = regions, checklist = checklist) + nanoparquet::write_parquet(x, nsr_table_path(nm, dir), compression = "gzip") + } + saveRDS(list(built = as.character(Sys.Date()), sources = srcs$source_name, + n_records = nrow(checklist), n_taxa = nrow(taxa), n_regions = nrow(regions)), + nsr_provenance_path("nsr", dir)) + + # spatial links between geographies, where GVS's GADM index is available + systems <- unique(regions$system) + if (length(systems) > 1 && file.exists(file.path(dir, "gadmindex-units-30s.tif"))) { + msg("Linking region systems (", paste(systems, collapse = ", "), ") ...") + try(nsr_build_region_links(dir = dir, quiet = quiet), silent = FALSE) + } else if (length(systems) > 1) { + warning("Sources use several geographies but GVS's GADM index is not in the cache, ", + "so they cannot be linked; queries by division name will only see their own system.", + call. = FALSE) + } + msg("Built: ", format(nrow(checklist), big.mark = ","), " checklist records, ", + format(nrow(taxa), big.mark = ","), " taxa, ", format(nrow(regions), big.mark = ","), + " regions.") + invisible(NSR_local_status(dir)) +} + +#' The GNRS political-division backbone, as NSR keys on it +#' @keywords internal +#' @noRd +nsr_gnrs_backbone <- function(dir) { + f <- function(x) file.path(dir, paste0("gnrs-", x, ".gz.parquet")) + if (!file.exists(f("country")) || !file.exists(f("state_province"))) { + stop("The GNRS backbone is not in ", dir, ".\n", + "Build it first with GNRS::GNRS_local_build(); NSR keys its divisions on GNRS ids.", + call. = FALSE) + } + list(country = as.data.frame(nanoparquet::read_parquet(f("country"))), + state = as.data.frame(nanoparquet::read_parquet(f("state_province")))) +} + +#' Political division ids for a set of country (and optionally state) names +#' +#' Internal. Exact, case- and accent-insensitive matching against the GNRS backbone; +#' anything unmatched comes back NA and is reported by the caller, rather than being +#' guessed at. +#' @keywords internal +#' @noRd +nsr_match_poldiv <- function(country, state = NULL, bb) { + norm <- function(x) { + x <- iconv(as.character(x), to = "ASCII//TRANSLIT") + tolower(trimws(gsub("[^A-Za-z0-9 ]", "", ifelse(is.na(x), "", x)))) + } + ci <- match(norm(country), norm(bb$country$country)) + out <- data.frame( + country_id = bb$country$country_id[ci], country = bb$country$country[ci], + state_province_id = NA_integer_, state_province = NA_character_, + stringsAsFactors = FALSE + ) + if (!is.null(state)) { + key <- paste(norm(country), norm(state)) + # the backbone spells some divisions bilingually ("New Brunswick/Nouveau-Brunswick"), + # so every variant is a valid target + variants <- function(col) { + parts <- strsplit(ifelse(is.na(col), "", col), "/", fixed = TRUE) + data.frame(row = rep(seq_along(parts), lengths(parts)), + name = norm(unlist(parts)), stringsAsFactors = FALSE) + } + v <- unique(rbind( + data.frame(row = seq_len(nrow(bb$state)), name = norm(bb$state$state_province)), + data.frame(row = seq_len(nrow(bb$state)), name = norm(bb$state$state_province_std)), + variants(bb$state$state_province), variants(bb$state$state_province_std))) + v$key <- paste(norm(bb$state$country)[v$row], v$name) + v <- v[!duplicated(v$key), , drop = FALSE] + si <- v$row[match(key, v$key)] + out$state_province_id <- bb$state$state_province_id[si] + out$state_province <- bb$state$state_province[si] + } + out$poldiv_level <- ifelse(!is.na(out$state_province_id), "state_province", + ifelse(!is.na(out$country_id), "country", NA_character_)) + out$poldiv_id <- ifelse(!is.na(out$state_province_id), paste0("S", out$state_province_id), + ifelse(!is.na(out$country_id), paste0("C", out$country_id), NA_character_)) + out +} + +#' Resolve checklist names to WCVP accepted species +#' +#' Internal. The id returned is the SPECIES' id, not an infraspecific one, so every +#' source keys on the same unit as POWO after its infraspecifics are rolled up. +#' @keywords internal +#' @noRd +nsr_resolve_names <- function(names_in, quiet = FALSE, species_index = NULL) { + r <- TNRS::TNRS_local(taxonomic_names = names_in, sources = "wcvp", matches = "best", + build_missing = FALSE, quiet = quiet) + acc <- r$Accepted_species + acc[is.na(acc) | !nzchar(acc)] <- NA_character_ + # prefer the species' own id where we know it, so subspecies land on their species + id <- ifelse(is.na(acc), NA_character_, r$Accepted_name_id) + if (!is.null(species_index)) { + sid <- species_index$taxon_id[match(acc, species_index$species_name)] + id <- ifelse(!is.na(sid), sid, id) + } + data.frame( + name = r$Name_submitted, + taxon_id = id, + species_name = acc, + family = r$Accepted_family, + genus = sub(" .*", "", acc), + rank = r$Accepted_name_rank, + stringsAsFactors = FALSE + ) +} diff --git a/R/local_cache.R b/R/local_cache.R new file mode 100644 index 0000000..215b91e --- /dev/null +++ b/R/local_cache.R @@ -0,0 +1,133 @@ +#' Cache directory for the offline NSR +#' +#' Internal. Shared with GNRS and GVS, so the three services build one cache: the +#' political divisions NSR keys on are GNRS's, and a user who has built either of the +#' others already has them. \code{options(NSR.cache_dir=)} overrides, then +#' \code{options(GNRS.cache_dir=)}, then \code{tools::R_user_dir("GNRS", "cache")}. +#' @keywords internal +#' @noRd +nsr_cache_dir <- function(create = FALSE) { + dir <- getOption("NSR.cache_dir", + getOption("GNRS.cache_dir", tools::R_user_dir("GNRS", which = "cache"))) + if (create && !dir.exists(dir)) dir.create(dir, recursive = TRUE, showWarnings = FALSE) + dir +} + +#' Paths of the local NSR tables +#' @keywords internal +#' @noRd +nsr_table_path <- function(table, dir = nsr_cache_dir()) { + file.path(dir, paste0("nsr-", table, ".gz.parquet")) +} + +#' @keywords internal +#' @noRd +nsr_provenance_path <- function(component, dir = nsr_cache_dir()) { + file.path(dir, paste0("nsr-", component, "-provenance.rds")) +} + +#' The tables the local NSR is made of +#' @keywords internal +#' @noRd +nsr_tables <- function() c("sources", "taxa", "regions", "checklist") + +#' Checklist sources the local NSR can build +#' +#' Internal. One entry per source: what it covers, whether it is comprehensive for its +#' area (so that absence there is informative), where it comes from and under what +#' licence. Sources are fetched on the user's machine at build time; nothing derived is +#' shipped with the package. +#' +#' \code{region_system} is the geography the source publishes against, kept as it is. +#' \code{is_comprehensive} says the source aims to list everything in its area; it does +#' not say which taxa it covers, which is measured from the data instead. +#' @keywords internal +#' @noRd +nsr_builtin_registry <- function() { + list( + powo = list( + source_name = "powo", + source_name_full = "Plants of the World Online / World Checklist of Vascular Plants (WCVP)", + version = "v15", + url = "https://sftp.kew.org/pub/data-repositories/WCVP/", + licence = "CC BY 4.0", + citation = "Govaerts, R. et al. The World Checklist of Vascular Plants (WCVP). Royal Botanic Gardens, Kew.", + is_comprehensive = TRUE, + region_system = "wgsrpd3", + coverage = "global, WGSRPD level 3 (the geography WCVP publishes against)", + local_file = "wcvp-v15.zip (from the TNRS cache, or supplied)" + ), + vascan = list( + source_name = "vascan", + source_name_full = "Database of Vascular Plants of Canada (VASCAN)", + version = NA_character_, + url = "https://data.canadensys.net/ipt/archive.do?r=vascan", + licence = "CC0 1.0", + citation = "Brouillet, L. et al. VASCAN, the Database of Vascular Plants of Canada.", + is_comprehensive = TRUE, + region_system = "gadm1", + coverage = "Canada, by province and territory", + local_file = "vascan_dwca.zip" + ), + # usda is deliberately absent: the USDA no longer serves the per-state native-status + # export the service was built from. The GBIF copy of PLANTS is names only, and the + # current API gives status by region (CAN / L48 / AK / HI), not by state, which adds + # nothing over POWO. See dev_notes/01-offline-nsr-design.md. + flbr = list( + source_name = "flbr", + source_name_full = "Flora e Funga do Brasil (Lista Oficial)", + version = NA_character_, + url = "https://api.checklistbank.org/dataset/2031/archive", + licence = "CC BY 4.0", + citation = "Flora e Funga do Brasil. Jardim Botanico do Rio de Janeiro.", + is_comprehensive = TRUE, + region_system = "gadm1", + coverage = "Brazil, by state; the level WGSRPD cannot reach (five level-3 regions)", + local_file = "flbr_dwca.zip" + ) + ) +} + +#' Which parts of the local NSR are built +#' +#' @param dir Cache directory. Defaults to the one shared with GNRS and GVS. +#' @return A data.frame with one row per source: whether it is built, how many checklist +#' records and taxa it contributed, and when it was built. +#' @export +NSR_local_status <- function(dir = nsr_cache_dir()) { + reg <- nsr_builtin_registry() + built_tables <- vapply(nsr_tables(), function(t) file.exists(nsr_table_path(t, dir)), logical(1)) + src <- if (built_tables[["sources"]]) { + as.data.frame(nanoparquet::read_parquet(nsr_table_path("sources", dir))) + } else NULL + out <- do.call(rbind, lapply(reg, function(s) { + row <- if (!is.null(src) && s$source_name %in% src$source_name) src[src$source_name == s$source_name, ] else NULL + data.frame( + source = s$source_name, coverage = s$coverage, licence = s$licence, + built = !is.null(row), + n_records = if (is.null(row)) NA_integer_ else as.integer(row$n_records), + n_taxa = if (is.null(row)) NA_integer_ else as.integer(row$n_taxa), + date_built = if (is.null(row)) NA_character_ else as.character(row$date_accessed), + stringsAsFactors = FALSE + ) + })) + rownames(out) <- NULL + attr(out, "cache_dir") <- dir + out +} + +#' Remove the local NSR tables +#' +#' @param dir Cache directory. +#' @param quiet Suppress messages? +#' @return Invisibly, the files removed. +#' @export +NSR_local_remove <- function(dir = nsr_cache_dir(), quiet = FALSE) { + f <- c(vapply(nsr_tables(), nsr_table_path, character(1), dir = dir), + nsr_table_path("region-links", dir), nsr_raster_path("wgsrpd3", dir), + nsr_raster_path("wgsrpd3_index", dir), nsr_provenance_path("nsr", dir)) + f <- f[file.exists(f)] + unlink(f) + if (!quiet) message("Removed ", length(f), " file(s) from ", dir) + invisible(f) +} diff --git a/R/local_import.R b/R/local_import.R new file mode 100644 index 0000000..9441fab --- /dev/null +++ b/R/local_import.R @@ -0,0 +1,219 @@ +# Source importers for the local NSR. +# +# Each returns a data.frame with one row per taxon x division x opinion: +# taxon_name the source's own name string (resolved later against WCVP) +# taxon_id the WCVP accepted id, where the source supplies one directly +# status "native", "introduced" or "present" +# is_cultivated 1 where the source says the taxon is cultivated in that division +# region_key ":", the source's OWN geography +# region_name the source's name for it +# +# Each source keeps the regions it publishes against: WCVP records WGSRPD level-3 areas +# (biogeographic, e.g. "Borneo"), VASCAN and Flora do Brasil record provinces and states +# (GADM level 1). Queries are resolved against those polygons, not against names; see +# local_regions.R. + +#' Read a delimited file quickly when data.table is installed +#' @keywords internal +#' @noRd +nsr_read_delim <- function(path_or_conn, sep = "\t", quote = "\"", select = NULL) { + if (requireNamespace("data.table", quietly = TRUE) && is.character(path_or_conn)) { + return(as.data.frame(data.table::fread(path_or_conn, sep = sep, quote = quote, + showProgress = FALSE, data.table = FALSE, + select = select, colClasses = "character"))) + } + x <- utils::read.delim(path_or_conn, sep = sep, quote = quote, colClasses = "character", + stringsAsFactors = FALSE, na.strings = "") + if (!is.null(select)) x <- x[, intersect(select, names(x)), drop = FALSE] + x +} + +#' Extract one member of a zip to a temporary file +#' @keywords internal +#' @noRd +nsr_unzip_member <- function(zip, member) { + td <- file.path(tempdir(), paste0("nsr-", basename(zip))) + dir.create(td, showWarnings = FALSE, recursive = TRUE) + utils::unzip(zip, files = member, exdir = td, overwrite = TRUE) + file.path(td, member) +} + +#' POWO / WCVP: global native and introduced ranges at WGSRPD level 3 +#' +#' Internal. Distributions are kept against the WGSRPD level-3 areas WCVP publishes +#' them against; the polygons do the rest. Matching those areas to political divisions +#' by NAME loses a quarter of them (Borneo, Sulawesi, Maluku, New Guinea, Lesser Sunda +#' Is. are not country or state names) and cannot express "Brazil Southeast" containing +#' four states, so it is not done. +#' +#' Names need no resolver: WCVP is the backbone, so \code{accepted_plant_name_id} is the +#' taxon id directly. +#' @keywords internal +#' @noRd +nsr_import_powo <- function(zip = NULL, bb, quiet = FALSE) { + if (is.null(zip)) zip <- nsr_find_wcvp_zip() + if (is.null(zip) || !file.exists(zip)) { + stop("Could not find the WCVP archive. Supply it as files = list(powo = \"wcvp-v15.zip\"), ", + "or build the TNRS 'wcvp' source, whose cache holds one.", call. = FALSE) + } + dist_f <- nsr_unzip_member(zip, "wcvp_distribution.csv") + names_f <- nsr_unzip_member(zip, "wcvp_names.csv") + if (!quiet) message(" reading WCVP ...") + dist <- nsr_read_delim(dist_f, sep = "|", + select = c("plant_name_id", "area_code_l3", "area", "introduced", + "extinct", "location_doubtful")) + nm <- nsr_read_delim(names_f, sep = "|", + select = c("plant_name_id", "accepted_plant_name_id", "parent_plant_name_id", + "taxon_name", "taxon_rank", "taxon_status", "family", "genus")) + + # distributions are recorded against a name; carry them to its accepted taxon + i <- match(dist$plant_name_id, nm$plant_name_id) + acc <- nm$accepted_plant_name_id[i] + acc[is.na(acc) | !nzchar(acc)] <- nm$plant_name_id[i][is.na(acc) | !nzchar(acc)] + # WCVP records distributions against accepted infraspecific taxa too (60k of them), so a + # subspecies' range hangs off the subspecies. Roll those up to the species: it is NSR's + # own rule (native propagates upward) and it is what a species-level query needs. + sp <- nsr_species_of(acc, nm) + j <- match(sp, nm$plant_name_id) + dist$taxon_id <- sp + dist$taxon_name <- nm$taxon_name[j] + dist$rank_published <- nm$taxon_rank[match(acc, nm$plant_name_id)] + dist <- dist[!is.na(dist$taxon_id) & !is.na(dist$taxon_name), , drop = FALSE] + + flag <- function(x) !is.na(x) & x %in% c("1", 1, TRUE, "TRUE") + out <- data.frame( + taxon_name = dist$taxon_name, taxon_id = dist$taxon_id, + species_name = dist$taxon_name, family = nm$family[j], genus = nm$genus[j], + rank = nm$taxon_rank[j], + status = ifelse(flag(dist$location_doubtful), "present", + ifelse(flag(dist$introduced), "introduced", "native")), + is_cultivated = 0L, + region_key = paste0("wgsrpd3:", dist$area_code_l3), region_name = dist$area, + stringsAsFactors = FALSE + ) + out[!is.na(out$region_key), , drop = FALSE] +} + +#' Walk accepted infraspecific taxa up to their species +#' +#' Internal. Up to three steps (a variety under a subspecies), stopping at the first +#' ancestor of species rank; anything already at species rank or above is left alone. +#' @keywords internal +#' @noRd +nsr_species_of <- function(ids, nm, max_steps = 3L) { + rank <- tolower(nm$taxon_rank) + is_infra <- !rank %in% c("species", "genus", "family", "") + out <- ids + for (step in seq_len(max_steps)) { + i <- match(out, nm$plant_name_id) + up <- !is.na(i) & is_infra[i] + if (!any(up)) break + par <- nm$parent_plant_name_id[i[up]] + par[is.na(par) | !nzchar(par)] <- out[up][is.na(par) | !nzchar(par)] + out[up] <- par + } + out +} + +#' Where TNRS keeps its copy of WCVP, if it has one +#' @keywords internal +#' @noRd +nsr_find_wcvp_zip <- function() { + if (!requireNamespace("TNRS", quietly = TRUE)) return(NULL) + d <- try(TNRS:::tnrs_cache_dir(), silent = TRUE) + if (inherits(d, "try-error") || !dir.exists(d)) return(NULL) + f <- list.files(d, pattern = "^wcvp.*\\.zip$", full.names = TRUE) + if (!length(f)) NULL else f[order(file.mtime(f), decreasing = TRUE)][1] +} + +#' VASCAN: Canada, by province and territory +#' @keywords internal +#' @noRd +nsr_import_vascan <- function(zip, bb, quiet = FALSE) { + if (is.null(zip) || !file.exists(zip)) { + stop("Supply the VASCAN archive as files = list(vascan = \"vascan_dwca.zip\"); ", + "download it from https://data.canadensys.net/ipt/archive.do?r=vascan", call. = FALSE) + } + tx <- nsr_read_delim(nsr_unzip_member(zip, "taxon.txt")) + di <- nsr_read_delim(nsr_unzip_member(zip, "distribution.txt")) + # carry each distribution to the accepted name of its taxon + acc <- tx$acceptedNameUsageID + acc[is.na(acc) | !nzchar(acc)] <- tx$id[is.na(acc) | !nzchar(acc)] + i <- match(di$id, tx$id) + j <- match(acc[i], tx$id) + di$taxon_name <- tx$scientificName[j] + # VASCAN splits the two questions: occurrenceStatus is present / excluded / doubtful / + # irregular / absent, establishmentMeans is native / introduced. Nativeness comes from + # the latter; the former only says whether the record counts at all. + occ <- tolower(ifelse(is.na(di$occurrenceStatus), "", di$occurrenceStatus)) + est <- tolower(ifelse(is.na(di$establishmentMeans), "", di$establishmentMeans)) + status <- ifelse(est %in% c("native", "endemic"), "native", + ifelse(est %in% c("introduced", "naturalized", "naturalised"), "introduced", + "present")) + status[occ %in% c("absent", "excluded")] <- NA_character_ + # locationID is "ISO3166-2:CA-NB"; the code beats the name, since VASCAN's spellings and + # GADM's differ and the list also covers places outside Canada + code <- sub("^ISO3166-2:", "", ifelse(is.na(di$locationID), "", di$locationID)) + gid1 <- bb$state$gid_1[match(ifelse(nzchar(code), sub("-", ".", toupper(code)), NA_character_), + bb$state$hasc_full)] + # the three VASCAN entries carrying no code: two are one GADM province, one another country + loc <- tolower(ifelse(is.na(di$locality), "", di$locality)) + gid1[is.na(gid1) & loc %in% c("newfoundland", "labrador")] <- + bb$state$gid_1[match("CA.NF", bb$state$hasc_full)] + gid0 <- rep(NA_character_, nrow(di)) + gid0[is.na(gid1) & loc == "greenland"] <- "GRL" + out <- data.frame(taxon_name = di$taxon_name, taxon_id = NA_character_, status = status, + is_cultivated = 0L, + region_key = ifelse(!is.na(gid1), paste0("gadm1:", gid1), + ifelse(!is.na(gid0), paste0("gadm0:", gid0), NA_character_)), + region_name = di$locality, stringsAsFactors = FALSE) + if (!quiet) { + miss <- unique(di$locality[is.na(out$region_key)]) + if (length(miss)) message(" ", length(miss), " VASCAN divisions unmatched: ", + paste(utils::head(miss, 6), collapse = ", ")) + } + out[!is.na(out$status) & !is.na(out$region_key) & !is.na(out$taxon_name), , drop = FALSE] +} + +#' Flora e Funga do Brasil: Brazil, by state +#' +#' Internal. \code{establishmentMeans} is Portuguese: NATIVA, NATURALIZADA, CULTIVADA. +#' Cultivated rows are kept as "present" and flagged, which is what the service's +#' \code{isCultivatedNSR} reports. Divisions come from the ISO 3166-2 code in +#' \code{locationID} (BR-SP ...), matched to the GNRS backbone's HASC code. +#' @keywords internal +#' @noRd +nsr_import_flbr <- function(zip, bb, quiet = FALSE) { + if (is.null(zip) || !file.exists(zip)) { + stop("Supply the Flora do Brasil archive as files = list(flbr = \"flbr_dwca.zip\"); ", + "download it from https://api.checklistbank.org/dataset/2031/archive", call. = FALSE) + } + tx <- nsr_read_delim(nsr_unzip_member(zip, "taxon.txt")) + di <- nsr_read_delim(nsr_unzip_member(zip, "distribution.txt")) + acc <- tx$acceptedNameUsageID + acc[is.na(acc) | !nzchar(acc)] <- tx$id[is.na(acc) | !nzchar(acc)] + i <- match(di$id, tx$id) + j <- match(acc[i], tx$id) + di$taxon_name <- tx$scientificName[j] + + em <- toupper(ifelse(is.na(di$establishmentMeans), "", di$establishmentMeans)) + status <- ifelse(em == "NATIVA", "native", + ifelse(em == "NATURALIZADA", "introduced", "present")) + cult <- as.integer(em == "CULTIVADA") + + # BR-SP -> GNRS hasc_full "BR.SP" -> GADM gid_1 + hasc <- sub("-", ".", toupper(ifelse(is.na(di$locationID), "", di$locationID))) + k <- match(hasc, bb$state$hasc_full) + out <- data.frame( + taxon_name = di$taxon_name, taxon_id = NA_character_, status = status, is_cultivated = cult, + region_key = ifelse(is.na(bb$state$gid_1[k]), NA_character_, paste0("gadm1:", bb$state$gid_1[k])), + region_name = bb$state$state_province[k], + stringsAsFactors = FALSE + ) + if (!quiet) { + miss <- unique(di$locationID[is.na(k)]) + if (length(miss)) message(" ", length(miss), " Brazilian location codes unmatched: ", + paste(utils::head(miss, 8), collapse = ", ")) + } + out[!is.na(out$region_key) & !is.na(out$taxon_name), , drop = FALSE] +} diff --git a/R/local_regions.R b/R/local_regions.R new file mode 100644 index 0000000..3c60385 --- /dev/null +++ b/R/local_regions.R @@ -0,0 +1,199 @@ +# Regions, and the spatial relations between them. +# +# Sources publish against different geographies: WCVP against WGSRPD level-3 areas +# (biogeographic), VASCAN and Flora do Brasil against GADM level-1 units (political). +# Both are kept as they are. A coordinate is located in each system directly; a query +# given only division names is carried between systems by a link table derived from the +# polygons themselves, so nothing is matched by name. + +#' Paths of the region rasters +#' @keywords internal +#' @noRd +nsr_raster_path <- function(what, dir = nsr_cache_dir()) { + switch(what, + wgsrpd3 = file.path(dir, "nsr-wgsrpd3-30s.tif"), + wgsrpd3_index = file.path(dir, "nsr-wgsrpd3-index.gz.parquet"), + stop("unknown raster: ", what) + ) +} + +#' Rasterize the WGSRPD level-3 polygons +#' +#' Internal. 30 arc-seconds, to match the GADM index GVS builds, so the two can be +#' cross-tabulated cell for cell. Polygons come from \code{rWCVPdata}, the published +#' TDWG shapes; nothing is downloaded here. +#' @keywords internal +#' @noRd +nsr_build_wgsrpd_raster <- function(dir = nsr_cache_dir(create = TRUE), resolution = 1 / 120, + quiet = FALSE) { + for (pkg in c("sf", "terra")) { + if (!requireNamespace(pkg, quietly = TRUE)) stop("Needs the '", pkg, "' package.", call. = FALSE) + } + pol <- nsr_wgsrpd_polygons() + idx <- data.frame(unit = seq_len(nrow(pol)), area_code_l3 = pol$LEVEL3_COD, + area_name = pol$LEVEL3_NAM, stringsAsFactors = FALSE) + nanoparquet::write_parquet(idx, nsr_raster_path("wgsrpd3_index", dir), compression = "gzip") + if (!quiet) message("Rasterizing ", nrow(pol), " WGSRPD level-3 areas at ", + round(resolution * 3600), " arc-seconds ...") + pol$unit <- idx$unit + tmpl <- terra::rast(terra::ext(-180, 180, -90, 90), resolution = resolution, crs = "EPSG:4326") + terra::rasterize(terra::vect(pol), tmpl, field = "unit", + filename = nsr_raster_path("wgsrpd3", dir), overwrite = TRUE, + wopt = list(datatype = "INT2U", gdal = c("COMPRESS=DEFLATE", "TILED=YES"))) + if (!quiet) message(" done") + invisible(nsr_raster_path("wgsrpd3", dir)) +} + +#' The TDWG level-3 polygons +#' @keywords internal +#' @noRd +nsr_wgsrpd_polygons <- function() { + if (requireNamespace("rWCVPdata", quietly = TRUE)) { + p <- try(rWCVPdata::wgsrpd3, silent = TRUE) + if (!inherits(p, "try-error")) return(sf::st_as_sf(p)) + } + if (requireNamespace("rWCVP", quietly = TRUE)) { + e <- new.env() + utils::data("wgsrpd3", package = "rWCVP", envir = e) + return(sf::st_as_sf(get("wgsrpd3", envir = e))) + } + stop("The WGSRPD level-3 polygons need the 'rWCVPdata' (or 'rWCVP') package.", call. = FALSE) +} + +#' Spatial links between region systems +#' +#' Internal. Cross-tabulates the WGSRPD level-3 raster against the GADM unit index GVS +#' builds, in blocks of rows, and records for each overlapping pair how much of each +#' region lies in the other. \code{fraction_from} is the share of the FROM region inside +#' the TO region, so \code{within} means the from-region is (almost) entirely inside. +#' +#' @return Invisibly, the link table. +#' @keywords internal +#' @noRd +nsr_build_region_links <- function(dir = nsr_cache_dir(), gadm_raster = NULL, block_rows = 500, + within_threshold = 0.95, quiet = FALSE) { + if (!requireNamespace("terra", quietly = TRUE)) stop("Needs the 'terra' package.", call. = FALSE) + if (is.null(gadm_raster)) gadm_raster <- file.path(dir, "gadmindex-units-30s.tif") + if (!file.exists(gadm_raster)) { + stop("The GADM unit index is not in ", dir, ".\n", + "Build it with GVS (GVS_local_build(\"index\")); NSR reuses it rather than making its own.", + call. = FALSE) + } + if (!file.exists(nsr_raster_path("wgsrpd3", dir))) nsr_build_wgsrpd_raster(dir, quiet = quiet) + + w <- terra::rast(nsr_raster_path("wgsrpd3", dir)) + g <- terra::rast(gadm_raster) + if (!isTRUE(all.equal(as.vector(terra::ext(w)), as.vector(terra::ext(g)))) || + terra::ncol(w) != terra::ncol(g) || terra::nrow(w) != terra::nrow(g)) { + stop("The WGSRPD and GADM rasters do not align; rebuild one at the other's resolution.", + call. = FALSE) + } + units <- as.data.frame(nanoparquet::read_parquet(file.path(dir, "gadmindex-units.gz.parquet"))) + widx <- as.data.frame(nanoparquet::read_parquet(nsr_raster_path("wgsrpd3_index", dir))) + + nr <- terra::nrow(w) + nc <- terra::ncol(w) + terra::readStart(w); terra::readStart(g) + on.exit({terra::readStop(w); terra::readStop(g)}, add = TRUE) + acc <- new.env(parent = emptyenv()) + if (!quiet) message("Cross-tabulating WGSRPD x GADM (", nr, " rows) ...") + for (s in seq(1L, nr, by = block_rows)) { + e <- min(nr, s + block_rows - 1L) + wv <- terra::readValues(w, row = s, nrows = e - s + 1L, col = 1L, ncols = nc) + gv <- terra::readValues(g, row = s, nrows = e - s + 1L, col = 1L, ncols = nc) + ok <- !is.na(wv) & !is.na(gv) & wv > 0 & gv > 0 + if (!any(ok)) next + key <- paste(wv[ok], gv[ok]) + tb <- table(key) + for (k in names(tb)) { + acc[[k]] <- (if (is.null(acc[[k]])) 0L else acc[[k]]) + as.integer(tb[[k]]) + } + rm(wv, gv, ok, key, tb) + } + keys <- ls(acc) + if (!length(keys)) stop("No overlap found between the two rasters.", call. = FALSE) + parts <- do.call(rbind, strsplit(keys, " ", fixed = TRUE)) + cells <- vapply(keys, function(k) acc[[k]], integer(1)) + d <- data.frame(w_unit = as.integer(parts[, 1]), g_unit = as.integer(parts[, 2]), + cells = as.integer(cells), stringsAsFactors = FALSE) + d$w_key <- paste0("wgsrpd3:", widx$area_code_l3[match(d$w_unit, widx$unit)]) + gi <- match(d$g_unit, units$unit) + d$g_key <- ifelse(!is.na(units$gid_1[gi]), paste0("gadm1:", units$gid_1[gi]), + paste0("gadm0:", units$gid_0[gi])) + # A query about a country must find links too, so the level-0 pairs are built as well - + # but each level needs its own denominators. Pooling them counts every cell twice (once + # for its state, once for its country), which halves every fraction and turns even + # coextensive regions into "overlaps". + d$g_key0 <- paste0("gadm0:", units$gid_0[gi]) + fracs <- function(gk) { + dd <- stats::aggregate(cells ~ w_key + g_key, + data = data.frame(w_key = d$w_key, g_key = gk, cells = d$cells, + stringsAsFactors = FALSE), FUN = sum) + wt <- stats::aggregate(cells ~ w_key, data = dd, FUN = sum) + gt <- stats::aggregate(cells ~ g_key, data = dd, FUN = sum) + dd$frac_w <- dd$cells / wt$cells[match(dd$w_key, wt$w_key)] + dd$frac_g <- dd$cells / gt$cells[match(dd$g_key, gt$g_key)] + dd + } + d <- rbind(fracs(d$g_key), fracs(d$g_key0)) + d <- d[!duplicated(paste(d$w_key, d$g_key)), , drop = FALSE] + + # Coextensive regions (California the GADM state and CAL the WGSRPD area) are the same + # place, and must not demote a native opinion to mere presence. The threshold is 0.95, + # not 0.99: these fractions come from 30 arc-second rasters of two independent + # coastlines, so Austria and WGSRPD's AUT score 0.990/0.988 rather than 1.0. + rel <- function(frac_from, frac_to) { + ifelse(frac_from >= within_threshold & frac_to >= within_threshold, "same", + ifelse(frac_from >= within_threshold, "within", + ifelse(frac_to >= within_threshold, "contains", "overlaps"))) + } + links <- rbind( + data.frame(from_region = d$w_key, to_region = d$g_key, + relation = rel(d$frac_w, d$frac_g), fraction = d$frac_w, stringsAsFactors = FALSE), + data.frame(from_region = d$g_key, to_region = d$w_key, + relation = rel(d$frac_g, d$frac_w), fraction = d$frac_g, stringsAsFactors = FALSE) + ) + nanoparquet::write_parquet(links, nsr_table_path("region-links", dir), compression = "gzip") + if (!quiet) message(" ", format(nrow(links), big.mark = ","), " links (", + sum(links$relation == "within"), " within, ", + sum(links$relation == "contains"), " contains, ", + sum(links$relation == "overlaps"), " overlaps)") + invisible(links) +} + +#' Locate coordinates in both region systems +#' +#' Internal. Raster lookup in the WGSRPD level-3 raster and in the GADM unit index, so a +#' record with coordinates is placed in each source's own geography without any crosswalk. +#' @return data.frame with \code{wgsrpd3} and \code{gadm} region keys (and the GADM +#' country key), one row per point. +#' @keywords internal +#' @noRd +nsr_locate_regions <- function(longitude, latitude, dir = nsr_cache_dir()) { + if (!requireNamespace("terra", quietly = TRUE)) stop("Needs the 'terra' package.", call. = FALSE) + n <- length(longitude) + out <- data.frame(wgsrpd3 = rep(NA_character_, n), gadm = NA_character_, gadm0 = NA_character_, + stringsAsFactors = FALSE) + ok <- is.finite(longitude) & is.finite(latitude) + if (!any(ok)) return(out) + xy <- cbind(longitude[ok], latitude[ok]) + + wf <- nsr_raster_path("wgsrpd3", dir) + if (file.exists(wf)) { + widx <- as.data.frame(nanoparquet::read_parquet(nsr_raster_path("wgsrpd3_index", dir))) + u <- as.integer(terra::extract(terra::rast(wf), xy)[, 1]) + out$wgsrpd3[ok] <- ifelse(is.na(u) | u == 0, NA_character_, + paste0("wgsrpd3:", widx$area_code_l3[match(u, widx$unit)])) + } + gf <- file.path(dir, "gadmindex-units-30s.tif") + if (file.exists(gf)) { + units <- as.data.frame(nanoparquet::read_parquet(file.path(dir, "gadmindex-units.gz.parquet"))) + u <- as.integer(terra::extract(terra::rast(gf), xy)[, 1]) + i <- match(u, units$unit) + out$gadm[ok] <- ifelse(is.na(i), NA_character_, + ifelse(!is.na(units$gid_1[i]), paste0("gadm1:", units$gid_1[i]), + paste0("gadm0:", units$gid_0[i]))) + out$gadm0[ok] <- ifelse(is.na(i), NA_character_, paste0("gadm0:", units$gid_0[i])) + } + out +} diff --git a/dev_notes/01-offline-nsr-design.md b/dev_notes/01-offline-nsr-design.md new file mode 100644 index 0000000..9d2a61b --- /dev/null +++ b/dev_notes/01-offline-nsr-design.md @@ -0,0 +1,494 @@ +# 05. Design: an offline NSR (native status) for the BIEN service stack + +Draft 2026-09-18, for BM's review before implementation. Sibling of `04-versioned-divisions-design.md`; +the political-division reference built there is the backbone here. Lives in `RGNRS/dev_notes` +only until an `RNSR` clone exists, then moves with the code. + +## Why this, and what already exists + +NSR is the last BIEN service with no offline implementation: + +| Service | Offline | Where | +|---|---|---| +| TNRS | yes | `RTNRS`: `TNRS_local()`, `local_build.R`, `local_backbone.R`, `tnrs_cache_dir()` | +| GNRS | yes | `RGNRS` branch `versioned-divisions`: `GNRS_local()`, history/cshapes components | +| GVS | yes | `RGVS` branch `local-implementation`: `GVS_local()`, raster index, centroids, CShapes | +| **NSR** | **no** | `NSR` 0.1.0 is API-only: `NSR()`, `NSR_simple()`, `NSR_sources()`, metadata helpers | + +The service's database is built by PHP under `ojalaquellueva/nsr/db` (the separate `nsr_db` repo is +retired), importing 12 checklists. There is no published build of the compiled database, and the +repository documents no licence for the source data, so, as with GADM and CShapes, **sources are +downloaded on the user's machine at build time and nothing derived is shipped**. + +## What the service does, and must be reproduced + +Output (probed live, `Pinus ponderosa` / United States / California), 18 columns: + + family, genus, species, country, state_province, county_parish, poldiv_full, poldiv_type, + native_status_country, native_status_state_province, native_status_county_parish, + native_status, native_status_reason, native_status_sources, isIntroduced, + isCultivatedNSR, is_cultivated_taxon, user_id + +Status codes: **N** native, **Ne** native and endemic, **I** introduced, **Ie** introduced inferred +from endemism elsewhere, **P** present without status, **A** absent from the checklists, **UNK** no +data. A status is returned per political level, plus an overall `native_status` taken at the lowest +division supplied. + +Propagation rules, from the service README: +- **Taxonomic:** native propagates *up* (a native variety makes the species, genus and family native + there); introduced propagates *down* (an introduced species makes its varieties introduced). With + no opinion for the taxon, higher ranks are consulted. +- **Political:** native propagates *up* (native in a county implies native in its state and country); + introduced propagates *down* (introduced in a country implies introduced in its states). +- **Absence:** absent from a checklist gives `A`; absent while endemic elsewhere gives `Ie`. +- `is_comprehensive` per source decides whether absence is informative at all: only a checklist that + claims to list everything in its area can turn "not listed" into evidence. + +## Architecture + +Mirrors the other two packages, so the three share one cache and one idiom. + +- Component `nsr` in the shared cache (`GNRS.cache_dir` / `R_user_dir("GNRS")`), tables prefixed + `nsr-`, provenance RDS per component, `nanoparquet` `.gz.parquet` files. +- `NSR_local_build(sources = c("powo", "vascan", "flbr"), dir, ...)`, one importer per + source behind a common interface (fetch -> standardise -> resolve names -> resolve divisions -> + append), so adding the remaining eight later is additive. +- `NSR_local(occurrence_dataframe, dir, ...)` with the service's input columns + (`species`/`genus`/`family`, `country`, `state_province`, `county_parish`) and its 18 output + columns, so existing code swaps one call for the other. +- **Names** resolve through `TNRS_local(sources = "wcvp")`, the resolver every other layer of this + project already uses; checklist names and query names go through the same call, so matching is + backbone-consistent by construction rather than by string equality. +- **Political divisions** resolve through `GNRS_local()`, which gives the country/state/county ids + and, new in the versioned-divisions work, **historical divisions**: a record labelled "USSR" or + "Zaire" resolves to the entity and its current successors, so native status can be looked up in + the successor whose polygon contains it. The live NSR cannot do this. + +### Cache tables + +| File | One row per | Columns | +|---|---|---| +| `nsr-sources` | source | source_name, full name, url, date_accessed, is_comprehensive, licence, citation, poldiv_level, divisions covered | +| `nsr-checklist` | taxon x division x source | taxon_id, rank, poldiv_id, poldiv_level, status (`native`/`introduced`/`present`), source_id | +| `nsr-taxa` | taxon | taxon_id (WCVP accepted id), family, genus, species, rank, is_cultivated_taxon | +| `nsr-regions` | region | region_key (`system:code`), system, code, name, level, GNRS ids and GADM gid where applicable | +| `nsr-region-links` | pair of regions | from_region, to_region, relation (`within`/`contains`/`overlaps`), fraction | +| `nsr-endemism` | taxon | taxon_id, endemic_poldiv_id, source (for the `Ie` rule) | + +Checklist rows key on resolved ids, never on strings, so a source's own spellings and division names +are a build-time problem only. + +## Sources, first release + +All four are public and redistributable-in-principle, and cover the places where WGSRPD level 3 (the +resolution WCVP gives us) is too coarse to say anything about sub-country status. + +| Source | What it gives | Format, size | Licence | +|---|---|---|---| +| **powo** | global native/introduced by WGSRPD level 3 | WCVP v15 **already in hand** (`data/wcvp_v15/wcvp_distribution.csv`), no download | CC BY 4.0 | +| ~~usda~~ | USA, by state | **not obtainable**: the per-state export is no longer served; the GBIF copy is names only and the API gives region-level status | CC0 | +| **vascan** | Canada, by province; native/introduced/ephemeral | DwC-A, `data.canadensys.net/ipt/archive.do?r=vascan`, ~10 MB, updated 2026-08-04 | CC0 | +| **flbr** | Brazil, by state; native/naturalised/cultivated | DwC-A, `ipt.jbrj.gov.br/jbrj/archive.do?r=lista_especies_flora_brasil`, ~50-100 MB | CC BY 4.0 | + +POWO is the service's global backbone and is the same Kew dataset as WCVP, at the same resolution; +we already apply its native filter in the garden pipeline (script 02 keeps `introduced == 0 & +extinct == 0 & location_doubtful == 0`). The other three add state-level status in the USA, Canada +and Brazil. Brazil matters most: WGSRPD divides it into five regions, so level 3 cannot distinguish +a species native to Amazonia from one introduced to São Paulo. + +**usda is dropped from the release (BM, 2026-09-19): the USDA no longer serves the file.** +NSR's own importer expects `usda_plants_native_status.csv`, an export of PLANTS' native status +by state. Checked 2026-09-18: the GBIF-hosted copy of the PLANTS database +(`hosted-datasets.gbif.org/datasets/usda.zip`, CC0, 2.1 MB) is a COLDP archive of names only - +`NameUsage.tsv` and `VernacularName.tsv`, no distribution file; the `csvdownload` endpoint the +old site used now returns the site shell; and the current API +(`plantsservices.sc.egov.usda.gov/api/PlantProfile?symbol=...`) returns native status only by +region - CAN, L48, AK, HI, PR - not by state, which is the resolution that would have added +anything over POWO. No public bulk route to per-state status was found. + +**Consequence, to state in the write-up:** in the USA the build has POWO alone, whose introduced +ranges are thin there - WCVP v15 holds no California record for *Quercus robur* at all, though it +is naturalised across the state. Five of the fifteen disagreements with the live service in the +API validation cited `usda`. US records will therefore skew towards `A` (absent) rather than `I` +where a species is introduced but unlisted, and `Ie` only rescues those whose native range is +confined elsewhere. Weakley's southeastern flora would cover part of the same ground if a route to +it appears. + +Deferred to a later release: conosur, fwi, ipane, mab, mexico, newguinea, tropicos, weakley. Each +needs its own acquisition route (several are web pages or an API with a key), and two (`mab` +endemic genera, `ipane` introduced-only) are small special-purpose lists rather than checklists. + +## Regions are polygons, not names (BM, 2026-09-18) + +**Each source keeps its own geography, and queries are resolved against it spatially.** +WCVP records distributions against WGSRPD level-3 areas, which are biogeographic, not political; +VASCAN and Flora do Brasil record provinces and states, which are GADM units. Crosswalking either +onto the other by *name* loses data and invents detail: + +- A first implementation matched level-3 unit names to GNRS divisions. Result: 79 of 370 units + matched a state, 194 fell back to country, and **97 (26%) matched nothing at all** - Borneo, + Sulawesi, Maluku, New Guinea, Lesser Sunda Is., Santa Cruz Is. Every POWO opinion for those areas + would have been dropped. +- The reverse fails too: "Brazil Southeast" is one level-3 unit covering four states, so a name + match can say nothing about Sao Paulo, although the containment is exact and knowable. + +Instead: + +| Source | Region system | Key | +|---|---|---| +| powo | WGSRPD level 3 (TDWG polygons, `rWCVPdata::wgsrpd3`) | `wgsrpd3:BZL` | +| vascan, flbr | GADM level 1, via the GNRS backbone's `gid_1` | `gadm1:BRA.25_1` | +| (later sources) | whichever they publish against | `:` | + +**A query with coordinates is resolved natively in both systems**, by raster lookup: the WGSRPD +level-3 raster (built for the garden pipeline) and the GADM unit index (built for GVS), both at 30 +arc-seconds with an exact fallback in boundary cells. No crosswalk is involved, which is the case +that matters for occurrence records. + +**A query with only division names** resolves through GNRS to GADM units, then to WGSRPD areas +through a **spatial link table** computed once at build time (`nsr-region-links`: from, to, relation, +fraction of the smaller unit's area inside the larger). The link table is derived by cross-tabulating +the two rasters, so it needs no new downloads and costs one pass over the cells. + +**Inheritance (revised 2026-09-18, BM): a place is judged by the polygons it lies IN.** +Evidence attaches to the polygon it was recorded for, and is not transferred to other polygons: + +- an opinion about a polygon **containing** the place applies to it (POWO's finest statement about + Guadeloupe is "native in the Leeward Islands"), and among those any native opinion wins; +- polygons **inside** the place describe only parts of it, so they are not its status. They answer + only when they all agree; otherwise the answer is `P`, the reason says the status varies, and + `n_subpolygons_native` / `n_subpolygons_introduced` give the split; +- polygons that merely **overlap** describe neither the place nor its parts, and answer `UNK`; +- `native_status_scope` records which of these produced the answer. + +Give coordinates and the question does not arise: the record is judged on the ground it sits on, +which is how the garden pipeline will use it. + +**Endemism is the deliberate exception** (BM): `Ne` and `Ie` are claims about the taxon's whole +range, not about one polygon, so they draw on evidence from elsewhere. A taxon confined to +California, found in Michigan, is introduced there whatever Michigan's checklists say. + +The earlier design let a sub-polygon set the whole place's status, which made *Polycarpon +tetraphyllum* native in the USA on one `powo:native` area against eleven `powo:introduced` ones. +That is the failure this rule removes. +while "introduced in Brazil Southeast" answers it as `I`. Nothing is invented and nothing is lost. + +## Resolution semantics + +For each row, after name and division resolution: + +1. Gather every checklist opinion for the taxon (or a higher rank, per the taxonomic rule) in the + division (or a containing/contained one, per the political rule), from all built sources. +2. Reduce to one status per political level. **Precedence (BM, 2026-09-18): if any source says + native, the taxon is native there.** Native is the hardest claim for a checklist to make by + accident, whereas "introduced" and "absent" are often an artefact of a list's scope or age. Where + sources disagree, the answer is native and the disagreement is recorded rather than hidden, in + `native_status_conflict` (whether opinions differed) and `native_status_opinions` (the per-source + verdicts, e.g. `powo:native; usda:introduced`). This mirrors `TNRS_local()`'s `Source_conflict`. + Failing a native opinion, prefer the finer political level, then explicit `introduced` over + `present`, then a comprehensive source over a non-comprehensive one. +3. Absence only yields `A`/`Ie` where a comprehensive source covers the division **and** the taxon + is evaluable, i.e. some source holds native-status information about it (see below). +4. `native_status` is the status at the lowest division supplied, and `isIntroduced` is 1 for `I` + and `Ie`. + +New columns beyond the service, marked as extensions so output stays a superset: +`poldiv_is_historical`, `resolved_in_successor` (which current division answered a historical one), +and `native_status_basis` (the source and division that decided it). + +## Absence: when "not listed" is evidence + +**Decisions (BM, 2026-09-18): POWO counts as comprehensive for the species it covers; and a taxon +about which no source holds any native-status information cannot be called Absent or Present at +all.** + +The live service does not make the second distinction, and it matters. Probed 2026-09-18: + +| Query | Live NSR | Comment | +|---|---|---| +| *Sphagnum palustre*, California | **A**, "Absent from all checklists for region", sources powo, usda | a moss; POWO covers vascular plants only | +| *Marchantia polymorpha*, California | **A**, same | a liverwort | +| *Quercus robur*, California | I, sources powo, usda, weakley | correct | +| *Pinus ponderosa*, Brazil | I, source flbr, `isCultivatedNSR` 1 | correct | +| *Zea mays*, Mexico | N, sources mexico, powo | correct | + +A bryophyte in California is reported Absent, when the truth is that no consulted checklist has an +opinion about bryophytes. Absence of evidence is returned as evidence of absence, and downstream +that becomes `isIntroduced` through the `Ie` rule. + +**The rule (BM's formulation): go by which taxa the sources hold native-status information about, +not by what the sources claim to cover.** A taxon is *evaluable* if it has at least one +native-status record in any built source, anywhere in the world. Absence is interpretable only for +evaluable taxa; for the rest the answer is `UNK`, with `native_status_reason` "no source holds +native-status information for this taxon". + +Why this rather than a declared taxonomic scope per source: +- **It does not need to know what a moss is.** *Sphagnum palustre* has no row in WCVP's + distributions, so it is not evaluable, and no `A` is emitted. The server never reasons about + taxonomic groups. +- **It survives new and user-supplied sources.** A curated scope field would have to be written, and + kept right, for every source anyone adds. Coverage here is a property of the data, computed at + build time as the distinct taxa in `nsr-checklist`. +- **It protects newly described species.** A species accepted in WCVP but with no distribution + records yet is not evaluable, so it returns `UNK` instead of being flagged Absent (and then + introduced) everywhere it is found. A declared-scope rule would get this wrong, since the species + is squarely inside POWO's stated scope. + +Implementation is one lookup: `evaluable(taxon_id)` is true when the taxon appears anywhere in +`nsr-checklist`. No scope column, no taxonomic-group inference, and the semantics stay the same as +sources are added: each new checklist can only widen the evaluable set. + +This is a deliberate divergence from the service, and one to report upstream: the same query that +returns `A` there returns `UNK` here. + +## Validation + +1. **Against the live service**, on a stratified sample of BIEN and GBIF records (target ~50,000 + rows spanning all four source regions and a tail of elsewhere): agreement on `native_status`, + a confusion matrix of codes, and disagreements broken down by source, division level and region. + The four-source build cannot match the service everywhere, so the honest target is: agreement + where only these four sources apply, and a measured, explained gap elsewhere. +2. **Against BIEN's stored `is_introduced`** on the 284M-record BIEN pull, the same check we ran for + geovalidity and centroids. +3. **Historical divisions**: records whose country no longer exists resolve and receive a status + through the successor, which the live service cannot do; report how many records that recovers. +4. **Garden pipeline**: the occurrence layer currently keeps a record only if it falls in the + species' native WGSRPD level-3 range (POWO). Re-run with NSR and report how many cells change, + which is the concrete answer to "does the finer native status matter for the paper?" + +## Built, 2026-09-18 (milestones 1-3) + +Branch `local-implementation`: `R/local_cache.R`, `R/local_build.R`, `R/local_import.R`, +`R/local_regions.R`, `R/NSR_local.R`, `tests/testthat/test-local-nsr.R` (24 expectations, all +passing). Built against the shared GNRS/GVS cache in 24 minutes: + +| Source | Records | Taxa | Geography | +|---|---|---|---| +| powo | 1,986,877 | 443,168 | 375 WGSRPD level-3 areas | +| flbr | 154,169 | 38,326 | 27 Brazilian states | +| vascan | 25,459 | 6,086 | 13 Canadian provinces + Greenland | + +Plus 415 regions and **13,742 spatial links** (3,330 within, 3,330 contains, 6,684 overlaps). +POWO needs no download (the WCVP archive is already in the TNRS cache) and no name resolution +(WCVP is the backbone, so `accepted_plant_name_id` is the taxon id). The other two resolve through +`TNRS_local()`: VASCAN matches 99.5% of names; FLBR 75%, the residue being the fungi and algae that +a vascular-plant backbone cannot match - which is exactly the evaluability rule's case. + +### Reference cases + +| Query | Local NSR | Live service | +|---|---|---| +| *Pinus ponderosa*, US/California | N (powo) | N | +| *Araucaria angustifolia*, Brazil/Sao Paulo | N (flbr + powo) | - | +| *Araucaria angustifolia*, Brazil/Amazonas | A | - | +| *Araucaria angustifolia*, by coordinates in Sao Paulo | N, resolved in both geographies | n/a | +| *Acer saccharum*, Canada/New Brunswick | N (vascan + powo) | - | +| *Zea mays*, Mexico | N (powo) | N | +| *Welwitschia mirabilis*, Namibia | N (powo) | N | +| *Sphagnum palustre* / *Marchantia polymorpha*, US/California | **UNK** | **A** | +| *Quercus robur*, US/California | **A** | **I** (powo, usda, weakley) | + +The bryophyte divergence is the intended one. The *Quercus robur* divergence is not: **WCVP holds no +California record for it at all** (22 introduced records elsewhere, none in CAL), so the service's +answer came from USDA and Weakley. POWO's introduced ranges are thinnest exactly where USDA would +cover, which is an argument for keeping `usda` in the first release rather than dropping it (see +Open questions). + +### Bugs worth remembering + +- **Pooling administrative levels in the link table.** Aggregating cells to both level 1 and level 0 + before computing fractions counted every cell twice, so California-the-state and CAL-the-WGSRPD-area + each looked half-inside the other and their relation degraded from `same` to `overlaps` - which + correctly demoted every native opinion to mere presence. Each level now has its own denominators. +- **A state query inheriting its country's answer.** Keys for both levels sat in one pool, so + Amazonas inherited Brazil's "native". The finest named place now answers alone; the country key + only fills the country column. +- **VASCAN's two status columns are the reverse of the obvious guess**: `occurrenceStatus` is + present/excluded/doubtful/irregular/absent, `establishmentMeans` is native/introduced. +- **The GNRS backbone spells some divisions bilingually** ("New Brunswick/Nouveau-Brunswick") and + Newfoundland's HASC is `CA.NF`, not `CA.NL`. + +### Validation against the live service (milestone 4) + +`scratchpad/nsr_validate.R`: a stratified sample of species x division (POWO natives, POWO +introduced, Brazilian states, Canadian provinces, and a tail of random species in random +countries), sent to both. The service drops roughly a third of API batches whatever their size, +so 150 of 300 rows were compared. + +**84% same decision** (78% identical code), against a service with four times as many sources: + +| local \ service | absent | introduced | native | present | unknown | +|---|---|---|---|---|---| +| absent | 6 | 4 | 0 | 0 | 0 | +| introduced | 5 | 49 | 1 | 1 | 1 | +| native | 9 | 1 | 71 | 0 | 1 | +| present | 1 | 0 | 0 | 0 | 0 | + +By stratum: powo 86.7% (n=120), vascan 86.7% (n=15), flbr 60% (n=5), no-source tail 60% (n=10). + +Disagreement classes, all explained: +- **local native, service absent (9).** Four cite no source at all (*Monanthotaxis schweinfurthii* + in DR Congo, *Galium spurium* in Syria, *Phleum* and *Trifolium spadiceum* in Russia); we cite + POWO v15. The service's POWO import is older, so these look like ours being more current. +- **local introduced, service absent (5).** US and Australian cases where POWO records the species + as introduced in an area inside the country and we inherit that upward; the service answers + absent. Ours is the more informative reading of the same data. +- **local absent, service introduced (4)** and **one each way on native/introduced**: sources we + lack (`usda`, `weakley`, `mab`, `fwi`, `conosur`) or POWO version differences. + +Two rule changes came out of reading the disagreements, both now in the code: +- **Inheritance is symmetric.** Any native opinion about a related region makes it native (BM's + precedence); failing that, any introduced opinion makes it introduced. Inferred answers say so in + `native_status_reason`, and `native_status_opinions` records the relation. +- **Links under 1% of a region are ignored**, so raster disagreement along a coastline carries + nothing. The `within` threshold is 0.95, not 0.99: these fractions come from 30 arc-second + rasters of two independently drawn coastlines, so coextensive regions score 0.987-0.990. + +Also added: **endemism** (`Ne`, `isEndemic`), read off the checklist - every region any source calls +the taxon native lies inside the queried place. Verified on single-area natives, and correctly +declining for *Araucaria angustifolia* (also Argentina, Paraguay) and *Welwitschia mirabilis* (also +Angola). Note it is a statement about the checklists, not a conservation claim: *Oryza sativa* in +China returns `Ne`, because POWO gives rice one native area. **Genus queries** are answered +directly (WCVP records distributions for 14,126 genera); what is refused is resolving a bare genus +onto some species of that genus. + +Bugs fixed this round, worth remembering: WCVP carries the same name at two ranks (a *variety* row +also called "Pinus ponderosa"), and only one holds the distributions, so names index to the +species-rank id that has opinions; and accepted infraspecific taxa (60,278 of them) keep their own +distributions, so they are rolled up to their species through `parent_plant_name_id`. + +### Status model completed: `Ie` (2026-09-18) + +Absence becomes introduction only where the taxon could not have been native: its whole native +range lies elsewhere **and is confined**, to one region or to one country. A widely native species +merely unrecorded here stays `A`, because absence alone is not evidence of introduction. + +| Query | Answer | Why | +|---|---|---| +| *Lepechinia calycina* (endemic to California) in Michigan | `Ie`, isIntroduced 1 | it could not be native there | +| *Lepechinia calycina* in California | `Ne`, isEndemic 1 | native and confined to it | +| *Homonoia* (widely native) in Michigan | `A` | absent, but could have been native | +| *Washingtonia* (California, Arizona, NW Mexico) in Michigan | `A` | native range not confined | + +Final validation: **85.3% same decision** (79.3% identical code) on 150 compared rows; by stratum +powo 86.7%, vascan 86.7%, no-source tail 80%, flbr 60% (n=5). `Ie` moved four records from absent +to introduced. + +The nine remaining "we native, service absent" rows were checked individually and are not artefacts +of loose inheritance: *Monanthotaxis schweinfurthii* in DR Congo is a **stated** `powo:native` while +the service returns absent citing nothing; *Galium spurium* in Syria is native within WGSRPD's +Lebanon-Syria unit, POWO's finest statement about Syria; *Trifolium spadiceum* is native in a +Russian sub-area. Over 1,500 country-level queries, inferred answers come through `within` and +`contains`, not weak overlaps. + +**One consequence of the any-native precedence worth stating in the write-up:** at country level a +single native area carries the whole country. *Polycarpon tetraphyllum* in the USA answers `N` from +one `powo:native(contains)` against eleven `powo:introduced(contains)`. The counts are visible in +`native_status_opinions`, and the alternative (majority, or requiring a stated opinion) would lose +the cases the rule exists for, but users filtering on `native_status` at country level should know +it. + +### Validation under the containment rule + +**82.0% same decision** (76.7% identical code), against 85.3% under the old inheritance. The whole +drop is the deliberate change: 9 rows where we now answer `P` (the status varies among the polygons +inside the queried place) while the service commits to one, and 2 where only partly overlapping +polygons have an opinion and we answer `UNK`. Excluding those eleven, agreement is ~89%. Nothing +accidental changed: the same `usda`/`weakley` gaps and POWO-version differences remain. + +By stratum: powo 82.5%, vascan 86.7%, no-source tail 80%, flbr 60% (n=5). + +This is a case where agreement with the service is the wrong target. The service propagates a +sub-region's status upward; we decline to, and say so. For records with coordinates - every record +in the garden pipeline - the two approaches coincide, because the point lies inside exactly one +polygon per system. + +## Milestones + +1. Cache scaffolding, `nsr-sources`, and the POWO importer from WCVP (no download); `NSR_local()` + returning service-shaped output for country-level queries. Validates against the live API on + country-level rows. +2. USDA, VASCAN and Flora do Brasil importers, with the state-level resolution and the political + propagation rules. +3. The taxonomic propagation rules, endemism (`Ie`) and conflict handling. (`is_cultivated_taxon` + was dropped - see open question 3.) +4. Validation 1-3; write-up. +5. Garden re-run (validation 4) and, separately, the remaining eight sources. + +## Open questions for BM + +1. **Conflict precedence.** Is the proposal above (finer division, explicit over present, + comprehensive over not, all opinions recorded) how the service behaves? If BM has the service's + SQL to hand it is faster than inferring it from probes. +2. **POWO as comprehensive.** Treating WCVP's absence as evidence of non-nativeness globally is what + makes `A`/`Ie` possible at all outside the three national checklists. Reasonable, or too strong? +3. Settled (BM, 2026-09-19): **`is_cultivated_taxon` is dropped.** It was already dead in + production - filled from table `cultspp`, loaded from `cultspp_staging` only + `if (exists_table(...))`, with no importer in the public repo populating it; probes return + `is_cultivated_taxon = 0` for both *Zea mays* and *Triticum aestivum*. The candidate + replacements were examined 2026-09-19 and none is worth carrying. + + **Why dropped rather than substituted (BM):** a use flag is not much use without a polygon + attached to it. Everything else NSR reports is an assertion about a taxon *in a place*, and + the containment semantics exist precisely so that polygon-level evidence is not spread to + other polygons. "Grown somewhere by someone" is not evidence that the record in hand was + planted, so a taxon-level flag would either sit unused or be read as if it were + place-specific. `isCultivatedNSR` already answers the version of the question that has a + geography, and is unaffected by this: it comes from checklists carrying a cultivated status + (FLBR flagged *Pinus ponderosa* in Brazil), reproduced from the same source field. + + What was examined, so it need not be examined again: + + - **Kew's World Checklist of Useful Plant Species** (Diazgranados et al. 2020, 40,292 species, + CC BY 4.0) is published only as an 11.3 MB PDF plus its EML (KNB, doi:10.5063/F1CV4G34) - + no data table. (The DOI in the earlier draft of this note, 10.34885/172, was wrong; that is + the State of the World's Plants and Fungi report.) **The PDF was parsed successfully on + 2026-09-19** - it encodes rank in the font, so the hierarchy is recoverable exactly - and the + result is `garden_variety_traits/data/nsr_sources/wcups_2020.csv`, produced by + `R_scripts/48_parse_wcups_pdf.py`. It reproduces every published control total (40,292 + records; 3 kingdoms / 6 phyla / 14 classes / 101 orders / 433 families / 6,737 genera; + 40,239 LSIDs; all ten use-category counts; the top five families and genera; 70 species with + all ten uses; 91 families and 2,790 genera with a single species). So availability is no + longer the obstacle - but the data are still ten use categories per taxon with **no + geography**, so the reason for dropping the flag is unchanged. + - **GRIN Taxonomy** (CC0, COLDP, `hosted-datasets.gbif.org/datasets/grin.zip`, 10.7 MB) does + carry economic uses in `TaxonProperty.tsv`: 17,111 taxa with a use, 9,836 in + cultivation-like classes (ornamental 9,155, plus food, forage, materials, fuel). Machine- + readable and usable - but taxon-level with **no geography at all**, which is the objection + above. Note too that economic use is not cultivation: the largest class after ornamental is + folklore medicine (4,812), largely wild-harvested. + - **Wiersema & Leon's World Economic Plants** has the same taxon-level-only shape. + + GRIN's region-level `cultivated` distribution status was checked as a way to give the flag a + geography, and is not fit for it: 1,437 taxa worldwide, 5,260 rows, 229 ISO3 countries and + nothing finer. Coverage is not the heavily cultivated combinations - *Triticum aestivum*, + *Oryza sativa*, *Glycine max* and *Malus domestica* have no distribution rows at all, and + *Zea mays*, *Manihot esculenta*, *Sorghum bicolor*, *Vitis vinifera* and *Theobroma cacao* + carry native ranges with zero cultivated entries, while *Hordeum vulgare* has 91. The reason + is that GRIN's distribution field documents germplasm provenance, not cultivation extent, so a + cultigen with no natural range often gets nothing; country totals (China 528, US 348, India + 299, Europe thin) track collecting effort rather than area under cultivation. The GBIF COLDP + export also flattens GRIN's sub-country geography onto ISO3 - visible as 71 identical + `iso:RUS` rows for a single taxon, one per Russian region - so even where it has content the + resolution is country-level. + + In the code: `NSR_local()` returns `is_cultivated_taxon = NA_integer_` rather than production's + `0`, so the column keeps its place in the service-shaped output without a never-populated + field being read as a real negative. `NSR_local_by_region()` omits it. + +4. Settled (BM, 2026-09-18): for the garden paper NSR is a **robustness check** on the existing + WCVP native-range filter, not a gate on the occurrence layer, to be revisited if validation 4 + shows it makes a material difference. +5. Settled: sources downloaded 2026-09-18 to `garden_variety_traits/data/nsr_sources/` (gitignored): + usda.zip (2.1 MB), VASCAN DwC-A, Flora do Brasil DwC-A. POWO comes from the WCVP v15 files + already in the repo's data directory. +6. **GRIN as a fourth source?** Open. Examining GRIN for the cultivated flag turned up + something more useful than the flag: `Distribution.tsv` holds 579,184 native/introduced rows + for 65,196 taxa across 229 countries, CC0, independent of POWO. It is country-level only (the + COLDP export flattens GRIN's sub-country geography onto ISO3), so it does not recover the + per-state US resolution that `usda` would have given, and the documented POWO-alone US skew + towards `A` stands either way. `iso:` maps directly onto the `gadm0:` keys the link table + already carries. Awaiting BM's call; not built. + +7. Settled: BM cloned `EnquistLab/RNSR` on 2026-09-19; work is on branch + `local-implementation` and this note lives in the clone. diff --git a/tests/testthat/test-local-nsr.R b/tests/testthat/test-local-nsr.R new file mode 100644 index 0000000..2662b97 --- /dev/null +++ b/tests/testthat/test-local-nsr.R @@ -0,0 +1,138 @@ +# Offline NSR. The pure helpers run anywhere; the rest need a built cache and are skipped +# unless the option NSR.test_cache names one. + +test_that("division names match through bilingual and standard spellings", { + bb <- list( + country = data.frame(country_id = 1L, country = "Canada", gid_0 = "CAN", + stringsAsFactors = FALSE), + state = data.frame(state_province_id = c(10L, 11L), country_id = 1L, country = "Canada", + state_province = c("New Brunswick/Nouveau-Brunswick", "Québec"), + state_province_std = c("New Brunswick/Nouveau-Brunswick", "Quebec"), + hasc_full = c("CA.NB", "CA.QC"), gid_1 = c("CAN.4_1", "CAN.11_1"), + stringsAsFactors = FALSE)) + r <- NSR:::nsr_match_poldiv(rep("Canada", 4), + c("New Brunswick", "Nouveau-Brunswick", "Quebec", "Atlantis"), bb) + expect_equal(r$state_province_id, c(10L, 10L, 11L, NA_integer_)) + expect_equal(r$poldiv_level, c("state_province", "state_province", "state_province", "country")) +}) + +test_that("a place is judged by the polygons containing it", { + db <- list(evaluable = "t1", covered = c("a", "b")) + op <- function(status, relation, source = "powo") { + data.frame(status = status, relation = relation, source_name = source, + is_cultivated = 0L, stringsAsFactors = FALSE) + } + # the polygon itself, and any polygon containing it, answer directly + expect_equal(NSR:::nsr_reduce(op("native", "same"), "a", TRUE, db)$code, "N") + expect_equal(NSR:::nsr_reduce(op("introduced", "within"), "a", TRUE, db)$code, "I") + r <- NSR:::nsr_reduce(op("native", "within"), "a", TRUE, db) + expect_equal(r$code, "N") + expect_equal(r$scope, "containing polygon") + # among containing polygons, any native opinion wins, and disagreement is recorded + r <- NSR:::nsr_reduce(rbind(op("native", "same"), op("introduced", "same", "usda")), + "a", TRUE, db) + expect_equal(r$code, "N") + expect_true(r$conflict) + # polygons INSIDE the place describe parts of it: they answer only when unanimous + r <- NSR:::nsr_reduce(rbind(op("native", "contains"), op("native", "contains")), + "a", TRUE, db) + expect_equal(r$code, "N") + expect_equal(r$scope, "sub-polygons agree") + r <- NSR:::nsr_reduce(rbind(op("native", "contains"), op("introduced", "contains")), + "a", TRUE, db) + expect_equal(r$code, "P") + expect_match(r$reason, "varies among the polygons") + expect_equal(r$n_sub_native, 1L) + expect_equal(r$n_sub_introduced, 1L) + # a containing polygon outranks the sub-polygons + r <- NSR:::nsr_reduce(rbind(op("introduced", "within"), op("native", "contains")), + "a", TRUE, db) + expect_equal(r$code, "I") + # merely overlapping polygons describe neither the place nor its parts + r <- NSR:::nsr_reduce(op("native", "overlaps"), "a", TRUE, db) + expect_equal(r$code, "UNK") + expect_match(r$reason, "partly overlapping") +}) + +test_that("absence is only read for taxa the sources know about", { + db <- list(evaluable = "t1", covered = "a") + # evaluable taxon, comprehensively listed region, no opinion -> absent + r <- NSR:::nsr_reduce(NULL, "a", TRUE, db) + expect_equal(r$code, "A") + # a taxon no source holds information about -> unknown, never absent + r <- NSR:::nsr_reduce(NULL, "a", FALSE, db) + expect_equal(r$code, "UNK") + expect_match(r$reason, "No source holds native status") + # evaluable taxon, but no comprehensive source covers the region -> unknown + expect_equal(NSR:::nsr_reduce(NULL, "z", TRUE, db)$code, "UNK") +}) + +cache <- getOption("NSR.test_cache", "") + +test_that("NSR_local answers the reference cases (needs a built cache)", { + skip_if(!nzchar(cache) || !file.exists(file.path(cache, "nsr-checklist.gz.parquet")), + "no built NSR cache (option NSR.test_cache)") + x <- data.frame( + species = c("Pinus ponderosa", "Araucaria angustifolia", "Araucaria angustifolia", + "Acer saccharum", "Sphagnum palustre", "Araucaria angustifolia"), + country = c("United States", "Brazil", "Brazil", "Canada", "United States", NA), + state_province = c("California", "Sao Paulo", "Amazonas", "New Brunswick", "California", NA), + county_parish = "", + latitude = c(NA, NA, NA, NA, NA, -23.5), + longitude = c(NA, NA, NA, NA, NA, -47.5), + stringsAsFactors = FALSE) + r <- NSR_local(x, dir = cache, quiet = TRUE) + expect_equal(r$native_status, c("N", "N", "A", "N", "UNK", "N")) + # Sao Paulo is answered by the Brazilian flora, which WGSRPD cannot reach + expect_match(r$native_status_sources[2], "flbr") + # a moss is unknown, not absent: POWO holds no information about it + expect_false(r$taxon_evaluable[5]) + # coordinates are resolved in each source's own geography + expect_equal(r$regions_matched[6], "coordinates") +}) + +test_that("region links relate the geographies sensibly (needs a built cache)", { + skip_if(!nzchar(cache) || !file.exists(file.path(cache, "nsr-region-links.gz.parquet")), + "no built NSR cache") + lk <- as.data.frame(nanoparquet::read_parquet(file.path(cache, "nsr-region-links.gz.parquet"))) + expect_true(all(lk$relation %in% c("same", "within", "contains", "overlaps"))) + expect_true(all(lk$fraction >= 0 & lk$fraction <= 1 + 1e-9)) + # Sao Paulo lies inside WGSRPD's Brazil Southeast + sp <- lk[lk$from_region == "gadm1:BRA.25_1" & lk$to_region == "wgsrpd3:BZL", ] + expect_equal(nrow(sp), 1) + expect_true(sp$relation %in% c("within", "same")) + # and the relation is recorded from both sides + rev <- lk[lk$from_region == "wgsrpd3:BZL" & lk$to_region == "gadm1:BRA.25_1", ] + expect_equal(nrow(rev), 1) + expect_true(rev$relation %in% c("contains", "same")) +}) + +test_that("absence becomes introduction only for taxa confined elsewhere", { + skip_if(!nzchar(cache) || !file.exists(file.path(cache, "nsr-checklist.gz.parquet")), + "no built NSR cache (option NSR.test_cache)") + db <- NSR:::nsr_local_db(cache) + nat <- db$checklist[db$checklist$status == "native", ] + per <- table(nat$taxon_id) + # a species POWO gives a single native area, queried far away + cal <- intersect(names(per)[per == 1], nat$taxon_id[nat$region_key == "wgsrpd3:CAL"]) + skip_if(!length(cal), "no single-area Californian native in this build") + sp <- db$taxa$species_name[match(cal[1], db$taxa$taxon_id)] + r <- NSR_local(data.frame(species = c(sp, sp), country = "United States", + state_province = c("Michigan", "California"), + county_parish = "", stringsAsFactors = FALSE), + dir = cache, quiet = TRUE) + expect_equal(r$native_status, c("Ie", "Ne")) + expect_equal(r$isIntroduced, c(1L, 0L)) + expect_equal(r$isEndemic, c(0L, 1L)) + expect_match(r$native_status_reason[1], "endemic to") + # a widely native species merely unrecorded there stays absent + wide <- names(per)[per > 20] + wide <- setdiff(wide, db$checklist$taxon_id[db$checklist$region_key == "wgsrpd3:MIC"]) + skip_if(!length(wide), "no widespread species absent from Michigan") + w <- db$taxa$species_name[match(wide[1], db$taxa$taxon_id)] + r2 <- NSR_local(data.frame(species = w, country = "United States", state_province = "Michigan", + county_parish = "", stringsAsFactors = FALSE), + dir = cache, quiet = TRUE) + expect_equal(r2$native_status, "A") + expect_equal(r2$isIntroduced, 0L) +}) From 88c5d0cb08527f59d058e4580bc4ad04d548bc6f Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Tue, 22 Sep 2026 16:03:25 -0400 Subject: [PATCH 2/8] Address PR review: packaging, resolver bugs, and the WCVP source name Packaging (R CMD check now returns Status: OK, from 1 NOTE and undeclared namespaces): - DESCRIPTION declares what the offline path actually calls: nanoparquet, terra, sf, TNRS, rWCVP and rWCVPdata in Suggests, stats/tools/utils in Imports. rWCVPdata is not on CRAN, so Additional_repositories names its drat. Everything offline stays a Suggest, behind a new nsr_need() guard, so an install with none of it keeps the API functions working. - NAMESPACE and man/ regenerated: NSR_local, NSR_local_build, NSR_local_status, NSR_local_remove and NSR_local_by_region were tagged @export but never exported, so the local API was unreachable from an installed package and the cached tests could not have run. - R/local_globals.R declares the data.table NSE column names. Resolver bugs: - The set path never populated country_code/state_code, so every NSR_local() call of 200 rows or more returned NA for native_status_country and native_status_state_province. The per-level codes are now a second and third pass of the same reducer, matching what the row path does per row; endemism is not applied to them, since Ne/Ie are claims about a whole range rather than about one level. - NSR_local_build() never built the WGSRPD raster for a single-system build, because the only call sat inside the link step. A wcvp-only cache had WGSRPD opinions and no raster, so coordinate queries silently matched nothing. It is now built whenever a source publishes against WGSRPD. - nsr_consulted_regions() ignored min_overlap while every neighbouring helper applied it, so a sliver link could turn a result into A or suppress an Ie. The threshold is now threaded through all four callers. - nsr_import_wcvp() filtered dist but not the parallel index j, so family, genus and rank came from unfiltered row positions while every other column came from the filtered frame - a recycling error, or silently wrong taxonomy where the lengths happened to divide. - Extinct WCVP records were read and then ignored, so a taxon that no longer occurs in a region came back native there. They are now dropped. Doubtful records stay as "present": too weak for native or introduced, but dropping them would read as a confident absence. - isCultivatedNSR was 0 for absence in the set path and NA in the row path, so the same query changed answer on batch size alone. Both now say NA. Likewise scope, conflict_type and the sub-polygon counts for rows with no taxon or no place. Contract: - NSR_local() carries user_id through, or generates sequential ids, on NSR()'s terms, so results join back to their input. - county_parish now warns that it is not resolved offline instead of returning a silent NA, and the help says so. - NSR_local_build() documented a download step it does not have, and defaulted to two sources that always stop without a local archive. It now defaults to the one source that can build unattended and the help states that the archives are not fetched for you. The Kew source is named wcvp rather than powo (BM): it is the WCVP archive that is read, and TNRS_local(sources = "wcvp") is what the sibling package calls the same data. "powo" is accepted as a synonym in sources= and files=. Known cost, recorded in the design note: native_status_sources no longer matches the live service's label for this source. tests/testthat/test-local-nsr-paths.R compares the row and set paths column for column on tables held in memory, so the 200-row switch cannot diverge again without a failure. nsr_index_db() was split out of nsr_local_db() to make that possible without a built cache. Not changed: Copilot's finding that terra::extract()[, 1] picks up an ID column. extract() adds ID for a SpatVector, not for the coordinate matrix passed here; verified against terra 1.9.50, where the call returns the raster value. Co-Authored-By: Claude Opus 5 --- DESCRIPTION | 14 +++- NAMESPACE | 5 ++ NEWS | 5 ++ R/NSR_local.R | 68 ++++++++++++++-- R/NSR_local_set.R | 48 ++++++++++-- R/local_build.R | 38 ++++++--- R/local_cache.R | 37 ++++++++- R/local_globals.R | 14 ++++ R/local_import.R | 14 +++- R/local_regions.R | 1 + dev_notes/01-offline-nsr-design.md | 8 +- man/NSR_local.Rd | 73 +++++++++++++++++ man/NSR_local_build.Rd | 55 +++++++++++++ man/NSR_local_by_region.Rd | 37 +++++++++ man/NSR_local_remove.Rd | 19 +++++ man/NSR_local_status.Rd | 18 +++++ tests/testthat/test-local-nsr-paths.R | 108 ++++++++++++++++++++++++++ 17 files changed, 529 insertions(+), 33 deletions(-) create mode 100644 R/local_globals.R create mode 100644 man/NSR_local.Rd create mode 100644 man/NSR_local_build.Rd create mode 100644 man/NSR_local_by_region.Rd create mode 100644 man/NSR_local_remove.Rd create mode 100644 man/NSR_local_status.Rd create mode 100644 tests/testthat/test-local-nsr-paths.R diff --git a/DESCRIPTION b/DESCRIPTION index 019ef2b..c895f59 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -14,14 +14,24 @@ Encoding: UTF-8 LazyData: true Imports: jsonlite, - httr + httr, + stats, + tools, + utils Suggests: knitr, data.table, + nanoparquet, rmarkdown, + rWCVP, + rWCVPdata, + sf, + terra, testthat, devtools, BIEN, + TNRS, vcr (>= 0.6.0) +Additional_repositories: https://matildabrown.github.io/drat VignetteBuilder: knitr -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 diff --git a/NAMESPACE b/NAMESPACE index 12632d5..0d3cc0f 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -4,6 +4,11 @@ export(NSR) export(NSR_citations) export(NSR_data_dictionary) export(NSR_from_coordinates) +export(NSR_local) +export(NSR_local_build) +export(NSR_local_by_region) +export(NSR_local_remove) +export(NSR_local_status) export(NSR_metadata) export(NSR_political_divisions) export(NSR_simple) diff --git a/NEWS b/NEWS index a985c9d..920d3f8 100644 --- a/NEWS +++ b/NEWS @@ -4,6 +4,11 @@ ## NEW FEATURES * Added new function to check native status using species name and coordinates (NSR_from_coordinates) +* Added an offline implementation: NSR_local() answers without an internet connection from a + local cache built by NSR_local_build(), with NSR_local_status(), NSR_local_remove() and + NSR_local_by_region() alongside it. The checklist sources are acquired on the user's + machine; nothing derived ships with the package. The offline path needs data.table, + nanoparquet, terra, sf and TNRS, all of which are Suggests. ## Bug fixes diff --git a/R/NSR_local.R b/R/NSR_local.R index c711254..54959f5 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -38,7 +38,13 @@ #' \code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none. #' @param min_overlap Ignore region links covering less than this share of a region. #' @param quiet Suppress progress messages? -#' @return A data.frame, one row per input row. +#' @return A data.frame, one row per input row, carrying \code{user_id} from the input +#' (or sequential ids where the input has none), as \code{\link{NSR}} does. +#' @section County and parish: +#' \code{county_parish} is accepted and echoed, but not resolved: the political-division +#' backbone stops at state and province, so \code{native_status_county_parish} is always +#' \code{NA} and the answer is given at the finest level that could be matched. A +#' warning says so. Coordinates are the way to ask a finer question. #' @export NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names = TRUE, min_overlap = 0.01, quiet = FALSE) { @@ -55,6 +61,23 @@ NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names lon <- num("longitude") lat <- num("latitude") n <- nrow(x) + + # user_id on NSR()'s terms, so a result can be joined back to its input either way + user_id <- if ("user_id" %in% names(x)) x[["user_id"]] else rep(NA, n) + if (all(is.na(user_id))) user_id <- seq_len(n) + if (any(duplicated(user_id))) { + stop("user_id should be either null or populated by unique values", call. = FALSE) + } + + # No county polygons exist offline, so county_parish cannot be answered. The service + # answers it, so say plainly that this one does not rather than returning a silent NA. + if (any(!is.na(county) & nzchar(county))) { + warning("county_parish is not resolved offline: the backbone carries no county ", + "polygons, so native_status_county_parish is NA and the answer is given at ", + "the finest level that could be matched. Give coordinates for a finer answer.", + call. = FALSE) + } + db <- nsr_local_db(dir) # ---- taxa ------------------------------------------------------------------------- @@ -109,6 +132,7 @@ NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names native_status_opinions = res$opinions, regions_matched = places$matched, taxon_evaluable = !is.na(taxon_id) & taxon_id %in% db$evaluable, + user_id = user_id, stringsAsFactors = FALSE ) rownames(out) <- NULL @@ -121,6 +145,7 @@ nsr_session <- new.env(parent = emptyenv()) #' @keywords internal #' @noRd nsr_local_db <- function(dir) { + nsr_need("nanoparquet") stamp <- paste(normalizePath(dir), file.mtime(nsr_table_path("checklist", dir)), file.mtime(nsr_table_path("region-links", dir))) hit <- nsr_session$db @@ -136,6 +161,20 @@ nsr_local_db <- function(dir) { db$links <- if (file.exists(lf)) as.data.frame(nanoparquet::read_parquet(lf)) else data.frame(from_region = character(0), to_region = character(0), relation = character(0), fraction = numeric(0), stringsAsFactors = FALSE) + db <- nsr_index_db(db) + nsr_session$db <- list(stamp = stamp, value = db) + db +} + +#' Derive every lookup index the resolvers use from the four cache tables +#' +#' Internal. Split out from \code{nsr_local_db()} so the resolvers can be exercised on +#' tables held in memory, without a built cache: the row path and the set path must give +#' the same answer for the same query, and that is only testable if a db can be made +#' without a multi-gigabyte build. +#' @keywords internal +#' @noRd +nsr_index_db <- function(db) { db$link_idx <- split(seq_len(nrow(db$links)), db$links$from_region) db$chk_idx <- split(seq_len(nrow(db$checklist)), paste(db$checklist$taxon_id, db$checklist$region_key)) @@ -170,7 +209,6 @@ nsr_local_db <- function(dir) { db$evaluable <- unique(db$checklist$taxon_id) db$comprehensive <- db$sources$source_name[db$sources$is_comprehensive %in% TRUE] db$covered <- unique(db$checklist$region_key[db$checklist$source_name %in% db$comprehensive]) - nsr_session$db <- list(stamp = stamp, value = db) db } @@ -244,6 +282,12 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { out$code[i] <- "UNK" out$reason[i] <- if (is.na(tid)) "Taxon not matched to a species in the backbone" else "Place not matched to any region" + # nothing was consulted, which is a definite "no disagreement, no sub-polygons" + # rather than an unknown; the set path says the same for these rows + out$scope[i] <- "none" + out$conflict_type[i] <- "none" + out$n_sub_native[i] <- 0L + out$n_sub_introduced[i] <- 0L next } op <- nsr_opinions_for(tid, keys, db, min_overlap) @@ -259,7 +303,7 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { op <- if (is.null(op)) add else Map(c, op, add) } } - consulted <- nsr_consulted_regions(keys, db) + consulted <- nsr_consulted_regions(keys, db, min_overlap) ev <- tid %in% db$evaluable r <- nsr_reduce(op, consulted, ev, db) if (identical(r$code, "N") && nsr_is_endemic(tid, keys, db, min_overlap)) { @@ -277,10 +321,12 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { for (f in names(r)) out[[f]][i] <- r[[f]] # per-level codes, for the service's columns out$country_code[i] <- if (length(ctry)) { - nsr_reduce(nsr_opinions_for(tid, ctry, db, min_overlap), nsr_consulted_regions(ctry, db), ev, db)$code + nsr_reduce(nsr_opinions_for(tid, ctry, db, min_overlap), + nsr_consulted_regions(ctry, db, min_overlap), ev, db)$code } else NA_character_ out$state_code[i] <- if (length(fine)) { - nsr_reduce(nsr_opinions_for(tid, fine, db, min_overlap), nsr_consulted_regions(fine, db), ev, db)$code + nsr_reduce(nsr_opinions_for(tid, fine, db, min_overlap), + nsr_consulted_regions(fine, db, min_overlap), ev, db)$code } else NA_character_ } out @@ -350,7 +396,7 @@ nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { if (is.null(rows)) return(NA_character_) nat <- unique(db$chk_region[rows][db$chk_status[rows] == "native"]) if (!length(nat)) return(NA_character_) - if (any(nat %in% nsr_consulted_regions(keys, db))) return(NA_character_) + if (any(nat %in% nsr_consulted_regions(keys, db, min_overlap))) return(NA_character_) nm <- function(k) { i <- match(k, db$regions$region_key) if (is.na(i)) k else db$regions$region_name[i] @@ -372,13 +418,19 @@ nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { #' Internal. Used to decide whether absence is interpretable: a query about California #' finds no source publishing against GADM California, but WGSRPD's CAL is the same place #' and is comprehensively listed, so absence there does mean something. +#' +#' \code{min_overlap} applies here for the same reason it applies to opinions: a link +#' covering a sliver of a region is not that region being listed, and letting one through +#' would turn "no opinion" into \code{A}, or suppress an \code{Ie}, on the strength of a +#' shared coastline. The set path filters its links once, when it builds them. #' @keywords internal #' @noRd -nsr_consulted_regions <- function(keys, db) { +nsr_consulted_regions <- function(keys, db, min_overlap = 0.01) { if (!length(keys)) return(character(0)) li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) if (!length(li)) return(keys) - unique(c(keys, db$link_to[li][db$link_rel[li] %in% c("same", "within", "contains")])) + keep <- db$link_rel[li] %in% c("same", "within", "contains") & db$link_frac[li] >= min_overlap + unique(c(keys, db$link_to[li][keep])) } #' What kind of disagreement is this? diff --git a/R/NSR_local_set.R b/R/NSR_local_set.R index f160563..4458670 100644 --- a/R/NSR_local_set.R +++ b/R/NSR_local_set.R @@ -86,10 +86,42 @@ NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, #' Resolve every query at once, with joins #' #' Internal. Needs \code{data.table}; \code{\link{NSR_local}} falls back to the row path -#' without it. Returns the same list of per-row vectors as \code{nsr_resolve_status()}. +#' without it. Returns the same list of per-row vectors as \code{nsr_resolve_status()}, +#' \emph{including} the per-level codes: the reduction is run once for the answer and once +#' more against each political level on its own, which is what the row path does per row. +#' Endemism is deliberately not applied to the per-level codes - \code{Ne} and \code{Ie} +#' are claims about the taxon's whole range rather than about one level - so that the two +#' implementations agree column for column. #' @keywords internal #' @noRd nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { + n <- length(taxon_id) + out <- nsr_resolve_status_set1(taxon_id, places, db, min_overlap, endemism = TRUE) + empty <- rep(list(character(0)), n) + # a level reports a code only when it was asked about AND there is a taxon to ask + # about: with no match in the backbone no lookup happened at any level, which is what + # the row path records + asked <- !is.na(taxon_id) + nc <- lengths(places$country) + nf <- lengths(places$fine) + if (any(nc > 0)) { + lv <- nsr_resolve_status_set1(taxon_id, list(fine = empty, country = places$country), + db, min_overlap, endemism = FALSE) + out$country_code <- ifelse(asked & nc > 0, lv$code, NA_character_) + } + if (any(nf > 0)) { + lv <- nsr_resolve_status_set1(taxon_id, list(fine = places$fine, country = empty), + db, min_overlap, endemism = FALSE) + out$state_code <- ifelse(asked & nf > 0, lv$code, NA_character_) + } + out +} + +#' One pass of the set resolver, against the polygons it is given +#' @keywords internal +#' @noRd +nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, + endemism = TRUE) { dt <- function(...) data.table::data.table(...) n <- length(taxon_id) blank <- rep(NA_character_, n) @@ -109,7 +141,7 @@ nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { dt(qid = rep(seq_len(n), nf), region_key = unlist(fine), relation = "same"), dt(qid = rep(seq_len(n), nc), region_key = unlist(ctry), relation = "same"))) # fixed below for rows that have a finer key - if (!nrow(Q)) return(nsr_set_finish(out, taxon_id, db, n)) + if (!nrow(Q)) return(nsr_set_finish(out, taxon_id, db, n, endemism = endemism)) Q[, relation := data.table::fifelse(has_fine[qid] & region_key %in% unlist(ctry), "within", relation)] Q[, taxon_id := taxon_id[qid]] Q <- Q[!is.na(taxon_id) & !is.na(region_key)] @@ -135,7 +167,7 @@ nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { # the opinions themselves O <- db$chk_dt[Q, on = c("taxon_id", "region_key"), nomatch = 0L, allow.cartesian = TRUE] - if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q)) + if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism)) O[, `:=`(inc = relation %in% c("same", "within"), sub = relation == "contains", ovl = relation == "overlaps")] @@ -222,13 +254,13 @@ nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { out$n_sub_native <- data.table::fifelse(is.na(res$n_sub_nat), 0L, as.integer(res$n_sub_nat)) out$n_sub_introduced <- data.table::fifelse(is.na(res$n_sub_int), 0L, as.integer(res$n_sub_int)) out$conflict <- !is.na(res$conflict_type) & res$conflict_type != "none" - nsr_set_finish(out, taxon_id, db, n, Q) + nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism) } #' Fill in the answers that need no opinions: absence, unknowns, and endemism #' @keywords internal #' @noRd -nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL) { +nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL, endemism = TRUE) { dt <- function(...) data.table::data.table(...) evaluable <- !is.na(taxon_id) & taxon_id %in% db$evaluable # the regions each query consulted, for the coverage test @@ -252,7 +284,9 @@ nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL) { out$conflict_type[no_answer] <- "none" out$n_sub_native[is.na(out$n_sub_native)] <- 0L out$n_sub_introduced[is.na(out$n_sub_introduced)] <- 0L - out$cultivated[out$code == "A"] <- 0L + # absence carries no opinion, so it carries no cultivation flag either: NA, matching + # nsr_reduce(). A 0 here would read as a checked negative and would differ from the + # row path for the same query, on nothing but batch size. # rows with no place at all, or no taxon none <- is.na(taxon_id) @@ -263,7 +297,7 @@ nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL) { } # --- endemism, the one rule allowed to look outside the queried polygon -------------- - if (!is.null(Q) && nrow(Q)) { + if (endemism && !is.null(Q) && nrow(Q)) { inside <- Q[relation %in% c("same", "contains"), .(qid, region_key)] nat <- db$native_dt tq <- dt(qid = seq_len(n), taxon_id = taxon_id)[!is.na(taxon_id)] diff --git a/R/local_build.R b/R/local_build.R index cd5e437..fc18c9f 100644 --- a/R/local_build.R +++ b/R/local_build.R @@ -1,10 +1,20 @@ #' Build the local native-status reference #' -#' Downloads (or reads) each checklist source, resolves its names against WCVP with +#' Reads each checklist source, resolves its names against WCVP with #' \code{TNRS::TNRS_local()} and its political divisions against the GNRS backbone, and -#' writes the shared cache tables. Sources are fetched on the user's machine; nothing +#' writes the shared cache tables. Sources are acquired on the user's machine; nothing #' derived ships with the package. #' +#' \strong{The archives are not downloaded for you.} \code{vascan} and \code{flbr} must +#' be supplied through \code{files}, and each stops with the URL to fetch if they are +#' not; \code{wcvp} is the exception, since it reads the WCVP archive already in the +#' TNRS cache when there is one. This is why \code{sources} defaults to \code{"wcvp"} +#' alone: it is the only source that can build unattended. +#' +#' \code{"powo"} is accepted as a synonym for \code{"wcvp"} - it is the name the live +#' service reports for the same Kew dataset - but \code{wcvp} is what gets written to +#' \code{native_status_sources}. +#' #' Tables written: \code{nsr-sources} (one per source, with licence and whether it is #' comprehensive), \code{nsr-taxa} (one per taxon, keyed on the WCVP accepted id), #' \code{nsr-regions} (one per region, in the geography its source publishes against) and @@ -14,14 +24,15 @@ #' #' @param sources Which sources to build. See \code{NSR_local_status()}. #' @param dir Cache directory, shared with GNRS and GVS. -#' @param files Named list of local archives to read instead of downloading, e.g. -#' \code{list(flbr = "flbr_dwca.zip")}. For \code{powo}, the WCVP zip; if absent, the -#' copy in the TNRS cache is used when there is one. +#' @param files Named list of local archives to read, e.g. +#' \code{list(flbr = "flbr_dwca.zip")}. Required for \code{vascan} and \code{flbr}. +#' For \code{wcvp}, the WCVP zip; if absent, the copy in the TNRS cache is used when +#' there is one. #' @param overwrite Rebuild sources that are already built? #' @param quiet Suppress progress messages? #' @return Invisibly, \code{NSR_local_status()}. #' @export -NSR_local_build <- function(sources = c("powo", "vascan", "flbr"), +NSR_local_build <- function(sources = "wcvp", dir = nsr_cache_dir(create = TRUE), files = list(), overwrite = FALSE, quiet = FALSE) { for (pkg in c("nanoparquet", "TNRS")) { @@ -30,7 +41,8 @@ NSR_local_build <- function(sources = c("powo", "vascan", "flbr"), } } reg <- nsr_builtin_registry() - sources <- match.arg(sources, names(reg), several.ok = TRUE) + sources <- match.arg(nsr_canonical_source(sources), names(reg), several.ok = TRUE) + if (length(files)) names(files) <- nsr_canonical_source(names(files)) msg <- function(...) if (!quiet) message(...) existing <- lapply(nsr_tables(), function(t) { @@ -51,7 +63,7 @@ NSR_local_build <- function(sources = c("powo", "vascan", "flbr"), msg("Importing ", s, " ...") t0 <- Sys.time() raw[[s]] <- switch(s, - powo = nsr_import_powo(files$powo, bb, quiet), + wcvp = nsr_import_wcvp(files$wcvp, bb, quiet), vascan = nsr_import_vascan(files$vascan, bb, quiet), flbr = nsr_import_flbr(files$flbr, bb, quiet) ) @@ -136,8 +148,16 @@ NSR_local_build <- function(sources = c("powo", "vascan", "flbr"), n_records = nrow(checklist), n_taxa = nrow(taxa), n_regions = nrow(regions)), nsr_provenance_path("nsr", dir)) - # spatial links between geographies, where GVS's GADM index is available + # The WGSRPD raster is what places coordinates in a WCVP region, so it is needed + # whenever a source publishes against WGSRPD - not only when there is a second + # geography to link it to. A wcvp-only build needs it just as much. systems <- unique(regions$system) + if ("wgsrpd3" %in% systems && !file.exists(nsr_raster_path("wgsrpd3", dir))) { + msg("Building the WGSRPD level-3 raster ...") + try(nsr_build_wgsrpd_raster(dir = dir, quiet = quiet), silent = FALSE) + } + + # spatial links between geographies, where GVS's GADM index is available if (length(systems) > 1 && file.exists(file.path(dir, "gadmindex-units-30s.tif"))) { msg("Linking region systems (", paste(systems, collapse = ", "), ") ...") try(nsr_build_region_links(dir = dir, quiet = quiet), silent = FALSE) diff --git a/R/local_cache.R b/R/local_cache.R index 215b91e..ce4ad0e 100644 --- a/R/local_cache.R +++ b/R/local_cache.R @@ -31,6 +31,36 @@ nsr_provenance_path <- function(component, dir = nsr_cache_dir()) { #' @noRd nsr_tables <- function() c("sources", "taxa", "regions", "checklist") +#' Resolve a source's accepted synonyms to its canonical name +#' +#' Internal. The Kew dataset is read as WCVP and reported as \code{wcvp}, which is what +#' it is and what \code{TNRS_local()} calls the same data. \code{powo}, the portal the +#' live service names it after, is accepted so existing calls and notes keep working. +#' @keywords internal +#' @noRd +nsr_canonical_source <- function(x) { + if (!length(x)) return(x) + alias <- c(powo = "wcvp") + i <- match(x, names(alias)) + unname(ifelse(is.na(i), x, alias[i])) +} + +#' Require the optional packages a local-NSR path needs +#' +#' Internal. The offline implementation is all Suggests: the API functions must keep +#' working on an install that has none of it, so every entry point says what is missing +#' rather than failing inside a \code{::} call. +#' @keywords internal +#' @noRd +nsr_need <- function(..., what = "The local NSR") { + for (pkg in c(...)) { + if (!requireNamespace(pkg, quietly = TRUE)) { + stop(what, " needs the '", pkg, "' package.", call. = FALSE) + } + } + invisible(TRUE) +} + #' Checklist sources the local NSR can build #' #' Internal. One entry per source: what it covers, whether it is comprehensive for its @@ -45,8 +75,8 @@ nsr_tables <- function() c("sources", "taxa", "regions", "checklist") #' @noRd nsr_builtin_registry <- function() { list( - powo = list( - source_name = "powo", + wcvp = list( + source_name = "wcvp", source_name_full = "Plants of the World Online / World Checklist of Vascular Plants (WCVP)", version = "v15", url = "https://sftp.kew.org/pub/data-repositories/WCVP/", @@ -72,7 +102,7 @@ nsr_builtin_registry <- function() { # usda is deliberately absent: the USDA no longer serves the per-state native-status # export the service was built from. The GBIF copy of PLANTS is names only, and the # current API gives status by region (CAN / L48 / AK / HI), not by state, which adds - # nothing over POWO. See dev_notes/01-offline-nsr-design.md. + # nothing over WCVP. See dev_notes/01-offline-nsr-design.md. flbr = list( source_name = "flbr", source_name_full = "Flora e Funga do Brasil (Lista Oficial)", @@ -95,6 +125,7 @@ nsr_builtin_registry <- function() { #' records and taxa it contributed, and when it was built. #' @export NSR_local_status <- function(dir = nsr_cache_dir()) { + nsr_need("nanoparquet") reg <- nsr_builtin_registry() built_tables <- vapply(nsr_tables(), function(t) file.exists(nsr_table_path(t, dir)), logical(1)) src <- if (built_tables[["sources"]]) { diff --git a/R/local_globals.R b/R/local_globals.R new file mode 100644 index 0000000..012f718 --- /dev/null +++ b/R/local_globals.R @@ -0,0 +1,14 @@ +# The offline resolver is written against data.table, which refers to columns by bare +# name. R CMD check cannot tell those from undefined globals, so they are declared here +# rather than silenced one call at a time. Keep this list in step with the set path. +utils::globalVariables(c( + ".", ".N", ".SD", ":=", + "code", "conflict_type", "consulted_here", "fraction", "from_region", + "has_int", "has_nat", "in_place", "is_cultivated", + "n_here", "n_in", "n_nat", "n_src_sets", "n_sub_int", "n_sub_nat", + "ovl", "ovl_srcs", "qid", "r_ord", "reason", "region_key", "relation", + "scope", "source_name", "srcs", "srcs_any", + "stated_int", "stated_nat", "stated_pre", "status", + "sub_sets", "sub_srcs", "sub_status", + "taxon_id", "to_region", "to_sys", "within_src" +)) diff --git a/R/local_import.R b/R/local_import.R index 9441fab..97d04e0 100644 --- a/R/local_import.R +++ b/R/local_import.R @@ -50,10 +50,10 @@ nsr_unzip_member <- function(zip, member) { #' taxon id directly. #' @keywords internal #' @noRd -nsr_import_powo <- function(zip = NULL, bb, quiet = FALSE) { +nsr_import_wcvp <- function(zip = NULL, bb, quiet = FALSE) { if (is.null(zip)) zip <- nsr_find_wcvp_zip() if (is.null(zip) || !file.exists(zip)) { - stop("Could not find the WCVP archive. Supply it as files = list(powo = \"wcvp-v15.zip\"), ", + stop("Could not find the WCVP archive. Supply it as files = list(wcvp = \"wcvp-v15.zip\"), ", "or build the TNRS 'wcvp' source, whose cache holds one.", call. = FALSE) } dist_f <- nsr_unzip_member(zip, "wcvp_distribution.csv") @@ -78,9 +78,17 @@ nsr_import_powo <- function(zip = NULL, bb, quiet = FALSE) { dist$taxon_id <- sp dist$taxon_name <- nm$taxon_name[j] dist$rank_published <- nm$taxon_rank[match(acc, nm$plant_name_id)] - dist <- dist[!is.na(dist$taxon_id) & !is.na(dist$taxon_name), , drop = FALSE] flag <- function(x) !is.na(x) & x %in% c("1", 1, TRUE, "TRUE") + # An extinct record is not a standing distribution. Left in, a taxon that no longer + # occurs in a region comes back native there, and counts as evaluable, which also moves + # absence and endemism answers elsewhere. A DOUBTFUL record is kept, but only as + # "present": it is too weak to support native or introduced, and dropping it would let + # the region read as a confident absence when the source in fact records a maybe. + keep <- !is.na(dist$taxon_id) & !is.na(dist$taxon_name) & !flag(dist$extinct) + dist <- dist[keep, , drop = FALSE] + j <- j[keep] # j indexes nm per dist row, so it is subset with it + out <- data.frame( taxon_name = dist$taxon_name, taxon_id = dist$taxon_id, species_name = dist$taxon_name, family = nm$family[j], genus = nm$genus[j], diff --git a/R/local_regions.R b/R/local_regions.R index 3c60385..962f4cf 100644 --- a/R/local_regions.R +++ b/R/local_regions.R @@ -48,6 +48,7 @@ nsr_build_wgsrpd_raster <- function(dir = nsr_cache_dir(create = TRUE), resoluti #' @keywords internal #' @noRd nsr_wgsrpd_polygons <- function() { + nsr_need("sf", what = "The WGSRPD polygons") if (requireNamespace("rWCVPdata", quietly = TRUE)) { p <- try(rWCVPdata::wgsrpd3, silent = TRUE) if (!inherits(p, "try-error")) return(sf::st_as_sf(p)) diff --git a/dev_notes/01-offline-nsr-design.md b/dev_notes/01-offline-nsr-design.md index 9d2a61b..e0be331 100644 --- a/dev_notes/01-offline-nsr-design.md +++ b/dev_notes/01-offline-nsr-design.md @@ -85,11 +85,17 @@ resolution WCVP gives us) is too coarse to say anything about sub-country status | Source | What it gives | Format, size | Licence | |---|---|---|---| -| **powo** | global native/introduced by WGSRPD level 3 | WCVP v15 **already in hand** (`data/wcvp_v15/wcvp_distribution.csv`), no download | CC BY 4.0 | +| **wcvp** | global native/introduced by WGSRPD level 3 | WCVP v15 **already in hand** (`data/wcvp_v15/wcvp_distribution.csv`), no download | CC BY 4.0 | | ~~usda~~ | USA, by state | **not obtainable**: the per-state export is no longer served; the GBIF copy is names only and the API gives region-level status | CC0 | | **vascan** | Canada, by province; native/introduced/ephemeral | DwC-A, `data.canadensys.net/ipt/archive.do?r=vascan`, ~10 MB, updated 2026-08-04 | CC0 | | **flbr** | Brazil, by state; native/naturalised/cultivated | DwC-A, `ipt.jbrj.gov.br/jbrj/archive.do?r=lista_especies_flora_brasil`, ~50-100 MB | CC BY 4.0 | +**The source is named `wcvp`, not `powo` (BM, 2026-09-22):** it is the WCVP archive that is +read, and it is what `TNRS_local(sources = "wcvp")` calls the same data, so the stack uses one +name for one dataset. `"powo"` is accepted as a synonym in `sources =` and in `files =`. The +cost is that `native_status_sources` no longer matches the live service's label for this source, +which is a known diff when validating `NSR_local()` against `NSR()`. + POWO is the service's global backbone and is the same Kew dataset as WCVP, at the same resolution; we already apply its native filter in the garden pipeline (script 02 keeps `introduced == 0 & extinct == 0 & location_doubtful == 0`). The other three add state-level status in the USA, Canada diff --git a/man/NSR_local.Rd b/man/NSR_local.Rd new file mode 100644 index 0000000..b22183f --- /dev/null +++ b/man/NSR_local.Rd @@ -0,0 +1,73 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/NSR_local.R +\name{NSR_local} +\alias{NSR_local} +\title{Determine native status without an internet connection} +\usage{ +NSR_local( + occurrence_dataframe, + dir = nsr_cache_dir(), + resolve_names = TRUE, + min_overlap = 0.01, + quiet = FALSE +) +} +\arguments{ +\item{occurrence_dataframe}{A data.frame with \code{species}, and either coordinates +or political division names (see Place).} + +\item{dir}{Cache directory, shared with GNRS and GVS.} + +\item{resolve_names}{Resolve submitted names against WCVP with +\code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none.} + +\item{min_overlap}{Ignore region links covering less than this share of a region.} + +\item{quiet}{Suppress progress messages?} +} +\value{ +A data.frame, one row per input row, carrying \code{user_id} from the input + (or sequential ids where the input has none), as \code{\link{NSR}} does. +} +\description{ +Offline Native Status Resolver. For each taxon and place, the status every built +checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s columns, +with four added. +} +\details{ +\strong{Place.} Give coordinates (\code{latitude}, \code{longitude}), political +division names (\code{country}, \code{state_province}, \code{county_parish}), or both. +Coordinates are the better input: each source is consulted in the geography it +publishes against (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil +against GADM states), with no crosswalk between them. Names are resolved to GADM +units through the GNRS backbone and then carried to other geographies by the spatial +link table. + +\strong{Precedence.} If any source says native, the answer is native: asserting +nativity is a positive claim, whereas "introduced" and "absent" are often artefacts of +a list's scope, age or purpose. Disagreement is recorded, not hidden, in +\code{native_status_conflict} and \code{native_status_opinions}. + +\strong{Which polygons answer.} A place is judged by the polygons it lies IN. An +opinion about a polygon containing the place applies to it (POWO's finest statement +about Guadeloupe is "native in the Leeward Islands"), and among those any native +opinion wins. Polygons INSIDE the place describe only parts of it, so they are not its +status: they answer only when they all agree, and otherwise the answer is \code{P} with +the reason saying the status varies and how many sub-polygons say what +(\code{n_subpolygons_native}, \code{n_subpolygons_introduced}). Give coordinates and +the question does not arise: the record is judged on the ground it sits on. +\code{native_status_scope} records which of these produced the answer. + +\strong{Endemism is the exception}, deliberately: \code{Ne} and \code{Ie} are claims +about the taxon's whole range rather than about one polygon, so they draw on evidence +from elsewhere. A record of a taxon confined to California, found in Michigan, is +introduced there however little Michigan's own checklists say. +} +\section{County and parish}{ + +\code{county_parish} is accepted and echoed, but not resolved: the political-division +backbone stops at state and province, so \code{native_status_county_parish} is always +\code{NA} and the answer is given at the finest level that could be matched. A +warning says so. Coordinates are the way to ask a finer question. +} + diff --git a/man/NSR_local_build.Rd b/man/NSR_local_build.Rd new file mode 100644 index 0000000..18e9cf7 --- /dev/null +++ b/man/NSR_local_build.Rd @@ -0,0 +1,55 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/local_build.R +\name{NSR_local_build} +\alias{NSR_local_build} +\title{Build the local native-status reference} +\usage{ +NSR_local_build( + sources = "wcvp", + dir = nsr_cache_dir(create = TRUE), + files = list(), + overwrite = FALSE, + quiet = FALSE +) +} +\arguments{ +\item{sources}{Which sources to build. See \code{NSR_local_status()}.} + +\item{dir}{Cache directory, shared with GNRS and GVS.} + +\item{files}{Named list of local archives to read, e.g. +\code{list(flbr = "flbr_dwca.zip")}. Required for \code{vascan} and \code{flbr}. +For \code{wcvp}, the WCVP zip; if absent, the copy in the TNRS cache is used when +there is one.} + +\item{overwrite}{Rebuild sources that are already built?} + +\item{quiet}{Suppress progress messages?} +} +\value{ +Invisibly, \code{NSR_local_status()}. +} +\description{ +Reads each checklist source, resolves its names against WCVP with +\code{TNRS::TNRS_local()} and its political divisions against the GNRS backbone, and +writes the shared cache tables. Sources are acquired on the user's machine; nothing +derived ships with the package. +} +\details{ +\strong{The archives are not downloaded for you.} \code{vascan} and \code{flbr} must +be supplied through \code{files}, and each stops with the URL to fetch if they are +not; \code{wcvp} is the exception, since it reads the WCVP archive already in the +TNRS cache when there is one. This is why \code{sources} defaults to \code{"wcvp"} +alone: it is the only source that can build unattended. + +\code{"powo"} is accepted as a synonym for \code{"wcvp"} - it is the name the live +service reports for the same Kew dataset - but \code{wcvp} is what gets written to +\code{native_status_sources}. + +Tables written: \code{nsr-sources} (one per source, with licence and whether it is +comprehensive), \code{nsr-taxa} (one per taxon, keyed on the WCVP accepted id), +\code{nsr-regions} (one per region, in the geography its source publishes against) and +\code{nsr-checklist} (one per taxon x region x source). \code{nsr-region-links}, the +spatial relations between geographies, is built too when the GADM index from GVS is in +the cache. +} diff --git a/man/NSR_local_by_region.Rd b/man/NSR_local_by_region.Rd new file mode 100644 index 0000000..fbfc6d1 --- /dev/null +++ b/man/NSR_local_by_region.Rd @@ -0,0 +1,37 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/NSR_local_set.R +\name{NSR_local_by_region} +\alias{NSR_local_by_region} +\title{Native status for taxa and polygons you have already resolved} +\usage{ +NSR_local_by_region( + taxon_id, + region_keys, + country_keys = NULL, + dir = nsr_cache_dir(), + min_overlap = 0.01 +) +} +\arguments{ +\item{taxon_id}{WCVP accepted ids (character), one per query.} + +\item{region_keys}{The polygons each record sits in, as a list of character vectors +(\code{"wgsrpd3:BZL"}, \code{"gadm1:BRA.25_1"}), or a single character vector for +one polygon each.} + +\item{country_keys}{Optional country polygons (\code{"gadm0:BRA"}), same shape.} + +\item{dir}{Cache directory.} + +\item{min_overlap}{Ignore region links covering less than this share of a region.} +} +\value{ +A data.frame with \code{native_status} and the same companion columns + \code{NSR_local()} returns. +} +\description{ +A pipeline that has located its records itself - in the WGSRPD raster and the GADM +index, as an occurrence workflow does once for all of its coordinates - already knows +the polygons each record sits in. This takes those directly, skipping name resolution +and point-in-polygon, and answers with the same rules as \code{\link{NSR_local}}. +} diff --git a/man/NSR_local_remove.Rd b/man/NSR_local_remove.Rd new file mode 100644 index 0000000..c3adef8 --- /dev/null +++ b/man/NSR_local_remove.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/local_cache.R +\name{NSR_local_remove} +\alias{NSR_local_remove} +\title{Remove the local NSR tables} +\usage{ +NSR_local_remove(dir = nsr_cache_dir(), quiet = FALSE) +} +\arguments{ +\item{dir}{Cache directory.} + +\item{quiet}{Suppress messages?} +} +\value{ +Invisibly, the files removed. +} +\description{ +Remove the local NSR tables +} diff --git a/man/NSR_local_status.Rd b/man/NSR_local_status.Rd new file mode 100644 index 0000000..9b63763 --- /dev/null +++ b/man/NSR_local_status.Rd @@ -0,0 +1,18 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/local_cache.R +\name{NSR_local_status} +\alias{NSR_local_status} +\title{Which parts of the local NSR are built} +\usage{ +NSR_local_status(dir = nsr_cache_dir()) +} +\arguments{ +\item{dir}{Cache directory. Defaults to the one shared with GNRS and GVS.} +} +\value{ +A data.frame with one row per source: whether it is built, how many checklist + records and taxa it contributed, and when it was built. +} +\description{ +Which parts of the local NSR are built +} diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R new file mode 100644 index 0000000..69008fa --- /dev/null +++ b/tests/testthat/test-local-nsr-paths.R @@ -0,0 +1,108 @@ +# NSR_local() switches from the row resolver to the data.table set resolver at 200 rows +# (option NSR.set_min). That switch must be invisible: the same query has to come back +# the same either way. It has not been, twice - the set path never filled the per-level +# columns, and it reported a cultivated flag of 0 where the row path reports NA - so this +# compares the two directly, on tables held in memory rather than a built cache. + +# a small world: Quebec in two geographies, inside Canada, plus a sliver of a region the +# query only just touches +fixture_db <- function() { + checklist <- data.frame( + taxon_id = c("t1", "t2", "t3", "t3", "t4", "t4"), + region_key = c("wgsrpd3:QUE", "gadm0:CAN", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:QUE", "gadm1:CAN.11_1"), + status = c("native", "introduced", "native", "native", "present", "native"), + source_name = c("wcvp", "vascan", "wcvp", "wcvp", "wcvp", "vascan"), + is_cultivated = c(0L, 0L, 0L, 0L, 1L, 0L), + stringsAsFactors = FALSE) + links <- data.frame( + from_region = c("gadm1:CAN.11_1", "wgsrpd3:QUE", "gadm1:CAN.11_1"), + to_region = c("wgsrpd3:QUE", "gadm1:CAN.11_1", "wgsrpd3:SLIVER"), + relation = c("same", "same", "within"), + fraction = c(0.99, 0.99, 0.001), + stringsAsFactors = FALSE) + taxa <- data.frame( + taxon_id = c("t1", "t2", "t3", "t4"), + species_name = c("Sp one", "Sp two", "Sp three", "Sp four"), + family = "Fam", genus = "Gen", rank = "species", stringsAsFactors = FALSE) + sources <- data.frame(source_name = c("wcvp", "vascan"), is_comprehensive = TRUE, + stringsAsFactors = FALSE) + regions <- data.frame( + region_key = c("wgsrpd3:QUE", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:SLIVER", + "gadm1:CAN.11_1", "gadm0:CAN"), + region_name = c("Quebec", "Ontario", "Brazil S", "Sliver", "Quebec", "Canada"), + level = c("state_province", "state_province", "state_province", "state_province", + "state_province", "country"), + stringsAsFactors = FALSE) + regions$system <- sub(":.*", "", regions$region_key) + NSR:::nsr_index_db(list(sources = sources, taxa = taxa, checklist = checklist, + regions = regions, links = links)) +} + +# every query names Quebec, inside Canada +fixture_places <- function(n_taxa) { + list(fine = rep(list("gadm1:CAN.11_1"), n_taxa), + country = rep(list("gadm0:CAN"), n_taxa)) +} + +test_that("the row path and the set path agree, column for column", { + skip_if_not_installed("data.table") + db <- fixture_db() + tid <- c("t1", "t2", "t3", "t4", NA_character_) + places <- fixture_places(length(tid)) + + row <- NSR:::nsr_resolve_status(tid, places, db) + set <- NSR:::nsr_resolve_status_set(tid, places, db) + + for (f in c("code", "country_code", "state_code", "scope", "conflict", "conflict_type", + "cultivated", "n_sub_native", "n_sub_introduced", "sources", "reason")) { + expect_equal(set[[f]], row[[f]], info = f) + } +}) + +test_that("the set path fills the per-level codes", { + skip_if_not_installed("data.table") + db <- fixture_db() + tid <- c("t1", "t2") + set <- NSR:::nsr_resolve_status_set(tid, fixture_places(2), db) + # native in Quebec via the WGSRPD twin; nothing said about Canada as a whole + expect_equal(set$state_code[1], "N") + expect_false(is.na(set$country_code[1])) + # introduced at country level propagates down to the state answer + expect_equal(set$country_code[2], "I") +}) + +test_that("a level that was not asked about has no code", { + skip_if_not_installed("data.table") + db <- fixture_db() + places <- list(fine = list(character(0)), country = list("gadm0:CAN")) + for (res in list(NSR:::nsr_resolve_status("t2", places, db), + NSR:::nsr_resolve_status_set("t2", places, db))) { + expect_true(is.na(res$state_code[1])) + expect_equal(res$country_code[1], "I") + } +}) + +test_that("absence leaves the cultivated flag unknown in both paths", { + skip_if_not_installed("data.table") + db <- fixture_db() + # t3 is evaluable (it is native in Ontario and in Brazil, so its range is not confined + # and the Ie rule does not fire) and Quebec is comprehensively listed, so a missing + # opinion here is a real absence rather than a gap + places <- fixture_places(1) + row <- NSR:::nsr_resolve_status("t3", places, db) + set <- NSR:::nsr_resolve_status_set("t3", places, db) + expect_equal(row$code, "A") + expect_equal(set$code, "A") + expect_true(is.na(row$cultivated)) + expect_true(is.na(set$cultivated)) +}) + +test_that("min_overlap keeps slivers out of the consulted set", { + db <- fixture_db() + keys <- "gadm1:CAN.11_1" + # the default threshold drops the 0.1% link but keeps the twin + expect_true("wgsrpd3:QUE" %in% NSR:::nsr_consulted_regions(keys, db)) + expect_false("wgsrpd3:SLIVER" %in% NSR:::nsr_consulted_regions(keys, db)) + # lower the threshold and it comes back, so it is the threshold doing the work + expect_true("wgsrpd3:SLIVER" %in% NSR:::nsr_consulted_regions(keys, db, min_overlap = 0)) +}) From 6f30047c083ee756386ca6155280c6f3091a9a82 Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Tue, 22 Sep 2026 16:38:22 -0400 Subject: [PATCH 3/8] Second review pass: link guard, GADM hierarchy, and build error handling The default build could not answer a single query by name. Names always resolve to GADM through the GNRS backbone, whatever geography a source publishes against, so a WGSRPD source is unreachable by name until the two are linked - but link construction was gated on length(systems) > 1, which a wcvp-only build never satisfies. Making wcvp the default in the previous commit turned that into the common case. The guard now asks whether linking is needed at all, not how many systems are present. Rasterization and link failures were caught by try() and then ignored, so a cache could be reported as built while being unable to place a coordinate or reach a WCVP opinion. Both now test for try-error and warn with the underlying message and the specific consequence, and the link step is skipped when the raster failed rather than failing again on the same cause. GADM parent/child containment is now resolved directly (BM's call). gid_1 carries its country in the key ("BRA.25_1" sits in "BRA"), so the relation is exact and is deliberately kept out of the link table: link rows carry overlap fractions and are filtered by min_overlap, and a state is a legitimately tiny share of its country - Distrito Federal is 0.07% of Brazil and would have been discarded as a sliver. Without this a country query could not see checklist rows published against that country's states, which is how VASCAN and Flora do Brasil publish, so native-up propagation, country coverage and confined-country Ie all failed at country level. Two things the row/set agreement test caught while wiring that up: - Only the downward direction belongs here. A first pass added parent edges too, which made a bare region query silently inherit a country it was never given; upward inheritance was already handled, where the caller asks for it, by the ancestor path and by country_keys. - nsr_confined_ranges() ran before the hierarchy was derived, so confined-country Ie quietly did nothing. The derivation moved above the data.table block. NSR_local_by_region() returns the per-level codes and native_status_conflict its help advertised, and the help now says plainly that it does not echo the input columns NSR_local() carries through, since it is given ids and polygon keys rather than a record. Three tests cover the hierarchy: a country query seeing the states inside it, min_overlap leaving the parent/child edges alone, and a range confined to one country's states resolving to that country. Not changed, and why: - "Regenerate NAMESPACE and exported documentation" is stale. All five exports and all five Rd files are on the branch already. - terra::extract()[, 1] is correct. extract() adds an ID column for a SpatVector, not for the coordinate matrix passed here. - Doubtful WCVP records stay as "present" rather than being dropped. Dropping them would turn a source's explicit "maybe here" into a confident absence. Extinct records are filtered, as of the last commit. Still open: coordinates and division names are unioned when both are given, so a record whose point and names disagree is answered from both. The service has no combined mode to follow - NSR() takes names only and NSR_from_coordinates() takes coordinates only - so this needs a decision rather than a precedent. Co-Authored-By: Claude Opus 5 --- R/NSR_local.R | 64 +++++++++++++++++++++++++-- R/NSR_local_set.R | 32 ++++++++++++-- R/local_build.R | 45 +++++++++++++++---- man/NSR_local_by_region.Rd | 7 ++- tests/testthat/test-local-nsr-paths.R | 63 ++++++++++++++++++++++---- 5 files changed, 185 insertions(+), 26 deletions(-) diff --git a/R/NSR_local.R b/R/NSR_local.R index 54959f5..fe19f3e 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -180,7 +180,26 @@ nsr_index_db <- function(db) { paste(db$checklist$taxon_id, db$checklist$region_key)) db$chk_env <- list2env(db$chk_idx, envir = new.env(hash = TRUE, parent = emptyenv())) db$link_env <- list2env(db$link_idx, envir = new.env(hash = TRUE, parent = emptyenv())) + # GADM's key carries its own hierarchy: gid_1 "BRA.25_1" sits in gid_0 "BRA". That + # containment is exact by construction, which is why it is kept OUT of the link table: + # link rows carry overlap fractions and are filtered by min_overlap, and a state is a + # legitimately tiny share of its country (Distrito Federal is 0.07% of Brazil), so as + # a fraction it would be discarded as a sliver. Only regions that actually carry + # checklist rows are listed, since those are the only ones an answer can draw on. + g1 <- unique(db$regions$region_key[db$regions$system == "gadm1"]) + g1 <- g1[g1 %in% db$checklist$region_key] + parent <- if (length(g1)) paste0("gadm0:", sub("\\..*$", "", sub("^gadm1:", "", g1))) else character(0) + db$gadm_parent <- stats::setNames(parent, g1) + db$gadm_children <- if (length(g1)) split(g1, parent) else list() if (requireNamespace("data.table", quietly = TRUE)) { + db$gadm_edges_dt <- if (length(g1)) { + data.table::data.table(from_region = c(g1, parent), to_region = c(parent, g1), + relation = rep(c("within", "contains"), each = length(g1))) + } else { + data.table::data.table(from_region = character(0), to_region = character(0), + relation = character(0)) + } + data.table::setkeyv(db$gadm_edges_dt, "from_region") db$chk_dt <- data.table::as.data.table(db$checklist) data.table::setkeyv(db$chk_dt, c("taxon_id", "region_key")) db$links_dt <- data.table::as.data.table(db$links) @@ -332,6 +351,29 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { out } +#' The GADM units directly below, or directly above, the regions given +#' +#' Internal. Exact containment read off the key, so no overlap fraction applies and +#' \code{min_overlap} does not filter it. See \code{nsr_index_db()}. +#' @keywords internal +#' @noRd +nsr_gadm_children <- function(keys, db) { + if (!length(keys) || !length(db$gadm_children)) return(character(0)) + k <- keys[startsWith(keys, "gadm0:")] + if (!length(k)) return(character(0)) + setdiff(as.character(unlist(db$gadm_children[k], use.names = FALSE)), keys) +} + +#' @keywords internal +#' @noRd +nsr_gadm_parents <- function(keys, db) { + if (!length(keys) || !length(db$gadm_parent)) return(character(0)) + k <- keys[startsWith(keys, "gadm1:")] + if (!length(k)) return(character(0)) + p <- as.character(unname(db$gadm_parent[k])) + setdiff(p[!is.na(p)], keys) +} + #' Opinions about one taxon bearing on a set of regions #' #' Internal. Direct opinions, plus opinions about regions the query's regions are inside @@ -358,6 +400,19 @@ nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01) { } } } + # The states inside a queried country: a source keyed on GADM states (VASCAN, Flora do + # Brasil) is otherwise invisible to a country query, which breaks native-up propagation. + # A state query does not come back here for its country - that is the ancestor path in + # nsr_resolve_status(), which takes direct rows only, so a country's OTHER states never + # reach it. + kids <- nsr_gadm_children(keys, db) + if (length(kids)) { + idx <- mget(paste(taxon_id, kids), db$chk_env, ifnotfound = list(NULL)) + if (sum(lengths(idx))) { + rows <- c(rows, unlist(idx, use.names = FALSE)) + rel <- c(rel, rep("contains", sum(lengths(idx)))) + } + } if (!length(rows)) return(NULL) list(status = db$chk_status[rows], source_name = db$chk_source[rows], is_cultivated = db$chk_cult[rows], relation = rel) @@ -378,6 +433,7 @@ nsr_is_endemic <- function(taxon_id, keys, db, min_overlap = 0.01) { li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) inside <- if (!length(li)) keys else c(keys, db$link_to[li][db$link_rel[li] %in% c("same", "contains") & db$link_frac[li] >= min_overlap]) + inside <- c(inside, nsr_gadm_children(keys, db)) all(nat %in% inside) } @@ -404,10 +460,10 @@ nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { if (length(nat) == 1) return(nm(nat)) containers <- lapply(nat, function(r) { li <- db$link_env[[r]] - if (is.null(li)) return(character(0)) + if (is.null(li)) return(nsr_gadm_parents(r, db)) keep <- db$link_rel[li] %in% c("same", "within") & db$link_frac[li] >= min_overlap & startsWith(db$link_to[li], "gadm0:") - db$link_to[li][keep] + c(db$link_to[li][keep], nsr_gadm_parents(r, db)) }) common <- Reduce(intersect, containers) if (length(common)) nm(common[1]) else NA_character_ @@ -428,9 +484,9 @@ nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { nsr_consulted_regions <- function(keys, db, min_overlap = 0.01) { if (!length(keys)) return(character(0)) li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) - if (!length(li)) return(keys) + if (!length(li)) return(unique(c(keys, nsr_gadm_children(keys, db)))) keep <- db$link_rel[li] %in% c("same", "within", "contains") & db$link_frac[li] >= min_overlap - unique(c(keys, db$link_to[li][keep])) + unique(c(keys, db$link_to[li][keep], nsr_gadm_children(keys, db))) } #' What kind of disagreement is this? diff --git a/R/NSR_local_set.R b/R/NSR_local_set.R index 4458670..2792adc 100644 --- a/R/NSR_local_set.R +++ b/R/NSR_local_set.R @@ -30,6 +30,13 @@ nsr_confined_ranges <- function(nat, db) { lk <- db$links_dt[relation %in% c("same", "within") & startsWith(to_region, "gadm0:"), list(region_key = from_region, country = to_region)] + if (!is.null(db$gadm_edges_dt) && nrow(db$gadm_edges_dt)) { + lk <- unique(data.table::rbindlist(list( + lk, + db$gadm_edges_dt[relation == "within", + list(region_key = from_region, country = to_region)]), + use.names = TRUE)) + } m <- merge(nat[taxon_id %in% multi$taxon_id], lk, by = "region_key", allow.cartesian = TRUE) per <- m[, list(n_in = data.table::uniqueN(region_key)), by = c("taxon_id", "country")] per <- merge(per, multi, by = "taxon_id") @@ -56,8 +63,11 @@ nsr_confined_ranges <- function(nat, db) { #' @param country_keys Optional country polygons (\code{"gadm0:BRA"}), same shape. #' @param dir Cache directory. #' @param min_overlap Ignore region links covering less than this share of a region. -#' @return A data.frame with \code{native_status} and the same companion columns -#' \code{NSR_local()} returns. +#' @return A data.frame with \code{native_status} and the status columns +#' \code{NSR_local()} returns, keyed on \code{taxon_id}. It does not echo the input +#' columns \code{NSR_local()} carries through (\code{species}, the division names, +#' \code{user_id}), since this entry point is given ids and polygon keys rather than +#' a record. #' @export NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, dir = nsr_cache_dir(), min_overlap = 0.01) { @@ -72,9 +82,13 @@ NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, res <- if (!is.null(db$chk_dt)) nsr_resolve_status_set(as.character(taxon_id), places, db, min_overlap) else nsr_resolve_status(as.character(taxon_id), places, db, min_overlap) data.frame(taxon_id = as.character(taxon_id), + native_status_country = res$country_code, + native_status_state_province = res$state_code, native_status = res$code, native_status_reason = res$reason, native_status_sources = res$sources, native_status_opinions = res$opinions, - native_status_scope = res$scope, native_status_conflict_type = res$conflict_type, + native_status_scope = res$scope, + native_status_conflict = res$conflict, + native_status_conflict_type = res$conflict_type, isIntroduced = as.integer(res$code %in% c("I", "Ie")), isEndemic = as.integer(res$code == "Ne"), isCultivatedNSR = res$cultivated, @@ -159,6 +173,18 @@ nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, .(qid, region_key = to_region, relation = relation, taxon_id)] Q <- data.table::rbindlist(list(Q, L), use.names = TRUE) } + # GADM's exact parent/child containment, joined on the query's OWN polygons only so a + # country's other states never reach a state-level query. Deliberately after the + # min_overlap filter: these edges carry no fraction because the containment is exact. + if (!is.null(db$gadm_edges_dt) && nrow(db$gadm_edges_dt)) { + H <- db$gadm_edges_dt[relation == "contains"][own, on = c(from_region = "region_key"), + nomatch = 0L, allow.cartesian = TRUE] + if (nrow(H)) { + Q <- data.table::rbindlist( + list(Q, H[to_region != from_region, .(qid, region_key = to_region, relation, taxon_id)]), + use.names = TRUE) + } + } # keep the strongest relation per (query, polygon) ord <- c(same = 1L, within = 2L, contains = 3L, overlaps = 4L) Q[, r_ord := ord[relation]] diff --git a/R/local_build.R b/R/local_build.R index fc18c9f..efbd99c 100644 --- a/R/local_build.R +++ b/R/local_build.R @@ -148,23 +148,50 @@ NSR_local_build <- function(sources = "wcvp", n_records = nrow(checklist), n_taxa = nrow(taxa), n_regions = nrow(regions)), nsr_provenance_path("nsr", dir)) + # ---- geography --------------------------------------------------------------------- # The WGSRPD raster is what places coordinates in a WCVP region, so it is needed # whenever a source publishes against WGSRPD - not only when there is a second # geography to link it to. A wcvp-only build needs it just as much. systems <- unique(regions$system) + raster_ok <- TRUE if ("wgsrpd3" %in% systems && !file.exists(nsr_raster_path("wgsrpd3", dir))) { msg("Building the WGSRPD level-3 raster ...") - try(nsr_build_wgsrpd_raster(dir = dir, quiet = quiet), silent = FALSE) + r <- try(nsr_build_wgsrpd_raster(dir = dir, quiet = quiet), silent = TRUE) + raster_ok <- !inherits(r, "try-error") + if (!raster_ok) { + warning("The WGSRPD level-3 raster could not be built: ", + conditionMessage(attr(r, "condition")), " +", + "The checklist tables are written, but until this raster exists no query ", + "can reach a WCVP opinion: coordinates cannot be placed in a WCVP region, ", + "and division names resolve to GADM, which needs the link table below. ", + "Install 'rWCVPdata' (or 'rWCVP'), then re-run NSR_local_build(overwrite = TRUE).", + call. = FALSE) + } } - # spatial links between geographies, where GVS's GADM index is available - if (length(systems) > 1 && file.exists(file.path(dir, "gadmindex-units-30s.tif"))) { - msg("Linking region systems (", paste(systems, collapse = ", "), ") ...") - try(nsr_build_region_links(dir = dir, quiet = quiet), silent = FALSE) - } else if (length(systems) > 1) { - warning("Sources use several geographies but GVS's GADM index is not in the cache, ", - "so they cannot be linked; queries by division name will only see their own system.", - call. = FALSE) + # Names always resolve to GADM through the GNRS backbone, whatever geography the + # sources use, so a WGSRPD source is unreachable by name until the two are linked. + # That holds for a wcvp-only build as much as a mixed one - hence no test on the + # number of systems, which is what used to leave the default build unanswerable. + needs_links <- "wgsrpd3" %in% systems || length(systems) > 1 + gadm_index <- file.path(dir, "gadmindex-units-30s.tif") + if (needs_links && raster_ok && file.exists(gadm_index)) { + msg("Linking region systems (", paste(c(systems, "gadm"), collapse = ", "), ") ...") + r <- try(nsr_build_region_links(dir = dir, quiet = quiet), silent = TRUE) + if (inherits(r, "try-error")) { + warning("The region link table could not be built: ", + conditionMessage(attr(r, "condition")), " +", + "The cache is written, but queries by division name will not reach sources ", + "published against another geography. Re-run NSR_local_build(overwrite = TRUE) ", + "once the cause is fixed.", call. = FALSE) + } + } else if (needs_links && raster_ok) { + warning("GVS's GADM index is not in ", dir, ", so the geographies cannot be linked. ", + "Queries by division name resolve to GADM and will not reach sources published ", + "against WGSRPD; coordinates still work. Build it with GVS, then re-run ", + "NSR_local_build(overwrite = TRUE).", call. = FALSE) } msg("Built: ", format(nrow(checklist), big.mark = ","), " checklist records, ", format(nrow(taxa), big.mark = ","), " taxa, ", format(nrow(regions), big.mark = ","), diff --git a/man/NSR_local_by_region.Rd b/man/NSR_local_by_region.Rd index fbfc6d1..806681a 100644 --- a/man/NSR_local_by_region.Rd +++ b/man/NSR_local_by_region.Rd @@ -26,8 +26,11 @@ one polygon each.} \item{min_overlap}{Ignore region links covering less than this share of a region.} } \value{ -A data.frame with \code{native_status} and the same companion columns - \code{NSR_local()} returns. +A data.frame with \code{native_status} and the status columns + \code{NSR_local()} returns, keyed on \code{taxon_id}. It does not echo the input + columns \code{NSR_local()} carries through (\code{species}, the division names, + \code{user_id}), since this entry point is given ids and polygon keys rather than + a record. } \description{ A pipeline that has located its records itself - in the WGSRPD raster and the GADM diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R index 69008fa..7d7dea9 100644 --- a/tests/testthat/test-local-nsr-paths.R +++ b/tests/testthat/test-local-nsr-paths.R @@ -8,11 +8,11 @@ # query only just touches fixture_db <- function() { checklist <- data.frame( - taxon_id = c("t1", "t2", "t3", "t3", "t4", "t4"), - region_key = c("wgsrpd3:QUE", "gadm0:CAN", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:QUE", "gadm1:CAN.11_1"), - status = c("native", "introduced", "native", "native", "present", "native"), - source_name = c("wcvp", "vascan", "wcvp", "wcvp", "wcvp", "vascan"), - is_cultivated = c(0L, 0L, 0L, 0L, 1L, 0L), + taxon_id = c("t1", "t2", "t3", "t3", "t4", "t4", "t5", "t5"), + region_key = c("wgsrpd3:QUE", "gadm0:CAN", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:QUE", "gadm1:CAN.11_1", "gadm1:CAN.11_1", "gadm1:CAN.4_1"), + status = c("native", "introduced", "native", "native", "present", "native", "native", "native"), + source_name = c("wcvp", "vascan", "wcvp", "wcvp", "wcvp", "vascan", "vascan", "vascan"), + is_cultivated = c(0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L), stringsAsFactors = FALSE) links <- data.frame( from_region = c("gadm1:CAN.11_1", "wgsrpd3:QUE", "gadm1:CAN.11_1"), @@ -28,10 +28,11 @@ fixture_db <- function() { stringsAsFactors = FALSE) regions <- data.frame( region_key = c("wgsrpd3:QUE", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:SLIVER", - "gadm1:CAN.11_1", "gadm0:CAN"), - region_name = c("Quebec", "Ontario", "Brazil S", "Sliver", "Quebec", "Canada"), + "gadm1:CAN.11_1", "gadm1:CAN.4_1", "gadm0:CAN"), + region_name = c("Quebec", "Ontario", "Brazil S", "Sliver", "Quebec", + "New Brunswick", "Canada"), level = c("state_province", "state_province", "state_province", "state_province", - "state_province", "country"), + "state_province", "state_province", "country"), stringsAsFactors = FALSE) regions$system <- sub(":.*", "", regions$region_key) NSR:::nsr_index_db(list(sources = sources, taxa = taxa, checklist = checklist, @@ -97,6 +98,52 @@ test_that("absence leaves the cultivated flag unknown in both paths", { expect_true(is.na(set$cultivated)) }) +# GADM keys carry their own hierarchy (BRA.25_1 sits in BRA). Without it, a country +# query cannot see checklist rows published against that country's states, which is how +# VASCAN and Flora do Brasil publish - so native-up propagation and confined-country Ie +# both fail at country level. +test_that("a country query sees the states inside it", { + skip_if_not_installed("data.table") + db <- fixture_db() + places <- list(fine = list(character(0)), country = list("gadm0:CAN")) + row <- NSR:::nsr_resolve_status("t4", places, db) + set <- NSR:::nsr_resolve_status_set("t4", places, db) + # t4's only Canadian opinion is VASCAN's, recorded against Quebec the GADM state. + # Without the hierarchy the country query cannot see it at all and answers A; with it + # the state answers, and since Quebec is t4's only native region it is endemic there. + expect_equal(row$code, "Ne") + expect_equal(set$code, "Ne") + expect_equal(row$scope, "sub-polygons agree") + expect_equal(set$scope, "sub-polygons agree") + expect_equal(set$sources, row$sources) +}) + +test_that("the hierarchy is exact, so min_overlap cannot filter it out", { + db <- fixture_db() + kids <- NSR:::nsr_gadm_children("gadm0:CAN", db) + expect_setequal(kids, c("gadm1:CAN.11_1", "gadm1:CAN.4_1")) + expect_equal(NSR:::nsr_gadm_parents("gadm1:CAN.4_1", db), "gadm0:CAN") + # a state is a tiny share of its country, so as an overlap fraction it would be + # dropped; the parent/child edges carry no fraction and are not filtered + expect_setequal(NSR:::nsr_gadm_children("gadm0:CAN", db), kids) + # regions with no checklist rows are not children of anything + expect_equal(NSR:::nsr_gadm_children("gadm0:BRA", db), character(0)) +}) + +test_that("a range confined to one country's states is confined to that country", { + skip_if_not_installed("data.table") + db <- fixture_db() + # t5 is native in Quebec and New Brunswick and nowhere else; asked about a Brazilian + # region that is comprehensively listed, it is absent there and endemic to Canada + places <- list(fine = list("wgsrpd3:BZL"), country = list(character(0))) + row <- NSR:::nsr_resolve_status("t5", places, db) + set <- NSR:::nsr_resolve_status_set("t5", places, db) + expect_equal(row$code, "Ie") + expect_equal(set$code, "Ie") + expect_match(row$reason, "Canada") + expect_match(set$reason, "Canada") +}) + test_that("min_overlap keeps slivers out of the consulted set", { db <- fixture_db() keys <- "gadm1:CAN.11_1" From 8dd05625658f6b67f598ccdd929733249bd250d2 Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Tue, 22 Sep 2026 16:55:48 -0400 Subject: [PATCH 4/8] Decide extinct records at query time, not at build time WCVP's `extinct` is a bare 0/1 with no date: checked against the v15 distribution table, whose only columns are plant_locality_id, plant_name_id, continent_code_l1, continent, region_code_l2, region, area_code_l3, area, introduced, extinct and location_doubtful. It means "considered no longer present as of this release" and cannot be compared against an occurrence's own date. Whether such a record counts is therefore a property of the question rather than of the data, so the previous commit's build-time filter was in the wrong place. The build now keeps the row with an is_extinct flag and NSR_local(exclude_extinct = TRUE) decides per call (BM): the default answers about the present day, FALSE models a past distribution. NSR_local_by_region() takes the same argument. 2,701 of 1,970,252 rows are affected, over 2,393 taxa, 906 of which have no surviving native record at all. Dropping the rows at build time also destroyed evidence the endemism rules need. A region a taxon has been lost from is still a region it was native to, so the Ie rule must not fire against it: the answer for a natively extirpated record is A (gone), never I or Ie (arrived). With the rows filtered out of the cache entirely, a taxon whose surviving range happened to be confined would have come back "introduced here (inferred)" for a record that was natively there - a worse answer than the N it replaced, and one this filter introduced. The endemism rules now read the taxon's whole native range, extinct records included, whatever exclude_extinct says; only the status opinions are filtered. The set path filters after its join rather than holding a second copy of the checklist, since the rows removed are 0.14% of it. Caches built before this default the column to 0, which is true of every row they hold. Two tests: that the switch moves the answer between A and N for a lost population, and that t7 (never present) resolves Ie while t6 (present and lost) resolves A - so the inference machinery is demonstrably live in the case it must not fire for. R CMD check: Status OK. Co-Authored-By: Claude Opus 5 --- NEWS | 5 +++ R/NSR_local.R | 46 +++++++++++++++------ R/NSR_local_set.R | 26 ++++++++---- R/local_build.R | 4 +- R/local_globals.R | 2 +- R/local_import.R | 22 +++++++---- dev_notes/01-offline-nsr-design.md | 17 ++++++++ man/NSR_local.Rd | 21 ++++++++++ man/NSR_local_by_region.Rd | 6 ++- tests/testthat/test-local-nsr-paths.R | 57 +++++++++++++++++++++++---- 10 files changed, 169 insertions(+), 37 deletions(-) diff --git a/NEWS b/NEWS index 920d3f8..0d947d1 100644 --- a/NEWS +++ b/NEWS @@ -9,6 +9,11 @@ NSR_local_by_region() alongside it. The checklist sources are acquired on the user's machine; nothing derived ships with the package. The offline path needs data.table, nanoparquet, terra, sf and TNRS, all of which are Suggests. +* NSR_local() gains exclude_extinct (default TRUE). WCVP flags a distribution extinct + without a date, so the cache keeps those records and the query decides: the default + answers about the present day, FALSE models a past distribution. Either way a region + a taxon has been lost from stays part of its native range, so such a record is + reported absent rather than introduced. ## Bug fixes diff --git a/R/NSR_local.R b/R/NSR_local.R index fe19f3e..42ae764 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -37,9 +37,25 @@ #' @param resolve_names Resolve submitted names against WCVP with #' \code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none. #' @param min_overlap Ignore region links covering less than this share of a region. +#' @param exclude_extinct Ignore distribution records the source marks extinct? Default +#' \code{TRUE}, which answers about the present day. Set \code{FALSE} to model a past +#' distribution, where a region the taxon has since been lost from still counts. See +#' Extinct records. #' @param quiet Suppress progress messages? #' @return A data.frame, one row per input row, carrying \code{user_id} from the input #' (or sequential ids where the input has none), as \code{\link{NSR}} does. +#' @section Extinct records: +#' WCVP marks a distribution extinct with a bare flag and no date, so the record says +#' only "considered no longer present as of this release". Whether it should count is +#' therefore a property of the question rather than of the data: a present-day status +#' wants it out, a historical distribution wants it in, and the cache keeps both so +#' \code{exclude_extinct} can decide per call. +#' +#' Excluding it never makes the taxon introduced there. The endemism rules read the +#' taxon's whole native range, extinct records included, whatever this argument says - +#' a region a taxon has been lost from is still a region it was native to, so \code{Ie} +#' ("introduced, inferred from endemism elsewhere") cannot fire against it. The answer +#' for such a record is \code{A}: gone, not foreign. #' @section County and parish: #' \code{county_parish} is accepted and echoed, but not resolved: the political-division #' backbone stops at state and province, so \code{native_status_county_parish} is always @@ -47,7 +63,7 @@ #' warning says so. Coordinates are the way to ask a finer question. #' @export NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names = TRUE, - min_overlap = 0.01, quiet = FALSE) { + min_overlap = 0.01, exclude_extinct = TRUE, quiet = FALSE) { if (!inherits(occurrence_dataframe, "data.frame")) { stop("occurrence_dataframe should be a data.frame", call. = FALSE) } @@ -99,8 +115,8 @@ NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names # ---- opinions --------------------------------------------------------------------- use_set <- !is.null(db$chk_dt) && nrow(occurrence_dataframe) >= getOption("NSR.set_min", 200) - res <- if (use_set) nsr_resolve_status_set(taxon_id, places, db, min_overlap) else - nsr_resolve_status(taxon_id, places, db, min_overlap) + res <- if (use_set) nsr_resolve_status_set(taxon_id, places, db, min_overlap, exclude_extinct) + else nsr_resolve_status(taxon_id, places, db, min_overlap, exclude_extinct) out <- data.frame( family = db$taxa$family[match(taxon_id, db$taxa$taxon_id)], @@ -175,6 +191,9 @@ nsr_local_db <- function(dir) { #' @keywords internal #' @noRd nsr_index_db <- function(db) { + # Caches built before extinct records were retained have no such column; they simply + # hold no extinct rows, so 0 is the truth for every row they do hold. + if (is.null(db$checklist$is_extinct)) db$checklist$is_extinct <- 0L db$link_idx <- split(seq_len(nrow(db$links)), db$links$from_region) db$chk_idx <- split(seq_len(nrow(db$checklist)), paste(db$checklist$taxon_id, db$checklist$region_key)) @@ -220,6 +239,7 @@ nsr_index_db <- function(db) { db$chk_status <- db$checklist$status db$chk_source <- db$checklist$source_name db$chk_cult <- db$checklist$is_cultivated + db$chk_extinct <- as.integer(db$checklist$is_extinct) %in% 1L db$chk_region <- db$checklist$region_key db$link_to <- db$links$to_region db$link_rel <- db$links$relation @@ -284,7 +304,8 @@ nsr_query_regions <- function(lon, lat, country, state, county, dir, db) { #' Gather and reduce every opinion bearing on each row #' @keywords internal #' @noRd -nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { +nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01, + exclude_extinct = TRUE) { n <- length(taxon_id) blank <- rep(NA_character_, n) out <- list(code = blank, reason = blank, sources = blank, opinions = blank, @@ -309,13 +330,14 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { out$n_sub_introduced[i] <- 0L next } - op <- nsr_opinions_for(tid, keys, db, min_overlap) + op <- nsr_opinions_for(tid, keys, db, min_overlap, exclude_extinct) # a coarser polygon the place sits in (its country) also contains it, so opinions # recorded ON that polygon apply; its OTHER sub-polygons do not, so only direct rows # are taken, never its links anc <- setdiff(ctry, keys) if (length(anc)) { rows <- unlist(mget(paste(tid, anc), db$chk_env, ifnotfound = list(NULL)), use.names = FALSE) + if (exclude_extinct && length(rows)) rows <- rows[!db$chk_extinct[rows]] if (length(rows)) { add <- list(status = db$chk_status[rows], source_name = db$chk_source[rows], is_cultivated = db$chk_cult[rows], relation = rep("within", length(rows))) @@ -340,11 +362,11 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01) { for (f in names(r)) out[[f]][i] <- r[[f]] # per-level codes, for the service's columns out$country_code[i] <- if (length(ctry)) { - nsr_reduce(nsr_opinions_for(tid, ctry, db, min_overlap), + nsr_reduce(nsr_opinions_for(tid, ctry, db, min_overlap, exclude_extinct), nsr_consulted_regions(ctry, db, min_overlap), ev, db)$code } else NA_character_ out$state_code[i] <- if (length(fine)) { - nsr_reduce(nsr_opinions_for(tid, fine, db, min_overlap), + nsr_reduce(nsr_opinions_for(tid, fine, db, min_overlap, exclude_extinct), nsr_consulted_regions(fine, db, min_overlap), ev, db)$code } else NA_character_ } @@ -381,9 +403,11 @@ nsr_gadm_parents <- function(keys, db) { #' inheritance rules. #' @keywords internal #' @noRd -nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01) { +nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01, exclude_extinct = TRUE) { if (!length(keys)) return(NULL) - rows <- unlist(mget(paste(taxon_id, keys), db$chk_env, ifnotfound = list(NULL)), use.names = FALSE) + live <- function(r) if (exclude_extinct && length(r)) r[!db$chk_extinct[r]] else r + rows <- live(unlist(mget(paste(taxon_id, keys), db$chk_env, ifnotfound = list(NULL)), + use.names = FALSE)) rel <- if (length(rows)) rep("same", length(rows)) else character(0) li <- unlist(mget(keys, db$link_env, ifnotfound = list(NULL)), use.names = FALSE) if (length(li)) { @@ -392,7 +416,7 @@ nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01) { good <- !(to %in% keys) & fr >= min_overlap & !(sub(":.*", "", to) %in% known) if (any(good)) { to <- to[good]; rl <- rl[good] - idx <- mget(paste(taxon_id, to), db$chk_env, ifnotfound = list(NULL)) + idx <- lapply(mget(paste(taxon_id, to), db$chk_env, ifnotfound = list(NULL)), live) n <- lengths(idx) if (sum(n)) { rows <- c(rows, unlist(idx, use.names = FALSE)) @@ -407,7 +431,7 @@ nsr_opinions_for <- function(taxon_id, keys, db, min_overlap = 0.01) { # reach it. kids <- nsr_gadm_children(keys, db) if (length(kids)) { - idx <- mget(paste(taxon_id, kids), db$chk_env, ifnotfound = list(NULL)) + idx <- lapply(mget(paste(taxon_id, kids), db$chk_env, ifnotfound = list(NULL)), live) if (sum(lengths(idx))) { rows <- c(rows, unlist(idx, use.names = FALSE)) rel <- c(rel, rep("contains", sum(lengths(idx)))) diff --git a/R/NSR_local_set.R b/R/NSR_local_set.R index 2792adc..8e9f059 100644 --- a/R/NSR_local_set.R +++ b/R/NSR_local_set.R @@ -63,6 +63,8 @@ nsr_confined_ranges <- function(nat, db) { #' @param country_keys Optional country polygons (\code{"gadm0:BRA"}), same shape. #' @param dir Cache directory. #' @param min_overlap Ignore region links covering less than this share of a region. +#' @param exclude_extinct Ignore distribution records the source marks extinct? See +#' \code{\link{NSR_local}}. #' @return A data.frame with \code{native_status} and the status columns #' \code{NSR_local()} returns, keyed on \code{taxon_id}. It does not echo the input #' columns \code{NSR_local()} carries through (\code{species}, the division names, @@ -70,7 +72,8 @@ nsr_confined_ranges <- function(nat, db) { #' a record. #' @export NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, - dir = nsr_cache_dir(), min_overlap = 0.01) { + dir = nsr_cache_dir(), min_overlap = 0.01, + exclude_extinct = TRUE) { if (!is.list(region_keys)) region_keys <- as.list(region_keys) if (is.null(country_keys)) country_keys <- vector("list", length(taxon_id)) if (!is.list(country_keys)) country_keys <- as.list(country_keys) @@ -79,8 +82,9 @@ NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, db <- nsr_local_db(dir) clean <- function(z) { z <- z[!is.na(z) & nzchar(z)]; if (!length(z)) character(0) else z } places <- list(fine = lapply(region_keys, clean), country = lapply(country_keys, clean)) - res <- if (!is.null(db$chk_dt)) nsr_resolve_status_set(as.character(taxon_id), places, db, min_overlap) - else nsr_resolve_status(as.character(taxon_id), places, db, min_overlap) + res <- if (!is.null(db$chk_dt)) + nsr_resolve_status_set(as.character(taxon_id), places, db, min_overlap, exclude_extinct) + else nsr_resolve_status(as.character(taxon_id), places, db, min_overlap, exclude_extinct) data.frame(taxon_id = as.character(taxon_id), native_status_country = res$country_code, native_status_state_province = res$state_code, @@ -108,9 +112,11 @@ NSR_local_by_region <- function(taxon_id, region_keys, country_keys = NULL, #' implementations agree column for column. #' @keywords internal #' @noRd -nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { +nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01, + exclude_extinct = TRUE) { n <- length(taxon_id) - out <- nsr_resolve_status_set1(taxon_id, places, db, min_overlap, endemism = TRUE) + out <- nsr_resolve_status_set1(taxon_id, places, db, min_overlap, endemism = TRUE, + exclude_extinct = exclude_extinct) empty <- rep(list(character(0)), n) # a level reports a code only when it was asked about AND there is a taxon to ask # about: with no match in the backbone no lookup happened at any level, which is what @@ -120,12 +126,14 @@ nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { nf <- lengths(places$fine) if (any(nc > 0)) { lv <- nsr_resolve_status_set1(taxon_id, list(fine = empty, country = places$country), - db, min_overlap, endemism = FALSE) + db, min_overlap, endemism = FALSE, + exclude_extinct = exclude_extinct) out$country_code <- ifelse(asked & nc > 0, lv$code, NA_character_) } if (any(nf > 0)) { lv <- nsr_resolve_status_set1(taxon_id, list(fine = places$fine, country = empty), - db, min_overlap, endemism = FALSE) + db, min_overlap, endemism = FALSE, + exclude_extinct = exclude_extinct) out$state_code <- ifelse(asked & nf > 0, lv$code, NA_character_) } out @@ -135,7 +143,7 @@ nsr_resolve_status_set <- function(taxon_id, places, db, min_overlap = 0.01) { #' @keywords internal #' @noRd nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, - endemism = TRUE) { + endemism = TRUE, exclude_extinct = TRUE) { dt <- function(...) data.table::data.table(...) n <- length(taxon_id) blank <- rep(NA_character_, n) @@ -193,6 +201,8 @@ nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, # the opinions themselves O <- db$chk_dt[Q, on = c("taxon_id", "region_key"), nomatch = 0L, allow.cartesian = TRUE] + # after the join, so the checklist is never copied just to drop 0.14% of its rows + if (exclude_extinct && nrow(O) && "is_extinct" %in% names(O)) O <- O[is_extinct == 0L] if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism)) O[, `:=`(inc = relation %in% c("same", "within"), sub = relation == "contains", diff --git a/R/local_build.R b/R/local_build.R index efbd99c..d46c090 100644 --- a/R/local_build.R +++ b/R/local_build.R @@ -105,7 +105,9 @@ NSR_local_build <- function(sources = "wcvp", d$taxon_id[fill] <- res$taxon_id[match(d$taxon_name[fill], res$name)] d <- d[!is.na(d$taxon_id) & !is.na(d$region_key), , drop = FALSE] data.frame(taxon_id = d$taxon_id, region_key = d$region_key, status = d$status, - is_cultivated = d$is_cultivated, source_name = s, stringsAsFactors = FALSE) + is_cultivated = d$is_cultivated, + is_extinct = if (is.null(d$is_extinct)) 0L else as.integer(d$is_extinct), + source_name = s, stringsAsFactors = FALSE) }) checklist <- unique(do.call(rbind, c(list(keep_old(existing$checklist, "source_name")), chk))) diff --git a/R/local_globals.R b/R/local_globals.R index 012f718..25288e3 100644 --- a/R/local_globals.R +++ b/R/local_globals.R @@ -4,7 +4,7 @@ utils::globalVariables(c( ".", ".N", ".SD", ":=", "code", "conflict_type", "consulted_here", "fraction", "from_region", - "has_int", "has_nat", "in_place", "is_cultivated", + "has_int", "has_nat", "in_place", "is_cultivated", "is_extinct", "n_here", "n_in", "n_nat", "n_src_sets", "n_sub_int", "n_sub_nat", "ovl", "ovl_srcs", "qid", "r_ord", "reason", "region_key", "relation", "scope", "source_name", "srcs", "srcs_any", diff --git a/R/local_import.R b/R/local_import.R index 97d04e0..fff73aa 100644 --- a/R/local_import.R +++ b/R/local_import.R @@ -80,12 +80,17 @@ nsr_import_wcvp <- function(zip = NULL, bb, quiet = FALSE) { dist$rank_published <- nm$taxon_rank[match(acc, nm$plant_name_id)] flag <- function(x) !is.na(x) & x %in% c("1", 1, TRUE, "TRUE") - # An extinct record is not a standing distribution. Left in, a taxon that no longer - # occurs in a region comes back native there, and counts as evaluable, which also moves - # absence and endemism answers elsewhere. A DOUBTFUL record is kept, but only as - # "present": it is too weak to support native or introduced, and dropping it would let - # the region read as a confident absence when the source in fact records a maybe. - keep <- !is.na(dist$taxon_id) & !is.na(dist$taxon_name) & !flag(dist$extinct) + # Extinct records are KEPT and flagged rather than dropped here. WCVP's `extinct` is a + # bare 0/1 with no date - it means "considered no longer present as of this release" - + # so whether it should count is a property of the question, not of the build: present-day + # status wants it out, a historical distribution wants it in. NSR_local(exclude_extinct) + # decides at query time, which is only possible if the build keeps the row. Dropping it + # here would also lose the fact that the taxon was ONCE native there, and the Ie rule + # needs that fact or it will call a natively extirpated record introduced. + # A DOUBTFUL record is kept as "present": too weak to support native or introduced, and + # dropping it would let the region read as a confident absence when the source in fact + # records a maybe. + keep <- !is.na(dist$taxon_id) & !is.na(dist$taxon_name) dist <- dist[keep, , drop = FALSE] j <- j[keep] # j indexes nm per dist row, so it is subset with it @@ -95,7 +100,7 @@ nsr_import_wcvp <- function(zip = NULL, bb, quiet = FALSE) { rank = nm$taxon_rank[j], status = ifelse(flag(dist$location_doubtful), "present", ifelse(flag(dist$introduced), "introduced", "native")), - is_cultivated = 0L, + is_cultivated = 0L, is_extinct = as.integer(flag(dist$extinct)), region_key = paste0("wgsrpd3:", dist$area_code_l3), region_name = dist$area, stringsAsFactors = FALSE ) @@ -171,7 +176,7 @@ nsr_import_vascan <- function(zip, bb, quiet = FALSE) { gid0 <- rep(NA_character_, nrow(di)) gid0[is.na(gid1) & loc == "greenland"] <- "GRL" out <- data.frame(taxon_name = di$taxon_name, taxon_id = NA_character_, status = status, - is_cultivated = 0L, + is_cultivated = 0L, is_extinct = 0L, region_key = ifelse(!is.na(gid1), paste0("gadm1:", gid1), ifelse(!is.na(gid0), paste0("gadm0:", gid0), NA_character_)), region_name = di$locality, stringsAsFactors = FALSE) @@ -214,6 +219,7 @@ nsr_import_flbr <- function(zip, bb, quiet = FALSE) { k <- match(hasc, bb$state$hasc_full) out <- data.frame( taxon_name = di$taxon_name, taxon_id = NA_character_, status = status, is_cultivated = cult, + is_extinct = 0L, region_key = ifelse(is.na(bb$state$gid_1[k]), NA_character_, paste0("gadm1:", bb$state$gid_1[k])), region_name = bb$state$state_province[k], stringsAsFactors = FALSE diff --git a/dev_notes/01-offline-nsr-design.md b/dev_notes/01-offline-nsr-design.md index e0be331..6c2aa69 100644 --- a/dev_notes/01-offline-nsr-design.md +++ b/dev_notes/01-offline-nsr-design.md @@ -90,6 +90,23 @@ resolution WCVP gives us) is too coarse to say anything about sub-country status | **vascan** | Canada, by province; native/introduced/ephemeral | DwC-A, `data.canadensys.net/ipt/archive.do?r=vascan`, ~10 MB, updated 2026-08-04 | CC0 | | **flbr** | Brazil, by state; native/naturalised/cultivated | DwC-A, `ipt.jbrj.gov.br/jbrj/archive.do?r=lista_especies_flora_brasil`, ~50-100 MB | CC BY 4.0 | +**Extinct records are kept and filtered at query time (BM, 2026-09-22):** WCVP's +`extinct` is a bare 0/1 with no date or year - checked against the v15 distribution table, +whose only columns are `plant_locality_id, plant_name_id, continent_code_l1, continent, +region_code_l2, region, area_code_l3, area, introduced, extinct, location_doubtful`. It +therefore means "considered no longer present as of this release" and cannot be compared +against an occurrence's own date. Whether it should count is a property of the question, +so the build keeps the row with an `is_extinct` flag and `NSR_local(exclude_extinct = TRUE)` +decides per call: the default answers about the present day, `FALSE` models a past +distribution. 2,701 of 1,970,252 rows are affected, over 2,393 taxa; 906 of those taxa have +no surviving native record at all. + +Independently of that switch, the endemism rules always read the taxon's whole native +range, extinct records included. A region a taxon has been lost from is still a region it +was native to, so `Ie` must not fire against it - the answer for a natively extirpated +record is `A` (gone), never `I`/`Ie` (arrived). Filtering at build time, as the garden +pipeline's script 02 does, would have destroyed the evidence needed for that distinction. + **The source is named `wcvp`, not `powo` (BM, 2026-09-22):** it is the WCVP archive that is read, and it is what `TNRS_local(sources = "wcvp")` calls the same data, so the stack uses one name for one dataset. `"powo"` is accepted as a synonym in `sources =` and in `files =`. The diff --git a/man/NSR_local.Rd b/man/NSR_local.Rd index b22183f..985ac5a 100644 --- a/man/NSR_local.Rd +++ b/man/NSR_local.Rd @@ -9,6 +9,7 @@ NSR_local( dir = nsr_cache_dir(), resolve_names = TRUE, min_overlap = 0.01, + exclude_extinct = TRUE, quiet = FALSE ) } @@ -23,6 +24,11 @@ or political division names (see Place).} \item{min_overlap}{Ignore region links covering less than this share of a region.} +\item{exclude_extinct}{Ignore distribution records the source marks extinct? Default +\code{TRUE}, which answers about the present day. Set \code{FALSE} to model a past +distribution, where a region the taxon has since been lost from still counts. See +Extinct records.} + \item{quiet}{Suppress progress messages?} } \value{ @@ -63,6 +69,21 @@ about the taxon's whole range rather than about one polygon, so they draw on evi from elsewhere. A record of a taxon confined to California, found in Michigan, is introduced there however little Michigan's own checklists say. } +\section{Extinct records}{ + +WCVP marks a distribution extinct with a bare flag and no date, so the record says +only "considered no longer present as of this release". Whether it should count is +therefore a property of the question rather than of the data: a present-day status +wants it out, a historical distribution wants it in, and the cache keeps both so +\code{exclude_extinct} can decide per call. + +Excluding it never makes the taxon introduced there. The endemism rules read the +taxon's whole native range, extinct records included, whatever this argument says - +a region a taxon has been lost from is still a region it was native to, so \code{Ie} +("introduced, inferred from endemism elsewhere") cannot fire against it. The answer +for such a record is \code{A}: gone, not foreign. +} + \section{County and parish}{ \code{county_parish} is accepted and echoed, but not resolved: the political-division diff --git a/man/NSR_local_by_region.Rd b/man/NSR_local_by_region.Rd index 806681a..58befab 100644 --- a/man/NSR_local_by_region.Rd +++ b/man/NSR_local_by_region.Rd @@ -9,7 +9,8 @@ NSR_local_by_region( region_keys, country_keys = NULL, dir = nsr_cache_dir(), - min_overlap = 0.01 + min_overlap = 0.01, + exclude_extinct = TRUE ) } \arguments{ @@ -24,6 +25,9 @@ one polygon each.} \item{dir}{Cache directory.} \item{min_overlap}{Ignore region links covering less than this share of a region.} + +\item{exclude_extinct}{Ignore distribution records the source marks extinct? See +\code{\link{NSR_local}}.} } \value{ A data.frame with \code{native_status} and the status columns diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R index 7d7dea9..d259bec 100644 --- a/tests/testthat/test-local-nsr-paths.R +++ b/tests/testthat/test-local-nsr-paths.R @@ -7,12 +7,25 @@ # a small world: Quebec in two geographies, inside Canada, plus a sliver of a region the # query only just touches fixture_db <- function() { + # one row per taxon x region x source; is_extinct only ever set by WCVP + spec <- c( + "t1 wgsrpd3:QUE native wcvp 0 0", + "t2 gadm0:CAN introduced vascan 0 0", + "t3 wgsrpd3:ONT native wcvp 0 0", + "t3 wgsrpd3:BZL native wcvp 0 0", + "t4 wgsrpd3:QUE present wcvp 1 0", + "t4 gadm1:CAN.11_1 native vascan 0 0", + "t5 gadm1:CAN.11_1 native vascan 0 0", + "t5 gadm1:CAN.4_1 native vascan 0 0", + # t6 was native in Quebec and has been lost from it; it survives only in Brazil + "t6 wgsrpd3:QUE native wcvp 0 1", + "t6 wgsrpd3:BZL native wcvp 0 0", + # t7 is the contrast: native ONLY in Brazil, never recorded in Quebec at all + "t7 wgsrpd3:BZL native wcvp 0 0") + f <- do.call(rbind, strsplit(trimws(spec), "[[:space:]]+")) checklist <- data.frame( - taxon_id = c("t1", "t2", "t3", "t3", "t4", "t4", "t5", "t5"), - region_key = c("wgsrpd3:QUE", "gadm0:CAN", "wgsrpd3:ONT", "wgsrpd3:BZL", "wgsrpd3:QUE", "gadm1:CAN.11_1", "gadm1:CAN.11_1", "gadm1:CAN.4_1"), - status = c("native", "introduced", "native", "native", "present", "native", "native", "native"), - source_name = c("wcvp", "vascan", "wcvp", "wcvp", "wcvp", "vascan", "vascan", "vascan"), - is_cultivated = c(0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L), + taxon_id = f[, 1], region_key = f[, 2], status = f[, 3], source_name = f[, 4], + is_cultivated = as.integer(f[, 5]), is_extinct = as.integer(f[, 6]), stringsAsFactors = FALSE) links <- data.frame( from_region = c("gadm1:CAN.11_1", "wgsrpd3:QUE", "gadm1:CAN.11_1"), @@ -21,8 +34,8 @@ fixture_db <- function() { fraction = c(0.99, 0.99, 0.001), stringsAsFactors = FALSE) taxa <- data.frame( - taxon_id = c("t1", "t2", "t3", "t4"), - species_name = c("Sp one", "Sp two", "Sp three", "Sp four"), + taxon_id = paste0("t", 1:7), + species_name = paste("Sp", 1:7), family = "Fam", genus = "Gen", rank = "species", stringsAsFactors = FALSE) sources <- data.frame(source_name = c("wcvp", "vascan"), is_comprehensive = TRUE, stringsAsFactors = FALSE) @@ -144,6 +157,36 @@ test_that("a range confined to one country's states is confined to that country" expect_match(set$reason, "Canada") }) +# WCVP's extinct flag carries no date, so whether it counts belongs to the question +# rather than to the build: the cache keeps the row and the query decides. +test_that("exclude_extinct decides whether a lost population still answers", { + skip_if_not_installed("data.table") + db <- fixture_db() + places <- list(fine = list("gadm1:CAN.11_1"), country = list(character(0))) + for (f in list(NSR:::nsr_resolve_status, NSR:::nsr_resolve_status_set)) { + gone <- f("t6", places, db, 0.01, exclude_extinct = TRUE) + past <- f("t6", places, db, 0.01, exclude_extinct = FALSE) + # present day: it is not there any more + expect_equal(gone$code, "A") + # as a historical distribution: it was native there + expect_equal(past$code, "N") + } +}) + +test_that("a lost population is absent, never introduced", { + skip_if_not_installed("data.table") + db <- fixture_db() + places <- list(fine = list("gadm1:CAN.11_1"), country = list(character(0))) + # t7 is native only in Brazil and was never in Quebec, so Quebec is an inferred + # introduction - this is the machinery that must NOT fire for t6 + expect_equal(NSR:::nsr_resolve_status("t7", places, db)$code, "Ie") + expect_equal(NSR:::nsr_resolve_status_set("t7", places, db)$code, "Ie") + # t6's surviving range is equally confined to Brazil, but Quebec is part of its native + # range whether or not the extinct record answers, so absence there is loss, not arrival + expect_equal(NSR:::nsr_resolve_status("t6", places, db)$code, "A") + expect_equal(NSR:::nsr_resolve_status_set("t6", places, db)$code, "A") +}) + test_that("min_overlap keeps slivers out of the consulted set", { db <- fixture_db() keys <- "gadm1:CAN.11_1" From 73a4ace187ba6e82a09591e8b11802eac03ace0a Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Tue, 22 Sep 2026 17:17:17 -0400 Subject: [PATCH 5/8] Answer from division names by default; compare coordinates, never pool them Two related changes to how NSR_local() decides where a record is. Coordinates are now opt-in. use_coordinates defaults to FALSE (BM), and any latitude and longitude present are echoed but not consulted unless it is set. Placing every record by point is a raster lookup per call, the most expensive thing in the query path, and it was being paid whether or not the caller wanted it; and turning coordinates into political divisions is GVS's job, so doing it here too means two implementations of one lookup drifting apart. The coordinate path stays available, because consulting each source in its own geography with no crosswalk is the finer question - but it is now a choice rather than a default. When both are given, they are compared rather than pooled. The region keys from the point and from the names used to be unioned, so a record whose two statements of place disagreed was answered from evidence for both and described neither. Coordinates are not made authoritative instead: a transposed longitude and a mistyped province are equally easy mistakes. The two are compared in the geography they share, GADM, at each level separately - so a right country with a wrong state is caught as readily as a wrong country - and a level only one of them speaks to is not a disagreement. Where they agree the wider key set answers, as before. Where they do not, place_conflict is TRUE, native_status is UNK, and native_status_coordinates and native_status_names carry the two answers for the caller to judge. The service has no combined mode to follow here: NSR() takes names and NSR_from_coordinates() takes coordinates, and the two are separate endpoints. The existing NSR_from_coordinates() example already renames country to country_declared so declared and inferred divisions can be compared, which is the same instinct. Those columns stay in the output whether or not coordinates are used, rather than the schema changing with an argument: place_conflict is FALSE and the two per-basis columns are NA. A message, suppressible with quiet, says when coordinates were supplied and not used, since occurrence data usually carries them and their being ignored should not be silent. nsr_places_agree() and nsr_split_keys() are split out and tested directly, since nsr_query_regions() needs a built cache and they do not. The cached reference test asks for coordinates explicitly and now says so; it also asserts the new default, that the row with coordinates and no division names has no place at all while the named rows are unchanged. R CMD check: Status OK. Co-Authored-By: Claude Opus 5 --- NEWS | 8 ++ R/NSR_local.R | 137 +++++++++++++++++++++----- dev_notes/01-offline-nsr-design.md | 25 +++++ man/NSR_local.Rd | 32 ++++-- tests/testthat/test-local-nsr-paths.R | 27 +++++ tests/testthat/test-local-nsr.R | 13 ++- 6 files changed, 210 insertions(+), 32 deletions(-) diff --git a/NEWS b/NEWS index 0d947d1..8015925 100644 --- a/NEWS +++ b/NEWS @@ -14,6 +14,14 @@ answers about the present day, FALSE models a past distribution. Either way a region a taxon has been lost from stays part of its native range, so such a record is reported absent rather than introduced. +* NSR_local() answers from political division names by default. Coordinates are + echoed but not used unless use_coordinates = TRUE: placing every record by point + costs a raster lookup per call, and turning coordinates into political divisions is + what GVS already does. +* Where a record gives both coordinates and division names and the two resolve to + different places, NSR_local() no longer pools them. Neither is treated as + authoritative: place_conflict flags the record, native_status is UNK, and the two + bases are answered separately in native_status_coordinates and native_status_names. ## Bug fixes diff --git a/R/NSR_local.R b/R/NSR_local.R index 42ae764..b28dbbf 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -4,13 +4,26 @@ #' checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s columns, #' with four added. #' -#' \strong{Place.} Give coordinates (\code{latitude}, \code{longitude}), political -#' division names (\code{country}, \code{state_province}, \code{county_parish}), or both. -#' Coordinates are the better input: each source is consulted in the geography it -#' publishes against (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil -#' against GADM states), with no crosswalk between them. Names are resolved to GADM -#' units through the GNRS backbone and then carried to other geographies by the spatial -#' link table. +#' \strong{Place.} By default the answer comes from political division names +#' (\code{country}, \code{state_province}, \code{county_parish}), which are resolved to +#' GADM units through the GNRS backbone and carried to each source's own geography by the +#' spatial link table. Any \code{latitude} and \code{longitude} present are echoed but +#' not used: placing every record by point costs a raster lookup per call, and turning +#' coordinates into political divisions is what GVS already does. Resolve them there and +#' pass the divisions here, or set \code{use_coordinates = TRUE} to have this function do +#' it. +#' +#' With \code{use_coordinates = TRUE}, a record is placed by its point in each source's +#' own geography (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil against +#' GADM states), with no crosswalk between them - the finer question, where the data +#' supports it. Given coordinates and names together, the two are independent claims +#' about where the record is and either can be wrong: a transposed longitude and a +#' mistyped province are equally easy mistakes, so neither is authoritative. Where they +#' agree, both are used. Where they disagree, \code{place_conflict} is \code{TRUE}, +#' \code{native_status} is \code{UNK}, and the two bases are answered separately in +#' \code{native_status_coordinates} and \code{native_status_names} for you to judge. The +#' service has no combined mode to follow here: \code{\link{NSR}} takes names and +#' \code{\link{NSR_from_coordinates}} takes coordinates. #' #' \strong{Precedence.} If any source says native, the answer is native: asserting #' nativity is a positive claim, whereas "introduced" and "absent" are often artefacts of @@ -36,6 +49,9 @@ #' @param dir Cache directory, shared with GNRS and GVS. #' @param resolve_names Resolve submitted names against WCVP with #' \code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none. +#' @param use_coordinates Place each record by its \code{latitude} and \code{longitude} +#' as well as by its division names? Default \code{FALSE}: coordinates cost a raster +#' lookup per call, and resolving them to political divisions is GVS's job. See Place. #' @param min_overlap Ignore region links covering less than this share of a region. #' @param exclude_extinct Ignore distribution records the source marks extinct? Default #' \code{TRUE}, which answers about the present day. Set \code{FALSE} to model a past @@ -63,7 +79,8 @@ #' warning says so. Coordinates are the way to ask a finer question. #' @export NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names = TRUE, - min_overlap = 0.01, exclude_extinct = TRUE, quiet = FALSE) { + use_coordinates = FALSE, min_overlap = 0.01, exclude_extinct = TRUE, + quiet = FALSE) { if (!inherits(occurrence_dataframe, "data.frame")) { stop("occurrence_dataframe should be a data.frame", call. = FALSE) } @@ -111,12 +128,47 @@ NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names } # ---- place ------------------------------------------------------------------------ - places <- nsr_query_regions(lon, lat, country, state, county, dir, db) + if (!use_coordinates && any(is.finite(lon) & is.finite(lat)) && !quiet) { + message("Coordinates present but not used: answering from division names. ", + "Set use_coordinates = TRUE to place each record by its point instead.") + } + places <- nsr_query_regions(lon, lat, country, state, county, dir, db, use_coordinates) # ---- opinions --------------------------------------------------------------------- use_set <- !is.null(db$chk_dt) && nrow(occurrence_dataframe) >= getOption("NSR.set_min", 200) - res <- if (use_set) nsr_resolve_status_set(taxon_id, places, db, min_overlap, exclude_extinct) - else nsr_resolve_status(taxon_id, places, db, min_overlap, exclude_extinct) + resolve <- function(tid, pl) if (use_set) + nsr_resolve_status_set(tid, pl, db, min_overlap, exclude_extinct) else + nsr_resolve_status(tid, pl, db, min_overlap, exclude_extinct) + res <- resolve(taxon_id, places) + + # ---- records whose two statements of place disagree --------------------------------- + # Neither basis is taken as authoritative, so the pooled answer - which drew on + # evidence for both places and describes neither - is withdrawn, and each basis is + # answered on its own. The two columns are populated only for these rows; elsewhere + # they are NA because there is nothing to compare, and native_status is the answer. + status_xy <- rep(NA_character_, n) + status_nm <- rep(NA_character_, n) + cf <- places$conflict + if (any(cf)) { + status_xy[cf] <- resolve(taxon_id[cf], nsr_split_keys(places$xy[cf]))$code + status_nm[cf] <- resolve(taxon_id[cf], nsr_split_keys(places$nm[cf]))$code + res$code[cf] <- "UNK" + res$reason[cf] <- paste0( + "Coordinates and division names resolve to different places, so no single answer ", + "is given; see native_status_coordinates and native_status_names") + res$scope[cf] <- "none" + res$sources[cf] <- NA_character_ + res$opinions[cf] <- NA_character_ + res$country_code[cf] <- NA_character_ + res$state_code[cf] <- NA_character_ + res$cultivated[cf] <- NA_integer_ + res$conflict[cf] <- FALSE + res$conflict_type[cf] <- "none" + res$n_sub_native[cf] <- 0L + res$n_sub_introduced[cf] <- 0L + if (!quiet) message(sum(cf), " record(s) whose coordinates and division names ", + "disagree; see place_conflict.") + } out <- data.frame( family = db$taxa$family[match(taxon_id, db$taxa$taxon_id)], @@ -146,6 +198,9 @@ NSR_local <- function(occurrence_dataframe, dir = nsr_cache_dir(), resolve_names n_subpolygons_native = res$n_sub_native, n_subpolygons_introduced = res$n_sub_introduced, native_status_opinions = res$opinions, + native_status_coordinates = status_xy, + native_status_names = status_nm, + place_conflict = places$conflict, regions_matched = places$matched, taxon_evaluable = !is.na(taxon_id) & taxon_id %in% db$evaluable, user_id = user_id, @@ -258,47 +313,85 @@ nsr_index_db <- function(db) { #' table. Returns, per row, the regions to consult and how each relates to the query. #' @keywords internal #' @noRd -nsr_query_regions <- function(lon, lat, country, state, county, dir, db) { +nsr_query_regions <- function(lon, lat, country, state, county, dir, db, + use_coordinates = FALSE) { n <- length(lon) direct <- vector("list", n) fine <- vector("list", n) ctry <- vector("list", n) + xy <- vector("list", n) + nm <- vector("list", n) + conflict <- rep(FALSE, n) label <- rep(NA_character_, n) level <- rep("country", n) matched <- rep("none", n) - loc <- if (any(is.finite(lon) & is.finite(lat))) nsr_locate_regions(lon, lat, dir) else NULL + loc <- if (use_coordinates && any(is.finite(lon) & is.finite(lat))) + nsr_locate_regions(lon, lat, dir) else NULL bb <- try(nsr_gnrs_backbone(dir), silent = TRUE) pd <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, state, bb) else NULL pd0 <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, NULL, bb) else NULL + at <- function(k, p) grep(p, k, value = TRUE) for (i in seq_len(n)) { - keys <- character(0) - if (!is.null(loc)) keys <- c(keys, stats::na.omit(c(loc$wgsrpd3[i], loc$gadm[i], loc$gadm0[i]))) + kx <- if (!is.null(loc)) + as.character(stats::na.omit(c(loc$wgsrpd3[i], loc$gadm[i], loc$gadm0[i]))) else character(0) + kn <- character(0) if (!is.null(pd)) { gid1 <- if (!is.na(pd$state_province_id[i])) bb$state$gid_1[match(pd$state_province_id[i], bb$state$state_province_id)] else NA_character_ gid0 <- if (!is.na(pd0$country_id[i])) bb$country$gid_0[match(pd0$country_id[i], bb$country$country_id)] else NA_character_ - if (!is.na(gid1)) keys <- c(keys, paste0("gadm1:", gid1)) - if (!is.na(gid0)) keys <- c(keys, paste0("gadm0:", gid0)) + if (!is.na(gid1)) kn <- c(kn, paste0("gadm1:", gid1)) + if (!is.na(gid0)) kn <- c(kn, paste0("gadm0:", gid0)) } - keys <- unique(keys[!is.na(keys)]) + xy[[i]] <- kx + nm[[i]] <- kn + conflict[i] <- !nsr_places_agree(kx, kn) + keys <- unique(c(kx, kn)) + keys <- keys[!is.na(keys)] # the finest place the query actually names, kept apart from its country: a question # about Amazonas must not inherit Brazil's answer + ctry[[i]] <- at(keys, "^gadm0:") fine[[i]] <- grep("^gadm0:", keys, value = TRUE, invert = TRUE) - ctry[[i]] <- grep("^gadm0:", keys, value = TRUE) direct[[i]] <- keys has_state <- any(grepl("^gadm1:", keys)) || (!is.null(pd) && !is.na(pd$state_province_id[i])) level[i] <- if (has_state) "state_province" else "country" - matched[i] <- if (!length(keys)) "none" else + matched[i] <- if (!length(keys)) "none" else if (conflict[i]) "conflict" else if (!is.null(loc) && !is.na(loc$gadm[i])) "coordinates" else "names" label[i] <- paste(stats::na.omit(c(country[i], if (!is.na(state[i]) && nzchar(state[i])) state[i])), collapse = ":") if (!nzchar(label[i]) && length(keys)) label[i] <- paste(keys, collapse = " + ") } - list(direct = direct, fine = fine, country = ctry, label = label, level = level, - matched = matched) + list(direct = direct, fine = fine, country = ctry, xy = xy, nm = nm, + conflict = conflict, label = label, level = level, matched = matched) +} + +#' Do coordinates and division names describe the same place? +#' +#' Internal. Coordinates and names are two independent claims about where a record is, +#' and either can be wrong: a transposed longitude and a mistyped province are equally +#' easy mistakes, so neither is authoritative. They are compared in the geography they +#' share, GADM, and at each level separately, so a right country with a wrong state is +#' caught as readily as a wrong country. A level only one of them speaks to is not a +#' disagreement. Where they agree, pooling the keys costs nothing and the wider set +#' answers; where they do not, \code{\link{NSR_local}} answers each on its own. +#' @keywords internal +#' @noRd +nsr_places_agree <- function(kx, kn) { + lvl <- function(p) { + a <- grep(p, kx, value = TRUE); b <- grep(p, kn, value = TRUE) + !length(a) || !length(b) || identical(sort(a), sort(b)) + } + lvl("^gadm0:") && lvl("^gadm1:") +} + +#' Split a set of region keys into the finest level named and its country +#' @keywords internal +#' @noRd +nsr_split_keys <- function(keys) { + list(fine = lapply(keys, function(z) grep("^gadm0:", z, value = TRUE, invert = TRUE)), + country = lapply(keys, function(z) grep("^gadm0:", z, value = TRUE))) } #' Gather and reduce every opinion bearing on each row diff --git a/dev_notes/01-offline-nsr-design.md b/dev_notes/01-offline-nsr-design.md index 6c2aa69..915a124 100644 --- a/dev_notes/01-offline-nsr-design.md +++ b/dev_notes/01-offline-nsr-design.md @@ -90,6 +90,31 @@ resolution WCVP gives us) is too coarse to say anything about sub-country status | **vascan** | Canada, by province; native/introduced/ephemeral | DwC-A, `data.canadensys.net/ipt/archive.do?r=vascan`, ~10 MB, updated 2026-08-04 | CC0 | | **flbr** | Brazil, by state; native/naturalised/cultivated | DwC-A, `ipt.jbrj.gov.br/jbrj/archive.do?r=lista_especies_flora_brasil`, ~50-100 MB | CC BY 4.0 | +**Division names are the default input; coordinates are opt-in (BM, 2026-09-22):** +`NSR_local()` answers from `country`/`state_province` unless `use_coordinates = TRUE`. +Two reasons. Placing every record by point costs a raster lookup per call, which is the +single most expensive thing the query path can do. And resolving coordinates to political +divisions is GVS's job: a pipeline that needs it already runs GVS, and duplicating that +step here would mean two implementations of the same lookup drifting apart. The +coordinate path stays available because it asks a finer question - it consults each +source in its own geography with no crosswalk - but it is a deliberate choice rather than +a default. + +**Coordinates and names are compared, not pooled (BM, 2026-09-22):** the first +implementation unioned the region keys from both, so a record whose point and whose +division names disagreed was answered from evidence for two places and described neither. +Coordinates are not taken as authoritative instead: a transposed longitude and a mistyped +province are equally easy mistakes, and the existing `NSR_from_coordinates()` example +already renames `country` to `country_declared` precisely so declared and inferred +divisions can be compared. The two are therefore compared in the geography they share, +GADM, at each level separately - a right country with a wrong state is caught as readily +as a wrong country - and a level only one of them speaks to is not a disagreement. Where +they agree the wider key set answers, as before. Where they do not, `place_conflict` is +`TRUE`, `native_status` is `UNK`, and `native_status_coordinates` and +`native_status_names` carry the two answers. There is no service behaviour to follow +here: `NSR()` takes names only and `NSR_from_coordinates()` takes coordinates only, so +the API never meets the case. + **Extinct records are kept and filtered at query time (BM, 2026-09-22):** WCVP's `extinct` is a bare 0/1 with no date or year - checked against the v15 distribution table, whose only columns are `plant_locality_id, plant_name_id, continent_code_l1, continent, diff --git a/man/NSR_local.Rd b/man/NSR_local.Rd index 985ac5a..b16f159 100644 --- a/man/NSR_local.Rd +++ b/man/NSR_local.Rd @@ -8,6 +8,7 @@ NSR_local( occurrence_dataframe, dir = nsr_cache_dir(), resolve_names = TRUE, + use_coordinates = FALSE, min_overlap = 0.01, exclude_extinct = TRUE, quiet = FALSE @@ -22,6 +23,10 @@ or political division names (see Place).} \item{resolve_names}{Resolve submitted names against WCVP with \code{TNRS::TNRS_local()}? Names already matching WCVP accepted names need none.} +\item{use_coordinates}{Place each record by its \code{latitude} and \code{longitude} +as well as by its division names? Default \code{FALSE}: coordinates cost a raster +lookup per call, and resolving them to political divisions is GVS's job. See Place.} + \item{min_overlap}{Ignore region links covering less than this share of a region.} \item{exclude_extinct}{Ignore distribution records the source marks extinct? Default @@ -41,13 +46,26 @@ checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s c with four added. } \details{ -\strong{Place.} Give coordinates (\code{latitude}, \code{longitude}), political -division names (\code{country}, \code{state_province}, \code{county_parish}), or both. -Coordinates are the better input: each source is consulted in the geography it -publishes against (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil -against GADM states), with no crosswalk between them. Names are resolved to GADM -units through the GNRS backbone and then carried to other geographies by the spatial -link table. +\strong{Place.} By default the answer comes from political division names +(\code{country}, \code{state_province}, \code{county_parish}), which are resolved to +GADM units through the GNRS backbone and carried to each source's own geography by the +spatial link table. Any \code{latitude} and \code{longitude} present are echoed but +not used: placing every record by point costs a raster lookup per call, and turning +coordinates into political divisions is what GVS already does. Resolve them there and +pass the divisions here, or set \code{use_coordinates = TRUE} to have this function do +it. + +With \code{use_coordinates = TRUE}, a record is placed by its point in each source's +own geography (WCVP against WGSRPD level-3 areas, VASCAN and Flora do Brasil against +GADM states), with no crosswalk between them - the finer question, where the data +supports it. Given coordinates and names together, the two are independent claims +about where the record is and either can be wrong: a transposed longitude and a +mistyped province are equally easy mistakes, so neither is authoritative. Where they +agree, both are used. Where they disagree, \code{place_conflict} is \code{TRUE}, +\code{native_status} is \code{UNK}, and the two bases are answered separately in +\code{native_status_coordinates} and \code{native_status_names} for you to judge. The +service has no combined mode to follow here: \code{\link{NSR}} takes names and +\code{\link{NSR_from_coordinates}} takes coordinates. \strong{Precedence.} If any source says native, the answer is native: asserting nativity is a positive claim, whereas "introduced" and "absent" are often artefacts of diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R index d259bec..ac8aed3 100644 --- a/tests/testthat/test-local-nsr-paths.R +++ b/tests/testthat/test-local-nsr-paths.R @@ -196,3 +196,30 @@ test_that("min_overlap keeps slivers out of the consulted set", { # lower the threshold and it comes back, so it is the threshold doing the work expect_true("wgsrpd3:SLIVER" %in% NSR:::nsr_consulted_regions(keys, db, min_overlap = 0)) }) + +# Coordinates and names are compared, not pooled, when they disagree: neither is +# authoritative, and a pooled answer would describe neither place. +test_that("place agreement is judged per level, in the geography the two share", { + xy <- c("wgsrpd3:QUE", "gadm1:CAN.11_1", "gadm0:CAN") + # the same place, said twice + expect_true(NSR:::nsr_places_agree(xy, c("gadm1:CAN.11_1", "gadm0:CAN"))) + # right country, wrong province - the case a country-only comparison would miss + expect_false(NSR:::nsr_places_agree(xy, c("gadm1:CAN.4_1", "gadm0:CAN"))) + # wrong country + expect_false(NSR:::nsr_places_agree(xy, c("gadm0:BRA"))) + # a level only one basis speaks to is not a disagreement + expect_true(NSR:::nsr_places_agree(xy, c("gadm0:CAN"))) + expect_true(NSR:::nsr_places_agree("wgsrpd3:QUE", c("gadm1:CAN.11_1", "gadm0:CAN"))) + # nothing to compare + expect_true(NSR:::nsr_places_agree(character(0), character(0))) +}) + +test_that("keys split into the finest level named and its country", { + s <- NSR:::nsr_split_keys(list(c("wgsrpd3:QUE", "gadm1:CAN.11_1", "gadm0:CAN"), + "gadm0:BRA", character(0))) + expect_equal(s$fine[[1]], c("wgsrpd3:QUE", "gadm1:CAN.11_1")) + expect_equal(s$country[[1]], "gadm0:CAN") + expect_equal(s$fine[[2]], character(0)) + expect_equal(s$country[[2]], "gadm0:BRA") + expect_equal(s$fine[[3]], character(0)) +}) diff --git a/tests/testthat/test-local-nsr.R b/tests/testthat/test-local-nsr.R index 2662b97..cb854a2 100644 --- a/tests/testthat/test-local-nsr.R +++ b/tests/testthat/test-local-nsr.R @@ -81,14 +81,21 @@ test_that("NSR_local answers the reference cases (needs a built cache)", { latitude = c(NA, NA, NA, NA, NA, -23.5), longitude = c(NA, NA, NA, NA, NA, -47.5), stringsAsFactors = FALSE) - r <- NSR_local(x, dir = cache, quiet = TRUE) + r <- NSR_local(x, dir = cache, quiet = TRUE, use_coordinates = TRUE) expect_equal(r$native_status, c("N", "N", "A", "N", "UNK", "N")) # Sao Paulo is answered by the Brazilian flora, which WGSRPD cannot reach expect_match(r$native_status_sources[2], "flbr") - # a moss is unknown, not absent: POWO holds no information about it + # a moss is unknown, not absent: WCVP holds no information about it expect_false(r$taxon_evaluable[5]) # coordinates are resolved in each source's own geography expect_equal(r$regions_matched[6], "coordinates") + + # by default coordinates are not consulted at all, so the row that has nothing else + # to go on has no place; the rows with division names are unaffected + d <- NSR_local(x, dir = cache, quiet = TRUE) + expect_equal(d$native_status[1:5], r$native_status[1:5]) + expect_equal(d$regions_matched[6], "none") + expect_equal(d$native_status[6], "UNK") }) test_that("region links relate the geographies sensibly (needs a built cache)", { @@ -113,7 +120,7 @@ test_that("absence becomes introduction only for taxa confined elsewhere", { db <- NSR:::nsr_local_db(cache) nat <- db$checklist[db$checklist$status == "native", ] per <- table(nat$taxon_id) - # a species POWO gives a single native area, queried far away + # a species WCVP gives a single native area, queried far away cal <- intersect(names(per)[per == 1], nat$taxon_id[nat$region_key == "wgsrpd3:CAL"]) skip_if(!length(cal), "no single-area Californian native in this build") sp <- db$taxa$species_name[match(cal[1], db$taxa$taxon_id)] From 531b8725aa11ef781b1923b003d84a022a949f00 Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Thu, 24 Sep 2026 11:41:27 -0400 Subject: [PATCH 6/8] addressing copilot stuff --- R/NSR_local.R | 7 +++ R/local_altdiv.R | 69 ++++++++++++++++++++++++++++++ tests/testthat/test-local-altdiv.R | 47 ++++++++++++++++++++ 3 files changed, 123 insertions(+) create mode 100644 R/local_altdiv.R create mode 100644 tests/testthat/test-local-altdiv.R diff --git a/R/NSR_local.R b/R/NSR_local.R index b28dbbf..6440bd4 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -331,6 +331,12 @@ nsr_query_regions <- function(lon, lat, country, state, county, dir, db, bb <- try(nsr_gnrs_backbone(dir), silent = TRUE) pd <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, state, bb) else NULL pd0 <- if (!inherits(bb, "try-error")) nsr_match_poldiv(country, NULL, bb) else NULL + # a declared state that belongs to another division system (a Norwegian county from + # after the 2018 or 2020 reform, a lan under the name records write) contributes the + # GADM units it covers, so the question is answered below country level instead of + # falling back to the country + iso0 <- if (!is.null(pd0)) bb$country$iso[match(pd0$country_id, bb$country$country_id)] else rep(NA_character_, n) + alt_keys <- nsr_altdiv_keys(iso0, state, dir) at <- function(k, p) grep(p, k, value = TRUE) for (i in seq_len(n)) { @@ -343,6 +349,7 @@ nsr_query_regions <- function(lon, lat, country, state, county, dir, db, gid0 <- if (!is.na(pd0$country_id[i])) bb$country$gid_0[match(pd0$country_id[i], bb$country$country_id)] else NA_character_ if (!is.na(gid1)) kn <- c(kn, paste0("gadm1:", gid1)) + if (is.na(gid1) && length(alt_keys[[i]])) kn <- c(kn, alt_keys[[i]]) if (!is.na(gid0)) kn <- c(kn, paste0("gadm0:", gid0)) } xy[[i]] <- kx diff --git a/R/local_altdiv.R b/R/local_altdiv.R new file mode 100644 index 0000000..530f474 --- /dev/null +++ b/R/local_altdiv.R @@ -0,0 +1,69 @@ +# =========================================================================== +# Divisions that belong to another division system +# +# A large part of the record stream declares a state that is not a GADM division: +# Swedish landskap and lappmarker, Watsonian vice-counties, Norwegian counties from +# after the 2018 and 2020 reforms. GNRS records these as units in their own right +# and, where it knows which GADM units they cover, writes that extent to the shared +# cache (altdiv-*). NSR reads those tables the same way it reads the GNRS reference: +# straight from the cache, without depending on the package. +# +# Only a unit whose extent is known is of any use here: it turns a declared name +# that matched nothing into the set of divisions the question is really about, so +# native status can be answered below country level. A unit with no extent yet adds +# no regions, which leaves the query where it was - answered at country level. +# =========================================================================== + +#' The alternative divisions whose extent is known, keyed by country and name +#' +#' Internal. NULL when the component has not been built. Cached per directory +#' for the session, like the other reference tables. +#' @keywords internal +#' @noRd +nsr_altdiv <- function(dir) { + key <- paste0("altdiv:", dir) + if (!is.null(nsr_session[[key]])) return(nsr_session[[key]]) + f <- function(x) file.path(dir, paste0("altdiv-", x, ".gz.parquet")) + if (!all(file.exists(f(c("units", "names", "extent"))))) return(NULL) + units <- as.data.frame(nanoparquet::read_parquet(f("units"))) + names_ <- as.data.frame(nanoparquet::read_parquet(f("names"))) + extent <- as.data.frame(nanoparquet::read_parquet(f("extent"))) + # only exact names, and only units that have an extent: a regex system (the + # vice-counties) has none yet, and a unit without one adds no regions + names_ <- names_[names_$match == "exact" & names_$entity_key %in% extent$entity_key, , drop = FALSE] + iso <- units$country_iso[match(names_$entity_key, units$entity_key)] + out <- list( + key = paste(toupper(iso), tolower(trimws(names_$name)), sep = "\u001f"), + entity = names_$entity_key, + extent = split(extent$gid, extent$entity_key), + level = vapply(split(extent$level, extent$entity_key), function(x) as.integer(x[1]), integer(1)) + ) + nsr_session[[key]] <- out + out +} + +#' Region keys for a declared division that belongs to another system +#' +#' Internal. Vectorised over rows; returns a list of character vectors, empty +#' where the declared name is not such a division or its extent is unknown. +#' @keywords internal +#' @noRd +nsr_altdiv_keys <- function(country_iso, name, dir) { + n <- length(name) + out <- vector("list", n) + for (i in seq_len(n)) out[[i]] <- character(0) + a <- nsr_altdiv(dir) + if (is.null(a)) return(out) + nm <- tolower(trimws(ifelse(is.na(name), "", name))) + cc <- toupper(ifelse(is.na(country_iso), "", country_iso)) + hit <- match(paste(cc, nm, sep = "\u001f"), a$key) + for (i in which(!is.na(hit))) { + k <- a$entity[hit[i]] + gid <- a$extent[[k]] + # NSR's regions are countries and states; a county-level extent has nothing to map + # onto yet, so it is skipped rather than answered at the wrong level + if (!length(gid) || a$level[[k]] != 1L) next + out[[i]] <- paste0("gadm1:", gid) + } + out +} diff --git a/tests/testthat/test-local-altdiv.R b/tests/testthat/test-local-altdiv.R new file mode 100644 index 0000000..3e21366 --- /dev/null +++ b/tests/testthat/test-local-altdiv.R @@ -0,0 +1,47 @@ +context("divisions that belong to another division system") + +# Self-contained: the tables GNRS writes to the shared cache, in a temp directory. + +dir <- file.path(tempdir(), "nsr-altdiv") +unlink(dir, recursive = TRUE) +dir.create(dir, recursive = TRUE, showWarnings = FALSE) +w <- function(x, table) nanoparquet::write_parquet(x, file.path(dir, paste0("altdiv-", table, ".gz.parquet")), + compression = "gzip") +w(data.frame(system = c("no-fylke", "se-landskap", "gb-vice-county"), + entity_key = c("NO-VIKEN", "SE-LS-Uppland", "GB-VC"), + country_iso = c("NO", "SE", "GB"), + name = c("Viken", "Uppland", "Watsonian vice-county"), + kind = c("superseded", "parallel", "parallel"), + extent_known = c(TRUE, FALSE, FALSE), stringsAsFactors = FALSE), "units") +w(data.frame(entity_key = c("NO-VIKEN", "SE-LS-Uppland", "GB-VC"), + name = c("Viken", "Uppland", "^VC ?[0-9]+"), + match = c("exact", "exact", "regex"), stringsAsFactors = FALSE), "names") +w(data.frame(entity_key = "NO-VIKEN", level = 1L, + gid = c("NOR.1_1", "NOR.4_1", "NOR.2_1"), stringsAsFactors = FALSE), "extent") + +test_that("a division with a known extent contributes the units it covers", { + k <- nsr_altdiv_keys("NO", "Viken", dir) + expect_equal(k[[1]], c("gadm1:NOR.1_1", "gadm1:NOR.4_1", "gadm1:NOR.2_1")) +}) + +test_that("a division whose extent is unknown contributes nothing", { + # recognised by GNRS, but there is no set of GADM units to ask about yet, so the + # query stays where it was rather than being answered at the wrong level + expect_equal(nsr_altdiv_keys("SE", "Uppland", dir)[[1]], character(0)) + expect_equal(nsr_altdiv_keys("GB", "VC57 Derbyshire", dir)[[1]], character(0)) +}) + +test_that("an ordinary GADM division is left alone, and matching is per country", { + expect_equal(nsr_altdiv_keys("NO", "Akershus", dir)[[1]], character(0)) + # "Viken" is a Norwegian entity: the same string under another country is not it + expect_equal(nsr_altdiv_keys("SE", "Viken", dir)[[1]], character(0)) + expect_equal(nsr_altdiv_keys(NA_character_, NA_character_, dir)[[1]], character(0)) +}) + +test_that("without the component built, nothing is claimed", { + empty <- file.path(tempdir(), "nsr-altdiv-empty") + unlink(empty, recursive = TRUE) + dir.create(empty, recursive = TRUE, showWarnings = FALSE) + expect_null(nsr_altdiv(empty)) + expect_equal(nsr_altdiv_keys("NO", "Viken", empty)[[1]], character(0)) +}) From 54b791c30c7e15a8384d6a1f5301459252c12c72 Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Thu, 24 Sep 2026 12:49:22 -0400 Subject: [PATCH 7/8] Third review pass: R 4.0 guard, cache schema, and four row/set divergences Packaging and build: - nsr_cache_dir() called tools::R_user_dir(), which arrived in R 4.0.0, while the package declares R (>= 3.5.0). On older R the local entry points failed in their default argument before doing any work. It now says what is missing and how to set it. Inventing a different default there would have put the cache where GNRS and GVS do not look and given each package its own copy, which is worse than an error; raising the package minimum would bind the API functions, which do not need it. - Adding a source to a cache built before extinct records were retained failed on the bind: the rows kept from the old checklist have no is_extinct column, and rbind() sees two schemas. nsr_index_db() already defaulted the column on read, but that is too late for the rebuild. The retained rows are normalised first. - Two warning strings held a real newline where the escape was meant. They parse, and R CMD check has been passing, but they were not what was written. Four more places the row and set paths disagreed, each found by asking the fixture rather than by reading: - A matched taxon with no place keys: the row path says "Place not matched to any region", the set path fell through to the generic absence logic and said "No comprehensive checklist covers this polygon". - conflict_type was NA in the set path for any result nothing disagreed about - a country-only query, or one whose opinions are all contains or overlaps - where the row path returns the contract value "none". - Absence was not read against the country a state was given with, although that country's opinions already feed the answer, so a state with no coverage of its own inside a comprehensively listed country was UNK in the row path and A in the set path. - A country row is now its own container for the confined-range test. A taxon native both to gadm0:CAN and to a Canadian state had no common container, so the Ie inference was missed in both paths: gadm0 has no parent to look up and its links only reach other geographies. Four tests cover them. NSR_local()'s help listed "four added columns" and now names what it actually returns. R CMD check: Status OK. Not changed, and why: - terra::extract()[, 1] is correct for the coordinate matrix passed here; the ID column is a SpatVector behaviour. - "Regenerate NAMESPACE" remains stale: the exports and Rd files have been on the branch since the first review pass. - The GADM hierarchy is deliberately not in the link table. Link rows carry overlap fractions and are filtered by min_overlap, and a state is a legitimately tiny share of its country; the containment is exact, so it is resolved directly instead. - Doubtful WCVP records stay "present": dropping them would turn a source's explicit "maybe here" into a confident absence. Co-Authored-By: Claude Opus 5 --- R/NSR_local.R | 40 ++++++++++++++++---- R/NSR_local_set.R | 32 ++++++++++++++-- R/local_build.R | 12 ++++-- R/local_cache.R | 16 +++++++- man/NSR_local.Rd | 14 ++++++- tests/testthat/test-local-nsr-paths.R | 54 +++++++++++++++++++++++++-- 6 files changed, 145 insertions(+), 23 deletions(-) diff --git a/R/NSR_local.R b/R/NSR_local.R index 6440bd4..2401a1f 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -1,8 +1,18 @@ #' Determine native status without an internet connection #' #' Offline Native Status Resolver. For each taxon and place, the status every built -#' checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s columns, -#' with four added. +#' checklist gives it, reduced to one answer. +#' +#' Output carries \code{\link{NSR}}'s columns, and adds: \code{isEndemic}; +#' \code{native_status_conflict} and \code{native_status_conflict_type}, saying whether +#' and how the sources disagreed; \code{native_status_scope}, which of the polygons +#' bearing on the place produced the answer; \code{n_subpolygons_native} and +#' \code{n_subpolygons_introduced}; \code{native_status_opinions}, every opinion consulted; +#' \code{native_status_coordinates}, \code{native_status_names} and \code{place_conflict} +#' (see Place); \code{regions_matched}, what the place was matched on; and +#' \code{taxon_evaluable}, whether any source holds native status for the taxon at all - +#' which is what separates \code{A} (absent) from \code{UNK} (no data). The coordinates +#' given are echoed as \code{latitude} and \code{longitude}. #' #' \strong{Place.} By default the answer comes from political division names #' (\code{country}, \code{state_province}, \code{county_parish}), which are resolved to @@ -444,7 +454,10 @@ nsr_resolve_status <- function(taxon_id, places, db, min_overlap = 0.01, op <- if (is.null(op)) add else Map(c, op, add) } } - consulted <- nsr_consulted_regions(keys, db, min_overlap) + # the ancestor country supplies opinions above, so absence has to be read against it + # too, or a state with no coverage of its own inside a comprehensively listed country + # answers UNK here and A in the set path, which reads Q's country row + consulted <- nsr_consulted_regions(c(keys, anc), db, min_overlap) ev <- tid %in% db$evaluable r <- nsr_reduce(op, consulted, ev, db) if (identical(r$code, "N") && nsr_is_endemic(tid, keys, db, min_overlap)) { @@ -489,13 +502,24 @@ nsr_gadm_children <- function(keys, db) { #' @keywords internal #' @noRd nsr_gadm_parents <- function(keys, db) { - if (!length(keys) || !length(db$gadm_parent)) return(character(0)) + if (!length(keys)) return(character(0)) k <- keys[startsWith(keys, "gadm1:")] - if (!length(k)) return(character(0)) - p <- as.character(unname(db$gadm_parent[k])) + p <- if (length(k) && length(db$gadm_parent)) + as.character(unname(db$gadm_parent[k])) else character(0) setdiff(p[!is.na(p)], keys) } +#' The countries a set of regions lies in, a country counting as its own +#' +#' Internal. For the confined-range test only. A taxon native both to a country row +#' and to one of that country's states is confined to it, so gadm0 has to answer for +#' itself here - it has no parent to look up, and its links only reach other geographies. +#' @keywords internal +#' @noRd +nsr_gadm_countries <- function(keys, db) { + unique(c(nsr_gadm_parents(keys, db), keys[startsWith(keys, "gadm0:")])) +} + #' Opinions about one taxon bearing on a set of regions #' #' Internal. Direct opinions, plus opinions about regions the query's regions are inside @@ -584,10 +608,10 @@ nsr_endemic_elsewhere <- function(taxon_id, keys, db, min_overlap = 0.01) { if (length(nat) == 1) return(nm(nat)) containers <- lapply(nat, function(r) { li <- db$link_env[[r]] - if (is.null(li)) return(nsr_gadm_parents(r, db)) + if (is.null(li)) return(nsr_gadm_countries(r, db)) keep <- db$link_rel[li] %in% c("same", "within") & db$link_frac[li] >= min_overlap & startsWith(db$link_to[li], "gadm0:") - c(db$link_to[li][keep], nsr_gadm_parents(r, db)) + c(db$link_to[li][keep], nsr_gadm_countries(r, db)) }) common <- Reduce(intersect, containers) if (length(common)) nm(common[1]) else NA_character_ diff --git a/R/NSR_local_set.R b/R/NSR_local_set.R index 8e9f059..582d3fc 100644 --- a/R/NSR_local_set.R +++ b/R/NSR_local_set.R @@ -37,6 +37,13 @@ nsr_confined_ranges <- function(nat, db) { list(region_key = from_region, country = to_region)]), use.names = TRUE)) } + # a country row is its own country: a taxon native both to gadm0:CAN and to one of + # Canada's states is confined to Canada, and without this it has no common container + own <- unique(nat$region_key[startsWith(nat$region_key, "gadm0:")]) + if (length(own)) { + lk <- unique(data.table::rbindlist( + list(lk, dt(region_key = own, country = own)), use.names = TRUE)) + } m <- merge(nat[taxon_id %in% multi$taxon_id], lk, by = "region_key", allow.cartesian = TRUE) per <- m[, list(n_in = data.table::uniqueN(region_key)), by = c("taxon_id", "country")] per <- merge(per, multi, by = "taxon_id") @@ -163,7 +170,9 @@ nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, dt(qid = rep(seq_len(n), nf), region_key = unlist(fine), relation = "same"), dt(qid = rep(seq_len(n), nc), region_key = unlist(ctry), relation = "same"))) # fixed below for rows that have a finer key - if (!nrow(Q)) return(nsr_set_finish(out, taxon_id, db, n, endemism = endemism)) + has_place <- nf > 0 | nc > 0 + if (!nrow(Q)) return(nsr_set_finish(out, taxon_id, db, n, endemism = endemism, + has_place = has_place)) Q[, relation := data.table::fifelse(has_fine[qid] & region_key %in% unlist(ctry), "within", relation)] Q[, taxon_id := taxon_id[qid]] Q <- Q[!is.na(taxon_id) & !is.na(region_key)] @@ -203,7 +212,8 @@ nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, O <- db$chk_dt[Q, on = c("taxon_id", "region_key"), nomatch = 0L, allow.cartesian = TRUE] # after the join, so the checklist is never copied just to drop 0.14% of its rows if (exclude_extinct && nrow(O) && "is_extinct" %in% names(O)) O <- O[is_extinct == 0L] - if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism)) + if (!nrow(O)) return(nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism, + has_place = has_place)) O[, `:=`(inc = relation %in% c("same", "within"), sub = relation == "contains", ovl = relation == "overlaps")] @@ -290,13 +300,17 @@ nsr_resolve_status_set1 <- function(taxon_id, places, db, min_overlap = 0.01, out$n_sub_native <- data.table::fifelse(is.na(res$n_sub_nat), 0L, as.integer(res$n_sub_nat)) out$n_sub_introduced <- data.table::fifelse(is.na(res$n_sub_int), 0L, as.integer(res$n_sub_int)) out$conflict <- !is.na(res$conflict_type) & res$conflict_type != "none" - nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism) + # "none" is the contract value for a result nothing disagreed about; NA here would be + # read as unknown, and the row path says "none" for the same query + out$conflict_type[!is.na(out$code) & is.na(out$conflict_type)] <- "none" + nsr_set_finish(out, taxon_id, db, n, Q, endemism = endemism, has_place = has_place) } #' Fill in the answers that need no opinions: absence, unknowns, and endemism #' @keywords internal #' @noRd -nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL, endemism = TRUE) { +nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL, endemism = TRUE, + has_place = TRUE) { dt <- function(...) data.table::data.table(...) evaluable <- !is.na(taxon_id) & taxon_id %in% db$evaluable # the regions each query consulted, for the coverage test @@ -324,6 +338,16 @@ nsr_set_finish <- function(out, taxon_id, db, n, Q = NULL, endemism = TRUE) { # nsr_reduce(). A 0 here would read as a checked negative and would differ from the # row path for the same query, on nothing but batch size. + # A taxon that was matched but has nowhere to look is not the same as a place no + # source covers, and the row path distinguishes them; the generic absence logic above + # would otherwise report this as uncovered. + noplace <- !is.na(taxon_id) & !rep_len(has_place, n) + if (any(noplace)) { + out$code[noplace] <- "UNK" + out$reason[noplace] <- "Place not matched to any region" + out$scope[noplace] <- "none" + } + # rows with no place at all, or no taxon none <- is.na(taxon_id) if (any(none)) { diff --git a/R/local_build.R b/R/local_build.R index d46c090..0de0ed0 100644 --- a/R/local_build.R +++ b/R/local_build.R @@ -99,6 +99,12 @@ NSR_local_build <- function(sources = "wcvp", # ---- assemble --------------------------------------------------------------------- keep_old <- function(x, col) if (is.null(x)) NULL else x[!x[[col]] %in% todo, , drop = FALSE] + # A cache built before extinct records were retained has no is_extinct column, and the + # rows kept from it are rbind()ed with new rows that do; give it the column first, or + # adding a source to an older cache fails on the bind. + if (!is.null(existing$checklist) && is.null(existing$checklist$is_extinct)) { + existing$checklist$is_extinct <- 0L + } chk <- lapply(names(raw), function(s) { d <- raw[[s]] fill <- is.na(d$taxon_id) @@ -162,8 +168,7 @@ NSR_local_build <- function(sources = "wcvp", raster_ok <- !inherits(r, "try-error") if (!raster_ok) { warning("The WGSRPD level-3 raster could not be built: ", - conditionMessage(attr(r, "condition")), " -", + conditionMessage(attr(r, "condition")), "\n", "The checklist tables are written, but until this raster exists no query ", "can reach a WCVP opinion: coordinates cannot be placed in a WCVP region, ", "and division names resolve to GADM, which needs the link table below. ", @@ -183,8 +188,7 @@ NSR_local_build <- function(sources = "wcvp", r <- try(nsr_build_region_links(dir = dir, quiet = quiet), silent = TRUE) if (inherits(r, "try-error")) { warning("The region link table could not be built: ", - conditionMessage(attr(r, "condition")), " -", + conditionMessage(attr(r, "condition")), "\n", "The cache is written, but queries by division name will not reach sources ", "published against another geography. Re-run NSR_local_build(overwrite = TRUE) ", "once the cause is fixed.", call. = FALSE) diff --git a/R/local_cache.R b/R/local_cache.R index ce4ad0e..b0762a4 100644 --- a/R/local_cache.R +++ b/R/local_cache.R @@ -4,11 +4,23 @@ #' political divisions NSR keys on are GNRS's, and a user who has built either of the #' others already has them. \code{options(NSR.cache_dir=)} overrides, then #' \code{options(GNRS.cache_dir=)}, then \code{tools::R_user_dir("GNRS", "cache")}. +#' +#' \code{tools::R_user_dir()} arrived in R 4.0.0, and the package supports R 3.5 for the +#' API functions. Rather than invent a different default on old R - which would put the +#' cache somewhere GNRS and GVS do not look, and silently give each package its own copy - +#' this says what is missing and how to set it. The API functions never reach here. #' @keywords internal #' @noRd nsr_cache_dir <- function(create = FALSE) { - dir <- getOption("NSR.cache_dir", - getOption("GNRS.cache_dir", tools::R_user_dir("GNRS", which = "cache"))) + fallback <- function() { + if (!is.null(getNamespace("tools")[["R_user_dir"]])) { + return(tools::R_user_dir("GNRS", which = "cache")) + } + stop("The local NSR needs R >= 4.0.0 for the shared cache location, or an explicit ", + "one: set options(GNRS.cache_dir = ) to the directory GNRS and GVS use, or pass ", + "dir = to this function.", call. = FALSE) + } + dir <- getOption("NSR.cache_dir", getOption("GNRS.cache_dir", fallback())) if (create && !dir.exists(dir)) dir.create(dir, recursive = TRUE, showWarnings = FALSE) dir } diff --git a/man/NSR_local.Rd b/man/NSR_local.Rd index b16f159..a985495 100644 --- a/man/NSR_local.Rd +++ b/man/NSR_local.Rd @@ -42,10 +42,20 @@ A data.frame, one row per input row, carrying \code{user_id} from the input } \description{ Offline Native Status Resolver. For each taxon and place, the status every built -checklist gives it, reduced to one answer. Output carries \code{\link{NSR}}'s columns, -with four added. +checklist gives it, reduced to one answer. } \details{ +Output carries \code{\link{NSR}}'s columns, and adds: \code{isEndemic}; +\code{native_status_conflict} and \code{native_status_conflict_type}, saying whether +and how the sources disagreed; \code{native_status_scope}, which of the polygons +bearing on the place produced the answer; \code{n_subpolygons_native} and +\code{n_subpolygons_introduced}; \code{native_status_opinions}, every opinion consulted; +\code{native_status_coordinates}, \code{native_status_names} and \code{place_conflict} +(see Place); \code{regions_matched}, what the place was matched on; and +\code{taxon_evaluable}, whether any source holds native status for the taxon at all - +which is what separates \code{A} (absent) from \code{UNK} (no data). The coordinates +given are echoed as \code{latitude} and \code{longitude}. + \strong{Place.} By default the answer comes from political division names (\code{country}, \code{state_province}, \code{county_parish}), which are resolved to GADM units through the GNRS backbone and carried to each source's own geography by the diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R index ac8aed3..6b211cd 100644 --- a/tests/testthat/test-local-nsr-paths.R +++ b/tests/testthat/test-local-nsr-paths.R @@ -21,7 +21,10 @@ fixture_db <- function() { "t6 wgsrpd3:QUE native wcvp 0 1", "t6 wgsrpd3:BZL native wcvp 0 0", # t7 is the contrast: native ONLY in Brazil, never recorded in Quebec at all - "t7 wgsrpd3:BZL native wcvp 0 0") + "t7 wgsrpd3:BZL native wcvp 0 0", + # t8 is native to a country row AND to one of that country's states + "t8 gadm0:CAN native vascan 0 0", + "t8 gadm1:CAN.11_1 native vascan 0 0") f <- do.call(rbind, strsplit(trimws(spec), "[[:space:]]+")) checklist <- data.frame( taxon_id = f[, 1], region_key = f[, 2], status = f[, 3], source_name = f[, 4], @@ -34,8 +37,8 @@ fixture_db <- function() { fraction = c(0.99, 0.99, 0.001), stringsAsFactors = FALSE) taxa <- data.frame( - taxon_id = paste0("t", 1:7), - species_name = paste("Sp", 1:7), + taxon_id = paste0("t", 1:8), + species_name = paste("Sp", 1:8), family = "Fam", genus = "Gen", rank = "species", stringsAsFactors = FALSE) sources <- data.frame(source_name = c("wcvp", "vascan"), is_comprehensive = TRUE, stringsAsFactors = FALSE) @@ -223,3 +226,48 @@ test_that("keys split into the finest level named and its country", { expect_equal(s$country[[2]], "gadm0:BRA") expect_equal(s$fine[[3]], character(0)) }) + +# The two paths must agree on the answers that need no opinions at all, not just on the +# ones that do. Each of these diverged. +test_that("a matched taxon with nowhere to look is not an uncovered place", { + skip_if_not_installed("data.table") + db <- fixture_db() + nowhere <- list(fine = list(character(0)), country = list(character(0))) + row <- NSR:::nsr_resolve_status("t1", nowhere, db) + set <- NSR:::nsr_resolve_status_set("t1", nowhere, db) + expect_equal(row$reason, "Place not matched to any region") + expect_equal(set$reason, row$reason) + expect_equal(set$code, row$code) + expect_equal(set$scope, row$scope) +}) + +test_that("nothing disagreeing is \"none\", not unknown", { + skip_if_not_installed("data.table") + db <- fixture_db() + # a country-only query is answered by the states inside it, so no opinion is recorded + # against the polygon itself and the conflict summary has nothing to aggregate + p <- list(fine = list(character(0)), country = list("gadm0:CAN")) + expect_equal(NSR:::nsr_resolve_status("t4", p, db)$conflict_type, "none") + expect_equal(NSR:::nsr_resolve_status_set("t4", p, db)$conflict_type, "none") +}) + +test_that("absence is read against the country a state was given with", { + skip_if_not_installed("data.table") + db <- fixture_db() + # CAN.9_1 carries no checklist rows and no links, but it was asked about as part of + # Canada, whose rows the answer already draws on - so absence there is interpretable + p <- list(fine = list("gadm1:CAN.9_1"), country = list("gadm0:CAN")) + expect_equal(NSR:::nsr_resolve_status("t3", p, db)$code, "A") + expect_equal(NSR:::nsr_resolve_status_set("t3", p, db)$code, "A") +}) + +test_that("a country counts as its own container for a confined range", { + skip_if_not_installed("data.table") + db <- fixture_db() + # t8 is native to Canada the country row and to Quebec the state; asked about Brazil + # it is absent, and its whole native range is confined to Canada + p <- list(fine = list("wgsrpd3:BZL"), country = list(character(0))) + expect_equal(NSR:::nsr_resolve_status("t8", p, db)$code, "Ie") + expect_equal(NSR:::nsr_resolve_status_set("t8", p, db)$code, "Ie") + expect_match(NSR:::nsr_resolve_status("t8", p, db)$reason, "Canada") +}) From 5cebff510aa997259b5405391ef7e9e117fb5ece Mon Sep 17 00:00:00 2001 From: Brian Maitner Date: Thu, 24 Sep 2026 13:01:13 -0400 Subject: [PATCH 8/8] Decline the service's higher-rank rule, and say why The design note took the service's taxonomic rule wholesale, including "with no opinion for the taxon, higher ranks are consulted". The resolver never did that, so the note described something the code does not do. The code is right and the note was wrong (BM). The rule's two halves are not the same claim. Upward is sound and is implemented, at build time in nsr_species_of(): a distribution recorded against an accepted infraspecific taxon is rolled up to its species, because a subspecies of X occurring natively somewhere means X occurs natively there - the subspecies is an X. About 60k WCVP rows reach their species this way. Downward from a higher rank is the ecological fallacy in its division form: knowing a genus is native to a region tells us nothing about any particular species in it. A genus native to Quebec means some species of it is native to Quebec, and the species in hand is as readily the introduced one. Native is the strongest claim the resolver makes, and making it at genus level would manufacture it at scale for precisely the taxa the sources are silent about. So a species with no opinion of its own is UNK, or A where a comprehensive source covers the place, never native-by-genus; taxon_evaluable separates them. A genus query is still answered directly from the distributions WCVP records against genera, which is answering what was asked rather than inferring downward from it. This is a permanent divergence from the live service and is recorded as one, alongside the wcvp/powo source-name difference, for validation against NSR(). The note now says so in the propagation rules, in the resolution semantics, in the milestone that claimed to implement it, and in a section of its own. NSR_local()'s help gains a Rank paragraph, since a caller could otherwise reasonably expect the service's behaviour. A test pins it: the genus answers when the genus is asked about, the species does not borrow that answer, and it is not evaluable on the strength of it. R CMD check: Status OK. Co-Authored-By: Claude Opus 5 --- R/NSR_local.R | 10 +++++++ dev_notes/01-offline-nsr-design.md | 41 +++++++++++++++++++++++---- man/NSR_local.Rd | 10 +++++++ tests/testthat/test-local-nsr-paths.R | 29 +++++++++++++++++++ 4 files changed, 85 insertions(+), 5 deletions(-) diff --git a/R/NSR_local.R b/R/NSR_local.R index 2401a1f..f830787 100644 --- a/R/NSR_local.R +++ b/R/NSR_local.R @@ -40,6 +40,16 @@ #' a list's scope, age or purpose. Disagreement is recorded, not hidden, in #' \code{native_status_conflict} and \code{native_status_opinions}. #' +#' \strong{Rank.} A distribution recorded against an accepted infraspecific taxon counts +#' for its species: a subspecies of \emph{X} native somewhere means \emph{X} is native +#' there. The reverse is not done. A species with no opinion of its own is not given its +#' genus's status, because a genus being native to a region says nothing about which of its +#' species are - it is the introduced one as readily as the native one. Such a species is +#' \code{UNK}, or \code{A} where a comprehensive source covers the place; +#' \code{taxon_evaluable} distinguishes them. A genus \emph{query} is still answered +#' directly, from the distributions WCVP records against genera. \code{\link{NSR}} does +#' consult higher ranks, so the two differ here by design. +#' #' \strong{Which polygons answer.} A place is judged by the polygons it lies IN. An #' opinion about a polygon containing the place applies to it (POWO's finest statement #' about Guadeloupe is "native in the Leeward Islands"), and among those any native diff --git a/dev_notes/01-offline-nsr-design.md b/dev_notes/01-offline-nsr-design.md index 915a124..fab93f7 100644 --- a/dev_notes/01-offline-nsr-design.md +++ b/dev_notes/01-offline-nsr-design.md @@ -37,7 +37,8 @@ division supplied. Propagation rules, from the service README: - **Taxonomic:** native propagates *up* (a native variety makes the species, genus and family native there); introduced propagates *down* (an introduced species makes its varieties introduced). With - no opinion for the taxon, higher ranks are consulted. + no opinion for the taxon, higher ranks are consulted. **We reproduce the first and decline the + last - see below.** - **Political:** native propagates *up* (native in a county implies native in its state and country); introduced propagates *down* (introduced in a country implies introduced in its states). - **Absence:** absent from a checklist gives `A`; absent while endemic elsewhere gives `Ie`. @@ -225,8 +226,10 @@ while "introduced in Brazil Southeast" answers it as `I`. Nothing is invented an For each row, after name and division resolution: -1. Gather every checklist opinion for the taxon (or a higher rank, per the taxonomic rule) in the - division (or a containing/contained one, per the political rule), from all built sources. +1. Gather every checklist opinion for the taxon in the division (or a containing/contained one, + per the political rule), from all built sources. Opinions recorded against a taxon's own + infraspecifics are already part of it, having been rolled up at build time; opinions recorded + against a *higher* rank are not consulted (see "Higher ranks are not consulted", below). 2. Reduce to one status per political level. **Precedence (BM, 2026-09-18): if any source says native, the taxon is native there.** Native is the hardest claim for a checklist to make by accident, whereas "introduced" and "absent" are often an artefact of a list's scope or age. Where @@ -451,6 +454,33 @@ sub-region's status upward; we decline to, and say so. For records with coordina in the garden pipeline - the two approaches coincide, because the point lies inside exactly one polygon per system. +### Higher ranks are not consulted (BM, 2026-09-24) + +The service's taxonomic rule has two halves and they are not the same claim. + +Upward is sound and **is** implemented, at build time in `nsr_species_of()`: a distribution recorded +against an accepted infraspecific taxon is rolled up to its species, because a subspecies of *X* +occurring natively in a region means *X* occurs natively in that region. The subspecies is an *X*. +Roughly 60k WCVP distribution rows hang off infraspecifics and reach the species this way. + +Downward from a higher rank is **not** implemented, and deliberately. "With no opinion for the +taxon, higher ranks are consulted" means answering about a species from its genus's range, and +knowing that a genus is native to a region tells us nothing about any particular species in it: +that is inferring a property of a part from a property of the whole - the ecological fallacy, in +its division form. A genus native to Quebec means *some* species of it is native to Quebec, and the +species in hand is as likely to be the introduced one. Native is the strongest claim the resolver +can make, and making it at genus level would manufacture it at scale for exactly the taxa about +which the sources are silent. + +So a species with no opinion of its own is `UNK` (or `A` where a comprehensive source covers the +place), never native-by-genus. `taxon_evaluable` says which it is. A genus *query* is still answered +directly, because WCVP records distributions against genera and those rows are kept keyed on the +genus - that is answering the question asked, not inferring downward from it. + +**This is a known and permanent divergence from the service**, and belongs with the `wcvp`/`powo` +source-name difference in any validation against `NSR()`: where the live service returns a +genus-derived status for an unlisted species, we return `UNK` or `A`. + ## Milestones 1. Cache scaffolding, `nsr-sources`, and the POWO importer from WCVP (no download); `NSR_local()` @@ -458,8 +488,9 @@ polygon per system. country-level rows. 2. USDA, VASCAN and Flora do Brasil importers, with the state-level resolution and the political propagation rules. -3. The taxonomic propagation rules, endemism (`Ie`) and conflict handling. (`is_cultivated_taxon` - was dropped - see open question 3.) +3. Taxonomic roll-up of infraspecifics, endemism (`Ie`) and conflict handling. (`is_cultivated_taxon` + was dropped - see open question 3; consulting higher ranks was declined - see "Higher ranks are + not consulted".) 4. Validation 1-3; write-up. 5. Garden re-run (validation 4) and, separately, the remaining eight sources. diff --git a/man/NSR_local.Rd b/man/NSR_local.Rd index a985495..20d5889 100644 --- a/man/NSR_local.Rd +++ b/man/NSR_local.Rd @@ -82,6 +82,16 @@ nativity is a positive claim, whereas "introduced" and "absent" are often artefa a list's scope, age or purpose. Disagreement is recorded, not hidden, in \code{native_status_conflict} and \code{native_status_opinions}. +\strong{Rank.} A distribution recorded against an accepted infraspecific taxon counts +for its species: a subspecies of \emph{X} native somewhere means \emph{X} is native +there. The reverse is not done. A species with no opinion of its own is not given its +genus's status, because a genus being native to a region says nothing about which of its +species are - it is the introduced one as readily as the native one. Such a species is +\code{UNK}, or \code{A} where a comprehensive source covers the place; +\code{taxon_evaluable} distinguishes them. A genus \emph{query} is still answered +directly, from the distributions WCVP records against genera. \code{\link{NSR}} does +consult higher ranks, so the two differ here by design. + \strong{Which polygons answer.} A place is judged by the polygons it lies IN. An opinion about a polygon containing the place applies to it (POWO's finest statement about Guadeloupe is "native in the Leeward Islands"), and among those any native diff --git a/tests/testthat/test-local-nsr-paths.R b/tests/testthat/test-local-nsr-paths.R index 6b211cd..44d12a6 100644 --- a/tests/testthat/test-local-nsr-paths.R +++ b/tests/testthat/test-local-nsr-paths.R @@ -271,3 +271,32 @@ test_that("a country counts as its own container for a confined range", { expect_equal(NSR:::nsr_resolve_status_set("t8", p, db)$code, "Ie") expect_match(NSR:::nsr_resolve_status("t8", p, db)$reason, "Canada") }) + +# Knowing a genus is native to a region says nothing about which of its species are, so +# a species with no opinion of its own is not given its genus's status. The live service +# does consult higher ranks; this is a deliberate divergence, pinned here so it cannot +# drift back in. See dev_notes/01-offline-nsr-design.md, "Higher ranks are not consulted". +test_that("a species does not inherit its genus's native status", { + chk <- data.frame(taxon_id = "g1", region_key = "wgsrpd3:QUE", status = "native", + source_name = "wcvp", is_cultivated = 0L, is_extinct = 0L, + stringsAsFactors = FALSE) + db <- NSR:::nsr_index_db(list( + sources = data.frame(source_name = "wcvp", is_comprehensive = TRUE, + stringsAsFactors = FALSE), + taxa = data.frame(taxon_id = c("g1", "s1"), species_name = c("Gen", "Gen spec"), + family = "Fam", genus = "Gen", rank = c("genus", "species"), + stringsAsFactors = FALSE), + checklist = chk, + regions = data.frame(region_key = "wgsrpd3:QUE", region_name = "Quebec", + level = "state_province", system = "wgsrpd3", + stringsAsFactors = FALSE), + links = data.frame(from_region = character(0), to_region = character(0), + relation = character(0), fraction = numeric(0), + stringsAsFactors = FALSE))) + p <- list(fine = list("wgsrpd3:QUE"), country = list(character(0))) + # the genus answers when the genus is what was asked about + expect_equal(NSR:::nsr_resolve_status("g1", p, db)$code, "Ne") + # the species does not borrow it, and is not evaluable on the strength of it + expect_equal(NSR:::nsr_resolve_status("s1", p, db)$code, "UNK") + expect_false("s1" %in% db$evaluable) +})