In forestgeo/fgeo.biomass: Calculate Biomass with the Best Available Allometries

knitr::opts_chunk$set(
  cache = TRUE,
  collapse = TRUE,
  comment = "#>",
  fig.path = "man/figures/README-"
)

set.seed(1)

Calculate biomass

The goal of fgeo.biomass is to calculate biomass using ForestGEO data and equations from either the BIOMASS package or the allodb package.

The BIOMASS package is applicable to tropical forests. It was first published on CRAN in 2016 and on Methods on Ecology and Evolution in 2017. fgeo.biomass provides the main features of BIOMASS with a simpler interface, consistent with all fgeo packages.
The allodb package is work in progress, and aims to provide expert-selected allometric equations, both for tropical and temperate forests. fgeo.biomass provides a simple interface to automate the process of finding the right equation(s) for each stem and computing biomass.

Installation

Install the development version of fgeo.biomass with:

# install.packages("devtools")
devtools::install_github("forestgeo/fgeo.biomass")

Setup

In addition to the fgeo.biomass package we will use dplyr and ggplot2 for data wrangling and plotting.

library(ggplot2)
library(dplyr)
library(fgeo.biomass)

fgeo.biomass wrapping BIOMASS

We'll use data from the Barro Colorado Island, Panama (BCI). We first pick alive trees and drop missing dbh values as we can't calculate biomass for them.

if (!requireNamespace("bciex", quietly = TRUE)) {
  stop(
    "For this example, you must first install the bciex package with:\n",
    "devtools::install_github('forestgeo/bciex')"
  )
}

bci_tree <- as_tibble(bciex::bci12t7mini) %>% 
  filter(status == "A", !is.na(dbh))
bci_tree

We also need species data.

bci_species <- as_tibble(bciex::bci_species)
bci_species

add_tropical_biomass() adds biomass to your census data.

biomass <- add_tropical_biomass(bci_tree, bci_species)
biomass

You may also provide a specific region or latitude and longitude.

biomass <- add_tropical_biomass(
  bci_tree, 
  bci_species,
  latitude = 9.154965, 
  longitude = -79.845884
)

biomass %>% 
  select(biomass, everything())

propagate_errors() allows you to propagate errors.

str(
  propagate_errors(biomass)
)

model_height() allows you to create a height model, which you can use to propagate height errors. This is what the entire pipeline looks like:

model <- model_height(bci_tree)

errors <- bci_tree %>% 
  add_tropical_biomass(bci_species) %>% 
  propagate_errors(height_model = model)

str(errors)

