Highlights

  • The OGF vegetation, with more biomass to lose, was more vulnerable.

  • Soil biogeochemical parameters were the least resistance to the hurricane, especially in SF.

  • Increase in intensity/frequency of hurricanes could lead this TDF toward a high nutrient limitation.

Introduction

At present, all ecosystems face the adverse effects associated with the inter- and intra-annual climate variability produced by Global Climate Change (GCC) (Allen and others 2017). Standing out from these climatic phenomena that produce the most damage to ecosystems are the Extreme Climatic Events (ECE) which include cold and heat waves, frosts, torrential rains, and intense windstorms (Jentsch and Beierkuhnlein 2008; IPCC 2012; Bojórquez and others 2019, 2021). Climate change has produced both an increase in variance and a displacement in the average value of temperature and precipitation, which has favored an increase in the intensity and frequency of ECEs in all ecosystems around the planet (IPCC 2012; Martínez-Yrízar and others 2018). The main problem with ECEs is that these events are highly stochastic and, when they occur, can modify the structure, the function, and the stability of any ecosystem (Jentsch and Beierkuhnlein 2008).

Among the tropical forest types, the Tropical Dry Forest (TDF), characterized by its elevated number of endemic species and its high beta diversity at different spatial scales (that is, continental, regional, and local scales), is considered one of the most threatened forests worldwide. TDFs are currently challenged by direct human activities that promote land-use change, extraction of resources, introduction of invasive species, and eventually alterations in fire patterns, and even frost-induced damage, as well as general shifts in climate regimes (Janzen 1988; Miles and others 2006; Linares-Palomino and others 2011; Pulla and others 2015, Álvarez-Yépiz and others 2018, Bojórquez and others 2021). Additionally, ECEs represent a strong stress for this ecosystem (Lugo 2008).

Cyclonic storms are the ECEs that most frequently affect many TDF and impose high stress on vegetation and soil (Lugo and others 2000; Cantrell and others 2014). For the year 2100, a significant increase in the intensity of tropical cyclones in the Eastern North Pacific is projected (IPCC 2014; Walsh and others 2016), and this is precisely the zone where most of Mexico’s TDF is located (Challenger and Soberón 2008; Villaseñor and Ortiz 2014).

The response of an ecosystem to a hurricane depends on factors such as species diversity, functional diversity and, quite probably, the development degree of the ecosystem (Jentsch and Beierkuhnlein 2008). Specifically, the effects of high-intensity hurricanes on vegetation include changes in the assembly of communities (Kreyling and others 2011), massive sudden tree mortality, delayed patterns of tree mortality, and alternative patterns of forest regeneration (Lugo and others 2000; Ibanez and others 2019). Also, immediately after the hurricane (days), the increase in canopy openness produces a change in canopy structure and function (Shiels and González 2014; Parker and others 2018). On the other hand, the effects of hurricanes on soil include the alteration of (1) the quantity and quality of the Soil Organic Matter (SOM), by way of unusual litterfall pulse, (2) the microclimatic soil conditions, (3) the activity of macro-, meso-, and microfauna, and (4) the SOM transformation, from fragmentation down to enzymatic mineralization (Lugo 2000; De Vries and Shade 2013; Cantrell and others 2014, Gavito and others 2018). In addition, excess water from precipitation produces rainfall–runoff erosion and nutrient loss by leaching (Nearing and others 2004; Lugo 2008; Plante and Conant 2014), and there is a reduction in the size of nutrient reservoirs such as carbon (C), nitrogen (N), and phosphorus (P), both in the vegetation and in the soil (Herbert and others 1999). All these effects can be increased by anthropogenic activities that change the composition, height, and density of the stand, and modify the inputs of SOM (Harper and Macdonald 2002).

In synthesis, hurricanes affect the ecosystem stability (the ability of the ecosystem to continue functioning), which is defined by two components: (1) resistance (that is, the amount of change caused by the disturbance) and (2) resilience (that is, the speed at which the system will return to its pre-disturbance level) (Pimm 1984; Orwin and Wardle 2004). In this sense, the response of ecosystems to hurricanes will depend on the duration and the intensity of the extreme event but also on the time that has passed since the last ECE (Jentsch and Beierkuhnlein 2008).

The assessment of resistance and resilience on TDF components after a hurricane may be helpful to elucidate the most vulnerable components. Unfortunately, most studies that have evaluated the impact of high-intensity hurricanes so far have focused on the effects over the floristic composition and structure and have mostly considered conserved forests (fragmented and/or continuous) (v.gr. Herbert and others 1999; Imbert and Portecop 2008; Lugo 2008; Lewis and Bannar-Martin 2012; Carrington and others 2015, but see McGroddy and others 2013). Several studies have focused on secondary forests (v. gr. Pascarella and others 2004; Flynn and others 2010; Hogan and others 2016), but they rarely attempt to explain the structural and functional damage over the soil–tree vegetation system (but see Gavito and others 2018). Even fewer report on the ecosystem conditions before and immediately after an ECE (Parker and others 2018).

In 2015, the category 4 Hurricane Patricia made landfall at the Chamela-Cuixmala TDF affecting the vegetation structure and depositing a large amount of wood debris (trunks and branches) and leaves on the soil (Jimenez-Rodriguez and others 2018). In the present study, we aimed to answer three questions: (1) How resistant is a Mexican tropical dry forest to a high-intensity hurricane? (2) Is this resistance associated with the degree of conservation of the forest? and (3) How the change in the occurrence of extreme climatic events in the last 30 years can affect the functioning of this ecosystem? For this purpose, we evaluated the resistance of the following features in two habitats (secondary forest and old-growth forest) of the Chamela-Cuixmala Tropical Dry Forest (TDF): (I) richness and diversity of the tree community, (II) vegetation structure and spatial distribution, and (III) soil biogeochemical variables. We hypothesize that the resistance of a TDF to a high-intensity hurricane will be modulated by its conservation degree and it will be higher in an old-growth forest than in a secondary forest, as a high diversity and structural complexity of the vegetation, as well as a conserved soil, will be associated to a lower vulnerability to disturbance.

