Wind Speed Controls Forest Structure in a Subtropical Forest Exposed to Cyclones: A Case Study Using an Individual-Based Model

Rau, E-Ping; Gardiner, Barry A.; Fischer, Fabian Jörg; Maréchaux, Isabelle; Joetzjer, Emilie; Sun, I-Fang; Chave, Jérôme

doi:10.3389/ffgc.2022.753100

ORIGINAL RESEARCH article

Front. For. Glob. Change, 03 May 2022

Sec. Forest Disturbance

Volume 5 - 2022 | https://doi.org/10.3389/ffgc.2022.753100

This article is part of the Research Topic Living with Tropical Storms in a Changing Climate View all 6 articles

Wind Speed Controls Forest Structure in a Subtropical Forest Exposed to Cyclones: A Case Study Using an Individual-Based Model

$\r\nE-Ping Rau$ E-Ping Rau¹

Barry A. Gardiner^2,3

Fabian Jörg Fischer¹

Isabelle Maréchaux⁴

Emilie Joetzjer⁵

I-Fang Sun⁶ $Jrme Chave*\r\n$ Jérôme Chave^1*

¹Laboratoire Evolution et Diversité Biologique (EDB, UMR 5174) CNRS, IRD, Université Toulouse III, Toulouse, France
²Laboratoire Interactions Sol Plante Atmosphère (UMR 1391 ISPA), Bâtiment C2, INRAE Centre de Bordeaux-Aquitaine, Villenave, France
³Institut Européen De La Forêt Cultivée, Cestas, France
⁴AMAP, Université de Montpellier, CIRAD, CNRS, INRAE, IRD, Montpellier, France
⁵UMR SILVA INRAE, AgroParisTech, Université de Lorraine, Nancy, France
⁶Department of Natural Resources and Environmental Studies, National Dong Hwa University, Hualien City, Taiwan

Extreme wind blowdown events can significantly modify the structure and composition of forests, and the predicted shift in tropical cyclone regimes due to climate change could strongly impact forests across the tropics. In this study, we coupled an individual-based and spatially-explicit forest dynamics model (TROLL) with a mechanistic model estimating wind damage as a function of tree size, traits, and allometry (ForestGALES). We assimilated floristic trait data and climate data from a subtropical forest site in Taiwan to explore the effect of wind regimes on forest properties. We found that the average canopy height and biomass stocks decreased as wind disturbance strength increased, but biomass stocks showed a nonlinear response. Above a wind intensity threshold, both canopy height and biomass drastically decreased to near-zero, exhibiting a transition to a non-forest state. Wind intensity strongly regulated wind impact, but varying wind frequency did not cause discernible effects. The implementation of within-stand topographic heterogeneity led to weak effects on within-stand forest structure heterogeneity at the study site. In conclusion, the intensity of wind disturbances can potentially greatly impact forest structure by modifying mortality. Individual-based modeling provides a framework in which to investigate the impact of wind regimes on mortality, other factors influencing wind-induced tree mortality, as well as interaction between wind and other forms of forest disturbance and human land use legacy.

Introduction

Natural disturbances shape forest structure, composition and dynamics, and play a substantial role in controlling the global carbon cycle (Reichstein et al., 2013; Pugh et al., 2019). Wind is a major source of natural disturbances and an important driver of tree mortality (Mitchell, 2013; Gardiner et al., 2016): while wind can cause partial tree damage such as branch snapping and defoliation, extreme wind events can lead to tree death, mainly through stem breakage and tree uprooting. In particular, extreme wind disturbances generated by tropical cyclones in regions such as the Caribbean (Lugo et al., 2000) and the northwestern Pacific (Lin et al., 2011) have been shown to control forest structure at a regional scale (Hogan et al., 2018; Ibanez et al., 2020), with one notable effect being the reduction of canopy height with increased cyclone frequency and intensity (Ibanez et al., 2019). Cyclone regimes are projected to change in the future (Lin et al., 2020), with a general increase in cyclone intensity and interregional differences in the trend of cyclone frequency (Knutson et al., 2020). In order to anticipate how such changes will affect community- and ecosystem-scale properties, it is essential to deepen our understanding of the mechanisms through which cyclonic winds regulate tree mortality, which remains one of the demographic processes that are less well understood in forest ecosystems (Bugmann et al., 2019; Koch et al., 2021).

However, it remains intrinsically difficult to reach general conclusions on factors driving wind-induced tree mortality based on empirical studies: these studies rely on observations of forest damage after individual cyclone events, which are both limited in their temporal extent, and driven by contingent factors such as cyclone trajectory and force (Lugo, 2008) and human land use legacy (Kulakowski et al., 2011; Schwartz et al., 2017). In various studies, a higher probability and level of damage has been found to be correlated to tree characteristics such as larger stem diameter (Ostertag et al., 2005; Halder et al., 2021), larger height (Dunham and Cameron, 2000; Vandecar et al., 2011), or lower wood density (Curran et al., 2008; Webb et al., 2014). However, the size dependence of damage can vary with tree species (Xi, 2015) and damage type (Everham and Brokaw, 1996), and some studies even did not observe any size effect on tree damage (Bellingham et al., 1996). In addition, wind impact on trees can vary at fine scales within a given forest and for a given species, as wind speed interacts in complex ways with terrain topography, often accelerating at hilltops and mountain ridges, and decelerating in leeward valleys (Ruel et al., 1998; Belcher et al., 2011). As a result, windward, exposed terrains can experience heavier wind damage than leeward regions and valleys (McEwan et al., 2011; Yap et al., 2016) (but see de Toledo et al., 2012).

Recently, a number of studies have contributed to a more mechanistic understanding of wind-tree interactions at a finer scale. To do so, they modelled the response of trees or tree-like structures when exposed to wind drag using the finite-element method (Moore and Maguire, 2008; Jackson et al., 2019b,c), and used continuously collected field data to relate tree motion to architectural properties (Jackson et al., 2019a,2021). These studies provide important insights into the mechanisms underlying wind-tree interaction. However, the response of trees to wind depends not only on individual characteristics, but also on the tree neighborhood, and more specifically on the presence and the relative size of neighboring trees, which influence wind sheltering: as an illustration, wind typically has a disproportionate impact on emergent, unsheltered trees (Hale et al., 2012; Seidl et al., 2014; Duperat et al., 2021). Hence, the impact of wind on a forest stand is more than the sum of the impacts on single trees: the heterogeneity of forest structure and composition must be accounted for when modeling the impacts of wind, especially for mixed sized forests.

