We can use FIA data to estimate how forest structure has changed over time (i.e., has density increased or declined) or to estimate the average rate at which individual trees grow. Here we show how we can use the vitalRates
function in rFIA to do just that.
For this specific example, we’ll focus on estimating average annual change in individual tree basal area and cumulative tree basal area per acre, by site productivity class (SITECLCD) and species across the states of Washington and Oregon.
Before we get started, you’ll need rFIA v0.3.2 to run the examples below. You can check which version you have installed with packageVersion('rFIA')
. If you have an older version, you can update with devtools::install_github('hunter-stanke/rFIA')
. Also, if you’d rather download the script that details the analysis presented below, you can get it here.
So first, we’ll need to download some data (you can skip this step if you already have FIA data for WA and OR). Here we use getFIA
to download and save the state subsets to our computer:
library(rFIA)
library(dplyr)
## Download some FIA data from the DataMart ------------------------------------
getFIA(states = c('OR', 'WA'),
dir = 'path/to/save/FIA/data',
load = FALSE) # Download, but don't load yet
Next we’ll read our data into R with readFIA
and take a most recent subset across states with clipFIA
:
## Setting up the database - "remote" in this case, but you can read all the data
## into RAM if you want to modify columns/ etc by setting inMemory = TRUE
db <- readFIA(dir = fiaPath,
states = c('OR', 'WA'), # If you keep all your data together
nCores = cores,
inMemory = FALSE) # Set to TRUE if your computer has enough RAM
## Take the most recent subset -------------------------------------------------
db <- clipFIA(db,
mostRecent = TRUE)
Now that we have our data loaded, we can start estimating tree growth rates! First, we’ll estimate net growth rates, i.e., we will include trees that have recruited or died in our estimates of tree growth. Note that net growth can be negative if mortality exceeds recruitment and growth on live trees. This type of growth estimate is most useful if we’re interested in characterizing population-level change. For example, if we’d like to know how the average rate at which cumulative basal area per acre has changed in our population of interest, we’d shoot for net growth. Here’s how we do it with rFIA, where net growth is indicated by the argument treeType = 'all'
:
## Net annual change
net <- vitalRates(db,
treeType = 'all', # "all" indicates net growth
grpBy = SITECLCD, # Grouping by site productivity class
bySpecies = TRUE, # also grouping by species
variance = TRUE)
If we are instead interested in characterizing the average annual growth rate of individual trees, we’d most likely want to exclude stems that died or recruited into the population between plot measurements. To do this with rFIA, simply set treeType = 'live'
:
## Annual growth on live trees that remained live
live <- vitalRates(db,
treeType = 'live', # "live" excludes mortality and recruitment
grpBy = SITECLCD, # Grouping by site productivity class
bySpecies = TRUE, # also grouping by species
variance = TRUE)
By default, vitalRates
will estimate average annual DBH, basal area, biomass, and net volume growth rates for individual stems, along with average annual basal area, biomass, and net volume growth per acre. Here we’ll focus in on basal area, so let’s simplify and pretty up our tables a bit:
## Net annual change
net <- net %>%
## Dropping unnecessary columns
select(SITECLCD, COMMON_NAME, SCIENTIFIC_NAME, SPCD,
BA_GROW, BA_GROW_AC, BA_GROW_VAR, BA_GROW_AC_VAR, nPlots_TREE, N) %>%
## Dropping rows with no live trees (growth was NA)
filter(!is.na(BA_GROW)) %>%
## Making SITECLCD more informative
mutate(site = case_when(SITECLCD == 1 ~ "225+ cubic feet/acre/year",
SITECLCD == 2 ~ "165-224 cubic feet/acre/year",
SITECLCD == 3 ~ "120-164 cubic feet/acre/year",
SITECLCD == 4 ~ "85-119 cubic feet/acre/year",
SITECLCD == 5 ~ "50-84 cubic feet/acre/year",
SITECLCD == 6 ~ "20-49 cubic feet/acre/year",
SITECLCD == 7 ~ "0-19 cubic feet/acre/year")) %>%
## Arrange it nicely
select(COMMON_NAME, SCIENTIFIC_NAME, SITECLCD, site, everything()) %>%
arrange(COMMON_NAME, SCIENTIFIC_NAME, SITECLCD, site)
## Annual growth on live trees that remained live
live <- live %>%
## Dropping unnecessary columns
select(SITECLCD, COMMON_NAME, SCIENTIFIC_NAME, SPCD,
BA_GROW, BA_GROW_AC, BA_GROW_VAR, BA_GROW_AC_VAR, nPlots_TREE, N) %>%
## Dropping rows with no live trees (growth was NA)
filter(!is.na(BA_GROW)) %>%
## Making SITECLCD more informative
mutate(siteProd = case_when(SITECLCD == 1 ~ "225+ cubic feet/acre/year",
SITECLCD == 2 ~ "165-224 cubic feet/acre/year",
SITECLCD == 3 ~ "120-164 cubic feet/acre/year",
SITECLCD == 4 ~ "85-119 cubic feet/acre/year",
SITECLCD == 5 ~ "50-84 cubic feet/acre/year",
SITECLCD == 6 ~ "20-49 cubic feet/acre/year",
SITECLCD == 7 ~ "0-19 cubic feet/acre/year")) %>%
## Arrange it nicely
select(COMMON_NAME, SCIENTIFIC_NAME, SITECLCD, siteProd, everything()) %>%
arrange(COMMON_NAME, SCIENTIFIC_NAME, SITECLCD, siteProd)
Here BA_GROW
gives us annual basal area growth per tree in square feet/year, and BA_GROW_AC
gives us average basal area growth per acre in square feet/acre/year. But maybe we’d prefer units to be square centimeters instead - just remember to multiply the variance by the square of the conversion factor!
## Net annual change
net <- net %>%
## Convert to square centimeters instead of square feet
mutate(BA_GROW = BA_GROW * 929.03,
BA_GROW_AC = BA_GROW_AC * 929.03,
BA_GROW_VAR = BA_GROW_VAR * (929.03^2),
BA_GROW_AC_VAR = BA_GROW_AC_VAR * (929.03^2))
## Annual growth on live trees that remained live
live <- live %>%
## Convert to square centimeters instead of square feet
mutate(BA_GROW = BA_GROW * 929.03,
BA_GROW_AC = BA_GROW_AC * 929.03,
BA_GROW_VAR = BA_GROW_VAR * (929.03^2),
BA_GROW_AC_VAR = BA_GROW_AC_VAR * (929.03^2))