Materials and Methods

Study Region and Sampling Sites

The study was carried out in and around the Chamela-Cuixmala Biosphere Reserve (CCBR; 19°22′–19°35′ N, 104°56′–105°03′ W), in the Pacific coast of Jalisco, Mexico (Figure 1). The CCBR and the surroundings comprise approximately 921.60 km2, with an elevation that does not surpass 200 masl. The climate is warm subhumid (Aw0i, García 2004) with a marked seasonal fluctuation in the precipitation regime and a dry season that lasts approximately seven months (from November to May). The average annual temperature is 24.6 °C, and the annual rainfall is 915 ± 311 mm (García-Oliva and others 2002). The vegetation is predominantly deciduous tropical forest and riparian forest (Lott and Atkinson 2002).

Figure 1
figure 1

Hurricane Patricia trajectory and distribution of study sites. The polygon represents the limits of the Chamela-Cuixmala Biosphere Reserve. Triangles (▲) represent Old-Growth Forest (OGF) sites, and circles (●) represent Secondary Forest (SF) sites.

We used eight sampling sites representing two successional stages, four Old-Growth Forest (OGF) sites (LIM2, TEJ1, TEJ2, and UNAM), and four Secondary Forest (SF) sites (HN2, LFH, LIM1, and NAC2). The sites were previously established for the study of the functionality and diversity of tropical plant communities, and the description of the criteria for establishing the study sites can be consulted in detail in Jimenez-Rodriguez and others (2018). Each site was defined by a plot with an area of 1000 m2 (50 × 20 m) and was separated from each other by a distance ranging from 1 to 32 km (Figure 1). When selecting the locations, the study plots were evaluated for the potential effect of the variation in elevation, slope, and aspect, which are tightly related to variation in environmental conditions and resources availability, such as insolation and soil water (Balvanera and others 2002). All sampling sites were located at an average mean height of 143 masl (range: 200 masl), with slopes averaging 10° (range 15°), mostly facing toward the south and the southeast. To make sites more comparable, we avoided north-facing slopes because they can show higher heterogeneity in insolation and plant community composition in relation to altitude (Balvanera and others 2002).

Hurricane Patricia made landfall on October 23, 2015, at the end of the rainy season. Patricia was classified as a category 4 hurricane (in the hurricane wind scales of Saffir–Simpson) and is currently considered one of the strongest recorded on the Pacific coast of Mexico. Coincidentally, our eight sampling sites were located in the same area of influence of Hurricane Patricia (11,062.67 km2 area of maximum influence), where the sustained winds averaged 119 km/h (Jimenez-Rodriguez and others 2018).

Vegetation and Soil Sampling

In each site, before and after the hurricane (September 2015, accumulated rainfall of 1394 mm and September 2016, accumulated rainfall of 727 mm, respectively), we sampled all woody plants with a diameter at breast height (DBH; 1.30 m) equal to or greater than 2.5 cm. The information recorded for each sampled individual was (1) the species name (as accepted in the WFO Plant List of 2023), (2) the diameter at breast height (DBH), (3) the total number of branches at breast height, and (4) the height. Additionally, we georeferenced (in UTM units, Zone 13N) all woody plants with a GPS Mobile Mapper 20 (Spectra Precision, 50-cm error after the post-processing differential correction).

During the vegetation sampling, the soil sampling was also carried out. At each sampling site, five soil samples were taken with a soil-core sampler from the first 15 cm of soil depth and were combined to produce one compound sample per site. A total of eight composite samples were collected for each year. The soil samples were stored in black plastic bags at 4 °C until laboratory analysis.

Biogeochemical and Soil Enzymatic Activity Analyses

For the determination of total nutrients, 50 g of fresh soil was dried and milled with a pestle and mortar. The Total Carbon (TC) was measured from a subsample process of 15 mg. To determine the total nitrogen (TN) and phosphorus (TP), subsamples of 5 g were digested in a mixture of concentrated H2SO4, H2O2 (30%), and K2SO4 plus CuSO4, the latter acting as a catalyst at 360 °C. One aliquot of the extract was used to determine the TP and TN.

The dissolved inorganic nitrogen forms (ammonium NH4+ and nitrate NO3) and phosphorus (orthophosphate, PO4+) were extracted from 20 g of fresh field soil samples with deionized water after shaking for 45 min and then filtering through a Whatman filter paper No. 42 and a 0.45 nitrocellulose membrane. One aliquot of the filtrate was used to determine the dissolved nitrogen and phosphorus (Jones and Willett 2006).

The total organic carbon (TOC) was determined by combustion and coulometric detection (Huffman 1977) with a Total Carbon Analyzer (UIC Mod. CM5012; Chicago, USA). The nitrogen (N) forms were determined by the macro-Kjeldahl method (Bremner 1996), and the phosphorus (P) was determined by the molybdate colorimetric method, following ascorbic acid reduction (Murphy and Riley 1962). Both N and P forms were measured with colorimetric analyses using a Bran Luebbe AutoAnalyzer III (Norderstedt, Germany).