If you pass latitude and longitude to add_tropical_biomass(), and then you pass aheight_modeltopropagate_errors()`, then you will need to ignore the coordinates. On an interactive session, you should see something like this:

if (interactive()) {
  errors <- bci_tree %>% 
    add_tropical_biomass(
      bci_species, 
      latitude = 9.154965, 
      longitude = -79.845884
    ) %>% 
    propagate_errors(height_model = model)

  str(errors)
}

add_wood_density() adds wood density to your census data. It is not limited to tropical forests, and has support for all of these regions: r glue::glue_collapse(fgeo.biomass:::wd_regions(), sep = ", ", last = ", and ").

wood_density <- add_wood_density(bci_tree, bci_species)

wood_density %>% 
  select(starts_with("wd_"), everything())

The BIOMASS package provides a tool to correct taxonomic names. fgeo.biomass does not include that feature. You may use BIOMASS directly or the more focused taxize package.

fgeo.biomass wrapping allodb

Warning

These features are not ready for research. We are now building a Minimum Viable Product, with just enough features to collect feedback from alpha users and redirect our effort. The resulting biomass is still meaningless.

We'll use the add_biomass() with these inputs:

A ForestGEO-like stem or tree table.
A species table (internally used to look up the Latin species names from the species codes in the sp column of the census table).

We'll use data from the Smithsonian Conservation Biology Institute, USA (SCBI). We first pick alive trees and drop missing dbh values as we can't calculate biomass for them.

census <- fgeo.biomass::scbi_tree1 %>% 
  filter(status == "A", !is.na(dbh))

census

We now use add_biomass() to add biomass to our census dataset.

species <- fgeo.biomass::scbi_species

with_biomass <- census %>% 
  add_biomass(species, site = "SCBI")

We are warned that we are using a tree-table (as opposed to a stem-table), and informed about how to interpret the resulting biomass values for trees and shrubs.

Some equations couldn't be found. There may be two reasons:

Some stems in the data belong to species with no matching species in allodb.
Some stems in the data belong to species that do match species in allodb but the available equations were designed for a dbh range that doesn't include actual dbh values in the data.

Here are the most interesting columns of the result:

with_biomass %>% 
  select(treeID, species, biomass)

Let's now visualize the relationship between dbh and bbiomass by species (black points), in comparison with agb (above ground biomass) values calculated with allometric equations for tropical trees (grey points).

with_biomass %>% 
  # Convert agb from [Mg] to [kg]
  mutate(agb_kg = agb * 1e3) %>% 
  ggplot(aes(x = dbh)) +
  geom_point(aes(y = agb_kg), size = 1.5, color = "grey") +
  geom_point(aes(y = biomass), size = 1, color = "black") +
  facet_wrap("species", ncol = 4) +
  ylab("Reference `agb` (grey) and calculated `biomass` (black) in [kg]") +
  xlab("dbh [mm]") +
  theme_bw()

Above, the species for which biomass couldn't be calculated show no black points, although they do show grey reference-points.

To better understand the distribution of biomass values for each species we can use a box-plot.

with_biomass %>% 
  ggplot(aes(species, biomass)) +
  geom_boxplot() +
  ylab("biomass [kg]") +
  coord_flip()

For some species the maximum dbh for which biomass was calculated is much lower than the maximum dbh value for which the reference agb was calculated. This is because most equations in allodb are defined for a specific range of dbh values. Eventually allodb might provide equations beyond the dbh limits currently available.

To explore this issue, here we use add_component_biomass() which allows us to see intermediary results that add_biomass() doesn't show.

detailed_biomass <- suppressWarnings(suppressMessages(
  add_component_biomass(census, species, site = "SCBI")
))

# Maximum `dbh` values by species
max_by_species <- detailed_biomass %>% 
  select(species, dbh_max_mm) %>% 
  group_by(species) %>% 
  arrange(desc(dbh_max_mm)) %>% 
  filter(row_number() == 1L) %>% 
  ungroup()

# `dbh` is above the maximum limit, so `biomass` is missing (agb has a value)
detailed_biomass %>% 
  filter(dbh > 1000) %>% 
  select(-dbh_max_mm) %>% 
  left_join(max_by_species) %>% 
  mutate(agb_kg = agb * 1e3) %>%
  select(species, biomass, agb, dbh, dbh_max_mm) %>% 
  arrange(species) %>%
  print(n = Inf)

Biomass via BIOMASS versus allodb

temperate_biomass <- add_biomass(census, species, site = "scbi")

# Warning: Aplying tropical equations to a temperate forest for comparison
tropical_biomass <- add_tropical_biomass(census, species)

dbh_biomsss <- tibble(
  dbh = temperate_biomass$dbh,
  species = temperate_biomass$species,
  temperate_biomass = temperate_biomass$biomass, 
  tropical_biomass = tropical_biomass$biomass
)

dbh_biomsss %>% 
  ggplot(aes(x = dbh)) +
  geom_point(aes(y = tropical_biomass), size = 1.5, color = "grey") +
  geom_point(aes(y = temperate_biomass), size = 1) +
  facet_wrap("species", ncol = 4) +
  ylab("Biomass [kg] (via the BIOMASS (grey) and allodb (black) packages)") +
  xlab("dbh [mm]") +
  theme_bw()