Appendix 2: Creating your own spatial data

Load the required pacakges

pacman::p_load(tidyverse, sf, rnaturalearth, patchwork, prioritizr)
source("utils-functions.R")

Create a dataset boundary

dat_bndry <- tibble(x = 100, y = seq(-50, 0, by = 1)) %>%
  bind_rows(tibble(x = seq(100, 160, by = 1), y = 0)) %>%
  bind_rows(tibble(x = 160, y = seq(0, -50, by = -1))) %>%
  bind_rows(tibble(x = seq(160, 100, by = -1), y = -50)) %>%
  as.matrix() %>%
  list() %>%
  st_polygon() %>%
  st_sfc(crs = "EPSG:4326") %>%
  st_sf() %>%
  st_make_valid()

Create planning units

dat_PUs <- st_make_grid(dat_bndry, cellsize = 2) %>%
  st_sf() %>%
  mutate(cellID = row_number()) # Add a cell ID reference

ggplot() +
  geom_sf(data = dat_PUs)

Create species data

dat_species_prob <- dat_PUs %>%
  st_sf() %>%
  mutate(
    Spp1 = runif(n = dim(.)[[1]]),
    Spp2 = runif(n = dim(.)[[1]]),
    Spp3 = runif(n = dim(.)[[1]]),
    Spp4 = runif(n = dim(.)[[1]]),
    Spp5 = runif(n = dim(.)[[1]])
  )
print(dat_species_prob)

Simple feature collection with 780 features and 6 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: 100 ymin: -50 xmax: 160 ymax: 2
Geodetic CRS:  WGS 84
First 10 features:
   cellID                       geometry        Spp1      Spp2       Spp3
1       1 POLYGON ((100 -50, 102 -50,... 0.073318975 0.9132673 0.64967977
2       2 POLYGON ((102 -50, 104 -50,... 0.232265136 0.5278280 0.71378633
3       3 POLYGON ((104 -50, 106 -50,... 0.005846191 0.1002031 0.70895919
4       4 POLYGON ((106 -50, 108 -50,... 0.049310110 0.9396838 0.96180804
5       5 POLYGON ((108 -50, 110 -50,... 0.132347638 0.3551796 0.14346341
6       6 POLYGON ((110 -50, 112 -50,... 0.927284935 0.5081837 0.40945422
7       7 POLYGON ((112 -50, 114 -50,... 0.021034951 0.6488039 0.01025041
8       8 POLYGON ((114 -50, 116 -50,... 0.990886714 0.4781844 0.93609304
9       9 POLYGON ((116 -50, 118 -50,... 0.012682000 0.4531465 0.80747452
10     10 POLYGON ((118 -50, 120 -50,... 0.742546350 0.0353536 0.39667087
         Spp4      Spp5
1  0.67805419 0.7889837
2  0.31436926 0.9522712
3  0.31235934 0.3874380
4  0.11227152 0.6601634
5  0.10585245 0.2963974
6  0.06569794 0.2242155
7  0.29624232 0.8025695
8  0.42655559 0.4251609
9  0.65161999 0.1906374
10 0.27303567 0.4910942

Convert to binary data

dat_species_bin <- dat_species_prob %>%
  as_tibble() %>%
  mutate(across(
     # Apply to all columns except geometry and cellID
    -any_of(c("cellID", "geometry")),
    ~ case_when(
      . >= 0.5 ~ 1,
      . < 0.5 ~ 0,
      is.na(.data) ~ 0
    )
  )) %>%
  st_as_sf()

Setup prioritization

Extract feature names for use in prioritization

col_name <- dat_species_bin %>%
  select(-"cellID") %>%
  st_drop_geometry() %>%
  colnames()

Add a cost layer

out_sf <- dat_species_bin %>%
  mutate(CostArea = rep(1, 780))

Create targets object

targets <- rep(0.3, length(col_name))

Create a conservation problem

dat_problem <- problem(out_sf,
                       features = col_name,
                       cost_column = "CostArea"
) %>%
  add_min_set_objective() %>%
  add_relative_targets(targets) %>%
  add_binary_decisions() %>%
  add_default_solver(verbose = FALSE)

Run prioritization

Solve and plot conservation problem

dat_soln <- dat_problem %>%
  solve.ConservationProblem()

(gg_sol <- splnr_plot_Solution(dat_soln))

Warning: Duplicated `override.aes` is ignored.