The use of GEDI canopy structure for explaining variation in tree species richness in natural forests

We used tree species richness data which were compiled from publicly available, previously published studies, by Keil and Chase (2019). Specifically, the dataset consisted of data on species richness from 1166 forest plots in tropical and temperate regions from across the world. These forest plots were originally selected by Keil and Chase (2019) to represent both temperate and tropical natural and semi-natural forests across as many biomes and biogeographic realms as possible given open access data availability. Nevertheless, as with most such analyses, there are considerable data gaps in large parts of Africa and Asia due to limited availability of forest plots from these regions. Only plots from natural or semi-natural forest (e.g. not from plantations, production forests, or severely disturbed/damaged forests, for example, as a result of logging, hurricanes, or bark beetle outbreaks) were used. These plots were taken from published database compilations (Phillips and Miller 2002, Phillips et al 2003, Ramesh et al 2010, Anderson-Teixeira et al 2015, Coelho de Souza et al 2016, Ricklefs and He 2016, Sullivan et al 2017), national forest inventory surveys (De Natale et al 2005, USDA 2016, IGN 2017), and extracted manually from primary sources (Beals 1965, Yamada 1976, Abbott 1984, Kohyama 1986, Cao and Zhang 1997, Linder et al 1997, Ansley and Battles 1998, Narayanan and Parthasarathy 1999, Maycock et al 2000, Adam 2001, Kohira et al 2001, Cairns et al 2003, Enoki 2003, Hirayama and Sakimoto 2003, Cheng-Yang et al 2004, Fashing et al 2004, Jing-Yun et al 2004, Sheil and Salim 2004, Shu-Qing et al 2004, Wu et al 2004, Bonino and Araujo 2005, Davis et al 2005, Sawada et al 2005, Splechtna et al 2005, Adekunle 2006, Graham 2006, Malizia and Grau 2006, Szwagrzyk and Gazda 2007, Lopes et al 2008, Yasuoka 2008, Addo-Fordjour et al 2009, Eichhorn 2010, Krishnamurthy et al 2010, Lalfakawma et al 2010, Nagel et al 2010, Namikawa et al 2010, Lü et al 2010, Eshete et al 2011, Round et al 2011, Sanchez et al 2013, Popradit et al 2015, van Do et al 2015, Wusheng et al 2015). For each plot we used the number of tree species (S) as our focal response variable. From the abovementioned studies we further extracted three covariates related to sampling design, which can all influence the observed S, and thus need to be controlled for in statistical analyses. These covariates were:

• The total area of the sampled plot (km²).
• The minimum diameter at breast height (minDBH in cm); only trees of equal or larger diameter than minDBH were measured and included in the analysis.
• The number of trees larger than the minDBH in the sampled plot.

For each plot, we considered six variables related to the climate, topography, and biogeographic history which were previously successfully used to explain tree species richness in Keil and Chase (2019). These were:

• WorldClim (Hijmans et al 2005) mean annual temperature (BIO1), mean isothermality (BIO3), and mean precipitation seasonality (BIO15), which are 1970–2000 spatially interpolated averages extracted from rasters with a resolution of 30 arc-seconds (1 km² around equator) available at www.worldclim.org/data/v1.4/worldclim14.html.
• Absolute difference between the lowest and highest elevation in a 30 arc-seconds raster cell within which a given plot is located, derived from Danielson and Gesch (2011).
• Island vs. mainland, a binary factor indicating if a plot is on an island or mainland (from Weigelt et al 2013), following a loose classification in which shelf islands (connected to mainland in the past) are classified as islands.
• Classification of biogeographic realms into Afrotropic, Australasian, Easter Palearctic, Indo-Malay, Nearctic, Neotropic, and Western Palearctic, following (Ricklefs and He 2016).

For further details on these environmental data, see supplementary material of Keil and Chase (2019).

