In this vignette, we demonstrate the use of the
BerkeleyForestsAnalytics package. We provide an example of
how you might use the BFA package to compile pre- and post-burn field
data. Through this example, we highlight some key components of the BFA
package:
The vignette is not a replacement for the README file, which covers the inputs and outputs of each function in detail. Additionally, the README gives detailed background information and references for the methods used in the package. We recommend that you review the README prior to or in conjunction with the vignette (Find README here).
To begin, we’ll load the required packages:
This vignette uses data from the Fire and Fire Surrogate (FFS) Study. In brief, FFS is an experimental study that was designed to evaluate the impacts of fire-only (prescribed fire), mechanical-only (mechanical thinning from below followed by mastication), and mechanical + fire (mechanical thinning from below followed by mastication followed by prescribed fire) treatments on forest structure, ecological function, and future fire behavior (Read more about the FFS Study here).
The data used in this vignette are from the fire-only (i.e., prescribed fire) stands. We used data from two time periods: before treatment (2001) and one year after treatment (2003). Note that the data were slightly modified for demonstrations purposes. Therefore, the outputs should not be taken to be actual findings from the FFS Study.
The tree data have the following columns:
id time:site combinedtime pre (pre-burn) or post (post-burn)site compartment (60, 340, or 400)plot plot in which the individual tree was
measuredexp_factor stems per hectarestatus live (1) or dead (0)decay decay class. 1-5 for standing dead trees. 0 for
live trees.species species of the individual tree, using
four-letter species codesdbh diameter at breast height in centimetersht total tree height in metersThe surface and ground fuels data have the following columns:
time pre (pre-burn) or post (post-burn)site compartment (60, 340, or 400)plot plot in which the individual transect was
measuredtransect azimuth of transect on which the fuel data
were collectedcount_1h count of 1-hour fuelscount_10h count of 10-hour fuelscount_100h count of 100-hour fuelslength_1h length of the sampling transect for 1-hour
fuels in meterslength_10h length of the sampling transect for 10-hour
fuels in meterslength_100h length of the sampling transect for
100-hour fuels in meterslength_1000h length of the sampling transect for
1000-hour fuels in metersssd_S sum-of-squared-diameters for sound 1000-hour
fuelsssd_R sum-of-squared-diameters for rotten 1000-hour
fuelslitter_depth litter depth in centimetersduff_depth duff depth in centimetersslope slope along the transect in percentFirst, we’ll use the SummaryBiomass() function to get
above-ground stem, bark, and branch tree biomass at the plot level.
Let’s investigate the input dataframe:
# Note that the example data used in this vignette is included with the package
# which is why we do not have to read in the data
head(vign_trees_1)## id time site plot exp_factor status decay species dbh ht
## 1 post_340 post 340 103 24.69 1 0 QUKE 85.3 35.2
## 2 post_340 post 340 103 24.69 1 0 QUKE 71.4 31.1
## 3 post_340 post 340 103 24.69 1 0 ABCO 34.3 22.5
## 4 post_340 post 340 103 24.69 1 0 ABCO 18.3 15.7
## 5 post_340 post 340 103 24.69 1 0 ABCO 52.8 33.0
## 6 post_340 post 340 103 24.69 1 0 ABCO 19.8 17.5Attempt 1: Now, let’s try using
SummaryBiomass(). We’ll keep the defaults for sp_codes (=
“4letter”), units (= “metric”), and results (= “by_plot):
tree_bio <- SummaryBiomass(data = vign_trees_1,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
decay_class = "decay",
species = "species",
dbh = "dbh",
ht = "ht")## Error in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : There are plots with a recorded expansion factor of 0, but with more than one row.
## Plots with no trees should be represented by a single row with site and plot filled in as appropriate and an exp_factor of 0.And we get an error message. It looks like there is an improper use of a 0 expansion factor. An expansion factor of 0 should only be used to represent a plot with no trees. Let’s look at where 0 expansion factors show up in the data:
## id time site plot exp_factor status decay species dbh ht
## 1 post_60 post 60 112 0 1 0 CADE 25.4 14.4
## 2 post_60 post 60 113 0 <NA> <NA> <NA> NA NA## id time site plot exp_factor status decay species dbh ht
## 1 post_60 post 60 112 24.69 0 3 QUKE 51.1 22.2
## 2 post_60 post 60 112 24.69 1 0 ABCO 69.9 30.9
## 3 post_60 post 60 112 24.69 1 0 ABCO 58.2 30.1
## 4 post_60 post 60 112 24.69 1 0 ABCO 23.6 17.6
## 5 post_60 post 60 112 24.69 1 0 ABCO 42.7 24.0
## 6 post_60 post 60 112 24.69 1 0 PSME 51.8 29.4
## 7 post_60 post 60 112 24.69 1 0 CADE 25.1 11.5
## 8 post_60 post 60 112 24.69 1 0 CADE 53.8 25.0
## 9 post_60 post 60 112 24.69 1 0 CADE 42.7 23.1
## 10 post_60 post 60 112 0.00 1 0 CADE 25.4 14.4It looks like a 0 expansion factor is properly used for post-60-113, but improperly used for post-60-112. We know what plot radius was used for larger trees, so we can confidently fill in the correct exp_factor here (24.69). The expansion factor will differ among studies. If a nested plot design was used, the expansion factor will differ among trees within the same plot.
Attempt 2: After correcting the expansion factor in the input data, let’s try again:
tree_bio <- SummaryBiomass(data = vign_trees_2,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
decay_class = "decay",
species = "species",
dbh = "dbh",
ht = "ht",
results = "by_plot")## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : There are dead trees with NA and/or zero decay class codes.
## The biomass of these dead trees will NOT be adjusted.
## Consider investigating these trees with mismatched status/decay class.
## ## Error in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : Not all species codes were recognized!
## Unrecognized codes: ABCCO SMEAnd we get another error message. It looks like there are some typos/transcription errors in the species codes. Looking at the list of unrecognized codes, we can tell that “ABCCO” should be “ABCO” (Abies concolor, commonly known as white fir) and “SME” should probably be “PSME” (Pseudotsuga menziesii, commonly known as Douglas-fir). Depending on the severity of the typo, you may want to go back to double check the species recorded on the original datasheet. In this case, the typos are fairly obvious. Let’s figure out where these typos occur in the data:
## id time site plot exp_factor status decay species dbh ht
## 1 post_340 post 340 108 24.69 1 0 ABCCO 12.7 11.3
## 2 post_340 post 340 109 24.69 1 0 SME 19.6 13.4Attempt 3: After correcting the species codes in the input data, let’s try again:
tree_bio <- SummaryBiomass(data = vign_trees_3,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
decay_class = "decay",
species = "species",
dbh = "dbh",
ht = "ht",
results = "by_plot")## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : There are dead trees with NA and/or zero decay class codes.
## The biomass of these dead trees will NOT be adjusted.
## Consider investigating these trees with mismatched status/decay class.
## ## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : There are missing DBH values in the provided dataframe.
## Trees with NA DBH will have NA biomass estimates.
## ## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : The allometric equations are for live trees with DBH >= 2.54cm and dead trees with DBH >= 12.7cm.
## You inputted live trees with DBH < 2.54cm. These trees will have NA biomass estimates.
## ## site plot live_Mg_ha dead_Mg_ha
## 1 post_340 103 470.15 10.08
## 2 post_340 104 241.31 0.00
## 3 post_340 108 156.67 0.00
## 4 post_340 109 169.09 0.00
## 5 post_340 111 189.39 4.31
## 6 post_340 112 823.50 14.57This time the function runs. However, we get some warning messages. Let’s look at the first warning, which tells us that there are trees with mismatched status/decay class. Recall that dead trees should have a decay class of 1-5 and live trees should have a decay class of NA or 0 (in this dataset we use 0 for live trees). Let’s figure out where mismatches occur in the data:
## id time site plot exp_factor status decay species dbh ht
## 1 pre_60 pre 60 27 24.69 0 0 CADE 42.2 21.6This CADE (Calocedrus decurrens, commonly known as incense cedar) was documented as dead (status = 0) but was assigned a decay class of 0 (which is reserved for live trees). We look back at the original datasheet (not shown here) and see that the tree was recorded as dead with a decay class of 0 in the field (which tells us this was not a transcription error). However, we also notice that a height to live crown base was recorded for the tree, indicating that the dead status was likely a recording error. Given that two pieces of information recorded for the tree point to a live status (i.e., a decay class of 0 and a height to live crown base), we can fairly confidently change the CADE’s status to live.
Attempt 4: After correcting for the mismatch in status/decay class, let’s try again:
tree_bio <- SummaryBiomass(data = vign_trees_4,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
decay_class = "decay",
species = "species",
dbh = "dbh",
ht = "ht",
results = "by_plot")## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : There are missing DBH values in the provided dataframe.
## Trees with NA DBH will have NA biomass estimates.
## ## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : The allometric equations are for live trees with DBH >= 2.54cm and dead trees with DBH >= 12.7cm.
## You inputted live trees with DBH < 2.54cm. These trees will have NA biomass estimates.
## ## site plot live_Mg_ha dead_Mg_ha
## 1 post_340 103 470.15 10.08
## 2 post_340 104 241.31 0.00
## 3 post_340 108 156.67 0.00
## 4 post_340 109 169.09 0.00
## 5 post_340 111 189.39 4.31
## 6 post_340 112 823.50 14.57The first warning message is gone (good!), but we still have two other warning messages. Let’s look at the next warning, which tells us there are missing DBH values in the dataframe. Let’s look at where NA DBH values show up in the data:
## id time site plot exp_factor status decay species dbh ht
## 1 pre_340 pre 340 115 24.69 1 0 PILA NA 24.6We look back at the original datasheet (not shown here) and see that DBH was recorded for the tree in the field. This was just a simple transcription error. However, in your own dataset you may have DBH values (or height values) that are truly missing. In the case of truly missing values, you may want to build a model that will allow you to predict DBH from total height (or total height from DBH). Such models may already exist for your study area.
Attempt 5: After filling in the missing DBH value, let’s try again:
tree_bio <- SummaryBiomass(data = vign_trees_5,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
decay_class = "decay",
species = "species",
dbh = "dbh",
ht = "ht",
results = "by_plot")## Warning in ValidateSumData(data_val = data, site_val = site, plot_val = plot, : The allometric equations are for live trees with DBH >= 2.54cm and dead trees with DBH >= 12.7cm.
## You inputted live trees with DBH < 2.54cm. These trees will have NA biomass estimates.
## ## site plot live_Mg_ha dead_Mg_ha
## 1 post_340 103 470.15 10.08
## 2 post_340 104 241.31 0.00
## 3 post_340 108 156.67 0.00
## 4 post_340 109 169.09 0.00
## 5 post_340 111 189.39 4.31
## 6 post_340 112 823.50 14.57There’s only one warning message left. The Forest Inventory and Analysis (FIA) Regional Biomass Equations are for live trees with a DBH >= 2.54 cm (1 in) and dead trees with a DBH >= 12.7 cm (5 in). It is not appropriate to use these equations to estimate biomass for trees with DBH below the specified cutoffs. We could filter these trees out of the dataset to avoid the warning message; however, that isn’t necessary. In this example, we’ll leave the small trees in and make note of the warning. If you need to calculate biomass for the smaller trees in your dataset, you will need to explore other options.
Next, we’ll use the ForestComp() and
ForestStr() functions to get forest composition and
structure at the plot level.
Let’s try using ForestComp() first:
for_comp <- ForestComp(data = vign_trees_5,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
species = "species",
dbh = "dbh",
relative = "ba",
units = "metric")## Error in ValidateCompData(data_val = data, site_val = site, plot_val = plot, : The "relative" parameter must be set to either "BA" or "density".We get an error message. It looks like we set the “relative” parameter to “ba” instead of “BA”. That’s easy to fix:
for_comp <- ForestComp(data = vign_trees_5,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
status = "status",
species = "species",
dbh = "dbh",
relative = "BA",
units = "metric")## The following species were present: ABCO CADE CONU NODE PILA PIPO PSME QUKE## site plot species dominance
## 1 post_340 103 QUKE 34.7
## 2 post_340 103 ABCO 39.2
## 3 post_340 103 CADE 26.1
## 4 post_340 103 PSME 0.0
## 5 post_340 103 NODE 0.0
## 6 post_340 103 PIPO 0.0
## 7 post_340 103 PILA 0.0
## 8 post_340 103 CONU 0.0
## 9 post_340 104 QUKE 16.3
## 10 post_340 104 ABCO 24.3
## 11 post_340 104 CADE 5.8
## 12 post_340 104 PSME 53.7
## 13 post_340 104 NODE 0.0
## 14 post_340 104 PIPO 0.0
## 15 post_340 104 PILA 0.0
## 16 post_340 104 CONU 0.0This time the function runs without any warning or error messages. However, there is a note with a list of all the species present in the input data. Notice in the output (we just show two plots here) that each species present is accounted for in each plot, even if the dominance for the specific species is 0. When you compile the data beyond the plot level (e.g. to the compartment level), it’s important that the 0 dominance values are captured.
Before moving on, let’s quickly look at the output for the plot with no trees:
## site plot species dominance
## 1 post_60 113 QUKE NA
## 2 post_60 113 ABCO NA
## 3 post_60 113 CADE NA
## 4 post_60 113 PSME NA
## 5 post_60 113 NODE NA
## 6 post_60 113 PIPO NA
## 7 post_60 113 PILA NA
## 8 post_60 113 CONU NANotice that the dominance is NA for all species. If there are no trees present, then % dominance is not applicable (i.e., NA). But why wouldn’t it make sense to say that there is 0% CADE? There’s no CADE on the plot, right? However, think about what % dominance actually is:
$\frac{BA_{sp,p}}{BA_{total,p}}$
where
If there are no trees present on plot p, then BAtotal, p will be 0. And anything divided by 0 is not defined. From both a mathematical and logical perspective, % dominance for a plot without trees is not applicable.
Now let’s try using ForestStr(). We’ll keep the default
for units (= “metric”):
for_str <- ForestStr(data = vign_trees_5,
site = "id",
plot = "plot",
exp_factor = "Exp_factor",
dbh = "dbh",
ht = "ht")## Error in ValidateStrData(data_val = data, site_val = site, plot_val = plot, : There is no column named "Exp_factor" in the provided dataframe.We get an error message. Right, we have a column named “exp_factor” but not “Exp_factor” (column names are case sensitive)! That’s easy to fix:
for_str <- ForestStr(data = vign_trees_5,
site = "id",
plot = "plot",
exp_factor = "exp_factor",
dbh = "dbh",
ht = "ht")
head(for_str)## site plot sph ba_m2_ha qmd_cm dbh_cm ht_m
## 1 post_340 103 691 74.04 36.9 32.1 19.1
## 2 post_340 104 296 41.38 42.2 36.2 17.5
## 3 post_340 108 963 41.17 23.3 21.4 11.5
## 4 post_340 109 938 35.22 21.9 16.2 9.3
## 5 post_340 111 296 40.52 41.7 38.2 18.9
## 6 post_340 112 667 103.76 44.5 38.5 23.7This time the function runs without any warning or error messages. But before moving on, let’s once again look at the output for the plot with no trees:
## site plot sph ba_m2_ha qmd_cm dbh_cm ht_m
## 1 post_60 113 0 0 NA NA NANotice that basal area (ba_m2_ha) and stems per hectare (sph) are both 0. This is correct for a plot with no trees. Then, average quadratic mean diameter (qmd_cm), average diameter at breast height (dbh_cm), and average height (ht_m) are all NA. These tree-level variables are not applicable for a plot without trees. For example, there is no “average diameter at breast height” if there are no diameters at breast height to take an average of. The average diameter at breast height is not 0 in this case! When you compile the data beyond the plot level (e.g., to the compartment level), it’s important that the 0s and NAs are accurately captured. A 0 qmd_cm here (at the plot level) would incorrectly pull down the compartment qmd_cm.
Finally, we’ll use the FineFuels(),
CoarseFuels(), and LitterDuff() functions to
get surface and ground fuel loads at the plot level.
Let’s investigate the input dataframe:
## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 post 340 103 17 4 1 1 1.83 1.83
## 2 post 340 103 89 11 2 1 1.83 1.83
## 3 post 340 104 13 4 0 1 1.83 1.83
## 4 post 340 104 60 3 0 0 1.83 1.83
## 5 post 340 108 7 11 1 3 1.83 1.83
## 6 post 340 108 80 14 4 3 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 11.34 0 64 0.00 0.0 7
## 2 3.05 11.34 0 0 0.00 0.0 15
## 3 3.05 11.34 0 0 0.00 0.0 6
## 4 3.05 11.34 0 313 0.50 0.5 8
## 5 3.05 11.34 0 289 0.75 0.5 7
## 6 3.05 11.34 0 0 0.00 0.0 2Attempt 1: Now, let’s try using
FineFuels(). We’ll keep the defaults for sp_codes (=
“4letter”) and units (= “metric”):
## Error in ValidateFWD(fuel_data_val = fuel_data, units_val = units): For fuel_data, there are repeat time:site:plot:transect observations.
## There should only be one observation/row for an individual transect at a specific time:site:plot.
## Investigate the following time:site:plot:transect combinations: post-400-9-159And we get an error message. It looks like there is a duplicate time:site:plot:transect observation for post-400-9-159. Let’s take a closer look:
## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 post 400 9 159 19 2 0 1.83 1.83
## 2 post 400 9 159 28 5 0 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 11.34 0 0 1.5 1.0 11
## 2 3.05 11.34 0 0 1.0 0.5 11## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 pre 400 9 159 32 8 0 1.83 1.83
## 2 pre 400 9 239 42 5 0 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 11.34 0 0 1 4.5 11
## 2 3.05 11.34 225 0 1 0.5 11Based on the field protocol, we know that there should be two transects per time:site:plot and that the two azimuths should be the same pre- and post-treatment. For the post-treatment observations, we can see that the counts, litter depth, and duff depth are not the same. Looking at the pre-treatment observations, we can see that the transect azimuths for site 400, plot 9 should be 159 and 239. It’s likely that one of the two post-400-9 transects should be changed to 239. We look back at the original datasheet (not shown here) and see that this is indeed the case.
Attempt 2: After correcting the transect azimuth in the input fuel data, let’s try again:
## Error in ValidateFWD(fuel_data_val = fuel_data, units_val = units): For fuel_data, count_100h must be a positive, whole number.We get another error message. It looks like there is an issue with one (or more) of the 100-hour counts. Let’s figure out where the issue occurs in the data:
vign_fuels_2 %>%
mutate(count_100h_check = abs(round(count_100h))) %>%
filter(count_100h != count_100h_check)## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 pre 340 13 215 35 13 1.1 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 11.34 81 100 3.5 3 7
## count_100h_check
## 1 1We can see that the count_100h value for pre-340-13-215 is 1.1, which is not a whole number. We look back at the original datasheet (not shown here) and see that this count should be 1. This was a transcription error.
Attempt 3: After correcting the count_100h value in the input fuel data, let’s try again:
## Error in ValidateMatches(fuel_match = step1, tree_match = step2): Tree and fuel data did not completely match!
## These time:site:plot combinations have tree data but no fuel data: pre-340-111
## These time:site:plot combinations have fuel data but no tree data: pre-340-11And we get yet another error message. Remember that there must be a one-to-one match between time:site:plot identities of tree and fuel data. See background information in the README file for more on why this one-to-one match is important. We known that there isn’t a plot 11 in the dataset, but that there is a plot 111. Pre_340_11 should probably be corrected to pre_340_111.
Attempt 4: After correcting the plot id in the input fuel data, let’s try again:
## Warning in ValidateFWD(fuel_data_val = fuel_data, units_val = units): For fuel_data, there are missing values in the count_10h column.
## For transects with NA 10h counts, 10h fuel load estimates will be NA.
## ## Warning in ValidateOverstory(tree_data_val = tree_data, sp_codes_val = sp_codes): Not all species codes were recognized! Unrecognized codes were converted to "UNTR" for unknown tree
## and will receive generic coefficients. Unrecognized codes: CONU NODE QUKE
## ## time site plot load_1h_Mg_ha load_10h_Mg_ha load_100h_Mg_ha load_fwd_Mg_ha
## 1 post 400 101 0.9029624 1.4033325 0.000000 2.306295
## 2 pre 400 101 2.2310300 4.7643799 5.382827 12.378237
## 3 post 60 101 0.3571794 0.0000000 4.563755 4.920934
## 4 pre 60 101 0.8695333 2.0983663 6.749952 9.717852
## 5 post 400 102 0.1562554 0.5359028 4.278554 4.970712
## 6 pre 400 102 0.2723569 3.2140721 2.850828 6.337257
## sc_length_1h sc_length_10h sc_length_100h
## 1 3.652542 3.652542 6.087570
## 2 3.652542 3.652542 6.087570
## 3 3.641661 3.641661 6.069435
## 4 3.641661 3.641661 6.069435
## 5 3.651810 3.651810 6.086350
## 6 3.651810 3.651810 6.086350This time the function runs. However, we get some warning messages. Let’s look at the first warning, which tells us that there are missing 10-hour counts. Let’s look at where NA count_10h values show up in the data:
## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 post 400 2 252 12 NA 1 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 11.34 81 81 0 0 5We look back the the original datasheet (not shown here) and see that this 10-hour count was not recorded in the field. It truly is missing, and there is little we can do about that at this point. Fortunately, there are two transects per time:site:plot. The other transect will still allow us to get a plot-level estimate of the 10-hour fuel load. The function will appropriately account for the missing 10-hour count.
The final warning message tells us that not all species codes were recognized. We look at the list of unrecognized codes: CONU (Cornus nuttallii, commonly known as pacific dogwood), NODE (Notholithocarpus densiflorus, commonly known as tanoak), and QUKE (Quercus kelloggii, commonly known as black oak). None of these are species code typos (but you should check for typos!). For the surface and ground fuel load calculations, we currently only have the necessary values for 19 Sierra Nevada conifer species (see background information in the README file for more detail on this topic). These three species are hardwoods and are, therefore, not included in the list of 19 Sierra Nevada conifers. This warning should not alarm us. Given the information currently available for the Sierra Nevada, the best we can do is assign generic “all species” values for these hardwoods.
Attempt 1: Now, let’s try using
CoarseFuels(). We’ll keep the defaults for sp_codes (=
“4letter”) and units (= “metric”):
## Error in ValidateCWD(fuel_data_val = fuel_data, units_val = units, sum_val = summed): For fuel_data, there are missing values in the length_1000h column.Another error message! It looks like there is a missing transect length. Let’s look at where NA length_1000h values show up in the data:
## time site plot transect count_1h count_10h count_100h length_1h length_10h
## 1 post 60 102 347 21 8 5 1.83 1.83
## length_100h length_1000h ssd_S ssd_R litter_depth duff_depth slope
## 1 3.05 NA 0 1771 0.5 1.25 9Based on the field protocol, we know that 1000-hour fuels were sampled for 11.34 meters along each transect. This is an easy fix.
Attempt 2: After filling in the missing transect length, let’s try again:
## Warning in ValidateOverstory(tree_data_val = tree_data, sp_codes_val = sp_codes): Not all species codes were recognized! Unrecognized codes were converted to "UNTR" for unknown tree
## and will receive generic coefficients. Unrecognized codes: CONU NODE QUKE
## ## time site plot load_1000s_Mg_ha load_1000r_Mg_ha load_cwd_Mg_ha
## 1 post 400 101 0.00000 0.00000 0.00000
## 2 pre 400 101 0.00000 8.02592 8.02592
## 3 post 60 101 0.00000 0.00000 0.00000
## 4 pre 60 101 18.77363 1.68116 20.45479
## 5 post 400 102 0.00000 0.00000 0.00000
## 6 pre 400 102 15.56532 0.00000 15.56532
## sc_length_1000s sc_length_1000r
## 1 22.63378 22.63378
## 2 22.63378 22.63378
## 3 22.56636 22.56636
## 4 22.56636 22.56636
## 5 22.62925 22.62925
## 6 22.62925 22.62925This time the function runs. We are already familiar with this
warning message (discussed in the FineFuels() function
section above).
Lastly, let’s try using LitterDuff(). We’ll keep the
defaults for sp_codes (= “4letter”), units (= “metric”), and measurement
(= “separate”):
## Warning in ValidateOverstory(tree_data_val = tree_data, sp_codes_val = sp_codes): Not all species codes were recognized! Unrecognized codes were converted to "UNTR" for unknown tree
## and will receive generic coefficients. Unrecognized codes: CONU NODE QUKE
## ## time site plot litter_Mg_ha duff_Mg_ha
## 1 post 400 101 8.207450 10.78962
## 2 pre 400 101 22.940513 82.86332
## 3 post 60 101 0.000000 0.00000
## 4 pre 60 101 37.114751 86.42521
## 5 post 400 102 9.304178 7.08900
## 6 pre 400 102 32.336705 31.92495The function runs. Once again, we are familiar with this particular
warning message (discussed in the FineFuels() function
section above).
Now that we have everything summarized at the plot level, let’s try compiling some of the data to the compartment level and then to the entire treatment (here the fire-only treatment) level.
Step 1: estimate tree biomass at the plot level
We already went through this process above, but let’s recall what the output dataframe looks like:
## site plot live_Mg_ha dead_Mg_ha
## 1 post_340 103 470.15 10.08
## 2 post_340 104 241.31 0.00
## 3 post_340 108 156.67 0.00
## 4 post_340 109 169.09 0.00
## 5 post_340 111 189.39 4.31
## 6 post_340 112 823.50 14.57Step 2: create all necessary columns for input into
the CompilePlots() function
tree_bio_2 <- tree_bio %>%
separate(site, c("time", "site")) %>% # separate id into time and site columns
mutate(trt_type = "fire") %>% # create a trt_type column
select(time, trt_type, site, plot, everything()) # organize columns as desired
head(tree_bio_2)## time trt_type site plot live_Mg_ha dead_Mg_ha
## 1 post fire 340 103 470.15 10.08
## 2 post fire 340 104 241.31 0.00
## 3 post fire 340 108 156.67 0.00
## 4 post fire 340 109 169.09 0.00
## 5 post fire 340 111 189.39 4.31
## 6 post fire 340 112 823.50 14.57Step 3: input data into
CompilePlots()
# keep the defaults for wt_data (= "not_needed")
tree_bio_sum <- CompilePlots(data = tree_bio_2,
design = "FFS")## time trt_type site avg_live_Mg_ha se_live_Mg_ha avg_dead_Mg_ha se_dead_Mg_ha
## 1 post fire 340 315.9475 37.45587 7.713500 3.2631872
## 2 post fire 400 290.6375 22.01452 1.537500 0.5935450
## 3 post fire 60 223.9916 26.24832 5.202632 2.9223176
## 4 pre fire 340 315.1080 35.60767 13.598500 5.2956273
## 5 pre fire 400 291.7030 22.18940 1.288000 0.5949558
## 6 pre fire 60 246.3984 23.62663 10.758421 9.8535120## time trt_type avg_live_Mg_ha se_live_Mg_ha avg_dead_Mg_ha se_dead_Mg_ha
## 1 post fire 276.8589 27.42481 4.817877 1.793207
## 2 pre fire 284.4031 20.16778 8.548307 3.721584Step 1: estimate fine fuel loads at the plot level
We already went through this process above, but let’s recall what the output dataframe looks like:
## time site plot load_1h_Mg_ha load_10h_Mg_ha load_100h_Mg_ha load_fwd_Mg_ha
## 1 post 400 101 0.9029624 1.4033325 0.000000 2.306295
## 2 pre 400 101 2.2310300 4.7643799 5.382827 12.378237
## 3 post 60 101 0.3571794 0.0000000 4.563755 4.920934
## 4 pre 60 101 0.8695333 2.0983663 6.749952 9.717852
## 5 post 400 102 0.1562554 0.5359028 4.278554 4.970712
## 6 pre 400 102 0.2723569 3.2140721 2.850828 6.337257
## sc_length_1h sc_length_10h sc_length_100h
## 1 3.652542 3.652542 6.087570
## 2 3.652542 3.652542 6.087570
## 3 3.641661 3.641661 6.069435
## 4 3.641661 3.641661 6.069435
## 5 3.651810 3.651810 6.086350
## 6 3.651810 3.651810 6.086350Step 2: create all necessary columns for input into
the CompileSurfaceFuels() function
FWD_2 <- FWD %>%
mutate(trt_type = "fire") %>% # create a trt_type column
select(time, trt_type, site, plot, everything()) # organize columns as desired
head(FWD_2)## time trt_type site plot load_1h_Mg_ha load_10h_Mg_ha load_100h_Mg_ha
## 1 post fire 400 101 0.9029624 1.4033325 0.000000
## 2 pre fire 400 101 2.2310300 4.7643799 5.382827
## 3 post fire 60 101 0.3571794 0.0000000 4.563755
## 4 pre fire 60 101 0.8695333 2.0983663 6.749952
## 5 post fire 400 102 0.1562554 0.5359028 4.278554
## 6 pre fire 400 102 0.2723569 3.2140721 2.850828
## load_fwd_Mg_ha sc_length_1h sc_length_10h sc_length_100h
## 1 2.306295 3.652542 3.652542 6.087570
## 2 12.378237 3.652542 3.652542 6.087570
## 3 4.920934 3.641661 3.641661 6.069435
## 4 9.717852 3.641661 3.641661 6.069435
## 5 4.970712 3.651810 3.651810 6.086350
## 6 6.337257 3.651810 3.651810 6.086350Step 3: input data into
CompileSurfaceFuels()
# keep the defaults for cwd_data (= "none), wt_data (= "not_needed"), and units (= "metric")
FWD_sum <- CompileSurfaceFuels(fwd_data = FWD_2,
design = "FFS")## time trt_type site avg_1h_Mg_ha se_1h_Mg_ha avg_10h_Mg_ha se_10h_Mg_ha
## 1 post fire 400 0.4969988 0.05976717 1.139383 0.2120926
## 2 pre fire 400 0.9968562 0.14076733 2.991177 0.3775630
## 3 post fire 60 0.4725502 0.07563463 1.735057 0.3326378
## 4 pre fire 60 1.2239823 0.20421903 4.355628 0.6468982
## 5 post fire 340 0.3601983 0.04450288 1.051188 0.2369040
## 6 pre fire 340 1.0927353 0.19285796 5.313776 0.9169993
## avg_100h_Mg_ha se_100h_Mg_ha
## 1 1.317892 0.3381811
## 2 4.438505 0.5900822
## 3 4.151810 1.3535638
## 4 9.035420 1.3362516
## 5 2.351752 0.5969664
## 6 6.336984 1.3793372## time trt_type avg_1h_Mg_ha se_1h_Mg_ha avg_10h_Mg_ha se_10h_Mg_ha
## 1 post fire 0.4432491 0.04212089 1.308543 0.2147717
## 2 pre fire 1.1045246 0.06583011 4.220194 0.6738877
## avg_100h_Mg_ha se_100h_Mg_ha
## 1 2.607151 0.8279885
## 2 6.603636 1.3336957