The activities of three enzymes (extracellular) involved in the cleavage of organic molecules with C, N, and P were measured: β-1,4-Glucosidase (BG), β-1,4-N-acetylglucosaminidase (NAG), and phosphomonoesterase (PME) with the assay techniques reported by Tabatabai and Bremner (1969), Eivazi and Tabatabai (1977), Eivazi and Tabatabai (1988), and Verchot and Borelli (2005). For all enzymes, we used 2 g of fresh soil and 30 ml of modified universal buffer (MUB) at pH 7 for enzyme extraction. All enzyme assays were incubated at 40 °C, the BG for 2 h, the NAG for 3 h, and the PME for 1.25 h. Following the incubation period, the tubes were centrifuged at 10,000 rpm for 2 min, and 750 μl of supernatant were recovered. We diluted the supernatant in 2 ml of deionized water with 75 μl of NaOH and measured the absorbance of p-nitrophenol (pNP) liberated at 410 nm on an Evolution 201 spectrophotometer (Thermo Scientific Inc.). Enzyme activities were expressed as nanomoles of pNP per gram of dry soil per hour (nmol pNP [g SDE]−1 h−1).

Data Analysis

We used two datasets, one from before and one from after the hurricane, to analyze changes in vegetation and soil attributes as a result of the hurricane. Trees that had died because of the hurricane were no longer represented in the post-hurricane database (Jimenez-Rodriguez and others 2018).

Characterization of Vegetation Response to the Hurricane

We characterized the response of the studied vegetation to the hurricane in relation to its species composition, diversity, and structure. We quantified the variation among tree community species composition and structure through Non-Metric Multidimensional Scaling (NMDS) ordinations (McCune and others 2002). For the ordinations, we consider two types of abundance matrices (sites × species), one using the number of individuals as a measure of the abundance of the species in each site and another considering their basal area. For the ordination based on the number of individuals, the Bray–Curtis coefficient was used as a measure of dissimilarity, while for the ordination based on the basal area, we used the Euclidian distance as a measure of dissimilarity. Prior to the calculation, we square-rooted the individual-based matrix and normalized the basal area-based matrix in order to avoid biases toward species with the highest abundance values and with the greatest variation in their abundances (Kindt and Coe 2005). In the analyses, we consider the set of species that, in each community, include 90% of the total number of individuals or comprise 90% of the total basal area, encompassing the species that have the greatest impact on ecosystem processes. The stress value, scaled from 0 to 100, was used as a measure of the degree of reliability of the resulting ordinations. From these ordinations, we used the sites’ coordinates (scores) as synthetic variables summarizing variations in community composition and structure before and after the hurricane.

To characterize the species diversity of the tree communities, we used the first three Hill numbers, of order q = 0 (species richness), q = 1 (exponential of Shannon–Wiener index, eH’), and q = 2 (inverse of Simpson index, 1/D). Along the Hill series, the parameter q quantifies how much it is discounted for the rare species when calculating the diversity indices—where a greater value of q corresponds to a greater discount of rare species—. This allowed us to generate a quantitative profile of species diversity for our tree communities considering both, species richness and evenness. We calculated the observed and estimated values for each of the indices considering the asymptotic estimators for the Hill numbers presented by Chao and others (2014). The values of these indices are reported in units of species, considering the Jost’s considerations to estimate the effective number of species and to use these indices as measures of true diversity (Jost 2006).

For the characterization of the vegetation structure, we used the following parameters: the total number of individuals (NI), the stand basal area (BA), the average tree height (ATH), the emergent tree layer height (ETH)—average tree height for the 10 tallest trees—, the number of branches (NB), the average nearest neighbor distance (NND), and the average density of individuals (intensity of points/trees, DEN), these last two parameters as measures of the individuals spatial distribution. Before quantifying the spatial pattern of trees, we checked for duplication of coordinates, to avoid violating the basic assumption behind the statistical methodology for analyzing spatial point processes, which is that the points of the study process can never be coincident (Baddeley and others 2016). Finally, due to the non-normality of the vegetation parameters, we used the Wilcoxon rank-sum test to evaluate for significant differences between OGF and SF sites.

Analysis of Soil Biogeochemical Variables

We used the analyses of variance (ANOVAs), followed by post-hoc Tukey tests, to evaluate for differences between successional stages (OGF vs. SF) and time (before vs. after the hurricane) on each of the biogeochemical parameters measured (p < 0.05). We also assessed the hurricane effect on biogeochemical parameters by performing a canonical discriminant function analysis using the program XLSTAT (Addinsoft 2019).

Evaluation of the Level of Resistance of Vegetation and Soil to the Hurricane

We calculated a resistance index for each of the vegetation and soil parameters, using the resistance index equation proposed by Orwin and Wardle (2004). The data obtained before the hurricane (September 2015) were considered as the control (C0), and the data obtained after the hurricane (September 2016) were considered as the disturbance treatment (P0). The resistance index (RS) was calculated as follows:

$$ {\text{RS}} = 1 - \left[ {\left( {2{*}\left| {D_{{\text{O}}} } \right|} \right)/\left( {C_{{\text{O}}} + \left| {D_{{\text{O}}} } \right|} \right)} \right] $$
(1)

where

$$ D_{O} = C_{O} - P_{O} $$
(2)

The RS is bounded by − 1 and + 1. When the resistance index value is + 1, it means that the disturbance had no effect on the measured variable, and therefore, it is considered to have maximum resistance. When the resistance index value is zero, the measured variable was 100% disturbed, which indicates low resistance. Finally, when the resistance index value is negative (as low as –1), it indicates that the measured variable was more than 100% disturbed which corresponds to the minimum resistance and also to a completely different behavior of the variable. To evaluate for statistical differences in the resistance index between successional stages (OGF vs. SF), we used the Mann–Whitney U Tests (p < 0.050).

Standardized Precipitation Index Calculation for the CCBR