In light of these challenges, dynamic forest models that are both spatially explicit and individual-based could be useful since they link external environmental conditions to community- and ecosystem-scale properties, through individual-level processes such as mortality, while explicitly accounting for spatiotemporal heterogeneity of forests (Seidl et al., 2014; Fischer et al., 2016). By including extreme wind disturbances in the climate forcing of these models, it is possible to represent mortality in a more mechanistic fashion, and generate stand-level predictions about forest response to wind disturbances, which could then be compared with forest inventory or remote sensing data. Several individual-based models (Ancelin et al., 2004; Schelhaas et al., 2007; Uriarte et al., 2009; Seidl et al., 2014; Kamimura et al., 2019) and land surface models (Chen et al., 2018) have incorporated windthrow or storm damage models. However, these models are heavily skewed toward temperate forests with relatively homogeneous stand structure. To the best of our knowledge, there is currently no forest dynamics model that simulates the effects of extreme wind disturbance at species-rich, structurally heterogeneous forests in the subtropics and the tropics, where tropical cyclones and downburst winds can be important disturbance factors.

In this study, we explored the community-wide consequences of wind disturbances, using a model-based approach. We coupled a wind disturbance module to an individual-based forest model, and created forcing conditions from tropical cyclone records and wind speed data. We parameterized this model for the Fushan Forest in Taiwan, located in a cyclone-prone region, and explored how extreme wind events impact tree mortality and forest structure and dynamics. Specifically, we addressed the following questions:

(1) How does wind-induced tree mortality impact forest structure and dynamics? We expect that the presence of wind-induced tree mortality would cause a reduction in average canopy height due to higher vulnerability of tall trees to wind damage. This could in turn increase light availability in the understory and facilitate tree establishment, leading to an increase in stem density. We also expect that total mortality as well as the fraction of mortality due to treefall would increase.

(2) How do changes in frequency and intensity of extreme wind events lead to shifts in forest dynamics and structure? We expect that higher wind frequency and intensity would result in more wind-induced mortality for large trees, which in turn would increase the extent of canopy height reduction and stem density increase.

(3) How does topographical heterogeneity within a forest stand influence wind damage patterns? We expect that the inclusion of topography would lead to a higher heterogeneity in forest structure, with lower canopy height and higher stem density in more exposed locations.

Materials and Methods

Description of the TROLL Model

The TROLL model is a spatially explicit individual-based model that divides the aboveground space of a forest stand into 3D cells of size 1 m³ (voxels). Solar irradiance inside each voxel is computed as the irradiance fraction transmitted from immediately above the focal voxel (Chave, 1999; Maréchaux and Chave, 2017). No more than one individual tree can establish in each 1 × 1 m pixel at any given time, and only self-standing stems ≥ 1 cm in DBH (diameter at breast height, measured at 1.3 m above ground, or 20 cm above deformities or buttresses if present) are explicitly modeled (herbaceous plants and lianas are not included). Seeds and seedlings < 1 cm in DBH are not directly modeled, but represented as part of a seed/seedling pool. Each modeled tree is characterized by its DBH, height, crown size (radius and depth), total leaf surface area, and age. Trees are assigned species-specific traits, which control photosynthesis, growth as well as other physiological and demographic processes. At each monthly time step, trees grow depending on the balance between carbon acquisition based on a C3 photosynthesis model (Farquhar et al., 1980) and respiration, with assimilated carbon allocated into wood or leaf production. Height-DBH allometric relations are prescribed at the species level, while the DBH-crown size allometry is assumed to be species-independent. Tree mortality in TROLL is represented by several processes: (i) stochastic tree death, here assumed to be negatively dependent on species-specific wood density (WD) (Wright et al., 2010): m_eff=m−α×WD, (α is positive, so m is the maximal possible value of the mortality rate), (ii) carbon starvation, when net assimilated carbon is negative over a consecutive period exceeding leaf lifespan, so that old leaves have all died while no new leaves could be produced, (iii) stochastic treefall, dependent on a species-specific tree height threshold related to maximum realized tree height, and (iv) secondary treefall, when a tree is affected by a falling neighboring tree. To more explicitly link tree mortality to exogenous factors such as wind disturbance, we updated the treefall module for this study to simulate treefall more mechanistically (cf. below).

Inputs for a typical run of TROLL (v3.0) include (1) monthly mean climate forcing data for the simulated location, (2) average trajectory of daily climatic variation, (3) species-specific plant traits for the simulated forest, and (4) species-independent general parameters, including those related to the stochastic tree death and treefall. The source code of TROLL (v3.0) is written in C++ and is available on https://github.com/troll-code/troll. On a computing cluster, each simulation of 200 × 200 m and 500 years uses around 15 min of computer processing unit (CPU) time.

Modeling Tree-Top Wind Speeds in a Mixed-Sized Forest

Horizontal wind speed decreases as one approaches the canopy top from above the canopy. This is modeled by the aerodynamic momentum transfer model above a vegetation canopy (Monteith and Unsworth, 2008, p310). The horizontal wind speed profile with height is represented with the following equation:

u (z) = \frac{u_{*}}{k} \ln (\frac{z - d}{z_{0}}) i f z \geq H (1a)

where u(z) (m s^–1) is the horizontal wind speed at height z (m) above ground or a displacement surface, k = 0.41 is the von Kármán constant, and u_* (m s^–1) is the friction velocity. Parameter d is called the zero-plane displacement height, and z₀the aerodynamic roughness. Both d and z₀ are dependent on canopy height H (m), and much research has explored the variation in both parameters d and z₀ over forest vegetation (Shaw and Pereira, 1982; Dorman and Sellers, 1989; Shuttleworth et al., 1989; Raupach et al., 1991). In general, d ranges from 0.7 H to 0.9 H, and z₀ from 0.04 H to 0.08 H for tropical forests. For this study, we chose the intermediate values of d = 0.8 H and z₀ = 0.06 H. These values are consistent with dense and tall evergreen forest vegetation as simulated in the SiB2 model (Sellers et al., 1996).

The aerodynamic momentum transfer model (1a) is only applicable above the canopy. In this study, we decided to apply an additional model to represent wind within the canopy (Inoue, 1963):

