+ "markdown": "---\ntitle: \"Bivariate maps with R\"\ndate: \"2025-2-11\"\nimage: assets/bivariate-map_bii.png\nimage-alt: bivariate raster map\ncategories: [R, cartography, spatial analysis, data visualization]\ndraft: false\n---\n\n\n\nMaps of raster (pixel) data can be very useful but typically you are limited to a displaying a single variable. Bivariate maps get around this by... \n\n\n\n::: {.cell}\n\n```{.r .cell-code}\n# vector required packages\npkgs <- c(\"terra\", \"ggplot2\", \"dplyr\", \"ggspatial\", \"rnaturalearth\",\n \"rnaturalearthdata\", \"tidyterra\", \"ggthemes\", \"extrafont\", \"magick\")\n\n#load packages\nlapply(pkgs, library, character.only = TRUE)\n\n# custom plot theme\ntheme_publication <- function(\n base_size=14,\n base_family=\"Fira Sans\"\n) {\n (ggthemes::theme_foundation(\n base_size=base_size,\n base_family=base_family)+\n theme(\n plot.title = element_text(face = \"bold\",size = rel(1.2),hjust = 0.5),\n text = element_text(),\n panel.background = element_rect(colour = NA),\n plot.background = element_rect(colour = NA),\n panel.grid.major = element_blank(),\n panel.grid.minor = element_blank(),\n panel.border = element_rect(colour = NA),\n axis.title = element_text(face = \"bold\",size = rel(0.8)),\n axis.title.y = element_text(angle=90,vjust =2),\n axis.title.x = element_text(vjust = -0.2),\n axis.text = element_text(),\n axis.line = element_line(colour=\"black\"),\n axis.ticks = element_line(),\n panel.grid.major = element_blank(),\n panel.grid.minor = element_blank(),\n legend.key = element_rect(colour = NA),\n legend.position = \"bottom\",\n legend.direction = \"horizontal\",\n legend.key.size= unit(0.5, \"cm\"),\n legend.title = element_text(face = \"bold\",\n hjust=0.5),\n plot.margin=unit(c(10,5,5,5),\"mm\"),\n strip.background=element_rect(colour=\"#f0f0f0\",fill=\"#f0f0f0\"),\n strip.text = element_text(face=\"bold\")\n )\n )\n}\n\n# load files of biodiversity intactness index in 2000 and 2020\nbio_2000 <- terra::rast(\"biodiversity_intactness_index/bii-2000_v2-1-1.tif\")\nbio_2020 <- terra::rast(\"biodiversity_intactness_index/bii-2020_v2-1-1.tif\")\n\n# subtract the 2020 layer from 2000\nbio_diff <- bio_2020 - bio_2000\n\n# get data of world map outline\nworld <- ne_countries(scale = \"medium\", returnclass = \"sv\")\n\n# remove antarctica\nworld <- world[world$admin != \"Antarctica\", ]\n\n# get extent of wdpa\next <- ext(wdpa)\n\n# reproject world outline to match crs\nworld <- project(world, bio_diff)\n\n# replace all values greater than 0 with 0\nbio_diff[bio_diff > 0] <- 0\n\n# Stack the difference layer and the 2000 layer\nbi_var <- c(bio_diff, bio_2000)\n\n# name the layers\nnames(bi_var) <- c(\"diff\", \"importance\")\n\n# convert to df\nbi_var_df <- as.data.frame(bi_var, xy = TRUE)\n\n# convert to absolute values\nbi_var_df$diff <- abs(bi_var_df$diff)\n\n# normalise the diff values to be between 0 and 1 such that values approaching\n# 0 are the lowest and the greatest negative values are the highest\nbi_var_df$diff <- scales::rescale(bi_var_df$diff, to = c(0, 1))\n\n# vector the number of classes to be used for breaks\nnum_class <- 3\n\n# classify both the difference and importance measures\nbi_var_class <- bi_class(bi_var_df,\n x = diff,\n y = importance,\n style = \"quantile\",\n dim = num_class)\n\n# Create quantile breaks\nbi_var_breaks <- bi_class_breaks(bi_var_df,\n x = diff, \n y = importance, \n style = \"quantile\",\n dim = num_class,\n dig_lab = 2,\n split = FALSE)\n\n# Set up a colour palette\npallet <- \"DkViolet2\"\n\n# Create bivariate map\nbi_var_map <- ggplot() +\n geom_raster(data = bi_var_class,\n aes(x = x,\n y = y,\n fill = bi_class),\n show.legend = FALSE)+\n bi_scale_fill(pal = pallet,\n dim = Full_num_class,\n flip_axes = FALSE,\n rotate_pal = FALSE,\n na.value = \"#F2F2F2\",\n aes(alpha = 0.5)) +\n theme_publication()+\n theme(axis.text = element_blank(),\n axis.title.x = element_blank(),\n axis.title.y = element_blank(),\n axis.ticks = element_blank(),\n panel.grid = element_blank(),\n axis.line = element_blank(),\n panel.grid.major = element_blank(),\n panel.grid.minor = element_blank(),\n panel.border = element_rect(colour = NA),\n panel.background = element_rect(fill='transparent'),\n plot.background = element_rect(fill='transparent', color=NA))\n\n# Create the legend for the bivariate map\nbi_var_legend <- bi_legend(pal = pallet, \n flip_axes = FALSE,\n rotate_pal = FALSE,\n dim = Full_num_class,\n breaks = bi_var_breaks,\n xlab = \"x\",\n ylab = \"y\",\n arrow = FALSE)+\n theme(text = element_text(size = 4, family = \"FiraSans\"),\n axis.text.x = element_text(angle = 0, hjust = 0),\n axis.text.y = element_text(angle = 0, hjust = 0),\n panel.background = element_rect(fill='transparent'),\n plot.background = element_rect(fill='transparent', color=NA))\n\n# combine using patchwork\nbi_var_map_legend <- bi_var_map +\n # legend inset on left\n inset_element(bi_var_legend,\n left = 0,\n bottom = 0,\n right = 0.2,\n top = 0.45,\n align_to = 'full')\n\n# save the map\nggsave(\"biodiversity_intactness_map.svg\",\n plot = bi_var_map_legend,\n width = 34,\n height = 19,\n dpi = 600, \n units = \"cm\",\n bg = \"transparent\")\n```\n:::",
0 commit comments