Chapter 4 Make

4.1 Load data

Location where data was saved from previous chapter (clip)

## default location:
load_location_default <- './clip_FIA/'
## custom location:
load_location_custom <- 'C:/Users/JJGross/Documents/R_projects/FIA_data/clip_FIA/' 
#set variable
load_location <- load_location_default

Load data and set as default argument variables

## Load most recent (MR) Appalachian trail (at) dataset:
at_MR <- read_rds(file = paste0(load_location, "at_FIA_MR.rds"))
## Load entire Appalachian trail (at) dataset:
at <- read_rds(file = paste0(load_location, "at_FIA.rds")) 
## Load eco
eco <- read_rds(file = paste0(load_location, "eco.rds")) 

## Set default variables for rFIA function:
database_variable <- at
polys_variable <- eco

4.2 Save location

Choose where final tables will be saved:

## default location:
save_default <- './summary_data/' 
## Alternatively, enter specific custom location:
save_custom <- 'C:/Users/JJGross/Documents/R_projects/FIA_data/summary_data/'

save_location <- save_default

## Create directory if it doesn't already exist
if (!dir.exists(save_location)) {dir.create(save_location)}

4.3 Metrics

Derive population estimates (at the plot- or the ecoregion-scale) for each forest health metric:

Live tree abundance

Species diversity of live trees

diversity

Tree vital rates

growth

Forest demographic rates

mortality

Regeneration abundance

Snags

Down woody material

downwoody

Invasive plant abundance

invasive

Stand structural stage distributions

Each of the following subsections illustrates how to use the variety of functions in rFIA to calculate these metrics.

4.3.1 Argument Variables

Are plot-level records desired within each SUBSECTION? If so, run the entire 03-make.Rmd with by_plot_variable set to TRUE. Otherwise, leave by_plot_variable <- FALSE to aggregate FIA plots at SUBSECTION level.

by_plot_variable <- FALSE means tables are summarized at SUBSECTION-level
by_plot_variable <- TRUE means tables are summarized at the FIA plot-level

by_plot_variable <- FALSE

most_recent <- FALSE means all FIA plots are included in the data tables
most_recent <- TRUE means only the most recent FIA plots are included

most_recent <- TRUE

The code chunk below automatically updates argument variables and folder save location based on above user input (blue boxes) to by_plot_variable and most_recent. Code chunk also specifies the number of processing cores to use during the execution of rFIA functions. See the rFIA webpage Tips for working with Big Data for more information.

## save tables to new folder if by_plot_variable is TRUE
if (by_plot_variable == TRUE) {
  # new save location folder
  ifelse(!dir.exists(file.path(save_location)), dir.create(file.path(save_location)), FALSE)
  save_location <- paste0(save_location, "by_plot_variable/")
  
}

## if 'most_recent == TRUE', than change database_variable to 'at_MR' database 
## and create new folder to save MR tables to. 
if (most_recent == TRUE) {
  database_variable <- at_MR
  # new save location folder
  ifelse(!dir.exists(file.path(save_location, "most_recent")), 
         dir.create(file.path(save_location, "most_recent")), FALSE)
  save_location <- paste0(save_location, "most_recent/")
}

## default number of processing cores to utilize in rFIA function arguments below.
cores_variable <- parallel::detectCores()-2

4.3.2 Population estimate method

Functions in rFIA can be used to derive population estimates of forest data using 5 unique estimators (i.e. different methods for panel combination) using the method argument in each function. For the purposes of forest population estimation for APPA, only the method = ANNUAL argument is used (i.e. no panel combination). For more information on panel combinations see the rFIA Alternative design-based estimators page.

## Set method variable for rFIA function arguments here:
method_variable <- "ANNUAL"

Note that when using method = "ANNUAL" all rFIA:: functions result with the following error message:

“Bad stratification, i.e., strata too small to compute variance of annual panels. If you are only interested in totals and/or ratio estimates, disregard this. However, if interested in variance (e.g., for confidence intervals) try using method =”TI”.”