u (z) = u (H) \exp (- α (1 - z / H)) i f z < H (1b)

where α=H/L_s≈3; empirical values of L_s are reported in Table 1 of Raupach et al. (1996). With this parameterization, horizontal wind speed u(z) at a height of d+z₀=0.86H and within the canopy at H/2 is 66% and 22% of u(H), respectively. To reduce computational burden, we assumed that trees whose height is lower than H/2 were not directly affected by wind. To account for canopy height heterogeneity and neighbor tree sheltering, we calculated canopy height H for every 20 × 20 m quadrat by taking the arithmetic mean of the top leaf-containing voxel layer of each pixel within the quadrat.

TABLE 1

Table 1. Slope estimates of linear models of each simulated quadrat-level forest structure metrics to the wind speed correction factor, with p-values in parentheses.

From the CRU-NCEP reanalysis dataset, we obtained empirical atmospheric wind speed values (cf. the “Climate forcing” section below), which are conventionally measured at 10 m above ground or above a displacement surface. In order to obtain the horizontal wind speed experienced by each tree, we applied the following equation, derived from eq. (1a), to convert horizontal atmospheric wind speed at d + 10, u(d + 10), into top-canopy wind speed u(H) within each 20 × 20 m quadrat:

u (H) = u (d + 10) \times \frac{\ln (\frac{H - d}{z_{0}})}{\ln (\frac{10}{z_{0}})}

= u (d + 10) \times \frac{1.204}{5.116 - ln (H)} (2)

Finally, for each tree with height h > H/2, we computed u(h) using either Eqs. (1a) or (1b).

Modeling Wind-Induced Tree Death: Critical Wind Speed

To simulate the impact of extreme wind disturbance, we assumed that wind-induced treefall is a stochastic process whose probability depends on the difference between the tree-top wind speed u(h) exerted on a tree of height h (m), and a tree-specific critical wind speed (CWS). The higher u(h) is relative to CWS, the more likely the tree is to fall and die. In this study, we used a logistic model to relate the wind-induced tree death probability to (u(h) – CWS) (Valinger and Fridman, 1999; Hale et al., 2015):

p = 1 / (1 + e^{- (u (h) - C W S)}) (3)

For the estimation of CWS, we followed the approach of ForestGALES, a wind damage risk model originally developed for even-sized and even-spaced plantation forests (Gardiner et al., 2000, 2008) but then adapted for use with individual trees in more complex forest structures (Hale et al., 2012, 2015; Duperat et al., 2021; Quine et al., 2021). ForestGALES models the biomechanical conditions under which the tree, simplified as an anchored vertical object, will break (causing stem breakage) or lose root anchorage (causing uprooting) when subject to a turning moment caused by wind loading.

Importantly, the formulation of Eq. (3) departs from that of ForestGALES, where the control of wind speed on mortality is defined to be u(H), the top-canopy wind speed. To account for the reduced impact of wind within forest canopies, the CWS in ForestGALES can be modulated by a competition index, so that trees within the canopy are less exposed to wind damage than top-canopy trees (Hale et al., 2012; Duperat et al., 2021; Quine et al., 2021). In the present formulation, CWS does not depend on the tree neighborhood, but the probability of wind-induced tree death depends on wind speed at tree height, as computed from Eqs. (1a) and (1b).

At our study site, field observations indicate that the proportion of tree uprooting is low compared to that of stem breakage (Supplementary Appendix A). We therefore focused on wind-induced stem breakage. The critical wind speed (CWS, m s^–1) inducing stem breakage and treefall is calculated by the following equation (Hale et al., 2015):

C W S = P \times {(\frac{M O R \times D B H^{3}}{T_{C}})}^{\frac{1}{2}} (4)

where DBH (m) is the tree diameter at breast height, and MOR (Pa, kg m^–1 s^–2) is the fresh-wood modulus of rupture. We estimated MOR from oven-dry wood density (WD, g cm^–3) through the following exponential relationship (Green et al., 1999, see Supplementary Appendix B):

M O R = 17.2 \times e^{(2.51 \times W D)} \times 10^{6} (5)

T_C is the turning moment coefficient (kg), which relates the square of mean wind speed (m² s^–2) to the maximum turning moment (kg m² s^–2) experienced by an individual tree. We used the following empirical equation for T_C (Hale et al., 2012):

T_{C} = τ \times D B H^{2} \times h (6)

where the constant τ (kg m^–3) = 111.7 (Hale et al., 2015). A higher T_C value means larger turning moment for a given wind speed.

The wind damage parameter P (unitless) is stand-specific and determines the overall susceptibility of a forest stand to wind-induced tree death: the smaller P value is, the lower the critical wind speed is for any given tree under the same condition, meaning that the forest is overall more susceptible to wind-induced tree death. It encapsulates multiple constants and corrective factors, as well as factors in the original equation (see Hale et al., 2015, Eqs. 11 and 12; also see Supplementary Appendix C) that are difficult to estimate, and whose interpretation is beyond the scope of this study. The value of P was thus tuned here by means of a sensitivity analysis (see below).

Climatic Forcing

The TROLL model requires the following monthly mean climatic variables: daytime and nighttime mean temperature, daytime mean irradiance, and daytime mean vapor pressure deficit (VPD). We used the CRU-NCEP dataset to provide the monthly climate forcing, a global gridded (0.5° × 0.5°) sub-daily (6-hourly) climate product spanning the 1901--2016 period (version 8; version 7 archived¹) (Viovy, 2018). It was constructed by combining observation-based CRU TS 3.2 data (Harris et al., 2014) and model-based NCEP-NCAR data (Kalnay et al., 1996). We used CRU-NCEP variables for the 1980–2016 period, for which the most observations are available (Kistler et al., 2001).

We estimated the monthly average frequency of cyclones occurring within a sufficiently close distance to the study site. For this, we used the IBTrACS dataset (International Best Track Archive for Climate Stewardship database; v04r00, archived²), which contains records of global tropical cyclones occurring since 1945 (Knapp et al., 2010). A common measure of the spatial extent of tropical cyclones is the mean radius of gale-force winds, by convention defined as 17.5 m s^–1. Based on the range of reported values in the literature, we assumed the mean radius of gale-force winds (R₁₇, km) to be 150 km (Chan and Chan, 2012; Weber et al., 2014; Lu et al., 2017). Thus, for a given site, we calculated the monthly mean frequency of recorded tropical cyclones occurring within a 150-km distance from the site over the period of 1987–2020. At Fushan, cyclone occurrence frequency was highest from July to September (∼ 0.5 cyclones per month); on average, total annual cyclone frequency was around 1.84 cyclones per year.