We used GEDI-derived variables describing canopy structure from the first two years of data collection: April 18, 2019–May 14, 2021. GEDI is a sampling instrument mounted to the International Space Station, collecting lidar data between 51.6° N and S. GEDI has three lasers that are split into four beams total. The beams are dithered across track, resulting in eight beams along-track. The footprint-spacing along track is 60 m and the across-track spacing is 600 m. At each sampling location, GEDI collects a full-waveform with a nominal footprint size of ∼25 m (Dubayah et al 2020). During the first two years of the mission, GEDI was programmed to collect data with the least spatial overlap, acquiring the densest possible sampling at the footprint level. Because of the sampling method, footprints do not necessarily overlap exactly with the sampled forest plots from Keil and Chase (2019), thus we used the information from GEDI footprints in square buffer zones around the original plot locations. This also ensures that we capture the canopy structure of the forest plots for those plots with a limited spatial location accuracy. GEDI footprints themselves have a geolocation accuracy of about 10 m. We included plots that have at least one GEDI footprint within the buffer zone around the plot centre. We evaluated three buffer zones (figure 1(b)):

• 1 × 1 km buffer zone (1 km², 888 plots).
• 2 × 2 km buffer zone (4 km², 978 plots).
• 4 × 4 km buffer zone (16 km², 1053 plots).

Zoom In Zoom Out Reset image size

There were some plots where no GEDI shots intersected any of the buffers because of its sampling method. Moreover, the number of plots available at the different buffer zones increases with the size of the buffer zone (varying from 888 to 1053 plots out of the original 1166 plots considered, figure 1(a)).

According to Marselis et al (2019, 2020), the most important GEDI metrics related to tree species richness were: canopy height, total plant area index and the plant area index along the vertical forest axis divided into 10 m vertical bins. We adopted this same set of metrics. However, to confirm that other canopy structure variables do not explain more of the variation in richness, we included the results of a principal component analysis (PCA) analysis performed on a large set of canopy structure metrics (180 metrics, appendix A). The final set of canopy metrics used here consists of the following variables:

• Canopy Height, measured as the relative height at which 98% of the energy was returned (RH98).
• Total plant area index (PAI) along the vertical forest axis
• The total PAI between 0 and 10 m vertically, m²/10 m³
• The total PAI between 10 and 20 m vertically, units in m²/10 m³
• The total PAI between 20 and 30 m vertically, units in m²/10 m³
• The total PAI between 30 and 40 m vertically, units in m²/10 m³
• The total PAI between 40 and 50 m vertically, units in m²/10 m³
• The total PAI above 50 m vertically, units in m²/10 m³
• The ground elevation

All these metrics were derived from the GEDI Level 2A and 2B datasets version 2 (Dubayah et al 2021a, 2021b). The L2A GEDI dataset contains the relative height (RH) products of the processed, geolocated waveform, including RH98 (essentially canopy height) and the ground elevation. The L2B dataset contains the vertical profile metrics; including the total PAI and the PAI along the vertical axis. GEDI L2A data contains information from four processing algorithms applied to the geolocated waveforms. We used the canopy height and ground elevation from the processing algorithm flagged as the 'best' algorithm for each shot in the GEDI dataset (Dubayah et al 2021a, 2021b). We also only included waveforms that were collected in the leaf on season, using the 'leaf off flag' in the GEDI data products. Lastly, only waveforms that were flagged '1' for the metric 'successful processing at the 2b data level' (vertical canopy structure, Dubayah et al 2021b) were considered reliable and used in this analysis.

For each of the 16 km², 4 km² and 1 km² areas, we calculated the average of the nine canopy structure metrics from all footprints within the buffer zone and joined them to the field and environmental data, creating three final datasets (one for each buffer zone).

2.2.1. Raw relationships

2.2.2. RFs

2.2.3. Collinearity and variation partitioning