4.4 Live tree abundance

The rFIA::tpa() function Produces tree per acre (TPA) and basal area per acre (BAA) estimates from FIA data, along with population totals for each variable. Options to group estimates by species, size class, and other variables defined in the FIADB.

Argument treeDomain = DIA \>= 5 results in estimates only utilizing tree diameters above 5 inches DBH.

4.4.1 tpa

## Trees per Acre (species and size classes lumped)
tpa <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA >= 5
) 
saveRDS(tpa, paste0(save_location, "tpa.rds"))

4.4.2 tpa_spp

## Trees per Acre by species (size-classes lumped)
tpa_spp <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySpecies = TRUE,
  bySizeClass = FALSE
) 
saveRDS(tpa_spp, paste0(save_location, "tpa_spp.rds"))

4.4.3 tpa_spp_sizecl

## Trees per Acre (by species and size-class)
tpa_spp_sizecl <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySpecies = TRUE,
  bySizeClass = TRUE
)
saveRDS(tpa_spp_sizecl, paste0(save_location, "tpa_spp_sizecl.rds"))

4.4.4 biomass

The rFIA::biomass() function produces estimates of volume (cu.ft./acre), biomass (tons/acre), and carbon (tons/acre) with options to group estimates by species, size class, and other variables defined in the FIADB.

## biomass (by species and size-class)
# biomass <- rFIA::biomass(db = database_variable,
#            byPlot = by_plot_variable,
#            polys = polys_variable, 
#            method = method_variable,
#            nCores = cores_variable,
#            returnSpatial = TRUE,
#            treeType = "live",
#            totals = TRUE,
#            treeDomain = DIA >= 5,
#            bySpecies = TRUE, 
#            bySizeClass = TRUE) 
# saveRDS(biomass, paste0(save_location, "biomass.rds"))

Note that the rFIA::biomass() function currently produces the following error and needs attention from package developers.

