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.

We also need species data.

add_tropical_biomass() adds biomass to your census data.

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
)
#> ✔ Guessing dbh in [mm].
#> ℹ You may provide the dbh unit manually via the argument`dbh_unit`.
#> ℹ Wood density given in [g/cm^3].
#> ✔ Using `latitude` and `longitude` (ignoring `region`).
#> ℹ Biomass is given in [kg].
#> ✔ Adding new columns:
#>   family, genus, species, wd_level, wd_mean, wd_sd, latitude, longitude, biomass

biomass %>% 
  select(biomass, everything())
#> # A tibble: 538 x 29
#>    biomass treeID stemID tag   StemTag sp    quadrat    gx    gy MeasureID
#>      <dbl>  <int>  <int> <chr> <chr>   <chr> <chr>   <dbl> <dbl>     <int>
#>  1 2397.      858      1 0008… ""      apei… 4402     899.  42         766
#>  2 1884.     1129      1 0011… ""      quar… 4308     867. 163.        995
#>  3  264.     2143      1 0021… ""      beil… 3715     744  305.       1829
#>  4  911.     2388     10 0023… 1       lueh… 3622     724. 447.       2007
#>  5  961.     4448      1 0044… ""      sima… 2321     477. 428.       3741
#>  6 2473.     5877      1 0059… ""      quar… 1303     280.  70.4      4800
#>  7  570.     6487      1 0065… ""      alse… 1108     221. 178.       5226
#>  8    2.12   8651      1 0105… ""      hyba… 4811     974. 228.       6832
#>  9   16.0    9480      1 0114… ""      fara… 4814     977. 290        7373
#> 10    2.49  10179     11 0121… <NA>    hyba… 4819     979. 395.       7898
#> # … with 528 more rows, and 19 more variables: CensusID <int>, dbh <dbl>,
#> #   pom <chr>, hom <dbl>, ExactDate <chr>, DFstatus <chr>, codes <chr>,
#> #   nostems <dbl>, date <dbl>, status <chr>, agb <dbl>, family <chr>,
#> #   genus <chr>, species <chr>, wd_level <chr>, wd_mean <dbl>,
#> #   wd_sd <dbl>, latitude <dbl>, longitude <dbl>

propagate_errors() allows you to propagate errors.

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:

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:

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: AfricaExtraTrop, AfricaTrop, Australia, AustraliaTrop, CentralAmericaTrop, China, Europe, India, Madagascar, Mexico, NorthAmerica, Oceania, SouthEastAsia, SouthEastAsiaTrop, SouthAmericaExtraTrop, SouthAmericaTrop, and World.

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:

  1. A ForestGEO-like stem or tree table.
  2. 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.

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

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:

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()
#> Warning: Removed 15028 rows containing missing values (geom_point).

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.

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")
))
#> ✔ Guessing dbh in [mm].
#> ℹ You may provide the dbh unit manually via the argument`dbh_unit`.
#> ℹ biomass values are given in [kg].
#> ✔ Guessing dbh in [mm].
#> ℹ You may provide the dbh unit manually via the argument`dbh_unit`.

# 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)
#> Joining, by = "species"
#> # A tibble: 23 x 5
#>    species                 biomass   agb   dbh dbh_max_mm
#>    <chr>                     <dbl> <dbl> <dbl>      <dbl>
#>  1 fagus grandifolia            NA 13.7  1030.        890
#>  2 fraxinus americana           NA 14.2  1053.        550
#>  3 liriodendron tulipifera      NA  8.24 1012.        650
#>  4 liriodendron tulipifera      NA 11.2  1159.        650
#>  5 liriodendron tulipifera      NA 10.3  1118.        650
#>  6 liriodendron tulipifera      NA 10.6  1135.        650
#>  7 liriodendron tulipifera      NA  8.48 1025.        650
#>  8 liriodendron tulipifera      NA 15.9  1365.        650
#>  9 liriodendron tulipifera      NA  8.12 1006.        650
#> 10 liriodendron tulipifera      NA 11.5  1173.        650
#> 11 liriodendron tulipifera      NA 11.5  1174.        650
#> 12 liriodendron tulipifera      NA  9.02 1054         650
#> 13 liriodendron tulipifera      NA 13.9  1280.        650
#> 14 quercus alba                 NA 15.0  1018.        890
#> 15 quercus rubra                NA 27.7  1418.        890
#> 16 quercus rubra                NA 28.2  1432.        890
#> 17 quercus rubra                NA 25.5  1366.        890
#> 18 quercus rubra                NA 17.3  1143.        890
#> 19 quercus rubra                NA 21.9  1272.        890
#> 20 quercus velutina             NA 16.1  1107         890
#> 21 quercus velutina             NA 26.6  1393.        890
#> 22 quercus velutina             NA 15.6  1092.        890
#> 23 quercus velutina             NA 31.6  1511.        890

Biomass via BIOMASS versus allodb

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()
#> Warning: Removed 15028 rows containing missing values (geom_point).