The most important result is that, in local forest plots at a global extent, we were unable to unequivocally attribute variation in species diversity to GEDI-derived forest structure, nor to other environmental predictors. More specifically, even though forest structure (together with covariates describing study design) explained up to 66% of variation in tree species diversity, environmental predictors related to climate, topography, and history can explain even more variation in diversity (up to 80%), and there is no increase in explained variation when GEDI and environmental predictors are used together in one model. This means that GEDI-derived forest structure can indeed be used to explain local tree species diversity, which is in line with the relationship between canopy height and tree species diversity found by Gatti et al (2017), and it underpins the hypothesis that canopy structure data may serve to predict variation in tree species richness (Marselis et al 2020). However, if other environmental predictors are available they may provide even better diversity estimates than the GEDI-derived variables, and we found no clear advantage in combining the strengths of the two groups of predictors for the datasets we examined.

Even though this is a negative result from a practical point of view–the GEDI-derived forest structure does not improve existing correlative near-global models of tree species diversity–there is a valuable conceptual inference. Independent effects of forest structure and environment on forest diversity (figure 5(a)) are unlikely, because otherwise, including forest structure should have unequivocally improved the overall models. Our results are more consistent with an alternative (causal) pathway in which the abiotic environment drives both forest structure and species diversity (figure 5(b)); with the possibility of a secondary effect of forest structure on species diversity (figure 5(c)) or the possibility of a two-way effect where tree species richness also drives forest structure itself (figure 5(d)). However, our correlative models cannot reject or confirm this secondary link, since the overlapping fraction of explained variation (figure 4) is too large. Thus, we conclude that the three pathways (figures 5(b)–(d)) are all plausible. This is in line with other studies (Frolking et al 2009, Pfeifer et al 2018, Senf et al 2020, Hakkenberg and Goetz 2021) which hypothesize that canopy structure should depend on environmental factors, soil and disturbance. Potential evidence for the secondary effect (figure 5(c)) was found by Hakkenberg and Goetz (2021) who studied the relation between canopy structure, climate and tree species richness in forests of the United States. Their study showed that a small portion of the remaining variance was explained by climate x structure interaction terms. Such a result could not be reproduced at this global scale using GEDI canopy structure, potentially due to the differences in canopy structure metrics available as Hakkenberg and Goetz (2021) used high resolution discrete return airborne lidar data at 20 × 20 m plot level, which points to the importance of spatial scale, discussed below.

Zoom In Zoom Out Reset image size

Our tree species richness dataset only included data on natural or semi-natural forests; there were no (severely) disturbed forests, orchards, production forests, or plantations/monocultures. Thus, our forest plots mostly have a compact and well-developed canopy, which potentially reduces the variation of forest structure across the dataset, limiting the scope for the GEDI-derived structure to explain diversity. In contrast, human-affected forests have distinctly 'unnatural', regular, or loose canopies (Tropek et al 2014). If such forests were also included in the analysis, we might expect that the GEDI-derived forest structure should explain more variation in richness, and could even outperform the other environmental predictors. Here we see a potential for joint use of environment and GEDI-derived variables to predictions of actual richness of all forests, including those modified by humans, as opposed to potential richness of natural plots as in this study and in Keil and Chase (2019).

Furthermore, GEDI footprints are not randomly distributed in our plots' buffer zones and they may not exactly overlap with the plots. For example, GEDI measurements may only cover a part of the buffer zone (e.g. figure 1(b)). More local, plot-specific forest structure information could provide additional/different information on the canopy structure than our buffered and spatially averaged metrics, but unfortunately the GEDI sampling density is not high enough to assert this. Encouragingly, however, we did not find considerable differences in model performance between the smallest buffer zone areas (1 km² and 4 km²) and a small improvement with the largest buffer zone area (16 km²). This is likely because structural variables derived from these different buffer areas are similar (appendix C), which indicates that averaging a handful of GEDI footprints in a 1 km² buffer zone is already representative of general forest structure in a larger area. This could be in line with previous findings suggesting that the canopy structure in tropical forests becomes uncorrelated over just hundreds of meters (Réjou-Méchain et al 2014), suggesting that all of our buffer zones describe a general forest structure instead of the plot-level forest structure. In this case, the larger number of canopy structure samples in 16 km² zones, compared to the smaller zones, would likely result in a more accurate measure of the general forest structure which may explain the improved model performance at this larger buffer zone.