Error in fcase(is.na(DIA), NA_real_, !is.na(DRYBIO_WDLD_SPP), DRYBIO_WDLD_SPP/(jTotal - : object ‘DRYBIO_WDLD_SPP’ not found

4.5 Species diversity of live trees

The rFIA::diversity() function produces estimates of diversity from FIA data. Returns Shannon’s Index (H), Shannon’s Equitability (Eh), and Richness (S) for alpha (mean/SE of stands), beta, and gamma diversity. Default behavior estimates species diversity, using TPA as a state variable and Species Code (SPCD) to groups of individuals.

4.5.1 diversity

## Diversity
diversity <- rFIA::diversity(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA >= 5
)
saveRDS(diversity, paste0(save_location, "diversity.rds"))

4.6 Tree growth rates

The rFIA::vitalRates() function computes estimates of average annual DBH (inches/ yr), basal area (sq. ft./ yr), biomass (short tons/ yr), and net volume (cu. ft./yr) growth rates for individual stems, along with average annual basal area and net volume growth per acre.

To calcuate by size class include argument bySizeClass =TRUE

To estimate net growth rates (include trees that have recruited or died in estimates), use treeType = 'all' (default)

To exclude stems that died or recruited into the population between plot measurements set treeType = 'live'

4.6.1 growth

## Tree DBH growth
growth <- rFIA::vitalRates(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySpecies = TRUE,
  variance = TRUE # note error message: "Bad stratification, i.e., strata too small to compute variance of annual panels." 
)
saveRDS(growth, paste0(save_location, "growth.rds"))

4.7 Forest demographic rates

The rFIA::growMort() function estimates of annual growth, recruitment, natural mortality, and harvest rates, along with population estimates for each variable.

Recruitment events are defined as when a live stem which is less than 5 inches DBH at time 1, grows to or beyond 5 inches DBH by time 2. This does NOT include stems which grow beyond the 5-inch diameter criteria and are then subject to mortality prior to remeasurement. Natural mortality is defined as when a live stem is subject to non-harvest mortality between successive measurement periods. Finally, harvest is defined as when a live stem is cut and removed between successive measurements.

4.7.1 mortality

To estimate mortality per species:

## mortality
mortality <- rFIA::growMort(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "all", #default, includes all stems, live and dead
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySpecies = TRUE,
  variance = TRUE # note error message: "Bad stratification, i.e., strata too small to compute variance of annual panels." 
)
saveRDS(mortality, paste0(save_location, "mortality.rds"))

Column mortality$MORT_PERC results in Inf for each row. Suspected error with rFIA::growMort() because none of the tables produced had any values in MORT_PREC

4.8 Regeneration abundance

The rFIA::seedlings() and rFIA::tpa() functions are used to estimate regeneration of seedlings and sapling, respectively. Note that querying saplings can be achieved within the rFIA::tpa() function by setting the tree domain argument treeDomain = DIA < 5. Stems below 1 inch DBH are not included in output of rFIA::tpa() because FIA seedlings (<1 inch DBH) are sampled differently and therefore queried using the distinct rFIA::seedlings() function.

4.8.1 seedlings

Seedling abundance (trees per acre) is computed using the distinct rFIA::seedlings() function, not the rFIA::tpa() function.

The FIA glossary defines seedlings as: “Live trees smaller than 1.0 inch (2.5 cm) d.b.h./d.r.c. that are at least 6 inches (15.2 cm) in height for softwoods and 12-inches (30.5 cm) in height for hardwoods.”

## Trees per acre of seedlings (<1 inch DBH)
seedlings <- rFIA::seedling(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  totals = TRUE
) 
saveRDS(seedlings, paste0(save_location, "seedlings.rds"))

4.8.2 seedlings_spp

## Trees per acre of seedlings by species
seedlings_spp <- rFIA::seedling(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  totals = TRUE,
  bySpecies = TRUE # by species
) 
saveRDS(seedlings_spp, paste0(save_location, "seedlings_spp.rds"))

4.8.3 saplings

Saplings must use rFIA::tpa() function with treeDomain = DIA < 5 argument To estimate sapling size class (1 to 4.9 inch DBH).

## Trees per acre of saplings (trees greater than 1 inch DBH, and less than 5 inch DBH)
saplings <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA < 5
) 
saveRDS(saplings, paste0(save_location, "saplings.rds"))

4.8.4 saplings_spp

## Trees per acre of saplings by species 
saplings_spp <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA < 5,
  bySpecies = TRUE, # by species
  bySizeClass = FALSE
) 
saveRDS(saplings_spp, paste0(save_location, "saplings_spp.rds"))

4.8.5 saplings_sizecl

## Trees per acre of saplings by size class 
saplings_sizecl <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "live",
  totals = TRUE,
  treeDomain = DIA < 5,
  bySpecies = FALSE,
  bySizeClass = TRUE # by size class
) 
saveRDS(saplings_sizecl, paste0(save_location, "saplings_sizecl.rds"))

4.9 Snags

rFIA::tpa() and rFIA::biomass() can be used to estimate snag abundance, volume, and percentage.

4.9.1 snag_abundance

treeType = "dead" bySizeClass = TRUE

## Snags per acre
snag_abundance <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "dead",
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySizeClass = TRUE
)

saveRDS(snag_abundance, paste0(save_location, "snag_abundance.rds"))

4.9.2 dead_and_live_tpa

treeType = "all" includes both live and dead trees in calculations. Use grpBy = STATUSCD to avoid lumping Live and Dead tree.

STATUSCD (From The FIA Database User Guide)
STATUSCD	Description
0	No status - Tree is not presently in the sample (remeasurement plots only). Tree was incorrectly tallied at the previous inventory, currently not tallied due to definition or procedural change, or is not tallied because it is located on a nonsampled condition (e.g., hazardous or denied). RECONCILECD = 5-9 required for remeasured annual inventory data but not for periodic inventory data.
1	Live tree
2	Dead tree
3	Removed - Cut and removed by direct human activity related to harvesting, silviculture or land clearing. This tree is assumed to be utilized.

