9 Derived metrics at the cloud level

9.1 Overview

The “cloud” level of regularization corresponds to the computation of derived metrics using all available points. As seen in section 8, calculating derived metrics for the whole point cloud is straightforward and users only need to provide a formula to calculate metric(s) of interest. For example, to calculate the average height (mean(Z)) of all points we can run the following:

LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
las <- readLAS(LASfile)
cloud_metrics(las, func = ~mean(Z)) # calculate mean height
#> [1] 13.27202

To calculate more than one metric a custom function needs to be used, or one of the pre-defined functions within lidR. To calculate the whole suite of 36 metrics defined in stdmetrics_z() we can use func = .stdmetrics_z. When several metrics are computed they are returned as a list.

metrics <- cloud_metrics(las, func = .stdmetrics_z)
str(head(metrics)) # output is a list
#> List of 6
#>  $ zmax    : num 30
#>  $ zmean   : num 13.3
#>  $ zsd     : num 7.45
#>  $ zskew   : num -0.476
#>  $ zkurt   : num 2.08
#>  $ zentropy: num 0.903

9.2 Applications

Point cloud metrics become interesting when computed for a set of plot inventories. In this case it can serves to compute a set of metrics for each plot, where known attributes have been measured in the field to construct a predictive model. Users could easily clip and loop through plot inventory files but we defined an all-in-one convenient function plot_metrics(). In the following example we load a .las and compute the metrics for each plot inventory using a shapefile of plot centers

LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
las <- readLAS(LASfile, filter = "-keep_random_fraction 0.5")
shpfile <- system.file("extdata", "efi_plot.shp", package="lidR")
inventory <- sf::st_read(shpfile, quiet = TRUE)
metrics <- plot_metrics(las, .stdmetrics_z, inventory, radius = 11.28)

Look at the content of inventory and metrics. inventory contains the plot IDs, their coordinates, and VOI a Value Of Interest. metrics contains 36 derived metrics for each plot combined with the inventory data

head(inventory)
#> Simple feature collection with 5 features and 2 fields
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: 684838.9 ymin: 5017796 xmax: 684976.6 ymax: 5017958
#> Projected CRS: NAD83 / UTM zone 17N
#>   IDPEP       VOI                 geometry
#> 1 PEPQ1 14.157140 POINT (684976.6 5017821)
#> 2 PEPQ2 12.720584 POINT (684923.9 5017958)
#> 3 PEPQ3 11.396656 POINT (684838.9 5017942)
#> 4 PEPQ4 11.597471   POINT (684855 5017891)
#> 5 PEPQ5  8.263425   POINT (684944 5017796)
head(metrics[,1:8])
#> Simple feature collection with 5 features and 8 fields
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: 684838.9 ymin: 5017796 xmax: 684976.6 ymax: 5017958
#> Projected CRS: NAD83 / UTM zone 17N
#>   IDPEP       VOI  zmax    zmean      zsd      zskew    zkurt  zentropy
#> 1 PEPQ1 14.157140 26.61 15.29461 7.710077 -0.9766295 2.641739 0.8644264
#> 2 PEPQ2 12.720584 25.21 16.39847 6.306690 -1.1466891 3.595908 0.8691998
#> 3 PEPQ3 11.396656 26.13 16.22342 6.390952 -0.7976664 2.829905 0.9139921
#> 4 PEPQ4 11.597471 25.56 12.77717 6.273880 -0.2144288 2.569128 0.9192199
#> 5 PEPQ5  8.263425 22.45  7.96093 5.984458  0.4195126 2.315890 0.9042922
#>                   geometry
#> 1 POINT (684976.6 5017821)
#> 2 POINT (684923.9 5017958)
#> 3 POINT (684838.9 5017942)
#> 4   POINT (684855 5017891)
#> 5   POINT (684944 5017796)

We have computed many metrics for each plot and we know the value of interest VOI. We can use that to build a linear model with some metrics. Here we have only 5 plots so it is not going to be big science

model <- lm(VOI~zsd+zmax, data = metrics)
summary(model)
#> 
#> Call:
#> lm(formula = VOI ~ zsd + zmax, data = metrics)
#> 
#> Residuals:
#>        1        2        3        4        5 
#> -0.02129  1.32581 -0.90887 -0.07268 -0.32297 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)
#> (Intercept) -17.9438     9.0365  -1.986    0.185
#> zsd           1.0961     1.1577   0.947    0.444
#> zmax          0.8896     0.4816   1.847    0.206
#> 
#> Residual standard error: 1.161 on 2 degrees of freedom
#> Multiple R-squared:  0.858,  Adjusted R-squared:  0.7159 
#> F-statistic:  6.04 on 2 and 2 DF,  p-value: 0.142

This example can be improved. In section 14 we will study how to extract a ground inventory and in section 16 we will study more in depth modelling presenting a complete workflow from the plot extraction to the mapping of the predictive model.