For the selected tropical cyclone records, we calculated on-site cyclonic wind speed (V_site, m s^–1) using the empirical function that relates the radial variation of the tangential wind speed beyond the radius of maximum sustained wind in mature tropical cyclones (Anthes, 1982; Hsu and Babin, 2005): $V_{s i t e} = 17.5 \times \sqrt{R_{17} / d}$ , where d (km) is the distance between the site and the cyclone center. Since d < R₁₇ by definition, it follows that V_site > 17.5 m s^–1. We then fitted the on-site cyclonic wind speed values of each month to a Weibull distribution, using the function fitdistr in the R package MASS (Venables and Ripley, 2002), and used the scale and shape parameters as input climate forcing variables. At Fushan, the year-round V_site distribution is right-skewed, with 1st quartile = 18.85 m s^–1, median = 20.60 m s^–1, and 3rd quartile = 25.18 m s^–1.

The coupling of TROLL to wind disturbances was performed as follows. At each time step, the occurrence of an extreme wind event was drawn from a cyclone occurrence probability, assuming that one extreme wind event at most can occur per time step. If an extreme wind event occurred, we drew a random wind speed value from the input wind speed distribution. We accounted for canopy heterogeneity and neighbor tree sheltering by only considering trees with height > H/2 as being exposed to wind disturbance, where H represents quadrat-wide average top canopy height, calculated by the mean value of the top leaf-containing voxel height of each pixel within the quadrat. For each exposed tree, we converted atmospheric wind speed, u(d + 10), to its tree-top wind speed, u(h) (Eqs. 1a, 1b, 2), and calculated its critical wind speed (Eq. 4). We then compared the two wind speed values, and determined if wind-induced tree death happens through a stochastic process, dependent on a logistic function of the difference between the two wind speeds (Eq. 3). Secondary treefall was modeled by assuming that when each tree dies, it falls in a random direction and increases the death rate in the impacted pixels.

Study Site and Parameterization

We parameterized the TROLL model with wind-induced tree mortality for a forest site in Taiwan. The Fushan Forest Dynamics Plot (FDP) is a 25-hectare (500 m × 500 m) plot in a moist broadleaf forest in the northeastern region of Taiwan (Su et al., 2007), part of ForestGEO (Forest Global Earth Observatory) (Condit, 1998; Anderson-Teixeira et al., 2015). Fushan is under strong influence of the northeasterly monsoon in winter and frequent typhoon visits in summer and autumn, with mean annual precipitation around 4200 mm, mean annual temperature around 18°C, and a mean relative humidity around 95%. Plot elevation ranges from 600 m to 733 m (Su et al., 2007). It was established in 2004 and censused every 5 years since then: all self-standing stems with a DBH ≥ 1 cm were identified, measured, tagged and mapped, with a total of 110 recorded tree species in the plot (Su et al., 2007).

Species-specific plant functional traits that TROLL requires as input parameters include leaf mass per area (LMA, g m^–2), nitrogen and phosphorus content (N_mass and P_mass, mg g^–1), wood density (WD, g cm^–3), maximum DBH (cm), DBH-height allometric parameters and regional relative abundance. These traits were measured at Fushan according to ForestGEO protocol (Iida et al., 2014) for 94 species, covering ca. 90% of the censused trees. Climatic variables were extracted from the CRU-NCEP dataset at the geographic coordinates closest to Fushan FDP (24° 45′ N, 121° 32′ E).

Stem density and aboveground biomass (AGB) estimations were available at Fushan from census data. In order to estimate on-site tree mortality, we used data from the annual tree mortality survey, which has been conducted at Fushan following ForestGEO protocol since 2017. The mortality survey records the number of tree deaths and cause of death (standing, uprooted, or broken) of a subset of censused trees (Arellano et al., 2021). Using the mortality data at Fushan spanning 2017 to 2020, we calculated the mean annual mortality rate for all trees DBH > 10 cm, as well as the proportion of mortality attributed to treefall. For the detailed protocol of attribution of mortality factors, see Supplementary Appendix A.

Unless otherwise specified, for all simulations, we simulated forest regeneration from bare soil for a reference plot area of 4 hectares (200 m × 200 m) for a duration of 500 years with a monthly timestep. Since we aimed to examine the effects of extreme wind on mature forest, the wind-induced tree mortality sub-model was activated after a 100-year spin up phase.

Sensitivity Analysis

The wind-induced tree mortality sub-model includes a single parameter P (unitless). To investigate model responses to its value, we conducted a sensitivity analysis by varying the value of P across a range of [0.005, 1], with a varying step of 0.005: for each value, we performed five replicates of TROLL simulations (a total of 1000 simulations).

We calculated the steady-state values (mean over the last 100 years of the simulation) of three structure metrics: stem density (DBH > 10 cm; N₁₀, tree ha^–1), Lorey’s height (basal area-weighted mean tree height, m) (Pourrahmati et al., 2018), and aboveground biomass (AGB, Mg ha^–1). For trees with DBH > 10 cm, we also calculated two mortality statistics: mean annual mortality and fraction of mortality due to treefall (%M_treefall). We reported mortality statistics at the onset of wind disturbance (year 101–200, i.e., first 100 years after wind submodule activation) and after reaching the steady state (year 401–500, i.e., last 100 years of simulation). We qualitatively described trends and sensitivity of these statistics in response to variation of parameter value.

Effects of Wind Frequency and Intensity

To examine how the frequency and intensity of extreme wind events influence forests, we performed two series of simulations: (1) we varied cyclone occurrence frequency from 0.1 to 2 times the empirical frequency (obtained from cyclone best-track data), with a varying step of 0.1, while maintaining wind intensity; (2) we varied the scale parameter of the wind speed distribution (parameter controlling mean and median of the Weibull distribution), from 0.1 to 10 times the empirical values (estimated from the cyclone best-track wind speed distribution), with a varying step of 0.1, while maintaining empirical frequency. Five replicates were performed for each condition (in total, 100 simulations for frequency and 500 simulations for intensity were performed). We set P = 0.7 based on the results of the sensitivity analysis, where simulation results are close to field observations and not near the forest tipping point (P < 0.3). As in previous steps, we calculated the steady-state values (mean over the last 100 years of the simulation) of stem density (DBH > 10 cm) (tree ha^–1), Lorey’s height (m), aboveground biomass (AGB, Mg ha^–1), as well as mean annual mortality and fraction of mortality due to treefall (%M_treefall) for trees with DBH > 10 cm. We qualitatively described trends and sensitivity of these statistics in response to variation of wind frequency and intensity.