## dead and live trees - to calculate percent snags
dead_and_live_tpa <- rFIA::tpa(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  treeType = "all",
  grpBy = STATUSCD,
  totals = TRUE,
  treeDomain = DIA >= 5,
  bySizeClass = TRUE
)

dead_and_live_tpa <- dead_and_live_tpa |>
  filter(STATUSCD == 1 | STATUSCD == 2) |> #only interested in live and dead 
  mutate(tree_status = case_when(STATUSCD == 1 ~ "Live",
                                 STATUSCD == 2 ~ "Dead",
                                 #anything else will result in NA
                                 TRUE ~ NA)) 
  
saveRDS(dead_and_live_tpa, paste0(save_location, "dead_and_live_tpa.rds"))

4.9.3 snag_volume

# ## Snag Volume 
# snag_volume <- rFIA::biomass(
#   db = database_variable,
#   byPlot = by_plot_variable,
#   polys = polys_variable,
#   method = method_variable,
#   nCores = cores_variable,
#   returnSpatial = TRUE,
#   treeType = "dead",
#   totals = TRUE,
#   treeDomain = DIA >= 5 
# )
# saveRDS(snag_volume, paste0(save_location, "snag_volume.rds"))

Note that the rFIA::biomass() function currently produces the following error and needs attention from package developers.

Error in fcase(is.na(DIA), NA_real_, !is.na(DRYBIO_WDLD_SPP), DRYBIO_WDLD_SPP/(jTotal - : object ‘DRYBIO_WDLD_SPP’ not found

4.10 Down woody material

rFIA::dwm() produces estimates of down woody material stocks. Estimates are returned by fuel class (duff, litter, 1HR, 10HR, 100HR, 1000HR, piles) for application in fuels management.

1HR fuels: small, fine woody debris 10HR fuels: medium, fine woody debris 100HR fuels: large, fine woody debris 1000HR fuels: coarse woody debris and slash piles duff: O horizon; all unidentifiable organic material above mineral soil, beneath litter litter: identifiable plant material which is downed and smaller than 10HR fuel class (1HR class includes standing herbaceous material).

4.10.1 downwoody

## Down woody material 
downwoody <- rFIA::dwm(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  totals = TRUE
)
saveRDS(downwoody, paste0(save_location, "downwoody.rds"))

4.11 Invasive plant abundance

The rFIA::invasive() function produces estimates of the areal coverage (%) of invasive species and frequency of plots invasive species were detected in.

4.11.1 invasive

## invasive
invasive <- rFIA::invasive(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  totals = TRUE,
  variance = TRUE # note error message: "Bad stratification, i.e., strata too small to compute variance of annual panels." 
)
saveRDS(invasive, paste0(save_location, "invasive.rds"))

4.12 Stand structural stage distributions

The rFIA::standStruct() function estimates forest structural stage distributions as percent forested land area in pole, mature, late, and mosaic stages and returns the stand structural stage distribution of an area of forest/timberland from FIA data.

Description: Estimates of forest structural stage distributions as percent forested land area in pole, mature, late, and mosaic stages
- Diameter Classes:
  - Pole: 11 - 25.9 cm
  - Mature: 26 - 45.9 cm
  - Large: 46+ cm
- Structural Stage Classification:
  - Pole Stage: > 67% BA in pole and mature classes, with more BA in pole than mature.
  - Mature Stage: > 67% BA in pole and mature classes, with more BA in mature than pole OR > 67% BA in mature and large classes, with more BA in mature.
  - Late-Successional Stage:: > 67% BA in mature and large classes, with more in large.
  - Mosiac:: Any plot not meeting above criteria.

4.12.1 ss

## stand structural stage
ss <- rFIA::standStruct(
  db = database_variable,
  byPlot = by_plot_variable,
  polys = polys_variable,
  method = method_variable,
  nCores = cores_variable,
  returnSpatial = TRUE,
  totals = TRUE
)
saveRDS(ss, paste0(save_location, "ss.rds"))