Finally, we assessed the effect of hurricanes on annual precipitation regimes, by calculating the Standardized Precipitation Index (SPI) as proposed by Zhang and others in 2009 (Table 1, Supplementary material S1). For this, we used a dataset of 38 years obtained from the manual meteorological station of the Chamela-Cuixmala Biosphere Reserve (CCBR), which provides the oldest meteorological record in the study region.

Table 1 Standardized Precipitation Index (SPI) Categories

Analyses were performed in R (R Core Team 2021), using different packages: “spatstat” for point pattern analysis—quantification of the spatial arrangement of trees in space—(Baddeley and Turner 2005; Baddeley and others 2016); “iNEXT” for the diversity indices calculations (Hsieh and others 2020); and “vegan” for NMDS ordination (Oksanen and others 2020).

Results

Resistance Regarding Tree Community Attributes

Our sampling allowed us to register a total of 1907 individuals representing 203 species in OGF sites and 1216 individuals representing 195 species in SF sites. After the passage of the hurricane, the numbers of individuals and species present in our sites decreased to 1764 individuals and 195 species in OGF and 1173 individuals and 189 species in SF. In general, the SF sites presented the greatest variation in terms of species richness (ranging from 10 to 68 species) and diversity (ranging from 4 to 43 for Shannon index, and from 3 to 29 for Simpson index), including both, the sites with the highest values of richness and diversity (LIM1 and HN2), and those with the lowest values (LFH) (Table S1). Conversely, the OGF sites were more similar to each other in terms of species richness and diversity. These patterns persisted after the passage of the hurricane. The high variation of SF probably precluded the detection of significant differences in diversity between successional stages, before and after the hurricane.

However, the OGF and SF sites differed significantly from each other in terms of the species composition, density of individuals, and basal area (Figure 2). These differences between successional stages also persisted after the passage of the hurricane. As in the case of species diversity, secondary forest also showed the greatest variation in terms of specific composition. The 10 most dominant species in the OGF sites were (in order of abundance) Croton roxanae Croizat, Croton suberosus Kunth, Apoplanesia paniculata C. Presl, Thouinia paucidentata Radlk, Lonchocarpus lanceolatus Benth, Heliocarpus pallidus Rose, Guapira petenensis (Lundell) Lundell, Caesalpinia eriostachys Benth, Serjania brachycarpa A.Gray ex Radlk, and Celosia monosperma Rose; while the 10 most dominant species in the SF sites were Haematoxylum brasiletto H. Karst, Acacia cochliacantha Willd, Caesalpinia coriaria (Jacq.) Willd, Croton alamosanus Rose, Psidium sartorianum (O.Berg) Nied, Lonchocarpus constrictus Pittier, Casearia corymbosa Kunth, Zanthoxylum fagara (L.) Sarg, Senna pallida (Vahl) H.S.Irwin & Barneby, and Guazuma ulmifolia Lam.

Figure 2
figure 2

Non-Metric Multidimensional Scaling Ordinations mapping plant communities’ dissimilarities regarding tree community composition and structure. Chart A shows the individual-based ordination (a) and chart B shows the basal area-based ordination (b). Ecosystem types: Old-Growth Forest (OGF) and the Secondary Forest (SF). Suffix A corresponds to the state of the communities prior to the arrival of the hurricane, while suffix B corresponds to the state of the communities after the passage of the hurricane. Results of the Wilcoxon tests contrasting OGF and SF site scores on the individual-based ordination: W = 15, P = 0.05 for axis 1 and W = 16, P = 0.03 for axis 2; and on the basal area-based ordination: W = 16, P = 0.03 for axis 1 and W = 12, P = 0.34 for axis 2. The same statistical results were obtained considering the scores corresponding to before and after the hurricane.

We did not observe differences between OGF and SF sites in their resistance index regarding species composition (NMDS scores), species richness, and species diversity (Shannon and Simpson indexes) (Figure 4a). For all these variables, the resistance index values were between 0.8 and 1.

Resistance Regarding Vegetation Structure

In general terms, the OGF, characterized by higher basal area and height of emergent trees, as well as by a shorter distance between individuals and a greater density in its distribution, presented a greater structural complexity than SF (Table 2). However, the differences between the two types of forests decreased after the passage of the hurricane, becoming indistinguishable in terms of their basal area and the spatial distribution of their individuals (NND and DEN) because, as a consequence of the high mortality of trees, OGF sites suffered a drastic reduction in its basal area, an increase in the distance between trees, and a decrease in their density. In fact, when we specifically contrasted the distribution of dead trees between OGF and SF sites, we found a higher density of these in the OGF (mean 0.0264) than in the SF (mean 0.0080) (Wilcoxon rank-sum test: Z = 2.02, P = 0.05; Figure 3).

Table 2 Vegetation Attributes at Old-Growth Forest (OGF) and Secondary Forest (SF) Site, Before and After the Hurricane
Figure 3
figure 3

Spatial distribution of trees at each study site, showing the trees that died after the hurricane event. Dots (●) represent living trees, and the plus sign (+) represents the dead trees. Charts A–D correspond to old-growth forest sites, and charts E–H correspond to secondary forest sites.

For the resistance index of structural attributes, we found that the total basal area and the average tree height showed a higher resistance in the SF than in the OGF (Figure 4b). In addition, the basal area in the OGF exhibited the lowest resistance value (0.63 ± 0.06) of all the structural attributes measured. Although we did not find significant differences between forest types for the height of emergent trees, the number of branches, and the number of individuals, we did observe a trend of lower resistance in the OGF (indices ranged from 0.63 to 0.95) than in the SF (indices ranged from 0.85 to 1, Figure 4b). Finally, while we did not detect significant differences between OGF and SF regarding the resistance for the density of individuals and the nearest neighbor distance, we did observe a trend toward a lower resistance in the OGF sites in relation to these spatial parameters (Figure 4b).