Effects of Topography on Wind Disturbances

In the original TROLL model, the effect of topography was not taken into account. Based on the knowledge that wind speed is altered over an uneven topography, we implemented quadrat-scale wind speed correction factors in the model that account for these topographical effects. For this, we used the Global Wind Atlas (GWA) data, which are produced through downscaling using the WAsP program (Mortensen et al., 2000), incorporating surface elevation and aerodynamic roughness lengths (Badger et al., 2015). We obtained 250 m × 250 m GWA pixels that fall in the area covered by the 1° × 1° CRU-NCEP pixel where the Fushan site is located: this represents a grid of 200 × 200 GWA pixels. We normalized the GWA wind speed values of the selected pixels so that the mean GWA wind speed is equal to the mean CRU-NCEP wind speed. We then resampled the GWA pixels to the 20m × 20m quadrat scale using bilinear interpolation with the resample function in the raster package (Hijmans, 2022), and selected the resampled pixels falling within the area of the Fushan plot: this represents a grid of 25 × 25 = 625 resampled pixels. We used the GWA wind speed values of these resampled pixels, normalized by the plot-wide mean, as the wind speed correction factor for each quadrat.

The wind speed correction factor ranged from 0.27 to 1.96, and was used as a proxy for topographic heterogeneity: if the topographic effect is implemented, the wind speed experienced at each quadrat is the plot-wide wind speed (randomly drawn from the input wind speed distribution) multiplied by this correction factor: a value above 1 means that the wind speed at that quadrat is considered to speed up (due to exposed terrain), and considered to slow down when the value is below 1.

We performed simulations with and without topographical effect at the Fushan site for a plot area of 25 hectares (500 m × 500 m), and examined the relationship between quadrat-level steady-state values (mean over the last 100 years of the simulation) of stem density (DBH > 10 cm) (tree ha^–1), Lorey’s height (m), aboveground biomass (AGB, Mg ha^–1) and the topographic effect by performing linear regressions for each statistics as a function of quadrat-level wind speed correction factor.

In order to compare simulated results with field observation, we also calculated field-measured quadrat-level stem density (DBH > 10 cm), Lorey’s height (m) and aboveground biomass. We plotted them against quadrat-level mean elevation (m), which ranged from 616 to 730 m, and performed linear regressions for each statistic as a function of the quadrat-level mean elevation.

Data Analysis

Data processing, statistical analysis and visualization were performed in R 3.3.0 (R Core Team, 2021). Apart from those already mentioned elsewhere, R packages ggplot2, cowplot, ncdf4, BIOMASS, geosphere, sp, and tidyr were used for this study (Pebesma and Bivand, 2005; Wickham, 2016; Réjou-Méchain et al., 2017; Wilke, 2020; Hijmans, 2021; Pierce, 2021; Wickham and Girlich, 2022).

Results

Sensitivity Analysis

As the P parameter value decreased (stronger wind-induced tree mortality), both average canopy height (Lorey’s height) and aboveground biomass (AGB) decreased, although AGB showed a hump-shaped pattern (Figure 1): at low P-values (P < 0.3), canopy height and biomass decreased to extremely low levels, suggesting a transition from forest to non-forest state. Stem density (N₁₀) only slightly decreased, but showed an abrupt increase around P = 0.15, where average canopy height and AGB reached extremely low levels, before decreasing again at lower P-values.

FIGURE 1

Figure 1. Summary statistics of the simulated forests, in relation to the critical wind speed parameter P (smaller P means stronger effect): (A) N₁₀, density of stems with DBH > 10 cm (trees ha^–1); (B) Lorey’s height, basal area-weighted mean tree height (m); (C) AGB, aboveground biomass (Mg ha^–1). Shaded areas represent the value range of five replicates for each simulation condition. Dashed lines represent simulation value with no wind disturbance, and Solid lines represent field observations.

Mean annual mortality and fraction of mortality due to treefalls (%M_treefall) increased as P decreased, although %M_treefall exhibited an erratic nonlinear response. Both statistics were markedly lower at the end of the simulation (after the steady state has been reached) than immediately after the onset of wind forcing (Figure 2). The fraction of mortality due to treefalls (%M_treefall) was overall around 75%, substantially higher than field observations (Figure 2, horizontal black line; Supplementary Appendix A).

FIGURE 2

Figure 2. (A) Mean annual mortality rate and (B) fraction of mortality attributed to treefalls (%M_treefall) for trees with DBH > 10 cm, in relation to the critical wind speed parameter P (smaller P means stronger effect). Shaded areas represent the value range of five replicates for each simulation condition. Black without border: mean values of the first 100 years after wind onset (year 101–200). Red with border: mean values of the last 100 years of simulation (year 401–500). Dashed lines represent simulation value with no wind disturbance, and solid lines represent field observations.

Effects of Wind Frequency and Intensity

None of the forest structure and mortality statistics showed significant change with varying wind frequency. On the other hand, as wind intensity increased, stem density (N₁₀), average canopy height (Lorey’s height) and aboveground biomass (AGB) all decreased (with a hump-shaped response for aboveground biomass), and mean mortality and the fraction of mortality due to treefalls (%M_treefall) increased, especially at higher wind intensity, where the near-zero level canopy height and biomass suggested a transition from forest to non-forest state as observed during the sensitivity test (Figure 3).

FIGURE 3

Figure 3. Simulated forest structure in relation to extreme wind frequency (A–E) and intensity (F–J), relative to empirical values from cyclone best-track data. For simulations with varying frequency, empirical wind intensity were used, and vice versa. Shaded areas represent the value ranges of five replicates for each simulation condition. (A,F) N₁₀, density of stems with DBH > 10 cm (trees ha^–1). (B,G) Lorey’s height, basal area-weighted mean tree height (m). (C,H) AGB, aboveground biomass (Mg ha^–1). (D,I) mean annual mortality rate. (E,J) %M_treefall, fraction of mortality due to treefalls. Dashed lines represent simulation value with no wind disturbance, and solid lines represent field observations.

