Nothing
Gap analysis workflow for multiple species"
In GapAnalysis: Conservation Indicators Using Spatial Information
knitr::opts_chunk$set(echo = TRUE)
# Check for pandoc and knitr status
has_pandoc <- rmarkdown::pandoc_available()
save_outputs <- !isTRUE(getOption('knitr.in.progress')) && has_pandoc
This workflow demonstrates how to import your own point data (occurrences) and
SDMs (rasters), check all inputs, and run the full gap analysis for multiple
species using the GapAnalysis R package. We provide a general approach using
a for-loop for processing and storing results for multiple species. Example
Cucurbita data is provided within the package library.
Note: This vignette recommends R version 4.5 or above. Windows users are also recommended to have Rtools 4.5 or above installed.
For each species, the workflow:
- Runs SRSex on all records.
- Checks inputs (occurrences, rasters, protected areas, ecoregions).
- Runs gap analysis metrics.
- Stores results as a table to an output excel file.
- Creates a conservation priority summary figure as an output png file.
- Generates interactive leaflet maps for both ex-situ and in-situ conservation as output HTML files.
At script end, the user is prompted to save the workflow by knitting to output HTML.
Download GapAnalysis repository (ALTERNATIVE)
Copy the GapAnalysis GitHub repository locally as an alternative method to
source functions and load example data if the GapAnalysis R package does not
install properly in the step above. See the
GitHub resource page
on copying a repository.
Easiest method: Download and extract the
GapAnalysis repository.
Select the green [Code] button on the main page to download as a ZIP.
You can also clone the repository if you have a GitHub account.
Data preparation
If you downloaded the repository locally, source functions and load data by
removing the hashes below. Your specific path will differ based on your working
directory and the relative path to the downloaded repository.
# ── Load libraries ─────────────────────────────────────────────────────────────
library(GapAnalysis)
library(dplyr)
library(tidyr)
library(stringr)
library(terra)
library(leaflet)
library(openxlsx)
library(ggplot2)
library(ggtext)
library(cowplot)
library(htmltools)
library(htmlwidgets)
library(kableExtra)
# ── Set outputs folder ────────────────────────────────────────────────────────
# NOTE: This path must match the path set in the knit: field in the YAML header
# All files will be saved here. Update both locations if changing the path, e.g.:
#outputs_folder <- "C:/Users/you/Desktop/GapAnalysis/outputs"
outputs_folder <- "outputs"
if (!dir.exists(outputs_folder)) dir.create(outputs_folder)
# ── Source functions from GapAnalysis repo (alternative method) ───────────────
#files <- list.files("R", pattern = "\\.[rR]$", full.names = TRUE)
#for (f in files) source(f)
# ── Load data from GapAnalysis repo (alternative method) ─────────────────────
#load("data/CucurbitaData.rda") # Occurrence data
#load("data/CucurbitaRasts.rda") # Raster list
#load("data/ProtectedAreas.rda") # Protected areas raster
#load("data/ecoregions.rda") # Ecoregions data
# ── Load your own data ────────────────────────────────────────────────────────
# NOTE: To use your own occurrence data and SDM rasters, remove the hashes
# below and update the file paths.
#occData <- read.csv("path/to/your_occurrences.csv") # Local occurrence data
#sdms <- list( # Local rasters read in as list
# rast("path/to/species1_sdm.tif"),
# rast("path/to/species2_sdm.tif"))
# NOTE: To download protected areas and ecoregions data locally, run getDatasets().
# This requires an active internet connection.
#getDatasets() # Local download of protected areas and ecoregions data
# ── Load example Cucurbita data (GapAnalysis library) ────────────────────────
# NOTE: To use your own data, hash out the two lines below.
occData <- CucurbitaData
sdms <- terra::unwrap(CucurbitaRasts)
# ── Prepare spatial layers ────────────────────────────────────────────────────
ecos <- terra::vect(ecoregions)
proAreas <- terra::unwrap(ProtectedAreas)
# Standardize CRS across all spatial layers to prevent mismatch warnings
target_crs <- "+proj=longlat +datum=WGS84"
ecos <- terra::project(ecos, target_crs)
proAreas <- terra::project(proAreas, target_crs)
# ── Prepare species list and SDM rasters ──────────────────────────────────────
taxa <- unique(occData$species)
sdmList <- lapply(seq_along(taxa), function(i) {
r <- sdms[[i]]
if (terra::crs(r, proj = TRUE) != terra::crs(ecos, proj = TRUE)) {
r <- terra::project(r, terra::crs(ecos))
}
r <- terra::subst(r, 0, NA)
r[!(is.na(r) | r == 1)] <- NA
r
})
names(sdmList) <- taxa
# Confirm species names
message("\nConfirm species for gap analysis:\n", paste(taxa, collapse = "\n"))
Run gap analysis
Example gap analysis workflow using a named list and for-loop to process and
store results for all species.
For details on each metric calculation, review the function documentation
in the GapAnalysis R folder.
results_list <- list()
for (i in seq_along(taxa)) {
taxon <- taxa[i]
occurrenceData <- occData[occData$species == taxon, ]
sdm <- sdms[[i]]
if (terra::crs(sdm, proj = TRUE) != terra::crs(ecos, proj = TRUE)) {
sdm <- terra::project(sdm, terra::crs(ecos))
}
sdm <- terra::subst(sdm, 0, NA)
sdm[!(is.na(sdm) | sdm == 1)] <- NA
# 1. Run SRSex first (on all records)
srsex <- SRSex(taxon = taxon, occurrenceData = occurrenceData)
# 2. Check inputs
occurrences <- checkOccurrences(csv = occurrenceData, taxon = taxon)
if (!inherits(sdm, "SpatRaster")) stop("sdm is not a SpatRaster")
sdm_checked <- checksdm(sdm)
proArea_cropped <- terra::crop(proAreas, terra::ext(sdm_checked))
proArea_checked <- checkProtectedAreas(protectedArea = proArea_cropped,
sdm = sdm_checked)
# Skip species if SDM raster has no values
if (all(is.na(terra::values(sdm_checked)))) {
message("Skipping ", taxon, ": SDM raster has no values.")
next
}
# 3. Generate G buffers
gBuffer <- generateGBuffers(
taxon = taxon,
occurrenceData = occurrences$data,
bufferDistM = 50000
)
# 4. Ex-situ analysis
grsex <- GRSex(taxon = taxon, sdm = sdm_checked, gBuffer = gBuffer)
ersex <- ERSex(
taxon = taxon,
sdm = sdm_checked,
occurrenceData = occurrences$data,
gBuffer = gBuffer,
ecoregions = ecos,
idColumn = "ECO_NAME",
limitByPoints = FALSE
)
fcsex <- FCSex(taxon = taxon, srsex = srsex, grsex = grsex, ersex = ersex)
# 5. In-situ analysis
srsin <- SRSin(
taxon = taxon,
sdm = sdm_checked,
occurrenceData = occurrences$data,
protectedAreas = proArea_checked
)
grsin <- GRSin(taxon = taxon, sdm = sdm_checked, protectedAreas = proArea_checked)
ersin <- ERSin(
taxon = taxon,
sdm = sdm_checked,
occurrenceData = occurrences$data,
protectedAreas = proArea_checked,
ecoregions = ecos,
idColumn = "ECO_NAME",
limitByPoints = FALSE
)
fcsin <- FCSin(taxon = taxon, srsin = srsin, grsin = grsin, ersin = ersin)
fcsmean <- FCSc_mean(taxon = taxon, fcsin = fcsin, fcsex = fcsex)
# 6. Store results
results_list[[taxon]] <- list(
srsex = srsex,
occurrences = occurrences,
sdm_checked = sdm_checked,
proArea_checked = proArea_checked,
eco_checked = ecoregions,
grsex = grsex,
ersex = ersex,
fcsex = fcsex,
srsin = srsin,
grsin = grsin,
ersin = ersin,
fcsin = fcsin,
fcsmean = fcsmean
)
}
Summary metrics tables
Results are organized as a named list, allowing review of individual species
metrics. The tables below summarize the main metrics for all species.
# ── Ex-situ metrics ───────────────────────────────────────────────────────────
exsitu_table <- do.call(rbind, lapply(names(results_list), function(taxon) {
r <- results_list[[taxon]]
data.frame(
Species = taxon,
SRS_exsitu = r$srsex$`SRS exsitu`,
GRS_exsitu = r$fcsex$`GRS exsitu`,
ERS_exsitu = r$fcsex$`ERS exsitu`,
FCS_exsitu = r$fcsex$`FCS exsitu`,
FCS_exsitu_score = r$fcsex$`FCS exsitu score`
)
}))
# ── In-situ metrics ───────────────────────────────────────────────────────────
insitu_table <- do.call(rbind, lapply(names(results_list), function(taxon) {
r <- results_list[[taxon]]
data.frame(
Species = taxon,
SRS_insitu = r$fcsin$`SRS insitu`,
GRS_insitu = r$fcsin$`GRS insitu`,
ERS_insitu = r$fcsin$`ERS insitu`,
FCS_insitu = r$fcsin$`FCS insitu`,
FCS_insitu_score = r$fcsin$`FCS insitu score`
)
}))
# ── Combined metrics ──────────────────────────────────────────────────────────
combined_table <- do.call(rbind, lapply(names(results_list), function(taxon) {
r <- results_list[[taxon]]
data.frame(
Species = taxon,
FCSc_min = r$fcsmean$FCSc_min,
FCSc_max = r$fcsmean$FCSc_max,
FCSc_mean = r$fcsmean$FCSc_mean,
FCSc_min_class = r$fcsmean$FCSc_min_class,
FCSc_max_class = r$fcsmean$FCSc_max_class,
FCSc_mean_class = r$fcsmean$FCSc_mean_class
)
}))
knitr::kable(exsitu_table, digits = 2, caption = "<span style='text-align:left; display:block;'>Ex-situ gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling()
knitr::kable(insitu_table, digits = 2, caption = "<span style='text-align:left; display:block;'>In-situ gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling()
knitr::kable(combined_table, digits = 2, caption = "<span style='text-align:left; display:block;'>Combined gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling()
plot_data <- do.call(rbind, lapply(names(results_list), function(taxon) {
r <- results_list[[taxon]]
data.frame(
ID = taxon,
SRSex = as.numeric(r$srsex[["SRS exsitu"]]),
GRSex = as.numeric(r$grsex$results[["GRS exsitu"]]),
ERSex = as.numeric(r$ersex$results[["ERS exsitu"]]),
FCSex = as.numeric(r$fcsex[["FCS exsitu"]]),
SRSin = as.numeric(r$srsin$results[["SRS insitu"]]),
GRSin = as.numeric(r$grsin$results[["GRS insitu"]]),
ERSin = as.numeric(r$ersin$results[["ERS insitu"]]),
FCSin = as.numeric(r$fcsin[["FCS insitu"]]),
FCSc_mean = as.numeric(r$fcsmean[["FCSc_mean"]])
)
}))
header_font_size <- 7
legend_label_font_mm <- 5
point_stroke <- 1.5
other_size <- 4.0
fcs_size <- 6.0
fcsc_size <- 10.5
all_species <- sort(unique(plot_data$ID))
format_id <- function(id) {
id |>
str_replace_all(" var\\. ", "</i> var. <i>") |>
str_replace_all(" x ", "</i> x <i>") |>
(\(x) paste0("<i>", x, "</i>"))()
}
formatted_species_levels <- format_id(all_species)
plot_data <- plot_data |>
mutate(
ID_italic = format_id(.data$ID),
ID_ordered = factor(.data$ID_italic, levels = rev(formatted_species_levels))
)
long_df <- plot_data |>
pivot_longer(
cols = c(SRSex, GRSex, ERSex, FCSex, SRSin, GRSin, ERSin, FCSin,
FCSc_mean),
names_to = "metric_raw",
values_to = "x"
) |>
filter(!is.na(.data$x)) |>
mutate(
metric_base = case_when(
str_detect(.data$metric_raw, "FCSc") ~ "FCSc_mean",
TRUE ~ str_remove(.data$metric_raw, "(ex|in)$")
),
origin_type = case_when(
str_detect(.data$metric_raw, "in$") ~ "in situ",
str_detect(.data$metric_raw, "ex$") ~ "ex situ",
TRUE ~ "combined"
),
point_size_dynamic = case_when(
.data$metric_base == "FCSc_mean" ~ fcsc_size,
.data$metric_base == "FCS" ~ fcs_size,
TRUE ~ other_size
),
plot_order = factor(.data$metric_base,
levels = c("FCSc_mean", "FCS", "SRS", "GRS", "ERS"))
) |>
arrange(.data$plot_order)
color_map <- c("FCSc_mean" = "#FF0000", "FCS" = "#000000",
"SRS" = "#0000FF", "GRS" = "#800080", "ERS" = "#228B22")
shape_map <- c("combined" = 18, "ex situ" = 16, "in situ" = 17)
panel_colors <- c("#ffb4b3", "#ffd380", "#ffff80", "#a8d2a8")
build_panel <- function(data, id_levels) {
num_levels <- length(id_levels)
header_y <- num_levels + 1.25
ymax <- num_levels + 2.0
ylim_top <- num_levels + 2.2
header_data <- tibble(
x = c(12.5, 37.5, 62.5, 87.5),
label = c("<b>Urgent<br>Priority</b>", "<b>High<br>Priority</b>",
"<b>Medium<br>Priority</b>", "<b>Low<br>Priority</b>")
)
ggplot(data, aes(x = .data$x, y = .data$ID_ordered)) +
annotate("rect",
xmin = c(0, 25, 50, 75), xmax = c(25, 50, 75, 100),
ymin = 0.5, ymax = ymax,
fill = panel_colors, alpha = 0.9) +
geom_segment(data = data.frame(y = seq_len(num_levels)),
aes(x = 0, xend = 100, y = y, yend = y),
color = "grey75", linewidth = 0.5, inherit.aes = FALSE) +
geom_segment(data = data.frame(x = c(0, 25, 50, 75, 100)),
aes(x = x, xend = x, y = 0.5, yend = ymax),
color = "black", linewidth = 1.5, inherit.aes = FALSE) +
geom_hline(yintercept = c(0.5, num_levels + 0.5, num_levels + 2.0),
color = "black", linewidth = 1.5) +
geom_jitter(aes(color = .data$metric_base, shape = .data$origin_type,
size = .data$point_size_dynamic),
stroke = point_stroke, width = 0.7, height = 0) +
geom_richtext(data = header_data,
aes(x = .data$x, y = header_y, label = .data$label),
size = header_font_size, color = "black",
fill = NA, label.color = NA,
lineheight = 0.9, hjust = 0.5, inherit.aes = FALSE) +
scale_color_manual(values = color_map, guide = "none") +
scale_shape_manual(values = shape_map, guide = "none") +
scale_size_identity() +
scale_x_continuous(breaks = c(0, 25, 50, 75, 100), expand = c(0, 0)) +
scale_y_discrete(drop = FALSE, limits = rev(id_levels)) +
coord_cartesian(xlim = c(0, 100), ylim = c(1, ylim_top), clip = "off") +
theme_minimal_grid() +
theme(
panel.grid = element_blank(),
axis.title = element_blank(),
axis.ticks.y = element_blank(),
axis.ticks.x = element_blank(),
axis.text.x = element_text(size = 16),
axis.text.y = element_markdown(size = 18, hjust = 1,
margin = margin(r = 15)),
plot.margin = margin(t = 30, r = 20, b = 20, l = 40)
)
}
main_plot <- build_panel(long_df, formatted_species_levels)
legend_data <- tribble(
~label, ~x, ~y, ~shape, ~color, ~symbol_size, ~label_size, ~vjust,
"<b>Conservation<br>Status</b>", 0.5, 1.8, NA_real_, "black", NA_real_, legend_label_font_mm, 1,
"FCSc-mean", 3.0, 1.7, 18, "#FF0000", fcsc_size, legend_label_font_mm, 0.5,
"FCS ex situ", 4.8, 1.7, 16, "black", fcs_size, legend_label_font_mm, 0.5,
"SRS ex situ", 6.6, 1.7, 16, "#0000FF", other_size, legend_label_font_mm, 0.5,
"GRS ex situ", 8.4, 1.7, 16, "#800080", other_size, legend_label_font_mm, 0.5,
"ERS ex situ", 10.2, 1.7, 16, "#228B22", other_size, legend_label_font_mm, 0.5,
"FCS in situ", 4.8, 1.3, 17, "black", fcs_size, legend_label_font_mm, 0.5,
"SRS in situ", 6.6, 1.3, 17, "#0000FF", other_size, legend_label_font_mm, 0.5,
"GRS in situ", 8.4, 1.3, 17, "#800080", other_size, legend_label_font_mm, 0.5,
"ERS in situ", 10.2, 1.3, 17, "#228B22", other_size, legend_label_font_mm, 0.5
)
legend_plot <- ggplot(legend_data, aes(x = .data$x, y = .data$y)) +
geom_point(aes(shape = .data$shape, color = .data$color,
size = .data$symbol_size), stroke = point_stroke) +
geom_richtext(aes(label = .data$label, x = .data$x + 0.15,
size = .data$label_size, vjust = .data$vjust),
hjust = 0, fill = NA, label.color = NA) +
scale_color_identity() +
scale_shape_identity() +
scale_size_identity() +
theme_void() +
coord_cartesian(xlim = c(0, 11.5), ylim = c(0.5, 2.2))
n_taxa <- length(formatted_species_levels)
plot_height <- max(10, n_taxa * 0.5)
final_plot <- plot_grid(main_plot, legend_plot,
ncol = 1, rel_heights = c(1, 0.1))
save_plot(
file.path(outputs_folder, "gapAnalysis_conservationPriority_figure.png"),
final_plot,
base_width = 12,
base_height = plot_height,
dpi = 300,
bg = "white")
message("Saved: ",
file.path(outputs_folder, "gapAnalysis_conservationPriority_figure.png"))
# ── Build table (1 row per species, all metrics as columns) ───────
all_metrics <- do.call(rbind, lapply(names(results_list), function(taxon) {
r <- results_list[[taxon]]
clean <- function(x) {
df <- as.data.frame(t(unlist(x)))
df <- df[, !grepl("(?i)^(taxon|species)$", names(df)), drop = FALSE]
return(df)
}
row <- cbind(
data.frame(Taxon = taxon),
clean(r$srsex),
clean(r$grsex$results),
clean(r$ersex$results),
clean(r$fcsex),
clean(r$srsin$results),
clean(r$grsin$results),
clean(r$ersin$results),
clean(r$fcsin),
clean(r$fcsmean)
)
row <- row[, !duplicated(names(row)), drop = FALSE]
return(row)
}))
# ── Round numeric columns ─────────────────────────────────────────────────────
for (col in names(all_metrics)) {
num <- suppressWarnings(as.numeric(as.character(all_metrics[[col]])))
if (all(is.na(num) == is.na(all_metrics[[col]]))) {
all_metrics[[col]] <- round(num, 2)
}
}
# ── Clean column names ────────────────────────────────────────────────────────
names(all_metrics) <- gsub("\\.", " ", names(all_metrics))
insert_pos <- which(names(all_metrics) == "Area in protected areas km2")
all_metrics <- cbind(
all_metrics[, seq_len(insert_pos), drop = FALSE],
"Area of model km2" = all_metrics[["Area of model km2"]],
all_metrics[, (insert_pos + 1):ncol(all_metrics), drop = FALSE]
)
names(all_metrics) <- dplyr::recode(names(all_metrics),
"Total records in SDM" = "Total records in model",
"FCSc_min" = "FCS minimum",
"FCSc_max" = "FCS maximum",
"FCSc_mean" = "FCSc mean",
"FCSc_min_class" = "FCS minimum priority category",
"FCSc_max_class" = "FCS maximum priority category",
"FCSc_mean_class" = "FCS mean priority category",
"Ecoregions with protected areas" = "Ecoregions within protected areas"
)
# ── Reorder columns ───────────────────────────────────────────────────────────
col_order <- c(
"Taxon", "Total records", "Total with coordinates",
"Total G records", "G records with coordinates",
"Total H records", "H records with coordinates",
"SRS exsitu", "G buffer areas in model km2", "Area of model km2",
"GRS exsitu", "Ecoregions within G buffer", "Ecoregions with records",
"ERS exsitu", "FCS exsitu", "FCS exsitu score",
"Total Observations", "Records in Protected areas", "Total records in model",
"SRS insitu", "Area in protected areas km2", "Area of model km2",
"GRS insitu", "Ecoregions within protected areas", "Ecoregions within model",
"ERS insitu", "FCS insitu", "FCS insitu score",
"FCS minimum", "FCS maximum", "FCSc mean",
"FCS minimum priority category", "FCS maximum priority category",
"FCS mean priority category"
)
names(all_metrics) <- make.unique(names(all_metrics), sep = "__TEMP__")
col_order_unique <- make.unique(col_order, sep = "__TEMP__")
all_metrics <- all_metrics[, col_order_unique, drop = FALSE]
names(all_metrics) <- gsub("__TEMP__\\d+$", "", names(all_metrics))
# ── Write to Excel ────────────────────────────────────────────────────────────
wb <- createWorkbook()
addWorksheet(wb, "gapAnalysis_results")
writeData(wb, "gapAnalysis_results", all_metrics, startRow = 1, startCol = 1)
addStyle(wb, "gapAnalysis_results",
style = createStyle(wrapText = TRUE, textDecoration = "bold"),
rows = 1,
cols = seq_len(ncol(all_metrics)),
gridExpand = TRUE
)
setRowHeights(wb, "gapAnalysis_results", rows = 1, heights = 60)
max_width <- max(nchar(as.character(all_metrics$Taxon))) + 2
setColWidths(wb, "gapAnalysis_results", cols = 1, widths = max_width)
saveWorkbook(wb,
file = file.path(outputs_folder, "gapAnalysis_results.xlsx"),
overwrite = TRUE
)
message("Saved: ", file.path(outputs_folder, "gapAnalysis_results.xlsx"))
clean_map <- function(map, occurrenceData = NULL, taxon = NULL) {
legend_call <- Filter(
function(c) !is.null(c$method) && c$method == "addLegend",
map$x$calls
)
labels <- if (length(legend_call) > 0) legend_call[[1]]$args[[1]][["labels"]] else NULL
is_insitu <- !is.null(labels) && any(grepl("Protected", labels))
if (is_insitu) {
layer_groups <- c("All ecos", "Eco gaps", "Distribution", "Protected areas")
} else {
layer_groups <- c("All ecos", "Eco gaps", "Distribution", "G buffer")
}
poly_idx <- 0
map$x$calls <- lapply(map$x$calls, function(call) {
if (!is.null(call$method) && call$method == "addPolygons") {
poly_idx <<- poly_idx + 1
if (poly_idx <= 2) call$args[[3]] <- layer_groups[poly_idx]
}
call
})
raster_idx <- 0
map$x$calls <- lapply(map$x$calls, function(call) {
if (!is.null(call$method) && call$method == "addRasterImage") {
raster_idx <<- raster_idx + 1
call$args[[4]] <- layer_groups[2 + raster_idx]
}
call
})
map$x$calls <- Filter(function(call) {
!(!is.null(call$method) && call$method == "addCircleMarkers")
}, map$x$calls)
map$x$calls <- Filter(function(call) {
!(!is.null(call$method) && call$method == "addLegend")
}, map$x$calls)
map$x$calls <- Filter(function(call) {
if (!is.null(call$method) && call$method == "addControl") {
html <- call$args[[1]]
if (is.character(html) && grepl("Ecoregions", html, ignore.case = TRUE)) {
return(FALSE)
}
}
TRUE
}, map$x$calls)
if (!is.null(occurrenceData) && !is.null(taxon)) {
occ <- occurrenceData[occurrenceData$species == taxon, ]
h_pts <- occ[occ$type == "H", ]
g_pts <- occ[occ$type == "G", ]
if (nrow(h_pts) > 0) {
map <- map |>
leaflet::addCircleMarkers(
lng = h_pts$longitude,
lat = h_pts$latitude,
radius = 2,
color = "#1184d4",
opacity = 1,
group = "H points (herbarium)"
)
}
if (nrow(g_pts) > 0) {
map <- map |>
leaflet::addCircleMarkers(
lng = g_pts$longitude,
lat = g_pts$latitude,
radius = 2,
color = "#6300f0",
opacity = 1,
group = "G points (genebank)"
)
}
}
all_groups <- c(layer_groups, "H points (herbarium)", "G points (genebank)")
map <- map |>
leaflet::addLayersControl(
overlayGroups = all_groups,
options = layersControlOptions(collapsed = FALSE)
)
map <- htmlwidgets::onRender(map, "
function(el, x) {
setTimeout(function() {
el.querySelectorAll('.leaflet-control-layers-overlays').forEach(function(s) {
var prev = s.previousElementSibling;
if (prev) prev.remove();
});
el.querySelectorAll('.leaflet-bottom.leaflet-left .leaflet-control')
.forEach(function(c) { c.remove(); });
}, 300);
}
")
map
}
for (taxon in names(results_list)) {
cat("###", taxon, "\n\n")
print(htmltools::tagList(
clean_map(
map = results_list[[taxon]]$ersex$map,
occurrenceData = occData,
taxon = taxon
)
))
cat("\n\n***\n\n")
print(htmltools::tagList(
clean_map(
map = results_list[[taxon]]$ersin$map,
occurrenceData = occData,
taxon = taxon
)
))
cat("\n\n***\n\n")
}
outputs_abs <- normalizePath(outputs_folder)
tmp_lib <- file.path(outputs_abs, "lib_tmp")
for (taxon in names(results_list)) {
clean_name <- tolower(gsub(" ", "_", taxon))
# Ex-situ map
ersex_map <- clean_map(
map = results_list[[taxon]]$ersex$map,
occurrenceData = occData,
taxon = taxon
)
if (!is.null(ersex_map)) {
file_name <- paste0(clean_name, "_exSitu_map.html")
htmlwidgets::saveWidget(
widget = ersex_map,
file = file.path(outputs_abs, file_name),
selfcontained = TRUE,
libdir = tmp_lib
)
if (dir.exists(tmp_lib)) unlink(tmp_lib, recursive = TRUE)
message("Saved: ", file.path(outputs_folder, file_name))
} else {
message("No ex-situ map available for: ", taxon)
}
# In-situ map
ersin_map <- clean_map(
map = results_list[[taxon]]$ersin$map,
occurrenceData = occData,
taxon = taxon
)
if (!is.null(ersin_map)) {
file_name <- paste0(clean_name, "_inSitu_map.html")
htmlwidgets::saveWidget(
widget = ersin_map,
file = file.path(outputs_abs, file_name),
selfcontained = TRUE,
libdir = tmp_lib
)
if (dir.exists(tmp_lib)) unlink(tmp_lib, recursive = TRUE)
message("Saved: ", file.path(outputs_folder, file_name))
} else {
message("No in-situ map available for: ", taxon)
}
}
message(
"\n",
"──────────────────────────────────────────────────────────────────────\n",
" WORKFLOW COMPLETE. WELL DONE!\n",
"──────────────────────────────────────────────────────────────────────\n",
"\n",
" >>> TO SAVE WORKFLOW AS HTML: Click the [Knit] button in RStudio <<<\n",
" at the top of the script.\n",
"\n",
" The HTML will save to the folder specified in the\n",
" knit: field of the YAML header.\n",
"\n",
" To change the save location, update the path in both:\n",
" - the knit: field in the YAML header\n",
" - the outputs_folder variable in Data Preparation\n",
"\n",
" Review your outputs folder for all saved files.\n",
"───────────────────────────────────────────────────────────────────────"
)
Try the GapAnalysis package in your browser
Any scripts or data that you put into this service are public.
GapAnalysis documentation built on July 24, 2026, 1:07 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
knitr::opts_chunk$set(echo = TRUE) # Check for pandoc and knitr status has_pandoc <- rmarkdown::pandoc_available() save_outputs <- !isTRUE(getOption('knitr.in.progress')) && has_pandoc
This workflow demonstrates how to import your own point data (occurrences) and SDMs (rasters), check all inputs, and run the full gap analysis for multiple species using the GapAnalysis R package. We provide a general approach using a for-loop for processing and storing results for multiple species. Example Cucurbita data is provided within the package library.
Note: This vignette recommends R version 4.5 or above. Windows users are also recommended to have Rtools 4.5 or above installed.
For each species, the workflow:
- Runs SRSex on all records.
- Checks inputs (occurrences, rasters, protected areas, ecoregions).
- Runs gap analysis metrics.
- Stores results as a table to an output excel file.
- Creates a conservation priority summary figure as an output png file.
- Generates interactive leaflet maps for both ex-situ and in-situ conservation as output HTML files.
At script end, the user is prompted to save the workflow by knitting to output HTML.
Download GapAnalysis repository (ALTERNATIVE)
Copy the GapAnalysis GitHub repository locally as an alternative method to source functions and load example data if the GapAnalysis R package does not install properly in the step above. See the GitHub resource page on copying a repository.
Easiest method: Download and extract the GapAnalysis repository. Select the green [Code] button on the main page to download as a ZIP.
You can also clone the repository if you have a GitHub account.
Data preparation
If you downloaded the repository locally, source functions and load data by removing the hashes below. Your specific path will differ based on your working directory and the relative path to the downloaded repository.
# ── Load libraries ───────────────────────────────────────────────────────────── library(GapAnalysis) library(dplyr) library(tidyr) library(stringr) library(terra) library(leaflet) library(openxlsx) library(ggplot2) library(ggtext) library(cowplot) library(htmltools) library(htmlwidgets) library(kableExtra) # ── Set outputs folder ──────────────────────────────────────────────────────── # NOTE: This path must match the path set in the knit: field in the YAML header # All files will be saved here. Update both locations if changing the path, e.g.: #outputs_folder <- "C:/Users/you/Desktop/GapAnalysis/outputs" outputs_folder <- "outputs" if (!dir.exists(outputs_folder)) dir.create(outputs_folder) # ── Source functions from GapAnalysis repo (alternative method) ─────────────── #files <- list.files("R", pattern = "\\.[rR]$", full.names = TRUE) #for (f in files) source(f) # ── Load data from GapAnalysis repo (alternative method) ───────────────────── #load("data/CucurbitaData.rda") # Occurrence data #load("data/CucurbitaRasts.rda") # Raster list #load("data/ProtectedAreas.rda") # Protected areas raster #load("data/ecoregions.rda") # Ecoregions data # ── Load your own data ──────────────────────────────────────────────────────── # NOTE: To use your own occurrence data and SDM rasters, remove the hashes # below and update the file paths. #occData <- read.csv("path/to/your_occurrences.csv") # Local occurrence data #sdms <- list( # Local rasters read in as list # rast("path/to/species1_sdm.tif"), # rast("path/to/species2_sdm.tif")) # NOTE: To download protected areas and ecoregions data locally, run getDatasets(). # This requires an active internet connection. #getDatasets() # Local download of protected areas and ecoregions data # ── Load example Cucurbita data (GapAnalysis library) ──────────────────────── # NOTE: To use your own data, hash out the two lines below. occData <- CucurbitaData sdms <- terra::unwrap(CucurbitaRasts) # ── Prepare spatial layers ──────────────────────────────────────────────────── ecos <- terra::vect(ecoregions) proAreas <- terra::unwrap(ProtectedAreas) # Standardize CRS across all spatial layers to prevent mismatch warnings target_crs <- "+proj=longlat +datum=WGS84" ecos <- terra::project(ecos, target_crs) proAreas <- terra::project(proAreas, target_crs) # ── Prepare species list and SDM rasters ────────────────────────────────────── taxa <- unique(occData$species) sdmList <- lapply(seq_along(taxa), function(i) { r <- sdms[[i]] if (terra::crs(r, proj = TRUE) != terra::crs(ecos, proj = TRUE)) { r <- terra::project(r, terra::crs(ecos)) } r <- terra::subst(r, 0, NA) r[!(is.na(r) | r == 1)] <- NA r }) names(sdmList) <- taxa # Confirm species names message("\nConfirm species for gap analysis:\n", paste(taxa, collapse = "\n"))
Run gap analysis
Example gap analysis workflow using a named list and for-loop to process and store results for all species.
For details on each metric calculation, review the function documentation in the GapAnalysis R folder.
results_list <- list() for (i in seq_along(taxa)) { taxon <- taxa[i] occurrenceData <- occData[occData$species == taxon, ] sdm <- sdms[[i]] if (terra::crs(sdm, proj = TRUE) != terra::crs(ecos, proj = TRUE)) { sdm <- terra::project(sdm, terra::crs(ecos)) } sdm <- terra::subst(sdm, 0, NA) sdm[!(is.na(sdm) | sdm == 1)] <- NA # 1. Run SRSex first (on all records) srsex <- SRSex(taxon = taxon, occurrenceData = occurrenceData) # 2. Check inputs occurrences <- checkOccurrences(csv = occurrenceData, taxon = taxon) if (!inherits(sdm, "SpatRaster")) stop("sdm is not a SpatRaster") sdm_checked <- checksdm(sdm) proArea_cropped <- terra::crop(proAreas, terra::ext(sdm_checked)) proArea_checked <- checkProtectedAreas(protectedArea = proArea_cropped, sdm = sdm_checked) # Skip species if SDM raster has no values if (all(is.na(terra::values(sdm_checked)))) { message("Skipping ", taxon, ": SDM raster has no values.") next } # 3. Generate G buffers gBuffer <- generateGBuffers( taxon = taxon, occurrenceData = occurrences$data, bufferDistM = 50000 ) # 4. Ex-situ analysis grsex <- GRSex(taxon = taxon, sdm = sdm_checked, gBuffer = gBuffer) ersex <- ERSex( taxon = taxon, sdm = sdm_checked, occurrenceData = occurrences$data, gBuffer = gBuffer, ecoregions = ecos, idColumn = "ECO_NAME", limitByPoints = FALSE ) fcsex <- FCSex(taxon = taxon, srsex = srsex, grsex = grsex, ersex = ersex) # 5. In-situ analysis srsin <- SRSin( taxon = taxon, sdm = sdm_checked, occurrenceData = occurrences$data, protectedAreas = proArea_checked ) grsin <- GRSin(taxon = taxon, sdm = sdm_checked, protectedAreas = proArea_checked) ersin <- ERSin( taxon = taxon, sdm = sdm_checked, occurrenceData = occurrences$data, protectedAreas = proArea_checked, ecoregions = ecos, idColumn = "ECO_NAME", limitByPoints = FALSE ) fcsin <- FCSin(taxon = taxon, srsin = srsin, grsin = grsin, ersin = ersin) fcsmean <- FCSc_mean(taxon = taxon, fcsin = fcsin, fcsex = fcsex) # 6. Store results results_list[[taxon]] <- list( srsex = srsex, occurrences = occurrences, sdm_checked = sdm_checked, proArea_checked = proArea_checked, eco_checked = ecoregions, grsex = grsex, ersex = ersex, fcsex = fcsex, srsin = srsin, grsin = grsin, ersin = ersin, fcsin = fcsin, fcsmean = fcsmean ) }
Summary metrics tables
Results are organized as a named list, allowing review of individual species metrics. The tables below summarize the main metrics for all species.
# ── Ex-situ metrics ─────────────────────────────────────────────────────────── exsitu_table <- do.call(rbind, lapply(names(results_list), function(taxon) { r <- results_list[[taxon]] data.frame( Species = taxon, SRS_exsitu = r$srsex$`SRS exsitu`, GRS_exsitu = r$fcsex$`GRS exsitu`, ERS_exsitu = r$fcsex$`ERS exsitu`, FCS_exsitu = r$fcsex$`FCS exsitu`, FCS_exsitu_score = r$fcsex$`FCS exsitu score` ) })) # ── In-situ metrics ─────────────────────────────────────────────────────────── insitu_table <- do.call(rbind, lapply(names(results_list), function(taxon) { r <- results_list[[taxon]] data.frame( Species = taxon, SRS_insitu = r$fcsin$`SRS insitu`, GRS_insitu = r$fcsin$`GRS insitu`, ERS_insitu = r$fcsin$`ERS insitu`, FCS_insitu = r$fcsin$`FCS insitu`, FCS_insitu_score = r$fcsin$`FCS insitu score` ) })) # ── Combined metrics ────────────────────────────────────────────────────────── combined_table <- do.call(rbind, lapply(names(results_list), function(taxon) { r <- results_list[[taxon]] data.frame( Species = taxon, FCSc_min = r$fcsmean$FCSc_min, FCSc_max = r$fcsmean$FCSc_max, FCSc_mean = r$fcsmean$FCSc_mean, FCSc_min_class = r$fcsmean$FCSc_min_class, FCSc_max_class = r$fcsmean$FCSc_max_class, FCSc_mean_class = r$fcsmean$FCSc_mean_class ) })) knitr::kable(exsitu_table, digits = 2, caption = "<span style='text-align:left; display:block;'>Ex-situ gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling() knitr::kable(insitu_table, digits = 2, caption = "<span style='text-align:left; display:block;'>In-situ gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling() knitr::kable(combined_table, digits = 2, caption = "<span style='text-align:left; display:block;'>Combined gap analysis metrics</span>", format = "html") %>% kableExtra::kable_styling()
plot_data <- do.call(rbind, lapply(names(results_list), function(taxon) { r <- results_list[[taxon]] data.frame( ID = taxon, SRSex = as.numeric(r$srsex[["SRS exsitu"]]), GRSex = as.numeric(r$grsex$results[["GRS exsitu"]]), ERSex = as.numeric(r$ersex$results[["ERS exsitu"]]), FCSex = as.numeric(r$fcsex[["FCS exsitu"]]), SRSin = as.numeric(r$srsin$results[["SRS insitu"]]), GRSin = as.numeric(r$grsin$results[["GRS insitu"]]), ERSin = as.numeric(r$ersin$results[["ERS insitu"]]), FCSin = as.numeric(r$fcsin[["FCS insitu"]]), FCSc_mean = as.numeric(r$fcsmean[["FCSc_mean"]]) ) })) header_font_size <- 7 legend_label_font_mm <- 5 point_stroke <- 1.5 other_size <- 4.0 fcs_size <- 6.0 fcsc_size <- 10.5 all_species <- sort(unique(plot_data$ID)) format_id <- function(id) { id |> str_replace_all(" var\\. ", "</i> var. <i>") |> str_replace_all(" x ", "</i> x <i>") |> (\(x) paste0("<i>", x, "</i>"))() } formatted_species_levels <- format_id(all_species) plot_data <- plot_data |> mutate( ID_italic = format_id(.data$ID), ID_ordered = factor(.data$ID_italic, levels = rev(formatted_species_levels)) ) long_df <- plot_data |> pivot_longer( cols = c(SRSex, GRSex, ERSex, FCSex, SRSin, GRSin, ERSin, FCSin, FCSc_mean), names_to = "metric_raw", values_to = "x" ) |> filter(!is.na(.data$x)) |> mutate( metric_base = case_when( str_detect(.data$metric_raw, "FCSc") ~ "FCSc_mean", TRUE ~ str_remove(.data$metric_raw, "(ex|in)$") ), origin_type = case_when( str_detect(.data$metric_raw, "in$") ~ "in situ", str_detect(.data$metric_raw, "ex$") ~ "ex situ", TRUE ~ "combined" ), point_size_dynamic = case_when( .data$metric_base == "FCSc_mean" ~ fcsc_size, .data$metric_base == "FCS" ~ fcs_size, TRUE ~ other_size ), plot_order = factor(.data$metric_base, levels = c("FCSc_mean", "FCS", "SRS", "GRS", "ERS")) ) |> arrange(.data$plot_order) color_map <- c("FCSc_mean" = "#FF0000", "FCS" = "#000000", "SRS" = "#0000FF", "GRS" = "#800080", "ERS" = "#228B22") shape_map <- c("combined" = 18, "ex situ" = 16, "in situ" = 17) panel_colors <- c("#ffb4b3", "#ffd380", "#ffff80", "#a8d2a8") build_panel <- function(data, id_levels) { num_levels <- length(id_levels) header_y <- num_levels + 1.25 ymax <- num_levels + 2.0 ylim_top <- num_levels + 2.2 header_data <- tibble( x = c(12.5, 37.5, 62.5, 87.5), label = c("<b>Urgent<br>Priority</b>", "<b>High<br>Priority</b>", "<b>Medium<br>Priority</b>", "<b>Low<br>Priority</b>") ) ggplot(data, aes(x = .data$x, y = .data$ID_ordered)) + annotate("rect", xmin = c(0, 25, 50, 75), xmax = c(25, 50, 75, 100), ymin = 0.5, ymax = ymax, fill = panel_colors, alpha = 0.9) + geom_segment(data = data.frame(y = seq_len(num_levels)), aes(x = 0, xend = 100, y = y, yend = y), color = "grey75", linewidth = 0.5, inherit.aes = FALSE) + geom_segment(data = data.frame(x = c(0, 25, 50, 75, 100)), aes(x = x, xend = x, y = 0.5, yend = ymax), color = "black", linewidth = 1.5, inherit.aes = FALSE) + geom_hline(yintercept = c(0.5, num_levels + 0.5, num_levels + 2.0), color = "black", linewidth = 1.5) + geom_jitter(aes(color = .data$metric_base, shape = .data$origin_type, size = .data$point_size_dynamic), stroke = point_stroke, width = 0.7, height = 0) + geom_richtext(data = header_data, aes(x = .data$x, y = header_y, label = .data$label), size = header_font_size, color = "black", fill = NA, label.color = NA, lineheight = 0.9, hjust = 0.5, inherit.aes = FALSE) + scale_color_manual(values = color_map, guide = "none") + scale_shape_manual(values = shape_map, guide = "none") + scale_size_identity() + scale_x_continuous(breaks = c(0, 25, 50, 75, 100), expand = c(0, 0)) + scale_y_discrete(drop = FALSE, limits = rev(id_levels)) + coord_cartesian(xlim = c(0, 100), ylim = c(1, ylim_top), clip = "off") + theme_minimal_grid() + theme( panel.grid = element_blank(), axis.title = element_blank(), axis.ticks.y = element_blank(), axis.ticks.x = element_blank(), axis.text.x = element_text(size = 16), axis.text.y = element_markdown(size = 18, hjust = 1, margin = margin(r = 15)), plot.margin = margin(t = 30, r = 20, b = 20, l = 40) ) } main_plot <- build_panel(long_df, formatted_species_levels) legend_data <- tribble( ~label, ~x, ~y, ~shape, ~color, ~symbol_size, ~label_size, ~vjust, "<b>Conservation<br>Status</b>", 0.5, 1.8, NA_real_, "black", NA_real_, legend_label_font_mm, 1, "FCSc-mean", 3.0, 1.7, 18, "#FF0000", fcsc_size, legend_label_font_mm, 0.5, "FCS ex situ", 4.8, 1.7, 16, "black", fcs_size, legend_label_font_mm, 0.5, "SRS ex situ", 6.6, 1.7, 16, "#0000FF", other_size, legend_label_font_mm, 0.5, "GRS ex situ", 8.4, 1.7, 16, "#800080", other_size, legend_label_font_mm, 0.5, "ERS ex situ", 10.2, 1.7, 16, "#228B22", other_size, legend_label_font_mm, 0.5, "FCS in situ", 4.8, 1.3, 17, "black", fcs_size, legend_label_font_mm, 0.5, "SRS in situ", 6.6, 1.3, 17, "#0000FF", other_size, legend_label_font_mm, 0.5, "GRS in situ", 8.4, 1.3, 17, "#800080", other_size, legend_label_font_mm, 0.5, "ERS in situ", 10.2, 1.3, 17, "#228B22", other_size, legend_label_font_mm, 0.5 ) legend_plot <- ggplot(legend_data, aes(x = .data$x, y = .data$y)) + geom_point(aes(shape = .data$shape, color = .data$color, size = .data$symbol_size), stroke = point_stroke) + geom_richtext(aes(label = .data$label, x = .data$x + 0.15, size = .data$label_size, vjust = .data$vjust), hjust = 0, fill = NA, label.color = NA) + scale_color_identity() + scale_shape_identity() + scale_size_identity() + theme_void() + coord_cartesian(xlim = c(0, 11.5), ylim = c(0.5, 2.2)) n_taxa <- length(formatted_species_levels) plot_height <- max(10, n_taxa * 0.5) final_plot <- plot_grid(main_plot, legend_plot, ncol = 1, rel_heights = c(1, 0.1)) save_plot( file.path(outputs_folder, "gapAnalysis_conservationPriority_figure.png"), final_plot, base_width = 12, base_height = plot_height, dpi = 300, bg = "white") message("Saved: ", file.path(outputs_folder, "gapAnalysis_conservationPriority_figure.png"))
# ── Build table (1 row per species, all metrics as columns) ─────── all_metrics <- do.call(rbind, lapply(names(results_list), function(taxon) { r <- results_list[[taxon]] clean <- function(x) { df <- as.data.frame(t(unlist(x))) df <- df[, !grepl("(?i)^(taxon|species)$", names(df)), drop = FALSE] return(df) } row <- cbind( data.frame(Taxon = taxon), clean(r$srsex), clean(r$grsex$results), clean(r$ersex$results), clean(r$fcsex), clean(r$srsin$results), clean(r$grsin$results), clean(r$ersin$results), clean(r$fcsin), clean(r$fcsmean) ) row <- row[, !duplicated(names(row)), drop = FALSE] return(row) })) # ── Round numeric columns ───────────────────────────────────────────────────── for (col in names(all_metrics)) { num <- suppressWarnings(as.numeric(as.character(all_metrics[[col]]))) if (all(is.na(num) == is.na(all_metrics[[col]]))) { all_metrics[[col]] <- round(num, 2) } } # ── Clean column names ──────────────────────────────────────────────────────── names(all_metrics) <- gsub("\\.", " ", names(all_metrics)) insert_pos <- which(names(all_metrics) == "Area in protected areas km2") all_metrics <- cbind( all_metrics[, seq_len(insert_pos), drop = FALSE], "Area of model km2" = all_metrics[["Area of model km2"]], all_metrics[, (insert_pos + 1):ncol(all_metrics), drop = FALSE] ) names(all_metrics) <- dplyr::recode(names(all_metrics), "Total records in SDM" = "Total records in model", "FCSc_min" = "FCS minimum", "FCSc_max" = "FCS maximum", "FCSc_mean" = "FCSc mean", "FCSc_min_class" = "FCS minimum priority category", "FCSc_max_class" = "FCS maximum priority category", "FCSc_mean_class" = "FCS mean priority category", "Ecoregions with protected areas" = "Ecoregions within protected areas" ) # ── Reorder columns ─────────────────────────────────────────────────────────── col_order <- c( "Taxon", "Total records", "Total with coordinates", "Total G records", "G records with coordinates", "Total H records", "H records with coordinates", "SRS exsitu", "G buffer areas in model km2", "Area of model km2", "GRS exsitu", "Ecoregions within G buffer", "Ecoregions with records", "ERS exsitu", "FCS exsitu", "FCS exsitu score", "Total Observations", "Records in Protected areas", "Total records in model", "SRS insitu", "Area in protected areas km2", "Area of model km2", "GRS insitu", "Ecoregions within protected areas", "Ecoregions within model", "ERS insitu", "FCS insitu", "FCS insitu score", "FCS minimum", "FCS maximum", "FCSc mean", "FCS minimum priority category", "FCS maximum priority category", "FCS mean priority category" ) names(all_metrics) <- make.unique(names(all_metrics), sep = "__TEMP__") col_order_unique <- make.unique(col_order, sep = "__TEMP__") all_metrics <- all_metrics[, col_order_unique, drop = FALSE] names(all_metrics) <- gsub("__TEMP__\\d+$", "", names(all_metrics)) # ── Write to Excel ──────────────────────────────────────────────────────────── wb <- createWorkbook() addWorksheet(wb, "gapAnalysis_results") writeData(wb, "gapAnalysis_results", all_metrics, startRow = 1, startCol = 1) addStyle(wb, "gapAnalysis_results", style = createStyle(wrapText = TRUE, textDecoration = "bold"), rows = 1, cols = seq_len(ncol(all_metrics)), gridExpand = TRUE ) setRowHeights(wb, "gapAnalysis_results", rows = 1, heights = 60) max_width <- max(nchar(as.character(all_metrics$Taxon))) + 2 setColWidths(wb, "gapAnalysis_results", cols = 1, widths = max_width) saveWorkbook(wb, file = file.path(outputs_folder, "gapAnalysis_results.xlsx"), overwrite = TRUE ) message("Saved: ", file.path(outputs_folder, "gapAnalysis_results.xlsx"))
clean_map <- function(map, occurrenceData = NULL, taxon = NULL) { legend_call <- Filter( function(c) !is.null(c$method) && c$method == "addLegend", map$x$calls ) labels <- if (length(legend_call) > 0) legend_call[[1]]$args[[1]][["labels"]] else NULL is_insitu <- !is.null(labels) && any(grepl("Protected", labels)) if (is_insitu) { layer_groups <- c("All ecos", "Eco gaps", "Distribution", "Protected areas") } else { layer_groups <- c("All ecos", "Eco gaps", "Distribution", "G buffer") } poly_idx <- 0 map$x$calls <- lapply(map$x$calls, function(call) { if (!is.null(call$method) && call$method == "addPolygons") { poly_idx <<- poly_idx + 1 if (poly_idx <= 2) call$args[[3]] <- layer_groups[poly_idx] } call }) raster_idx <- 0 map$x$calls <- lapply(map$x$calls, function(call) { if (!is.null(call$method) && call$method == "addRasterImage") { raster_idx <<- raster_idx + 1 call$args[[4]] <- layer_groups[2 + raster_idx] } call }) map$x$calls <- Filter(function(call) { !(!is.null(call$method) && call$method == "addCircleMarkers") }, map$x$calls) map$x$calls <- Filter(function(call) { !(!is.null(call$method) && call$method == "addLegend") }, map$x$calls) map$x$calls <- Filter(function(call) { if (!is.null(call$method) && call$method == "addControl") { html <- call$args[[1]] if (is.character(html) && grepl("Ecoregions", html, ignore.case = TRUE)) { return(FALSE) } } TRUE }, map$x$calls) if (!is.null(occurrenceData) && !is.null(taxon)) { occ <- occurrenceData[occurrenceData$species == taxon, ] h_pts <- occ[occ$type == "H", ] g_pts <- occ[occ$type == "G", ] if (nrow(h_pts) > 0) { map <- map |> leaflet::addCircleMarkers( lng = h_pts$longitude, lat = h_pts$latitude, radius = 2, color = "#1184d4", opacity = 1, group = "H points (herbarium)" ) } if (nrow(g_pts) > 0) { map <- map |> leaflet::addCircleMarkers( lng = g_pts$longitude, lat = g_pts$latitude, radius = 2, color = "#6300f0", opacity = 1, group = "G points (genebank)" ) } } all_groups <- c(layer_groups, "H points (herbarium)", "G points (genebank)") map <- map |> leaflet::addLayersControl( overlayGroups = all_groups, options = layersControlOptions(collapsed = FALSE) ) map <- htmlwidgets::onRender(map, " function(el, x) { setTimeout(function() { el.querySelectorAll('.leaflet-control-layers-overlays').forEach(function(s) { var prev = s.previousElementSibling; if (prev) prev.remove(); }); el.querySelectorAll('.leaflet-bottom.leaflet-left .leaflet-control') .forEach(function(c) { c.remove(); }); }, 300); } ") map }
for (taxon in names(results_list)) { cat("###", taxon, "\n\n") print(htmltools::tagList( clean_map( map = results_list[[taxon]]$ersex$map, occurrenceData = occData, taxon = taxon ) )) cat("\n\n***\n\n") print(htmltools::tagList( clean_map( map = results_list[[taxon]]$ersin$map, occurrenceData = occData, taxon = taxon ) )) cat("\n\n***\n\n") }
outputs_abs <- normalizePath(outputs_folder) tmp_lib <- file.path(outputs_abs, "lib_tmp") for (taxon in names(results_list)) { clean_name <- tolower(gsub(" ", "_", taxon)) # Ex-situ map ersex_map <- clean_map( map = results_list[[taxon]]$ersex$map, occurrenceData = occData, taxon = taxon ) if (!is.null(ersex_map)) { file_name <- paste0(clean_name, "_exSitu_map.html") htmlwidgets::saveWidget( widget = ersex_map, file = file.path(outputs_abs, file_name), selfcontained = TRUE, libdir = tmp_lib ) if (dir.exists(tmp_lib)) unlink(tmp_lib, recursive = TRUE) message("Saved: ", file.path(outputs_folder, file_name)) } else { message("No ex-situ map available for: ", taxon) } # In-situ map ersin_map <- clean_map( map = results_list[[taxon]]$ersin$map, occurrenceData = occData, taxon = taxon ) if (!is.null(ersin_map)) { file_name <- paste0(clean_name, "_inSitu_map.html") htmlwidgets::saveWidget( widget = ersin_map, file = file.path(outputs_abs, file_name), selfcontained = TRUE, libdir = tmp_lib ) if (dir.exists(tmp_lib)) unlink(tmp_lib, recursive = TRUE) message("Saved: ", file.path(outputs_folder, file_name)) } else { message("No in-situ map available for: ", taxon) } }
message( "\n", "──────────────────────────────────────────────────────────────────────\n", " WORKFLOW COMPLETE. WELL DONE!\n", "──────────────────────────────────────────────────────────────────────\n", "\n", " >>> TO SAVE WORKFLOW AS HTML: Click the [Knit] button in RStudio <<<\n", " at the top of the script.\n", "\n", " The HTML will save to the folder specified in the\n", " knit: field of the YAML header.\n", "\n", " To change the save location, update the path in both:\n", " - the knit: field in the YAML header\n", " - the outputs_folder variable in Data Preparation\n", "\n", " Review your outputs folder for all saved files.\n", "───────────────────────────────────────────────────────────────────────" )
Try the GapAnalysis package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.