For this example we’ll use forest inventory data, calculate volume and other variables for each plot.
library(forestmangr)
library(dplyr)
data(exfm9)
data_ex <- exfm9
data_ex
#> # A tibble: 900 × 14
#> MAP PROJECT STRATA GENCODE STRATA_AREA PLANTING_DATE SPACING PLOT
#> <chr> <int> <int> <fct> <int> <date> <fct> <int>
#> 1 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 2 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 3 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 4 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 5 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 6 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> # ℹ 894 more rows
#> # ℹ 6 more variables: MEASUREMENT_DATE <date>, PLOT_AREA <int>, PIT <int>,
#> # DBH <dbl>, TH <dbl>, OBS <fct>
The first step is to estimate the height of non measured trees. We’ll evaluate two hypsometric models. Henricksen’s: \[ Ln(H) = \beta_0 + \beta_1*Ln(H) \]
And Campos & Leite’s model, which uses dominant height: \[ Ln(H) = \beta_0 + \beta_1*\frac{1}{dbh} +
\beta_2*Ln(DH) \] In order to use this model, first we’ll need to
calculate the dominant height for each plot. To do this we’ll use the
dom_height
function. In it we’ll input the dataset, and
height, dbh, plot and observation variables. The observation variable
refers to the quality or classification of the tree, i.e., if it’s
normal, forked, dead or dominant. In addition to that, we’ll also input
the code used to distinguish dominant trees. In this dataset, the code
is "dom"
.
dom_height(df=data_ex,th="TH",dbh="DBH",plot="PLOT",obs="OBS",dom="D")
#> # A tibble: 10 × 2
#> PLOT DH
#> <int> <dbl>
#> 1 1 25.0
#> 2 2 23.8
#> 3 3 20.7
#> 4 4 20.0
#> 5 5 19.6
#> 6 7 24.4
#> # ℹ 4 more rows
Now that we’ve seen the dominant height of each plot, we can run the
function again, but this time set the merge_data
argument
as TRUE
, to bind the variable to the data directly:
data_ex <- dom_height(data_ex,"TH","DBH","PLOT","OBS","D",merge_data = TRUE)
head(as.data.frame(data_ex))
#> MAP PROJECT STRATA GENCODE STRATA_AREA PLANTING_DATE SPACING PLOT
#> 1 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 2 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 3 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 4 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 5 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 6 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> MEASUREMENT_DATE PLOT_AREA PIT DBH TH OBS DH
#> 1 2012-09-04 810 1 15.0 23.8 N 25.04
#> 2 2012-09-04 810 2 13.0 23.8 N 25.04
#> 3 2012-09-04 810 3 15.0 24.7 N 25.04
#> 4 2012-09-04 810 4 13.5 23.3 N 25.04
#> 5 2012-09-04 810 5 15.0 24.3 N 25.04
#> 6 2012-09-04 810 6 14.0 22.4 N 25.04
Now we can fit the hypsometric models. We’ll fit them using the
lm_table
function. the function
forestmangr::inv
will allow us to fit the Campos &
Leite model without creating a new variable for 1/dbh:
data_ex <- data_ex %>%
lm_table(log(TH) ~ inv(DBH) + log(DH),output="merge_est",est.name="CL") %>%
lm_table(log(TH) ~ log(DBH), output="merge_est",est.name="Henricksen")
head(data_ex)
#> MAP PROJECT STRATA GENCODE STRATA_AREA PLANTING_DATE SPACING PLOT
#> 1 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 2 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 3 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 4 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 5 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 6 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> MEASUREMENT_DATE PLOT_AREA PIT DBH TH OBS DH CL Henricksen
#> 1 2012-09-04 810 1 15.0 23.8 N 25.04 24.60879 22.98647
#> 2 2012-09-04 810 2 13.0 23.8 N 25.04 23.42624 20.75027
#> 3 2012-09-04 810 3 15.0 24.7 N 25.04 24.60879 22.98647
#> 4 2012-09-04 810 4 13.5 23.3 N 25.04 23.74891 21.31799
#> 5 2012-09-04 810 5 15.0 24.3 N 25.04 24.60879 22.98647
#> 6 2012-09-04 810 6 14.0 22.4 N 25.04 24.05250 21.87975
Ps: The lm_table
function checks if the model has log in
the y variable, and if it does, it removes it automatically when
estimating variables. Because of that, there’s no need to calculate the
exponential for the estimated variables.
We’ll evaluate the quality of the fitted values using
resid_plot
. Non measured trees will be removed
automatically:
Campos & Leite’s model was superior, thus, it will be used.
Now we can estimate the height of non measured trees, using
dplyr::mutate
and ifelse
:
data_ex <- data_ex %>%
mutate( TH_EST = ifelse(is.na(TH), CL, TH ), CL=NULL,Henricksen=NULL )
head(data_ex)
#> MAP PROJECT STRATA GENCODE STRATA_AREA PLANTING_DATE SPACING PLOT
#> 1 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 2 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 3 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 4 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 5 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 6 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> MEASUREMENT_DATE PLOT_AREA PIT DBH TH OBS DH TH_EST
#> 1 2012-09-04 810 1 15.0 23.8 N 25.04 23.8
#> 2 2012-09-04 810 2 13.0 23.8 N 25.04 23.8
#> 3 2012-09-04 810 3 15.0 24.7 N 25.04 24.7
#> 4 2012-09-04 810 4 13.5 23.3 N 25.04 23.3
#> 5 2012-09-04 810 5 15.0 24.3 N 25.04 24.3
#> 6 2012-09-04 810 6 14.0 22.4 N 25.04 22.4
To estimate the volume with bark, we’ll take a previously fitted equation, and save it in a data frame:
tabcoef_vwb <- data.frame(b0=-9.595863,b1=1.889372,b2=0.9071631)
tabcoef_vwb
#> b0 b1 b2
#> 1 -9.595863 1.889372 0.9071631
And do the same for volume without bark:
tabcoef_vwob <- data.frame(b0=-9.808975,b1=1.918007,b2=0.908154)
tabcoef_vwob
#> b0 b1 b2
#> 1 -9.808975 1.918007 0.908154
Now we’ll estimate volume, basal area and age:
data_ex <- data_ex %>%
mutate(CSA = pi*DBH^2/40000,
AGE = as.numeric(MEASUREMENT_DATE - PLANTING_DATE) / 30,
VWB = exp(tabcoef_vwb$b0 + tabcoef_vwb$b1*log(DBH) + tabcoef_vwb$b2*log(TH_EST) ),
VWOB = exp(tabcoef_vwob$b0 + tabcoef_vwob$b1*log(DBH) + tabcoef_vwob$b2*log(TH_EST) ) )
head(data_ex)
#> MAP PROJECT STRATA GENCODE STRATA_AREA PLANTING_DATE SPACING PLOT
#> 1 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 2 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 3 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 4 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 5 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> 6 FOREST1 1 2 FLO014 45 2007-04-01 3x3 1
#> MEASUREMENT_DATE PLOT_AREA PIT DBH TH OBS DH TH_EST CSA AGE
#> 1 2012-09-04 810 1 15.0 23.8 N 25.04 23.8 0.01767146 66.1
#> 2 2012-09-04 810 2 13.0 23.8 N 25.04 23.8 0.01327323 66.1
#> 3 2012-09-04 810 3 15.0 24.7 N 25.04 24.7 0.01767146 66.1
#> 4 2012-09-04 810 4 13.5 23.3 N 25.04 23.3 0.01431388 66.1
#> 5 2012-09-04 810 5 15.0 24.3 N 25.04 24.3 0.01767146 66.1
#> 6 2012-09-04 810 6 14.0 22.4 N 25.04 22.4 0.01539380 66.1
#> VWB VWOB
#> 1 0.2011052 0.1761617
#> 2 0.1534627 0.1338786
#> 3 0.2079921 0.1822010
#> 4 0.1616611 0.1411803
#> 5 0.2049342 0.1795194
#> 6 0.1670810 0.1460599
To summarise the plots, we’ll use plot_summarise
. We’ll
input the dataset, and variables such as plot, plot area, dbh, height,
total area, volume with bark, volume without bark, dominant height and
age. With this, the function will calculate for each plot the mean
diameter, quadratic diameter, mean height and mean dominant height. It
will also calculate to total number of individuals, basal area and
volume with and without bark for each plot, and extrapolate it to
hectares.
tab_invt <- plot_summarise(df=data_ex,plot="PLOT",plot_area="PLOT_AREA",dbh="DBH", th="TH_EST",total_area="STRATA_AREA",vwb="VWB",vwob="VWOB",dh="DH",age="AGE")
head(as.data.frame(tab_invt))
#> PLOT AGE STRATA_AREA PLOT_AREA DBH q TH_EST DH Indv Indvha
#> 1 1 66 45 810 14.1722 14.2196 24.0004 25.04 90 1111.111
#> 2 2 66 45 810 14.4494 14.5235 23.3723 23.76 90 1111.111
#> 3 3 63 45 810 12.8539 12.8902 19.9767 20.70 90 1111.111
#> 4 4 63 51 810 12.2833 12.3285 18.9652 19.98 90 1111.111
#> 5 5 63 51 810 12.6389 12.6939 18.8658 19.60 90 1111.111
#> 6 7 66 45 810 14.7584 14.8135 23.9202 24.40 90 1111.111
#> G G_ha VWB VWB_ha VWOB VWOB_ha
#> 1 1.4292 17.6450 16.5503 204.3245 14.4789 178.7514
#> 2 1.4744 18.2028 16.6473 205.5216 14.5737 179.9222
#> 3 1.1614 14.3388 11.5018 141.9973 10.0318 123.8493
#> 4 1.0744 13.2638 10.2121 126.0759 8.8958 109.8248
#> 5 1.1390 14.0618 10.7609 132.8508 9.3818 115.8243
#> 6 1.5339 18.9371 17.6205 217.5366 15.4333 190.5350