Figure 4
figure 4

Resistance index corresponding to plant community attributes, to vegetation structural attributes, and to soil biogeochemistry parameters for old-growth forest and secondary forest. Triangles (▲) represent Old-Growth Forest (OGF), and circles (●) represent Secondary Forest (SF). Plant community attributes (chart A): scores of tree communities in the axes 1 and 2 resulting from the ordinations (NMDS), considering as measures of species abundance, their number of individuals (Cind1 and Cind2), and their basal area (Cba1 and Cba2). These scores constitute a synthetic variable that reflects the variation between communities in terms of their composition and structure; observed and estimated species richness (SO and SE, respectively); observed and estimated Shannon diversity index (HO and HE, respectively); and observed and estimated Simpson diversity index (DO and DE, respectively). Vegetation structural attributes (chart B): basal area (BA); average tree height (ATH); emergent trees height (ETH); number of branches (NB); number of individuals (NI); average density of individuals (intensity, INT); and average nearest neighbor distance among trees (NND). Soil biogeochemistry parameters (chart C): percentage of water (%W); potential hydrogen (pH); β-1,4-glucosidase (BG); β-1,4-N-acetylglucosaminidase (NAG); phosphomonoesterase (PME); total organic carbon (TOC); total nitrogen (TN); total phosphorus (TP); the stoichiometry of C:N, C:P, and N:P; dissolved ammonium (NH4+); dissolved nitrate (NO3); orthophosphate (PO4); and the stoichiometry of NH4+:NO3. The asterisk (*) indicates significant differences (p < 0.05) among OGF and SF.

Resistance Regarding Soil Biogeochemical Variables

The discriminant analysis using the soil variables showed that two factors explain 98.96% of the variance and that sites can be separated in four groups (Figure 5a), in the first axis (F1) by time (before or after the hurricane) and in the second axis (F2) by successional stage (OGF and SF). The first discriminant function explained 89.78% of the between-group variance. The variables that contributed most to the first function were HPO4, %W, pH, TOC, BG, and NAG. The second discriminant function explained 9.71% of the between-group variance, and the variables that contributed most to the second function were NH4, NH4:NO3, N:P, PME, C:N, and C:P (Figure 5b).

Figure 5
figure 5

Discriminant analysis of ecosystem types, before and after the hurricane, regarding soil biogeochemistry parameters. Ecosystem types: Old-Growth Forest (OGF) and the Secondary Forest (SF). Suffix A corresponds to the state of the communities prior to the arrival of the hurricane, while suffix B corresponds to the state of the communities after the passage of the hurricane. Chart A shows the distribution of the study sites in the resulting two-dimensional space, where factors F1 and F2 explain 98.96% of the variance. Chart B shows the coordinates of the soil biogeochemistry parameters correlated with factors F1 and F2. Soil biogeochemistry parameters: percentage of water (%W); potential hydrogen (pH); β-1,4-glucosidase (BG); β-1,4-N-acetylglucosaminidase (NAG); phosphomonoesterase (PME); total organic carbon (TOC); total nitrogen (TN); total phosphorus (TP); the stoichiometry of C:N, C:P, and N:P; dissolved ammonium (NH4+); dissolved nitrate (NO3); orthophosphate (PO4); and the stoichiometry of NH4+:NO3.

The soil moisture and pH were higher before than one year after the hurricane (p < 0.01) for the OGF and SF (Tables 3, 4, and 5), while the resistance indices of these variables were not different between the stages. The resistance values for soil moisture ranged between 0.8 and 0.1 and for pH between 0.79 and 0.86, respectively (Figure 4c).

Table 3 Biogeochemical Parameters at Each Site, Before and After the Hurricane
Table 4 Results of the t test Contrasting the Old-growth Forest versus the Secondary Forest in Relation to the Biogeochemical Parameters
Table 5 Results of the t test Contrasting the Old-Growth and Secondary Forests Before versus After the Passage of the Hurricane in Relation to the Biogeochemical Parameters

For carbon, the mean values of total organic carbon (TOC) and the β-1,4-glucosidase activity (BG) changed significantly one year after the passage of the hurricane only at OGF (Tables 3, 4, and 5). In this stage, the TOC concentration increased after the hurricane and the BG activity decreased. The resistance index of the TOC was similar between the OGF and the SF with mean values of 0.42 and 0.48, respectively (Figure 4c), while for BG activity, the index was higher in the OGF (0.41 ± 0.1) than in the SF (− 0.05 ± 1).

For nitrogen, we detected differences for the concentrations of dissolved inorganic nitrogen (as ammonium DNH4) and for the DNH4+: DNO3 ratio (Tables 3, 4, and 5). Before the hurricane, the DNH4+ concentration and the DNH4+: DNO3 ratio were higher in the OGF than in the SF. Likewise, at the OGF, both parameters were higher before than after the hurricane (Tables 3, 4, and 5). In this sense, the initial difference between the stages in the DNH4+ concentration and the DNH4+: DNO3 ratio disappeared after the hurricane. For the resistance index of DNH4 and DNO3, we observed an inverse pattern (Figure 4c). While the DNH4 exhibited a higher resistance index in OGF (0.21 vs. − 0.38), the DNO3 showed a higher resistance in the SF (0.22 vs. − 0.32). Although the total N concentration and the NAG activity showed no significant differences, we did observe a trend of lower hurricane resistance in the OGF (0.2 and 0.27, respectively) than in the SF (0.4 and 0.36, respectively). For the NH4+: NO3 ratio, we did observe a trend of lower hurricane resistance in the SF (− 0.26), while the resistance in the OGF is null (− 0.01) (Figure 4c).