Effects of Topography on Wind Disturbances

Implementing the topographic effect, we found weak but significant relationships between topographic wind speed correction and forest structure. Results of linear regressions showed that stem density (DBH > 10 cm, N₁₀) was not significantly related to wind speed correction (p = 0.09), but average canopy height (Lorey’s height) decreased and aboveground biomass (AGB) increased significantly at quadrats with higher wind speed correction factors (representing more exposed terrain) (Figure 4 and Table 1). In contrast, the field-observed relationship between quadrat-scale forest structure statistics and elevation at Fushan was stronger and exhibited a different pattern: higher stem density and lower average canopy height were observed at more exposed quadrats, but biomass was not found to vary as a function of quadrat elevation (Figure 5 and Table 2).

FIGURE 4

Figure 4. Simulated quadrat-level forest structure statistics as a function of the topography-related wind speed correction factor of each quadrat. A correction factor > 1 means that the quadrat-level wind speed is accelerated,, and a correction factor < 1 means that the quadrat-level wind speed is decelerated. Shaded areas represent interquartile ranges (IQR), calculated within a moving window frame (0.025) across the whole x-axis value range and then linearly interpolated to the x-axis value at each quadrat. Black without border: without topography correction. Red with border: with topography correction. Solid lines represent the linear regression curve. (A) N₁₀, density of stems with DBH > 10 cm (trees ha^–1). (B) Lorey’s height, basal area-weighted mean tree height (m). (C) AGB, aboveground biomass (Mg ha^–1). For all structure statistics, the quadrat-wide values are converted to the corresponding hectare-wide values.

FIGURE 5

Figure 5. Quadrat-wide (20 m × 20 m) forest summary statistics at the Fushan site, as a function of the quadrat elevation above sea level. Shaded areas represent interquartile ranges (IQR), calculated within a moving window frame (2 m) across the whole x-axis value range and then linearly interpolated to the x-axis value at each quadrat. Solid lines represent the linear regression curve. (A) N₁₀, density of stems with DBH > 10 cm (trees ha^–1). (B) Lorey’s height, basal area-weighted mean tree height (m). (C) AGB, aboveground biomass (Mg ha^–1). The quadrat-wide values were converted to corresponding hectare-wide values.

TABLE 2

Table 2. Slope estimates of linear models of each quadrat-level forest structure metrics to the wind speed correction factor, with p-Values in parentheses.

Discussion

In this study, we included a model of wind-induced tree mortality in an individual-based forest model, and simulated how wind-induced tree mortality affects the structure and dynamics of a subtropical forest. We found that wind-induced tree mortality had a large negative impact on canopy height and a more complex influence on stem density and biomass. The comparison of tree mortality at the onset of wind disturbance and at the steady state reveals forest adaptation to repeated wind disturbances. Wind intensity was found to exert a strong control on wind impact, while wind frequency did not. Implementation of topographic heterogeneity showed a weak but significant effect on within-stand canopy height and biomass at the Fushan site: this implementation could serve as a basis to incorporate more complex wind-terrain interactions in individual-based models.

Effects of Wind-Induced Tree Death

In response to increasing susceptibility to wind-induced tree death, canopy height was found to decrease sharply. According to the wind damage model, T_C is proportional to DBH² × H (Eq. 6), and CWS_stembreak is proportional to (DBH³ / T_C)^1/2 (Eq. 4): it thus follows that CWS_stembreak is proportional to (DBH / H)^1/2. Hence, wind damage risk can be expected to be higher for trees with a smaller diameter-to-height ratio (Cremer et al., 1982), and the observed reduction of canopy height may be caused by the fact that taller trees are in general more exposed and less protected by the sheltering of neighboring trees, meaning that they experience higher tree-top wind speed and higher wind damage risk (Hale et al., 2012; Duperat et al., 2021). A corollary of this hypothesis is that tree height in the canopy may become more homogenous as the effect of wind-induced tree mortality strengthens (Van Bloem et al., 2007; Chi et al., 2015). The hump-shaped response of biomass, although unexpected, might be the result of the joint effect of the selection of trees with denser wood and the reduction of larger trees, which contribute the most to total forest biomass.

A transition from forest to non-forest state around P = 0.3 was apparent, where stem density increased abruptly as canopy height and biomass decreased to near-zero. This increase in stem density likely reflects a light condition that favors more smaller trees to establish as the forest canopy opens up. As the effect of wind disturbances grew even stronger, even these small-height trees started to be affected due to the absence of sheltering from taller trees, causing stem density to decrease as well.

Steady-state total mortality increased with disturbance intensity, but only starting from around the transition point. Steady-state fraction of mortality due to treefall decreased as disturbance intensity increased until after the transition point, where it increased drastically. In addition, both mortality statistics were markedly lower at the end of the simulation than at the onset of wind disturbance. These observations suggest that long-term exposure to wind disturbances results in modification to the forest structure through the higher survival of more resistant species, such as those with larger diameter-to-height ratio or higher wood density. This forest structural change could account for the low immediate mortality after a cyclone passage observed in some empirical studies (Walker, 1991; Bellingham et al., 1992).

In addition to changes in structure and species composition, plastic acclimation could also regulate forest response to wind disturbances, as shown by a modeling study (Kamimura et al., 2019). In this study, plastic acclimation was not taken into account as the model did not allow traits to vary plastically over the lifetime of an individual, and trait variations were not inheritable. In the future, it would be highly interesting to explore the eco-evolutionary dynamics that results from temporally plastic trait variability and the inheritance of more wind-resistant traits (Telewski and Moore, 2016).

Effects of Wind Frequency and Intensity

Increasing wind intensity negatively affects forest stature and biomass stocks, consistent with results of global-scale studies that showed a clear empirical relationship between higher cyclone intensity and lower overall forest stature (Ibanez et al., 2019); however, varying wind frequency did not cause any effect on forest structure. Forest response to disturbance must be considered by placing the return period of the disturbance events in perspective with the life span of the affected organism (Turner and Dale, 1998): as trees can live for decades or even centuries in a forest, they can be exposed to centennial disturbance events. Therefore, a forest could be capable of adapting to frequent cyclone visits without suffering catastrophic loss, as long as their intensity is within the tolerance range of the forest. On the other hand, extremely intense cyclones, even when occurring on a centennial basis, could be sufficient to cause severe consequences on forest structure and potentially even threaten the very persistence of the forest ecosystem.

