Andrew Plowright 2024-01-19
This vignette will take you through the process of generating tree
statistics over polygonal areas using the outputs from ForestTools
.
Start by loading in some test data. These are sample outputs from the
ForestTools
functions. kootenayTrees
are individual tree points
produced using vwf
. kootenayCrowns
are polygonal crowns generated by
mcws
. kootenayBlocks
are some polygonal areas of interest that we’ll
use for generating statistics.
# Load libraries
library(terra)
library(sf)
library(ForestTools)
library(dplyr)
library(magrittr)
# Load sample data
data("kootenayTrees", "kootenayBlocks", "kootenayCrowns", "kootenayCHM")
# Plot areas of interest and add trees
plot(unwrap(kootenayCHM))
plot(kootenayBlocks$geometry, add=T)
plot(kootenayTrees$geometry, pch=19, col="blue", add=T)
Calculating global statistics for the trees is fairly straightforward.
# Height
mean(kootenayTrees$height)
## [1] 5.251959
median(kootenayTrees$height)
## [1] 4.256785
max(kootenayTrees$height)
## [1] 13.49121
# Crown area
mean(kootenayCrowns$crownArea)
## [1] 9.006453
median(kootenayCrowns$crownArea)
## [1] 7
max(kootenayCrowns$crownArea)
## [1] 54.25
To calculate statistics by polygonal areas of interest, we’ll first use
the st_intersects
function to create subsets of trees for each
polygon. We’ll then calculate statistics for each of those subsets, and
return the result to polygons. Note: this can duplicate trees that are
contained in overlapping polygons.
# Create subsets of the trees for each polygonal area of interest
trees_by_poly <- kootenayBlocks %>%
st_intersects(kootenayTrees) %>%
lapply(function(x) kootenayTrees[x,])
# Calculate statistics for each polygonal area of interest
kootenayBlocks[["mean_height"]] <- sapply(trees_by_poly, function(trees) mean(trees$height))
kootenayBlocks[["max_height"]] <- sapply(trees_by_poly, function(trees) max(trees$height))
The same operation can be repeated for tree crowns. However, given that
the crowns are polygonal, we don’t want a single crown to be counted
twice between two adjoining areas of interest. The only variation
required here is that we use st_centroid
to compute centroids for each
crown before intersecting them with the areas of interest.
# Create subsets of the crowns for each polygonal area of interest
crowns_by_poly <- kootenayBlocks %>%
st_intersects(st_centroid(kootenayCrowns)) %>%
lapply(function(x) kootenayCrowns[x,])
# Calculate statistics for each polygonal area of interest
kootenayBlocks[["mean_crownArea"]] <- sapply(crowns_by_poly, function(crowns) mean(crowns$crownArea))
kootenayBlocks[["max_crownArea"]] <- sapply(crowns_by_poly, function(crowns) max(crowns$crownArea))
# View results
kootenayBlocks[,c("BlockID", "mean_height", "max_height", "mean_crownArea", "max_crownArea")]
## Simple feature collection with 3 features and 5 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 439689.1 ymin: 5526454 xmax: 439832.5 ymax: 5526563
## Projected CRS: WGS 84 / UTM zone 11N
## BlockID mean_height max_height mean_crownArea max_crownArea
## 0 101 7.338588 13.491207 13.683476 51.25
## 1 3308 3.208017 7.125149 4.847403 23.50
## 2 113 7.685308 12.583441 13.430990 54.25
## geometry
## 0 MULTIPOLYGON (((439707.9 55...
## 1 MULTIPOLYGON (((439832.5 55...
## 2 MULTIPOLYGON (((439832.4 55...
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.