Canopy Density and Roughness Differentiate Resistance of a Tropical Dry Forest to Major Hurricane Damage

Abstract

Tropical dry forest is vulnerable to increased climate variability with more frequent and severe storms. Studies of hurricane impact on tropical dry forest often focused on individual tree traits. How trees in tropical dry forests work together to combat wind damage is still unclear. To address this, we integrated ground-observed ecosystem structure from National Ecological Observation Network (NEON) with airborne-LiDAR images and analyzed resistance in forest structure of Guánica dry forest in Puerto Rico to major hurricanes in 2017 at the forest-stand level. Using each plot instead of the individual tree as the base unit, we regressed mean changes in stem height and fractions of lost or damaged stems at 15 plots on mean stem diameter, mean and standard deviation of stem height, stem density, and topography. Meanwhile, using the LiDAR-derived canopy heights, we compared the changes in canopy height before and after the hurricanes and regressed spatially the canopy height change on prior-hurricane tree cover, canopy height, and rugosity. We found that the damage was small in places with high stem density or high tree cover. Ground-observed damage in terms of height reduction significantly increased with the standard deviation of stem height, an index of roughness, but decreased with the mean stem diameter of the plots. LiDAR-detected damage in terms of reduction in canopy height was also found to decrease with tree cover and mean canopy height when the canopy height was small or moderate but increase with the rugosity. The fraction of lost stems significantly decreased with the stem density, and the fraction of damaged stems significantly increased with the roughness and the plot elevation. The collective parameters of forest stand quantified from ground-observation and LiDAR, such as stem density, tree cover, and canopy roughness or rugosity, highlighted mutual supports of trees and played important roles in resisting damages to the tropical dry forest during major hurricanes.

1. Introduction

Tropical and subtropical dry forest accounts for 42% of the global tropical and subtropical forest cover [1] and spans climate regions with annual rainfall ranging from 250 to 2000 mm, annual mean bio-temperature greater than 17 °C, and the annual ratio of potential evapotranspiration to precipitation greater than 1 [2]. The C stock of tropical dry forest accounts for 8–9% of the global estimate, yet the net primary production was 620 g C m⁻², which is only lower than that of the tropical wet forest [3,4]. There are always 2–6 months of dry season, and the water-limited nature makes tropical dry forests especially sensitive to climate variability and extreme climate events such as drought and storms [5,6,7,8].

The structure of the tropical dry forest is substantially different from that of the tropical wet and moist forests. The dry forest has, in general, lower tree richness than wet forest, and the low diversity is especially found in island forests [1]. According to the negative exponential relationship between root to shoot ratio and annual rainfall [9], the tropical dry forest has a larger root:shoot ratio (0.57–0.71) than moist and wet forests (0.17–0.35) to facilitate water and nutrient uptake [10], which might help resist root–soil failure during strong winds [11,12]. Limited by the hydraulic transporting capacity, the tropical dry forest also has a smaller stem height to diameter ratio than wet and moist forests have. For diameters in the ranges of 10–20 cm and 21–40 cm, the mean heights were 11.8 and 18.1 m for the tropical dry forest, 14.7 and 23.3 m for the moist forest, and 16.7 and 25.6 m for the wet forest, respectively [13]. Lower diversity, larger root:shoot ratio, and lower height:diameter (taper, lower slenderness) ratio make tropical dry forests differ from wet or moist forests in response and resistance to tropical storms.

Resistance to tropical storms and recovery of the tropical dry forest after disturbance also depend on other structural and functional traits. For example, wood density correlates to rupture modulus and elasticity modulus, and high wood density, associating with mechanical strength, might indicate less likely the wind snap (e.g., trunk breakage) but more likely the windthrown (e.g., uprooting) [14]. A tree species with a greater specific leaf area (SLA) may intercept more wind drag forces to rupture the stems and branches; however, thinner leaves are also easily torn by the wind and may be more advantageous to avoid damage to branches and trunks. On the other hand, a species with a smaller SLA has a longer leaf life span, stronger leaves, stronger ties with branches and twigs, and less likely to shed its leaves [15]. Since wind-induced load can be increased 2–3 times with foliage, compared to that without foliage, mechanically strong leaves make them hard to shed and may lead to severe damage to branches and trunks during storms [16]. After hurricane disturbance, tree species with greater SLA have a higher relative growth rate [17]. SLA and maximum height are positively associated with resprouting and recovery, whereas wood density and slenderness are negatively associated with resprouting and recovery [7].