Some global-scale comparative studies have demonstrated a relationship between structural or functional composition of forests and frequency of cyclone events (Hogan et al., 2018; Ibanez et al., 2019). Since in real life the frequency of cyclone events and the average experienced cyclone intensity are highly correlated, this observed relationship could be partially due to variation in intensity as well. On the other hand, if we take into account the phenomenon of delayed mortality due to wind-induced partial tree damage, there could potentially be a genuine effect of cyclone frequency. Finally, in this study, we tested the model at a forest site which has experienced high cyclone frequency in the past: it is possible that a forest site which has experienced low cyclone frequency in the past could respond more strongly to an increase in cyclone frequency. In light of predictions on future increase in tropical cyclone intensity (Knutson et al., 2020), this result suggests that forest ecosystems in cyclone-prone areas may risk undergoing profound structural modifications or even state shifts, and more intensive and preemptive monitoring would be needed to anticipate the consequences (Newman, 2019). Modeling scenarios of simultaneous modifications of wind frequency and intensity, as well as gradual changes of the wind regime over the course of the simulation (non-stationary extreme winds) should also be explored, in order to improve the realism of disturbance forcing.

Effect of Topography

Topography alters the impact of extreme wind on forests in complex ways. Many studies have found that exposed locations tended to be subject to more severe wind damage than sheltered locations (Ruel et al., 1998; Magnabosco Marra et al., 2014). In this study, we introduced a simple correction factor to model the effect of quadrat-wide topographic heterogeneity, based on the assumption that wind speed accelerates locally on exposed terrain such as hilltop and mountain ridges, and decelerates locally in sheltered terrain such as valleys (Belcher et al., 2011; Mitchell, 2013). The simulation results showed a weak but significant effect of topographic heterogeneity on canopy height and biomass stocks within the forest plot at Fushan, with higher average canopy height and aboveground biomass at more exposed quadrats. The increase of aboveground biomass at more exposed quadrats was consistent with the hump-shaped response of biomass as wind disturbance strength increased (Figure 1). However, field observations were not consistent with the prediction, showing higher stem density and lower average canopy height at higher quadrats, and no significant relationship between biomass and quadrat elevation (Figure 5 and Table 2).

The discrepancy between the predicted pattern and field-observed pattern strongly suggests that topographic effect on wind-induced mortality is not a simple function of elevation, and that it will be necessary to refine the model representation of wind-topography interaction. This could be improved, for example, by taking into account canopy surface model and wind direction, as well as tree acclimation to the local wind environment through resource reallocation and modification of allometry (Eloy et al., 2017; Dèfossez et al., 2021).

In the current study, wind speed correction was based on Global Wind Atlas data generated by the WasP model, which is the most applied model for predicting topographic effects on the local wind climate, in particular for the wind energy industry (Finnigan et al., 2020). The model is based on an analytical solution of wind flow over hills originally formulated by Jackson and Hunt (1975). It works best in regions with moderate slopes (≤1:10), which is the case at the Fushan site. Future work will require evaluation of the suitability of GWA data to represent spatial variation of wind speed at forest sites of various steepness and aerodynamic roughness.

Nevertheless, this preliminary exploration demonstrated how spatial variation of wind speed could be coupled to wind-induced tree mortality in an individual-based model, and could serve as a framework under which to further investigate how topography mediates the effect of wind disturbances.

Challenges of Model Representation of Wind-Induced Tree Mortality

The first and foremost challenge when simulating wind disturbance in a mixed-sized forest is the description of the wind profile: since the individual-based model TROLL does not prescribe fixed stand-level characteristics (notably top canopy height), many factors that control wind profile dynamically change across time and space. Wind speeds above the tree canopy are commonly modeled by a logarithmic profile (Raupach et al., 1991), which is the approach taken in the ForestGALES model. The aerodynamic parameters used in the logarithmic profile have been shown to depend on the plant area index, and as such, they are expected to vary seasonally. A detailed parameterization has recently been proposed based on remote-sensing LAI and canopy height products, which could help account for this variability and enhance realism of wind profile modeling (Liu et al., 2021).

In this study, we further assumed that the wind speed experienced by trees within the canopy follows an exponentially decreasing profile (Inoue, 1963), and that trees well below the displacement height in effect do not experience wind disturbance; other equations have been derived to describe the decreasing wind speed profile within the canopy (Raupach and Thom, 1981; Ancelin et al., 2004), and multiple measurements have also been made in different canopies (Raupach et al., 1996). Significant variations in wind profile could exist due to the heterogeneity of forest structural characteristics, such as leaf area index, stem density and spacing (Lalic et al., 2003; Moon et al., 2016; de Santana et al., 2017). In addition, the interaction of horizontal winds with the canopy structure creates turbulent eddy structures, whose effects have not yet been fully explored in the current model (Raupach et al., 1996; de Langre, 2008). For example, the gust structure of the wind could be accounted for by including a gust factor, as in the original ForestGALES model (Quine et al., 2021).

The original ForestGALES model estimated critical wind speed for both types of wind-induced damage, tree uprooting and stem breakage. In a preliminary test, we simulated both processes, but found that few tree uprooting occurred compared to stem breakage (Supplementary Appendix C). This is consistent with the empirical mortality survey at Fushan, where many more dead trees exhibit stem breakage than uprooting (Supplementary Appendix C). This has motivated our decision to simplify model representation and focus on stem breakage. One further motivation for ignoring tree uprooting is that it is largely controlled by root anchoring, which is currently still difficult to represent mechanistically and to parameterize for all tree species, due to our limited understanding on root traits and the physics underlying root-soil interaction. In the future, efforts should be made to devote more attention to the process of root anchoring.

Lastly, although we parameterized species-specific wind susceptibility using plant traits to the extent possible, and assumed that the free parameter P was identical for all the trees in a forest stand, P may in reality still be species-specific due to factors such as wood deformities, stem taper, relative allocation of total biomass to stem weight, as well as the capacity for defoliation that reduces wind loading. These factors could all contribute to model uncertainty, and further investigation is needed to better constrain them in future model developments.

Field Mortality Survey Data