For the phosphorous, we observed differences before and after the hurricane in the concentrations of total phosphorus (TP) and dissolved inorganic phosphorus (as phosphate, PO4+). The TP only showed differences at the OGF, and its concentration increased after the hurricane (Tables 3, 4, and 5). In contrast, the phosphate, whose concentration before the hurricane was higher in SF than in OGF, decreased dramatically in both stages, and the difference between them disappeared after the hurricane as PO4+ concentration dropped to zero in all the sites (Tables 3 and 4). Consequently, the resistance values of PO4+ were equal to zero for both forest types (Figure 4c). The phosphomonoesterase (PME) activity did not show significant differences between time or forest type (Tables 3, 4, and 5); however, it showed greater resistance in OGF than in SF, with values of 0.7 vs. 0.3, respectively (Figure 4c).

For the ratios C:P and N:P, no differences between sites or time were detected (Tables 3, 4, and 5), but both parameters showed higher resistance in the OGF than in the SF (for CP, 0.4 vs. − 0.01 and for NP, 0.4 vs. − 0.17, respectively; Figure 4c). For C:N, no differences between sites and time were detected. The resistance value did not differ between the stages either.

Standardized Precipitation Index (SPI) for the CCBR

The results of the classification of each year evaluated by the SPI and the hurricanes that occurred in the severely to extremely humid years (NOAA 2019; SMN 2019) are shown in Figure 6. We observed that in a 38-year period, there was no year as moist as 2015, the year that Hurricane Patricia occurred. We also observed that since the year 1998, there was an increase in the moderately humid years interspersed with years with moderate to extreme drought (that is, 373 mm in the 2001). Finally, we identified an increase in the frequency of moderately wet years in the second decade of this century.

Figure 6
figure 6

Standardized Precipitation Index calculated from 1980 to 2018 for the region of the Chamela-Cuixmala Biosphere Reserve.

Discussion

We found that the disturbance level mediates the differential resistance of the vegetation and soil attributes to high-intensity hurricanes. In general, vegetation attributes such as richness and diversity showed very high resistance one year after the hurricane, but several structural attributes did not, especially in the OGF. In this sense, SF vegetation appears to be more resistant, whereas the OGF, with more biomass to lose, appears to be more vulnerable. Conversely, most soil attributes showed low resistance in both stages but especially in SF. Then, we found that the short-term response of TDF to high-intensity hurricanes, in terms of above- and belowground processes, is in part dependent on their disturbance level. Furthermore, the low resistance of several vegetation and soil attributes (v.gr., trees basal area, PO4, NO3, and NH4 soil concentrations, NH4+: NO3 ratio, BG activity, among others) suggests that this TDF is facing a drastic change in functioning and that under a scenario of increasing intensity and frequency of hurricanes, its recovery to pre-hurricane conditions could not be possible.

Resistance in the Composition and Structure of the Tree Community

First, tree communities of OGF and SF were clearly different in terms of species composition of living trees before and after the hurricane (Figure 2). Remarkably, after a highly intense hurricane, this difference between habitats persisted as well as the dominance of families such as Leguminosae and Euphorbiaceae, although the abundance of the most dominant species notably decreased in OGF, suggesting higher levels of change in species composition of the mature forests and that their dominant species, especially those of high basal area, were strongly and differentially affected by the hurricane. Then, even though the species composition did not show great changes in both habitats one year after the hurricane, OGFs appear to be more vulnerable than SFs to high-intensity hurricanes, as these events could be purging the most massive species. This idea is supported by the higher frequency of dead trees found in OGF (Figure 3). However, the delayed mortality of several individuals and species could mask the real and long-term changes in species composition, species richness, and diversity (Xi and others 2008; Uriarte and others 2019).

We also observed a greater variation in species composition, richness, and alpha diversity in SF (Figure 2). This variation could arise from the stochastic nature of the early successional process in this ecosystem—driven by the diverse land-use histories/anthropogenic activities in the region—coupled with the high beta diversity and environmental heterogeneity characteristic of the TDF on the Pacific coast of Mexico (Alvarez-Añorve and others 2012; Jimenez-Rodriguez and others 2018). Likely due to these qualities of the system, the pronounced difference between habitats in inter-site variation persisted even after the hurricane.

Regarding the changes in vegetation richness and diversity, although some studies in TDF (Imbert and Portecop 2008; Ibanez and others 2019) have reported declines in the richness and diversity of the vegetation in response to a hurricane, we did not detect immediate differences after this hurricane, even in spite of the reduction, almost in half, for the percentage of individuals of OGF that included the most dominant species. As the decrease in the abundance of dominant species did not significatively affect the species diversity, this parameter showed a high level of resistance in OGF, at least in the short term. Indeed, we observed high resistance values in the richness and diversity of the tree community in both successional stages (Figure 4a). The latter implies that in a single high-intensity hurricane event, short-term changes in the diversity attributes of the communities of both forest types (OGF and SF) were not of great magnitude. A similar result was observed in a temperate forest affected by a category 3 hurricane, where the immediate damage effect was not significant and was assumed to be delayed as compositional changes in tree species may have a time lag (Xi and others 2008). On the other hand, five years after Hurricane Hugo (a category 3 storm) struck Puerto Rico in 1989, 30% of the stems that were broken and 55% of those that were uprooted during the storm had died, likely causing delayed changes in both the size structure and species composition of the Luquillo forest (Uriarte and others 2019).