Recent studies emphasized the important role of stand-level parameters in wind damages to forests [11,18]. Hart et al. [18] found that the random forest models on wind damage to forests provided high accuracy and the discriminatory power remained high with the removal of individual-level characteristics; however, the accuracy and the discriminatory power would be lost if the stand-level information is removed. Compared to a single-tree situation, trees in a closed forest stand provide sheltering effect and/or mechanically mutual supports during wind disturbance [11]. Moore and Lin [19] reported that stand density and height:diameter ratio are the most important factors associated with wind damage in terms of basal area loss. Compared to dense stands, winds inside a stand with widely spaced trees are higher [11]. Forest structure also affects its vulnerability to wind damage [12]. By analyzing the impacts of Typhoon Songda in northern Japan using the random forest algorithm, Morimoto et al. [20] reported that the natural mixed forest stand with a complex structure was less vulnerable to the catastrophic wind compared to the Abies plantation with a simple structure.

Monitoring wind damage in forests is mostly based on small-scale ground observations [7,21,22] and large-scale optical remote-sensed images [23,24]. Based on the fused vegetation index from Sentinel-2 and Landsat-8, the dynamic of forest greenness, indicating ecosystem functions, was assessed to explore how topography, drainage capability, and legacy of drought differentiated tropical ecosystem response to and recovery from major hurricanes [24]. Compared to ecosystem functions, ecosystem structures showed much slower recovery after hurricanes [22] and divergence [25]; thus, structural changes are better indicators of hurricane impacts in the short to long term. However, structural changes especially at scales larger than plots are still rare in existing studies. The recent rapid advance in remote sensing technology, particularly LiDAR, provides us necessary data and metrics such as canopy height, cover, or rugosity [26,27] to monitor structural changes at a larger scale and to explore underlying mechanisms of forests in response to storms.

The Caribbean region is prone to tropical storms. In September 2017, two successive major hurricanes impacted the region [24,28]. Recent research mostly focused on the impact of hurricanes on windward wet and moist forests [21], while the study of the impact on the leeward tropical dry forest is largely left untouched. On the other hand, compared to the trait-based individual-level studies mentioned above, the quantitative roles of stand-level parameters of tropical dry forest, such as stem density and canopy roughness, in the ecosystem resistance to tropical storms remain unclear.

This study aims to explore how the two major hurricanes impact a tropical dry forest at stand level by integrating damages to trees based on the National Ecological Observation Network (NEON) ground observations with changes in canopy height within the forest via the airborne LiDAR observations. We associated the impact of the hurricanes to the ground-measured stem density and stem diameter, unevenness of the stem height, and topography such as elevation, aspect, or topographic position, as well as to the LiDAR-sensed canopy height, tree cover, and rugosity. We hypothesized that more unevenness of the canopy or rugosity and higher plot elevation tend to bring more damage, whereas higher tree cover, higher stem density, and larger stem diameter tend to resist the damage. On the other hand, the rain brought by hurricanes may alleviate the drought stress during the years before the storms.

2. Materials and Methods

2.1. Study Area

NEON Guánica site (−66.86°W, 17.97°N, Figure 1) is located within the dry-forest life zone in southwest Puerto Rico near the coast of the Caribbean Sea. The site has a mean annual precipitation of 860 mm and a mean annual temperature of 25 °C [1]. Guánica dry forest was designated as a state forest in 1919 and a UNESCO International Biosphere Reserve in 1981. In 2017, Hurricane Irma passed by the tropical island Puerto Rico on 6 September with the sustaining wind speed of 298 km h⁻¹, and Hurricane Maria traversed the island on 20 September with 249 km h⁻¹ at its landfall, incurring enormous damage to the tropical ecosystem over the island (Figure 1). Prior to the two hurricanes, Puerto Rico experienced a severe drought in 2014–2016 with the island-wide water rationing lasting for months [29]. Historically, Hurricane Georges in 1998 and Hurricane San Felipe in 1928 affected the structures of the Guánica dry forest [30].