In this study, we retrieved data from the annual mortality surveys conducted in ForestGEO sites, in the hope of calibrating the model using empirical mortality data as a complement to inventory data. The survey data could inform us on the number of tree deaths and treefall that occur within the forest plot, and simulated mortality rates showed good correspondence with observed values calculated from the survey data, but inferring the cause of tree death from the observed damage modes of the dead trees has proven to be difficult. A tree may experience multiple damage successively or simultaneously, which would all contribute to causing treefall, or it may experience damage after its death: these different scenarios could not be distinguished by post hoc records of observed damage on dead trees. In addition, a large proportion of stem breakage and uprooting in the understory are likely to be caused by treefall events of neighboring trees, and not directly by wind.

As a consequence, it is possible that the field observed numbers of tree deaths with uprooting or stem breakage were an overestimation of the wind-induced uprooting or stem breakage events in reality. However, the observed fraction of tree deaths due to treefalls were considerably lower than the simulated values, suggesting that there might be additional acclimation and adaptation phenomena limiting treefalls, which were not accounted for in the current model. Nevertheless, the annual mortality survey data still contain the most detailed fine-scale records on tree mortality in species-rich natural forests available to date, and the observed patterns of damage modes, while not sufficient to calibrate the model, still served as a heuristic basis that provide information on patterns of tree death.

Perspectives

In a mixed-sized natural forest, natural thinning and gap dynamics frequently modify tree density and local stand height, causing changes in the level of physical sheltering from neighboring trees: the wind loading and damage risk of a tree may change considerably as its sheltering status changes even when its size does not (Schelhaas et al., 2007; Seidl et al., 2014; Quine et al., 2021). In the current implementation of the wind damage model, local sheltering is accounted for by considering a tree’s height with respect to local canopy top height, simulating within-canopy wind attenuation, and assigning a height threshold under which a tree is not considered to experience wind disturbances. A finer representation of the effect of local sheltering could be achieved by including a “competition index,” calculated based on a tree’s relative size to neighboring trees, in the formulation of the tree’s turning moment (T_C), so that local sheltering is represented as a continuous variation and not as a single cutoff point (Hale et al., 2012; Duperat et al., 2021), even though the choice of competition index is not trivial and requires careful consideration. In addition, the T_C equation is currently empirically fitted to a limited number of temperate plantation forest sites, and its transferability to the tropics and to mixed-sized natural forests needs to be examined in more detail: recent works of fine-scale measurement of tree movement in response to wind in natural forests represent the first step to overcoming this challenge (Jackson et al., 2021).

The current model could be further refined by including other aspects of wind impact on trees, such as successive damage that does not immediately cause death but increases delayed mortality (Walker, 1995; Tanner et al., 2014), coping mechanisms such as defoliation (reduction of wind drag at the cost of temporary lower productivity) and re-sprouting (allowing survival even after stem breakage). Other factors influencing tree mortality, such as preexisting stem rot or deformities, interactions with other forms of disturbance such as insect attacks, fire and drought (Seidl et al., 2011; Reichstein et al., 2013; Newman, 2019). Human land use legacy and fragmentation (Laurance and Curran, 2008; Uriarte et al., 2009; McGroddy et al., 2013; Schwartz et al., 2017) could also be considered conjointly with the wind disturbance model: this could potentially be useful for forest conservation and management, by evaluating how forest response to disturbance is modulated by different resource use and management strategies (Hogan et al., 2020). Parameter calibration and model validation could also be improved with the help of high-resolution satellite data that monitor wind gap formation and dynamics (Hayashi et al., 2015; Kislov and Korznikov, 2020; Ballère et al., 2021).

Finally, although the current study focuses on the impacts of wind on a cyclone-prone forest, the model we developed could also be applied to explore wind vulnerability of forests that are less accustomed to wind disturbances, as well as the effects of localized wind blowdown events, which are thought to shape the structure and dynamics of Amazonian forest (Magnabosco Marra et al., 2018; Negrón-Juárez et al., 2018; Peterson et al., 2019). The exploration of more forest settings and wind regimes could help us go beyond single forest plots and take into account landscape- and regional-level heterogeneity in the model (Seidl et al., 2014; Peereman et al., 2020), in order to explore the consequences of wind disturbance at a landscape or even regional level.

Conclusion

In this study, we explored the long-term effects of wind-induced tree mortality to a subtropical forest. The results indicate that wind disturbance could have strong negative effects on forest structure when intensity is high, which has important implications given the projected increase in tropical cyclone intensity. This modeling framework of wind-induced tree mortality, including the preliminary implementation of topographic effects, could serve as the basis for improving representation of the mortality processes in vegetation models, deepen our mechanistic understanding of how wind disturbance acts on forests over a large spatial scale in conjunction with other types of disturbance, as well as generate predictions on the future of natural forests in response to the changing wind regime.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

Author Contributions

EPR, BG, FF, IM, and JC contributed to conception and design of the study. EJ processed and provided climate forcing datasets. I-FS provided Fushan field census data. EPR performed the model simulations and wrote the first draft of the manuscript. EPR, BG, and JC wrote sections of the manuscript. All authors contributed to manuscript revision, read, and approved the submitted version.

Funding

This work was supported by the “Investissement d’Avenir” grants managed by the Agence Nationale de la Recherche (CEBA, ref. ANR-10-LABX-25-01; TULIP, ref. ANR-10-LABX-0041; ANAEE-France: ANR-11-INBS-0001). This work was also funded by the administrative region of Occitanie (France) and Centre National d’Études Spatiales (CNES).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

Simulations were run on the OLYMPE cluster at CALMIP, Toulouse. We thank Dr. Chia-Hao Chang-Yang at the National Sun Yat-sen University for helping to provide us with the mortality survey data at Fushan FDP. The data at Fushan FDP results from long-term fieldwork and research supported by the Taiwan Forestry Bureau, the Taiwan Forestry Research Institute and the Ministry of Science and Technology of Taiwan. We would like to express our gratitude to all field technicians and students who helped with the implementation and recensus of Fushan Forest Dynamics plots. We also thank the staff at Fushan Research Center for providing logistic support.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/ffgc.2022.753100/full#supplementary-material

Footnotes

References

Ancelin, P., Courbaud, B., and Fourcaud, T. (2004). Development of an individual tree-based mechanical model to predict wind damage within forest stands. Ecol. Manage. 203, 101–121. doi: 10.1016/j.foreco.2004.07.067