In terms of vegetation structure, the OGF was the least resistant to the hurricane. Before the hurricane, the basal area was higher in the OGF than in the SF, but after the hurricane, the difference between successional stages was greatly diminished (Table 2). Consequently, among all the structural variables evaluated, the basal area had the lowest resistance value, and this value was significantly lower in the OGF than in the SF. Additionally, the average height (AH) showed lower resistance in the OGF than in the SF, suggesting the existence of a relationship between the BA, the AH, and the level of damage caused by a high-intensity hurricane. If we also consider that in OGF, the density (DEN) of trees was higher (Figure 3) and that both descriptors of spatial distribution, NDD and DEN, showed a trend toward a low resistance to the hurricane (Figure 4b), such relationship can also be highly associated to the spatial distribution of the trees. In this sense, the greater damage and mortality associated with hurricane winds in sites with a high density of massive and tall trees, such as the OGF sites, can be due to a domino effect on the surrounding trees, as the fall of such tall and massive trees causes the damage, fall, and mortality of the adjacent ones. Certainly, the tallest trees are more vulnerable to mechanical instability under strong winds, as previously reported by Jimenez-Rodriguez and others (2018) in the study system.

The above confirms that the structure and spatial distribution of the OGF tree community play a central role in its lower resistance to a high-intensity hurricane. These results are consistent with the findings in a Caribbean forest where a high-intensity hurricane produced a setback in the successional stage of the old forest and favored the dominance of SF and its simpler structure (Lugo 2000; Imbert and Portecop 2008). Likewise, the canopy of Yucatán OGF in Mexico, remained significantly reduced two years after the Hurricane Wilma, whereas the canopy of SF quickly recovered to pre-hurricane conditions (Hasselquist and others 2010). Our result also confirms the idea that intense tropical cyclones may alter the TDF vegetation mainly through changes in density, height, basal area, and spatial distribution of trees (Jimenez-Rodriguez and others 2018; Ibanez and others 2019). These alterations would diminish the characteristic differences between OGF and earlier successional stages, favoring the dominance of SF like vegetation. Indeed, previous plant functional evaluations performed in TDFs found that the proportion of plant species combining the capacity to resist damage and harsh environmental conditions with a high resprouting ability (that is, TDF pioneer species) could increase in scenarios of increasing frequency and/or intensity of hurricanes (Van Bloem and others 2006; Jimenez-Rodriguez and others 2018). Furthermore, studies examining the effects of successive hurricanes on other tropical forests have found that the accelerated pace of change in storm disturbance regimes may exceed the capacity for adaptation of the forest communities, yielding a depauperate forest dominated by pioneer species (Uriarte and others 2019).

Resistance in the Soil Biogeochemical Variables

The biogeochemical variables showed a clear differentiation between habitats and behaved differently before and after the hurricane (Figure 5a). Overall, all the studied biogeochemical variables were affected in their resistance, indicating a low resistance of the soil biogeochemical processes to this hurricane (Figure 4c). Contrasting with the findings related to the vegetation, the low resistance values (even less than zero) found in these parameters suggest that SF soil could be highly vulnerable to ECEs, even more than OGF soil, probably because the loss of nutrients caused by storms is more accentuated in soils covered by a thinner layer of vegetation and consequently with lower retention of nutrients, as those that result from anthropogenic disturbances.

The increased concentration of TOC, TN, and TP in the OGF after the hurricane is probably due to a higher deposit of leaf litter and woody debris on the forest floor produced by the strong winds of hurricanes and transferring nutrients to the soil above normal levels. This finding is consistent with the increase in the nutrients flux detected after a nearly complete canopy defoliation occurred in this region after the Hurricanes Jova and Patricia (Jaramillo and others 2018, 2022; Parker and others 2018). Hurricane Jova (category 2) hit the Chamela TDF in the year 2011, and it immediately caused high amounts of dissolved organic C, N, and P, which were transferred by throughfall from vegetation to soil (Jaramillo and others 2018). Noticeably, in our secondary forests, this did not occur.

Prior to the Hurricane Patricia, the DNH4 + was significantly lower in the SF compared to the OGF (Table 3), indicating that this form may not be in excess in the disturbed forests, although the higher proportion of legumes with atmospheric nitrogen-fixing bacteria could help to explain this observation (Hasselquist and others 2010). However, one year after Patricia, the potentially available inorganic N and the BG activity were severely affected in the OGF soil (Table 3). Besides, the DNH4+, DNO3, BG activity, and NAG activity all demonstrated low resistance to the hurricane in both habitats (Figure 4c). We attribute this post-hurricane N disturbance, especially in the OGF soil, to a reduction of energy obtained from low BG activity, suggesting a shift in the processes regulated by the soil microbial community that results in a reduction of nutrient availability, which is a known response in water-limited ecosystems to fluctuation in precipitation events (Montiel and others 2017). Furthermore, it is likely that a nutrient loss by leaching also occurred, as when Hurricane Jova (category 2) hit the Chamela TDF in the year 2011 and immediately caused high amounts of organic dissolved C, N, and particularly P which were dispersed by runoff from vegetation to soil (Jaramillo and others 2018). Indeed, the low levels of inorganic N that we found one year after Patricia could reflect the effect of past hurricanes, such as Jova, that affected the nitrogen cycling resistance to these events (Jaramillo and others 2022).

Although phosphatases activity decreased slightly and presented medium resistance (a little higher in OGF), the available P (PO4) had zero resistance and reduced its concentration by 100% after the hurricane (Table 3), showing that the P was the most vulnerable soil nutrient to the hurricane in both OGF and SF. Again, the decrease in the bacterial community as well as the effect on the soil fungal community one year after the event could affect the dynamic of P that in this TDF is mainly controlled by biological processes (Campo and others 2001, Alvarez-Manjarrez and others 2021). Moreover, this could also reflect the action of past hurricanes, for example, Jaramillo and others (2022) found that up to 44 months after Patricia, the reductions in seasonal and annual N and P fluxes (and consequently in the productivity of this TDF) were still evident, differently to the relatively rapid recovery of such fluxes after Jova. They also found that the exceptionally high annual P flux attributed to Jova may have represented a threshold to the P cycle in this TDF. This confirms that resistance of biogeochemical processes after Hurricane Patricia was low implying a potential reduction in their resilience, at least in the short term.