2.2. NEON Plot Ground Data Analysis

NEON started the forest structure observation in 2015 at scattered plots (Figure 1). The 20 NEON tower plots are distributed close to the instrument tower, and each plot is 40 m × 40 m with two subplots sampled, each of which is 20 m × 20 m. The 20 distributed plots of 20 m × 20 m each are scattered in the farther distance but around the tower. We acquired the forest structure observations from the NEON sites [31], and the data showed that among the 40 plots, only 10 distributed and 5 tower plots had measurements both before and after the hurricanes. The distribution of these plots almost covered all the elevation ranges from lowland, mid slope, to top of the hill (Figure 1,

Table 1. The 15 plots used for vegetation structural change assessment before and after the major hurricanes in September 2017.

Plot ID	Plot Type	Time Before	Time After	Altitude (m)
GUAN_002	Distribute	August 2015	January 2018	89.1
GUAN_003	Distribute	April 2016	December 2017	182.9
GUAN_004	Distribute	February 2016	January 2018	133.5
GUAN_008	Distribute	March–April 2016	January 2018	65.8
GUAN_009	Distribute	May 2016	January 2018	116.5
GUAN_010	Distribute	May 2016	January 2018	177.6
GUAN_011	Distribute	May 2016	February 2018	201.2
GUAN_015	Distribute	April 2016	December 2017	207.8
GUAN_018	Distribute	February 2016	January 2018	162.6
GUAN_019	Distribute	May–June 2016	February 2018	131.7
GUAN_042	Tower	April 2017	February 2018	146.9
GUAN_043	Tower	April 2017	February 2018	115.9
GUAN_044	Tower	April 2017	February 2018	147.5
GUAN_045	Tower	April 2017	February 2018	145.3
GUAN_046	Tower	March 2017	February 2018	168.0

).