The temporal mismatches between collection of the field data, the environmental data, and the canopy structure data should also be kept in mind. Forest growth and (small) natural disturbances alter the forest over the shorter and longer term (Senf et al 2020), resulting in a mismatch of the forest structure and representation of the tree species diversity. This is a limitation that is difficult to overcome because of the time-consuming and expensive field data collection. Here we assume that the data used here (with the GEDI data from 2019 to 2021, climatic data from 1970 to 2000, and forest plot data collected between 1962 and 2016) represent the long-term situation in these natural and semi-natural forests. We assume that changes of climate and species temporal turnover in these forests are slow enough to still allow a correlative analysis involving the more recent GEDI data. The high explanatory power of both the climate and GEDI-derived variables support this assumption, although the real consequences of this temporal mismatch are unknown and should be investigated in the future. Additionally, it is yet unknown how the temporal dynamics of GEDI observations acquired from different times within the same leaf-on season, as well as across leaf-on seasons, may affect our results.

Future research should further explore the more realistic hypotheses (i.e. figures 5(b)–(d). We need to step back to properly examine the relationship between abiotic conditions and forest structure (such as recently performed in the US by Hakkenberg and Goetz (2021)). Now that the GEDI data are available, a global analysis of this relationship using its massive amount of remotely sensed structure data is possible. This should enable us to further untangle the (causal) relationships between climate, forest structure, and disturbance, improve our understanding of ecosystem functioning, and enable the analysis of the real usefulness of canopy structure for mapping tree species richness in both natural and disturbed forests.

This also points to the importance of spatial scale in general. The question of the role of structure vs. environment in driving patterns of tree species richness should be explored as a function of scale. This could, for example, help to determine whether there are characteristic scales of analysis where structure dominates beyond environmental variables or vice versa. Such an analysis would require high quality species richness data (e.g. at 25 m²), across a landscape, along with vertical canopy structure data at the same resolution in spatially continuous grids. Answering this question is exceptionally difficult because these kinds of continuous field data are rare. Furthermore, it has been shown that grain of disturbances is spatially very small, with mean disturbance sizes of 10–20 ha across continents, but with a median strongly skewed towards smaller patches (∼0.1 ha) (Taubert et al 2018), and not inconsistent with the decorrelation of structure over hundreds of meters noted above. Given the importance of small-scale disturbance to tree diversity, this suggests the need for fine-grained information on canopy structure at global scales. GEDI provides the basis for obtaining such observations when used in combination with other continuous remote sensing data at fine (30–100 m) resolutions, such as may be obtained from passive optical sensors (e.g. Landsat; Potapov et al 2021) or from radar sensors that measure 3D structure (e.g. TanDEM-X; Qi et al 2019). GEDI data can provide the basis for machine learning and other fusion methods that predict not only height and biomass at fine scales, but potentially the imputation of entire waveforms, opening the door to pervasive canopy structure predictions, at fine scales, across the landscape.

Some of the complex interactions can perhaps also be teased apart by focusing on fully mapped forest plots such as those in the CTFS/ForestGeo network (Anderson-Teixeira et al 2015), which offer rich data on exact location and growth form of every individual tree, collected over multiple time periods, with growth form and species of each known individual. If these were connected with remotely sensed forest structure using fusion with GEDI, an analysis teasing apart the interplay between forest structure, diversity, and abiotic conditions across a variety of scales may be possible.

Lastly, alternative analytical approaches that accommodate complex causal pathways, for example, structural equation modelling (Iriondo et al 2003, Lam and Maguire 2012), can be adopted. A similar analysis has been done to evaluate the direct and indirect effects of forest taxonomical and structural diversity, soil fertility and climate on aboveground carbon storage in neotropical forests (Poorter et al 2015).