knitr::opts_chunk$set(echo=T, message=F, warning=F)
library(prioritizr) # devtools::load_all("~/github/prioritizr")
library(bbnj)       # devtools::load_all() # devtools::install_local(force=T)
library(raster)
library(sf)
library(dplyr)
library(stringr)
library(glue)
library(here)
library(fs)
library(knitr)
library(formattable)
area   = raster::area
select = dplyr::select

rmd <- knitr::current_input(dir = T)
#rmd <- file.path(getwd(), "s01a.bio.10pct.gl.now.Rmd")
pfx <- rmd %>% path_ext_remove()

r basename(pfx)

# variables: target 30% ----
rel_target <- 0.3
redo       <- F

# planning unit: inverse of fishing profitability ----
r_fish_profits <- raster(s_fish_gfw, "mean_scaled_profits_with_subsidies") # plot(r_fish_profits)
r_pu <- r_fish_profits %>%
  gap_fill_raster() # %>%
  # rescale_raster(inverse=T) # plot(r_pu)

# biodiversity: now + 2100 ----
lyrs_bio_alltime <- names(s_bio_gmbi) %>% 
  setdiff(c("nspp_all", "nspp_na", "rls_all", "rls_na")) %>% 
  setdiff(c("nspp_all_2100", "nspp_na_2100", "rls_all_2100", "rls_na_2100"))
s_bio_alltime <- subset(s_bio_gmbi, lyrs_bio_alltime)

# features ----
s_features <- stack(
  r_vgpm,
  s_bio_alltime,
  #raster(s_fish_gfw, "mean_scaled_profits_with_subsidies") %>%
  #    gap_fill_raster() %>%
  #    rescale_raster(inverse=T),
  #raster(s_fish_ubc, "mcp_2004"),
  r_phys_seamounts,
  r_phys_vents,
  s_phys_scapes)
names(s_features) <- c(
  "bio_vgpm",
  names(s_bio_alltime),
  #"fish_profit.subs"
  #"fish_mcp.2004",
  "phys_seamounts",
  "phys_vents",
  glue("phys_scape.{1:11}"))

# problem ----
p <- problem(r_pu, s_features) %>%
  add_min_set_objective() %>% 
  add_relative_targets(rel_target)

# solve ----
s <- solve_log(p, pfx, redo=redo)
r_pu_area <- area(r_pu_id) %>%  # in km2
  mask(r_pu_id)
#tif <- "/Users/bbest/github/bbnj/inst/app/www/scenarios/s03a.biofish.30pct.gl.now_sol.tif"
tif <- glue("{pfx}_sol.tif")
r_sol <- raster(tif)
area_km2 <- cellStats(r_pu_area * r_sol, "sum")

Area ($km^2$): r format(area_km2, big.mark=",")

tbl_target_representation(glue("{pfx}_rep.csv"))


ecoquants/bbnj documentation built on April 4, 2023, 9:08 p.m.