NEON woody vegetation structural measurement includes stem height, diameter (mostly DBH ≥ 10 cm, diameter at breast height, https://data.neonscience.org/data-products/DP1.10098.001, accessed on 8 June 2021), crown diameter, etc. For this particular site with small trees and dense canopy of the dry forest, i.e., greater than 3000 stems per hectare, stem diameter ≥ 0.8 cm were measured. Among these metrics, stem height and diameter were mostly measured.

For each plot, we selected those individual stems labeled with “Live” as the plant status and stem height measured before the hurricanes and compared with their plant status and stem height after the hurricanes. We then analyzed these changes in plant status and stem height for each of the 15 plots and explored the underlying mechanisms that lead to the changes. The plant status after the hurricanes can be grouped into three categories: “Live”, “Lost”, and “Damaged”, and the latter two reflect a change in status. The “Lost” category included plant status marked as “Standing dead”, “Removed”, “Lost, fate unknown”, “Dead, broken bole”, and “Lost, presumed dead”; and the “Damaged” category included those marked as “Live, broken bole”, “Live, physically damaged”, and “Live, other damage”. The five tower plots were measured in early 2017 and early 2018 (Table 1). Nine out of the ten distributed plots were measured in 2016 and the end of 2017 to early 2018, and the remaining one (GUAN_0002) was measured in 2015 and early 2018. We assume the comparison between the measurements before and after is mostly caused by the major hurricanes in September 2017, but to a less extent, by the growth after the pulse rainfall brought by the hurricanes.

Multiple and logistic (quasibinomial) regressions were applied to explore the possible factors that contribute to the loss and/or damage of stems and the change in stem height. Quasibinomial is often used to eliminate the overdispersion of a generalized linear model. The explanatory variables are plot elevation and topographic position, the density of stems, the mean and standard deviation of stem height, and mean stem diameter. The topographic position index, representing the topographic position, was calculated as the elevation minus the average elevation of a circular neighborhood with a radius of five cells.

2.3. Canopy Height Change by Airborne Light Detection and Ranging (LiDAR)

Changes in canopy height were examined with the LiDAR-derived canopy height model (CHM). The NASA G-LiHT mission collected the LiDAR images along a few transects across the Guánica dry forest on 8 March 2017 [32], and the NEON Airborne Observation Platform (AOP) collected the LiDAR images covering the whole site on 28 May 2018. We selected the transect of PR_8March_2017_45_46, which crosses the forest and spans all the elevation range (Figure 2). The section of the transect within the NEON Guánica site was chosen for further analyses. The vertical and horizontal accuracies of both NEON [33,34] and G-LiHT images (glihtdata.gsfc.nasa.gov/, accessed on 8 June 2021) are within 1 m. Aside from the CHM images, we also used the digital surface model (DSM), digital terrain model (DTM), topographical slope and aspect, and rugosity based on the DSM. All images have a spatial resolution of 1 m. We obtained the data from the NEON data portal with the “neonUtilities” R Package and from the G-LiHT site (glihtdata.gsfc.nasa.gov/, accessed on 8 June 2021).

Both the NEON and NASA CHM were calibrated with the ground-measured stem height of the NEON plots located within the transect before performing the statistical comparison. Stem heights measured in early 2017 and early 2018 (from 1 January to 30 May) were used to calibrate the CHM of G-LiHT before the hurricanes and the CHM of NEON after the hurricanes, respectively. Changes in canopy height were obtained by subtracting the calibrated CHM before the hurricanes from the calibrated CHM after the hurricanes. We also calculated the tree cover as the proportion of pixels with canopy height above 2 m. Changes in canopy height, together with the prior-hurricane CHM, DSM, slope, aspect, tree cover, and rugosity, were averaged for a grid of 50 m within the transect.

A spatial error model was then applied to find whether and how the canopy height change depends on the prior-hurricane covariables including topography. The spatial error model [35] considers the spatial autocorrelation of error of a linear regression model as follows:

$$y=x\beta +u\hspace{1em}\hspace{1em}u=\lambda Wu+ϵ$$

where y is the dependent variable represented as an n × 1 matrix (n, number of cases in the dataset), x is one or multiple (m) independent variables represented by an n × p matrix with p = m + 1, β is a p × 1 coefficient matrix to be estimated, and u is an n × 1 spatial error matrix. W is an n × n weight matrix determined by the neighborhood structure, λ is the spatial autocorrelation coefficient, and ε is the normal i.i.d. (independent and identically distributed) residual.

Before proceeding, we averaged all these variables over a 50 m × 50 m grid. We then divided canopy height change by the canopy height prior to the hurricanes to obtain the relative height change ($\mathsf{\Delta }{h}_{cr}$), calculated a function of aspect as

$$f\left(\alpha \right)=-\cos \left(\alpha -\pi /4\right)$$

with α as the aspect angle, and normalized rugosity (r), elevation (Z), and canopy height (h) to [0,1] range. Tree cover is already in the range, and cosine of aspect is in the range from −1 to 1. We then regressed spatially the relative canopy height change on the normalized covariates with the spatial error model. The spatial correlation of errors was defined within eight nearest neighboring grid cells. Variables were selected with the ANOVA model comparison. Data preparation and analyses were performed in ArcGIS (ESRI, Redlands, CA, USA) and R [36]. Error assessment was performed using R² for the multivariate regression and the spatial error model and using the residual deviance versus the null-model (model with intercept only) deviance for the quasibinomial regression.

3. Results

3.1. Mean Change in Stem Height at Plot Level

The mean height change in the 15 plots showed that the reduction (negative change) prevailed (Figure 3) over the stimulated growth by the pulse rainfall brought by the hurricanes, 11 out of 15 plots exhibiting a reduction in mean stem height and 7 of them being significant (t-test at α = 0.05). The average mean height change at plot level was −0.14 m with a standard deviation of 0.22 m. The regression of the mean height change on topography, stem density, mean stem diameter, mean stem height, and standard deviation of the stem height measured before the hurricanes resulted in the equation via a stepwise selection of independent variables as follows:

$$\mathsf{\Delta }\overline{h}=-0.27\left(\pm 0.112\right)-0.264\left(\pm 0.1\right)s_h+0.094\left(\pm 0.022\right)\overline{D}$$

(1)

where $\mathsf{\Delta }\overline{h}$ is the mean stem height change with negative values for damage, $s_h$ is the standard deviation of stem height within the plot, and $\overline{D}$ is the mean stem diameter. Standard errors of the estimated parameters are shown in the parentheses. The regression yielded an $R^2$ of 0.59, and all the coefficients and the F-statistic were significant (p-value < 0.05). $s_h$ is a measurement of canopy unevenness connecting to the surface roughness. The equation implies that the roughness of the canopy tends to bring more reduction in stem height, whereas plot with trees having greater stem diameters tends to resist the reduction in stem height. Pooling the individual stems from all the 15 plots, we found a significant correlation coefficient of 0.157 between the height change and the stem diameter before the hurricanes, which means the larger the stem diameter, the less reduction was observed in height. Positive changes in plots 04, 08, 43, and 44 showed that the post-hurricane growth benefited from moisture brought by the storms exceeded the minor hurricane damage in these plots (Figure 3).

3.2. Lost and Damaged Stems

The number of affected stems (Figure 4), including the categories of “Lost” and “Damaged”, showed a wide range of variation with the mean ± standard deviation as 4 ± 3.7 for the number of lost stems and as 4.6 ± 4.1 for the number of damaged stems. The mean of the lost plus damaged stems was 8.6 (±4.9). The proportion of lost or damaged stems, that is, the number of lost or damaged stems divided by the total number of stems of the plot (Figure 4 lower panel), showed that the mean ± standard deviation of the fraction of the lost stem and that of the damaged stem were 2.2% ± 1.7% and 2.3% ± 1.9%, respectively. On average, 4.5% of plant stems were affected by the hurricanes.

After regressing the fraction of lost stems in the plots on pertinent covariates including plot topography, the density of stems, the mean and standard deviation of stem height, and mean stem diameter with assumed quasibinomial (logit) distribution, we obtained the following equation after stepwise selection of the covariates:

$$p_{Lost}=-2.87\left(\pm 0.21\right)-1.38\left(\pm 0.45\right)d$$

(2)

where $p_{Lost}$ is the fraction of lost stems, and $d$ is the stem density. The regression yielded residual deviance Equation (2) of 8.56, in comparison with the null deviance (the null model with intercept only) of 16.88. All the coefficients were significant (p-value < 0.01). Therefore, the proportion of lost stems decreased significantly with the stem density. Similarly, the regression of the fraction of damaged stems gave,

$$p_{Damage}=-6.79\left(\pm 1.36\right)+1.12\left(\pm 0.33\right)s_h+0.0082\left(\pm 0.0062\right)Elev$$

(3)

where $p_{Damage}$ is the fraction of damaged stems and $Elev$ is the elevation at the plot centroid. The regression had residual deviance of 22.68, in comparison with the null deviance of 46.4. All coefficients were significant at α = 0.01, except the coefficient in front of Elev which had a p-value of 0.21. Therefore, greater canopy unevenness and higher elevation were more likely to result in more damage.

3.3. Impact at Species Level

Table 2. The affected species (encountering lost or/and damaged stem), number of stems, number of affected, and proportion of affected in the 15 plots at NEON Guánica dry forest site.

Scientific Name	# Stems	# Affected	% Affected
Gymnanthes lucida Sw.	1714	18	1.1
Croton lucidus L.	537	15	2.8
Bucida buceras L.	610	13	2.1
Eugenia foetida Pers.	997	8	0.8
Thouinia striata Radlk. var. portoricensis (Radlk.) Votava & Alain	490	8	1.6
Leucaena leucocephala (Lam.) de Wit	487	6	1.2
Pisonia albida (Heimerl) Britton ex Standl.	541	6	1.1
Pithecellobium unguis-cati (L.) Benth.	171	6	3.5
Bursera simaruba (L.) Sarg.	305	5	1.6
Guaiacum sanctum L.	181	5	2.8
Coccoloba microstachya Willd.	223	3	1.3
Acacia farnesiana (L.) Willd.	82	2	2.4
Bunchosia glandulosa (Cav.) DC.	39	2	5.1
Capparis flexuosa (L.) L.	54	2	3.7
Colubrina arborescens (Mill.) Sarg.	29	2	6.9
Croton betulinus Vahl	6	2	33.3
Guaiacum officinale L.	262	2	0.8
Pictetia aculeata (Vahl) Urb.	173	2	1.2
Poitea florida (Vahl) Lavin	42	2	4.8
Rochefortia acanthophora (DC.) Griseb.	25	2	8.0
Schaefferia frutescens Jacq.	58	2	3.4
Swietenia mahagoni (L.) Jacq.	562	2	0.4
Capparis hastata Jacq.	112	1	0.9
Capparis indica (L.) Druce	22	1	4.5
Coccoloba diversifolia Jacq.	207	1	0.5
Comocladia dodonaea (L.) Urb.	35	1	2.9
Croton flavens L.	7	1	14.3
Croton sp.	40	1	2.5
Eugenia xerophytica Britton	30	1	3.3
Exostema caribaeum (Jacq.) Schult.	66	1	1.5
Guapira obtusata (Jacq.) Little	10	1	10.0
Guettarda krugii Urb.	122	1	0.8
Guettarda sp.	32	1	3.1
Lantana exarata Urb. & Ekman	81	1	1.2
Reynosia uncinata Urb.	39	1	2.6
Unknown plant	526	1	0.2

lists the species names with stems either lost or damaged (affected), the number of stems for each affected species in the site, the number of affected, and the percent of affected calculated as the number of affected divided by the number of stems. One of the species is marked as an unknown plant, and two of them have only their genus known. Due to the lack of complete information about the structural and functional traits, we could not analyze the relationship between impact and traits. However, the most affected species such as Gymnanthes lucida Sw., Eugenia foetida Pers., Bucida buceras L., Croton lucidus L., Thouinia striata, and Guaiacum sanctum L. are recorded as hardwood species (http://tropical.theferns.info/, accessed on 8 June 2021). The six species constituted more than 52% of all the affected stems. From the global wood density database [37], we found the wood density of 10 major species including Bursera simaruba (L.) Sarg., Cocoloba diversifolia Jacq., Exostema caribaeum (Jacq.) Schult., Guaiacum officinale L., Guaiacum sanctum L., Gymnathes lucida Sw., Leucaena leucocephala (Lam.) de Wit, Pictetia aculeata (Vahl) Urb., Swietenia mahagonia (L.) Jacq., and Thouinia striata. However, the percent of affected (fourth column in Table 2) did not have a significant correlation with the wood density of the 10 species. In the study of the impact of Hurricane Georges on the Guánica dry forest in 1998, van Bloem and Murphy [38] also found wood density did not play a significant role in resisting hurricane damage. However, they found the pioneer species in the forest were damaged more.

3.4. LiDAR-Derived Canopy Height Change along the Transect

Averaging the pixels across the width direction of the transect, we found that canopy height change was correlated with both prior-hurricane canopy height and tree cover. Canopy height change was mostly negative in the lower left of the transect with relatively lower tree cover and average canopy height but positive in part of the upper-right section with greater tree cover and average canopy height (Figure 5a,b). On average, the canopy height decreased by 0.23 m, in comparison with a 0.14 m decrease in stem height from plot measurement, probably because ground stem height only measures the top height of the stems, whereas LiDAR-derived canopy height includes both treetop and canopy foliage. At 50 m × 50 m grid scale, spatial regression of the relative canopy height change ($\mathsf{\Delta }{h}_{cr}$), i.e., canopy height change divided by the canopy height prior to the hurricanes, on the normalized covariates resulted in the following Equation:

$$\mathsf{\Delta }{h}_{cr}=-0.62+0.39 C+0.89h_c-0.62h_c{}^2-0.019f\left(\alpha \right)-0.047 r$$

(4)

where C and h_c are prior-hurricane tree cover and normalized canopy height, respectively, $f\left(\alpha \right)=-\cos \left(\alpha -\pi /4\right)$ with α as the aspect angle, and r is normalized rugosity. The model had a regression standard error of 0.03, in comparison with the range of relative canopy height as [–0.6, 0.17]. All the coefficients had p-values less than 0.01. The spatial error model did not provide a statistic similar to that of R²; however, we calculated a similar value of 0.81 based on the definition of the ordinary linear model (Figure 5d). The spatial autocorrelation coefficient, λ, was estimated as 0.79. The regression showed that $\mathsf{\Delta }{h}_{cr}$ increased with C and h_c in a second-order (squared) fashion; however, decreases with r derived from the digital surface model. As negative $\mathsf{\Delta }{h}_{cr}$ indicates decreased canopy height or damaged canopy, the result implied that canopy height reduction was resisted by the tree cover, canopy height when the canopy was low or moderate, but promoted by the rugosity, canopy height when the canopy was high. The aspect function $f\left(\alpha \right)$ made southwest as +1 and northeast as −1, so that the result showed that slopes facing southwest experienced more damage than slopes facing northeast, using Equation (4).

4. Discussion

Instead of using the individual-level tree traits to explain the hurricane damage, we integrated the NEON ground observations of vegetation structure and the airborne LiDAR-derived continuous structural information, which allowed us to use forest stands as base units to explore the underlying mechanisms of ecosystem resistance to hurricanes. Our analyses of the response of the tropical dry forest in Guánica to major hurricanes highlighted the importance of stand-level parameters such as stem density and rugosity. At the stand level, trees function interactively to resist wind damage. The results are aligned with recent modeling studies reporting that stand-level characteristics play more important roles in wind damage to forests, compared to the individual-level characteristics [11,18].

At the stand level, stem density is one of the primary parameters as seen in Equation (2). A forest can be viewed as a windbreak system. In general, the denser the stems are, the denser the canopy, and the more support among trees with the interaction of their branches and leaves in the wind. The bending of stems and branches and friction among leaves dissipate the wind energy. Stem density is especially important for tropical dry forests with small trees to resist winds. The NEON Guánica dry forest has a mean density of 3990 ± 2250 stems ha⁻¹ in the 15 plots, whereas the estimated mean density for the tropical moist forest is around 800 globally [39] and 2278 ha⁻¹ for a moist forest in the same island [40]. Even for a forest with large trees and low apparent density, two or multiple trees growing together can support each other [41]. For a stand with low density, wind speed inside the canopy will be higher than that of a denser stand [11], thus leading to more severe wind damage.

Mean stem diameter is another essential parameter to resist forest height reduction as seen in Equation (1). Wind force at the canopy constitutes a turning moment to snap or uproot trees. In the Amazon forest, the strength to sustain a large turning moment is positively correlated to the DBH of the experimental trees [11,42]. Trees with larger DBH are more resistant to snapping or uprooting. A recent study of tree adaptation and acclimation to wind based on long-term windthrow datasets also revealed that taper or lower height:diameter ratio was related to reduced wind damage [12,19]. The forest risk model ForestGALES [43] indicates that critical wind speed to rupture stems is proportional to the order of 3/2 of stem diameter and to the stem density, aligned with our findings in the tropical dry forest as seen in Equations (1) and (2).

The standard deviation of stem heights is an index for unevenness of the canopy height reflecting surface roughness or rugosity as derived from the digital surface model (Figure 5). We found that the rugosity is important in explaining the reduction of stem height or canopy height (Equations (1) and (4)) and the proportion of damaged stems (Equation (3)) during major hurricanes. Storm intensity increases with surface roughness [44], which implies that higher surface roughness tends to bring more damage. More surface roughness makes a boundary layer thinner [45] so that storm winds can traverse the forest with a faster speed at the lower height to bring more damage. Hence, a forest with a uniform stem height is resistant to hurricane damage. Surface roughness can be enhanced by forests with emergent high trees or exceedingly low stem density with heterogeneous height arrangements, which may create gaps. During strong winds, emergent high trees might be knocked down first due to the negative correlation between critical wind speed to topple trees and tree height, and then greater wind penetration into the stand through the newly created gaps might induce higher wind load to the remaining trees [11]. Therefore, more damages are expected as gaps propagate. As to a specific case of plantations, in general, they have uniform stem height and simple structure. When the storm wind is lower than the critical wind speed to topple trees, plantations are expected to resist wind damage due to low surface roughness. However, due to their simple structure, plantation trees have similar critical wind speed to topple trees. Therefore, when the wind exceeds the critical wind speed to topple plantation trees, few trees could resist the winds, and the structure might collapse [20].

Topography-related factors, such as elevation, topographic position, and aspect, reflect the exposure of forests to different wind speeds during storms. The higher the elevation, the stronger the hurricane winds encountered, thus more damage, e.g., Equation (3) [46]. The Guánica dry forest, located at the left side of hurricanes Irma and Maria, which spun counterclockwise and moved toward the northwest, might encounter the strongest winds coming from the west to southwest, and our results showed that more damages (Equation (4)) were associated with areas facing southwest.

Our model for LiDAR-detected change in canopy height along a transect Equation (4) basically confirms the finding with the ground plot measurements in Equations (1)–(3). Tree cover may be associated with stem density if the crown size is relatively uniform, which is approximately the case in this dry forest with dense stems of small trees. Higher tree cover might imply higher stem density and fewer gaps so that trees help each other to resist damage. Higher tree cover might also indicate a higher leaf area index (LAI), whereas gust wind is negatively related to LAI [11,47]. Therefore, a larger tree cover might weaken gust wind and reduce wind damage. Canopy rugosity derived from DSM Equation (4) and stem height unevenness Equations (1) and (3) both tend to enlarge wind damage. We found a significant Pearson’s correlation coefficient of 0.65 between stem diameter and stem height in the NEON ground measurements. Hence, a greater mean canopy height in the low and moderate range may indicate a greater mean stem diameter. We also found a significant correlation of 0.63 between tree cover and mean canopy height derived from LiDAR. Therefore, a higher canopy may also mean greater tree cover to resist wind damage. Mean canopy height is negatively correlated to the gust wind profile [11,47] and thus tends to protect the emergent trees, which agrees with the role of mean canopy height in Equation (4).

On average, 4.5% of stems in the Guánica forest were affected by the hurricanes in 2017, similar to the impact of Hurricane Gilbert on the subtropical dry forest in Jamaica in 1988 with 5% structural damage [48] and that of Hurricane David on a dry forest in Dominica in 1979 with 4% [49]. Hurricane Georges in 1998 caused 12.4% structural damage to the Guánica forest [38], probably because the path of Georges was much closer to the forest than the hurricanes in 2017; however, the 2% of lost stems in this study is close to the 2% mortality nine months after Hurricane Georges [30]. In the Guánica dry forest, the conspecific percentage of affected stems was on average 2.7% when only the known species having at least 10 stems were counted (Table 2).

Puerto Rico experienced a severe drought in 2014–2016, and the large pulse rain brought by the major hurricanes in 2017 might release the drought stress (e.g., positive height change in Figure 3) and help the recovery of the dry forest [24]. Climate models predict intensified climate variability [50,51], and tropical dry forests might face increasing extreme climate events [52]. While tropical dry forests have adapted to moisture-deficit environments and might show resistance to further drought, a recent review still pointed out that altered species distribution and ecosystem processes are likely in the altered rainfall regime [8]. A great challenge to be addressed by future studies would be the resistance and resilience of tropical dry forests in response to alternating severe drought and major storms.

5. Conclusions

Tropical dry forests are vulnerable to hurricane disturbances and function differently from tropical moist or wet forests in their response, resistance, and recovery after the disturbances. As individuals’ functional traits prevail in the literature on hurricanes’ impact, we hereby addressed it at forest stand level and tried to answer how trees in tropical dry forest work together to combat the wind damage. To achieve this, we integrated the NEON ground observations of ecosystem structure with the airborne-LiDAR-derived continuous structural information, which allowed us to use forest stand as a base unit to explore the underlying mechanisms of ecosystem resistance of the NEON Guánica dry forest to the two major hurricanes of Irma and Maria in 2017. Our results highlighted the mutual supports among trees during storms and the important roles of stand-level collective parameters such as stem density and rugosity in resisting wind damages to the tropical dry forest. At the stand level, trees function interactively to resist wind damage, especially in the case of tropical dry forests with low profile but dense canopy. The conclusion on the tropical dry forest reemphasizes the finding from recent modeling studies, namely, that stand-level characteristics are more essential in explaining wind damage to forests, compared to the individual-level characteristics. Coupling remote sensing with ground observations provides an effective and promising approach to advancing our mechanistic understanding of ecosystem responses to disturbances.