Campos & Leite’s model was superior, thus, it will be used.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 x 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
#> # ... with 894 more rows, and 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 x 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
#> # ... with 4 more rowsNow 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.04Now 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.87975Ps: 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:
resid_plot(data_ex, "TH", "CL","Henricksen") 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.4To 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.9071631And 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.908154Now 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.1460599To 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