Long-Term Effects on TDF Functionality

For the Chamela TDF, it has been found that the SF can recover a structure and floristic composition similar to that of the OGF after a period of between 13 and 15 years (Alvarez-Añorve and others 2012). Likewise, the rapid regeneration of biomass and structure has also been reported for other TDFs (Norden and others 2009; Chazdon and others 2016). However, it has also been reported that the recovery of plant functional attributes (as net photosynthesis or leaf thickness, among others) can be slower than the recovery of the structure (Alvarez-Añorve and others 2012). Furthermore, even the same species can function significantly different in an OGF and in a disturbed forest that has recovered its structure and composition.

Regarding the available soil nutrients that showed dramatic changes after the hurricane, we are uncertain of the time required for them to achieve a full recovery, but it will take longer after Patricia than after Jova, highlighting the importance of considering not only the consequences of high-intensity storms on forest resistance and recovery but also the legacy effects from previous hurricane disturbances (Jaramillo and others 2022). In the same sense, a study in another ecosystem that is also limited by P showed that 30 years after the disturbance of a grassland soil, the PO4 did not show signs of resilience, while the NH4 had a low resilience, and these conditions led the system toward a slow recovery of their soil functionality (Hernández-Becerra and others 2016). We consider that under the predicted increase in intensity/frequency of hurricanes in the Eastern North Pacific, in this TDF, the resilience of available nutrients will be compromised and soil functionality can be severely altered, but this needs to be investigated further.

In this TDF, the main problem lies in the recovery time or resilience of the functional attributes of both vegetation and soil, as feedbacks between aboveground and belowground processes have the potential to influence forest recovery. For example, it has been reported that loss of nitrogen can limit the aboveground productivity, which in turn may reduce fine root production. This reduction could have pronounced effects on the ability of trees to absorb nutrients from aboveground hurricane litterfall, thereby inhibiting the regeneration of the canopy (Hasselquist and others 2010). This is particularly true for the regeneration of the most massive elements. Consequently, the changes and low resistance in belowground processes may exert extensive feedback on the recovery of aboveground forest canopy. In this scenario, the feedback on the recovery process might be more pronounced for the SF due to the lower resistance of its belowground processes (Figure 7).

Figure 7
figure 7

Model exploring the impact of increased hurricane intensity and/or frequency on the tropical dry forests of the Mexican Pacific Coast. It illustrates the loss of soil function and decline in vegetation structural complexity in both old-growth and secondary forests. Secondary forests, which start with lower soil functionality and vegetation structural complexity, are likely to evolve into more degraded communities compared to old-growth forests under these conditions, due to the feedback of belowground processes on the recovery of aboveground forest canopy. However, in the long term, both forest types might transition into scrubbier communities.

Previous studies on Chamela TDF have assessed resilience in both functional and structural variables over extended periods. They have also reported the behavior of these variables seven years after Hurricane Jova. As a result, these studies were able to evaluate resilience under two conditions: pre-perturbation and several years post-perturbation (Álvarez-Yépiz and others 2018; Gavito and others 2018; Parker and others 2018; Tapia-Palacio and others 2018). We emphasize the importance of long-term assessments of the variables presented in this study. This would provide data on post-perturbation conditions and facilitate the evaluation of resilience in both structure and functionality of soil and vegetation following Hurricane Patricia in future research.

We observed with the SPI that, at least for 40 years (since 1980), an ECE like the one in 2015 has not occurred. However, not only the intensity but also the sum with the previous ECEs of lower intensity could have gradually degraded both structure and function of the ecosystem (Beard and others 2005). We also observed that since the year 1998, there was an increase in severe humid years interspersed with moderate to severe droughts and is reported for a forest of Puerto Rico that the variability in rainfall regimes and drought influenced both above- and belowground ecosystem processes directly (Beard and others 2005). Finally, we observed that there has been a tendency since 2007 to increase precipitation. From 2012 to 2014, in addition to the incidence of hurricanes during the cyclone season, the increase in annual humidity has also been explained by abnormal rains in the dry season from December to May (Maass and others 2018; Martínez-Yrízar and others 2018).

In synthesis, we consider that short-term changes in abiotic and biotic factors due to a high-intensity hurricane will likely lead to substantial long-term impacts in both structure and function of vegetation and soil. Furthermore, if the climatic disturbances return at intervals shorter than the recovery time (resilience) of the vegetation and soil of this ecosystem (which is currently happening), in the long term, we could expect changes in the structure and function toward a more humid but less diverse forest, lacking several species typical of the OGF and therefore with a simpler structure (that is, a lower height canopy or scrubbier appearance) and more legumes dominance, as this clade can contend with the limitation of nitrogen in the soil, favoring the aboveground productivity and forest recovery. Moreover, while in the long term, the OFGs could change their functioning, vegetation structure, and species composition, and these changes could be more drastic in the SFs which, even in the short term, were more heavily affected in their soil biogeochemical processes and could face greater nutrient limitations (Figure 7). Then, the increase in the intensity and frequency of hurricanes over the pacific coast of Mexico could lead to this ecosystem toward a high nutrient limitation (especially by phosphorus) for the plants and consequently toward a loss of functioning, especially in the SF. This eventually could produce a severe degradation in fundamental attributes and functions of the ecosystem. Considering this, it is crucial to halt or significantly reduce the transformation of this tropical dry forest to maintain, at least, its most basic functions and ecosystem services under the future scenarios of climatic change.