Light penetration and topography shape juvenile tree species assemblies in the understory of the tropical Andean cloud forest

Abstract

The floristic composition of the understory plays a fundamental role in the long-term conservation of the diversity, structure, and function of mountain cloud forests in the Andes. We evaluated the relationship between the understory tree floristic composition of four types of predefined cloud forests and the canopy structure, the light transmitted to the understory, and the effect of topography. Through multivariate analysis, we found an environmental gradient correlated with light penetration into the understory and a gradient associated with the slope and, to a lesser extent, with the elevation. Then, we identified floristically well-differentiated ecological groups in response to environmental conditions; however, the groups only partially coincided with the understory composition of the predefined forests. We found environmental response species groups such as Roupala obovata and Beilschmiedia sulcata that are indicator species of sites with lower light penetration into the understory but with steeper slopes and higher elevation. In comparison, Clusia multiflora and Zanthoxylum quinduense to be the main indicator species from sites with greater light penetration into the understory and lower slope and elevation. These findings support appropriate species selection when implementing restoration strategies in forest landscape restoration plans.

Introduction

Tropical humid forests represent the most diverse ecosystems on the planet, and amidst these, tropical mountain cloud forests (TMCFs) have the greatest biological diversity and fragility (Brown & Kappelle 2001); however, they are subject to high deforestation rates and damaging human interventions and the effects of climate change (Aide et al. 2019, Fadrique et al. 2018, Rodríguez-Morales et al. 2009). The TMCFs vegetation develops typically on rugged topography under conditions of persistent cloudiness, high rainfall, relatively low temperatures, and high humidity (Hogan & Machado 2002, Aparecido et al. 2018, Fahey et al. 2016, Luna et al. 2001, Oliveira et al. 2014). Low solar radiation reaching the canopy is perhaps the most limiting resource for plant growth and reproduction in these ecosystems (Pearcy 2007, Malhi et al. 2017), mainly due to high cloudiness and terrain inclination. A complex combination of topographical variables such as elevation, slope, and aspect, influence local climate variables such as air temperature, precipitation, wind, and solar radiation (Aparecido et al. 2018, Cavelier 1996, Fahey et al. 2016, Liu et al. 2023, Oliveira et al. 2014). In addition, the likelihood of natural disturbances increases with slope, substantially changing the local light regime. The interactions between the light environment and vegetation in the forest understory play a fundamental role in shaping forest ecosystems’ structure, diversity, and function, influencing the future forest community structure (Huo et al. 2014). The understory of TMCFs typically shows very low photosynthetic photon fluxes, which can suddenly increase after tree falls and landslides (Clark et al. 2016, Quevedo-Rojas et al. 2018, Schwarzkopf et al. 2011). Disturbances also determine mountain ecosystems’ structure, composition, and functioning (Dar & Parthasarathy 2022). Cloud immersion affects the amount of solar radiation that reaches TMCFs affecting temperature and humidity. The dynamics of structure and composition of these forests are ultimately the result of the relationships of climate, topographical variables (influencing soil structure, nutrients, and water relations), and floristic composition, especially the understory juvenile tree species. Even so, the role of juvenile tree communities in forest development is one of the least studied areas in the forest ecology of cloud forests (Dar & Parthasarathy 2022, Homeier et al. 2010, Rahman et al. 2017). This study evaluated the floristic composition of the juvenile tree species in the understory of four cloud forest types in the Venezuelan Andes. These forests were classified from aerial photographs (Rangel 2005) according to combinations of canopy height (tall canopy, 25–30 m height; medium, 20–25 m; low, 15–20 m); upper canopy density (dense to sparse), and the number of strata observed from the ground as indicators of structural complexity. We hypothesize that the understories of these forests differ in their juvenile tree floristic composition due to differences in canopy structure/floristic composition, light penetration, and topographic variables; therefore, we anticipate the existence of species or groups of indicator species in the understory of environmental conditions characteristic of each forest. Our goal was to answer the following questions: (i) Does the understory of these forests differ in the floristic composition of juvenile trees? (ii) Do canopy structure, light, and topographical variables constitute environmental gradients associated with variations in the juvenile tree floristic composition across the understory of the forest types? and (iii) Are there species or groups of indicator species characterizing the understory of each forest?

Materials and methods

Study area

The study was carried out in the San Eusebio University Forest (SEUF) which contains representative vegetation communities of the Andean cloud forest ecosystem. The SEUF (ca. 369 ha) is located in the Andes mountains of Mérida (8°37′00″ N–71°21′00″ W), Venezuela (Figure 1). The forest has an irregular structure, with broadleaf and conifer species, being the latter (Podocarpaceae) a remarkable component of the canopy. It is rich in tree species colonized by abundant and diverse epiphytes (Lamprecht & Veillon 1967). Besides the Podocarpaceae; Lauraceae, Melastomataceae, Euphorbiaceae, and Myrtaceae are the predominant tree families (Ramos & Plonczak 2007). Elevation ranges from 2220 to 2440 meters above sea level (masl) (Schwarzkopf et al. 2011) and the climate is characterized by persistent cloudiness occurring year-round, an average annual precipitation of 1,500 mm with a short dry season between January and March, an average year temperature of 14.9°C, and relative humidity above 90%. The general topography goes from rounded hills with gentle slopes to steep and abrupt hills and ravines (Márquez 1990).

Figure 1. Location of the SEUF study area in Mérida, Venezuela. In the figure at the right is the legend for the total area and the chosen forests: tall dense forest (TDF), low dense forest (LDF), medium-tall medium-dense forest (MMDF), and low sparse forest (LSF). The yellow parallel lines are the virtual transects, and the green points represent the points where plots are located (adapted from Rangel 2005).

The area is a pristine primary forest (Rollet 1984), although some anthropic disturbance occurred in the 1950s (Schwarzkopf et al. 2011). The four forest types identified by Rangel (2005) were classified as tall dense forest (TDF), medium-tall medium-dense forest (MMDF), low dense forest (LDF), and low sparsely-dense forest (LSF) (Figure 1). The TDF is characterized by a tall continuous canopy (25–30 m) and a complex structure with three well-differentiated strata: an upper dense stratum between 25 and 30 m high, with sparse emergent trees that can reach 40 m in height; an intermediate stratum between 20 and 25 m tall, formed by medium-size tree species and smaller trees of emergent species; and a lower stratum 10–15 m tall formed by small trees, shrubs, and regeneration of large trees. The MMDF is characterized by an upper stratum 20–25 m tall formed by an almost continuous canopy with medium-sized emergent trees. A second stratum is 15–18 m tall, and a third stratum is around 10 m tall composed of small trees, young medium-sized trees, and shrubs. The LDF has two strata, the upper dense one forming a continuous canopy 18–20 m tall and a second stratum 10–15 m tall. The LSF has a single open stratum of 8–10 m, with sparse emergent trees rarely reaching 15–20 m tall (Figure S1). Gaps in these forests are usually small with areas of 50–400 m² (Quevedo et al. 2016). Larger gaps formed by disturbances such as landslides are quickly covered by bamboo (Chusquea sp.) and vines that often limit the regeneration of tree species.

Sampling

We used the vegetation map of Rangel (2005) for delimiting the areas to sample in the selected forests, excluding sites previously identified as being close to the land-farm frontier, sites logged in the past, or with signs of previous silvicultural intervention (e.g., presence of planted trees). Further, we established a stratified sampling design with proportional allocation, i.e., the sample size was determined according to each forest’s relative weight (area). On the map, a starting point was randomly selected in each forest and identified by its geographic coordinates. From these points, ‘virtual’ sampling grids were laid out using MapSetToolkit v.1.77 (Anon 2022). The ‘virtual’ grids consisted of parallel transects of variable length (100–1,000 m), traced 100 m apart and perpendicularly oriented in the predominant slope direction. Also, the starting point of each transect was randomly chosen, and successive points were located systematically every 30 m (Figure 1). The base map, grids, and points were transferred to a Garmin 60 CSX GPS, so in the field, the sample points could be accessed through the shortest/easiest path, minimizing the aperture of trails in the forest. The geographic coordinates (Universal Transverse Mercator -UTM-) and elevation (masl) were stored for each sample point. The position error (±5 m) was minimized with the “waypoint averaging” option. At each point, a 1-m radius circular plot (3.14 m²) was laid out, and the juveniles of tree species 30–150 cm tall were identified and measured, excluding resprouts. Individuals below 30 cm tall were discarded since they have a large chance of death from droughts, floods, and pathogens. For each juvenile tree, we recorded the species name and total height. In addition, large plots (200 m², r = 7.92 m) were established surrounding a fraction of the 3.14 m² plots (one in five). In these large plots, all trees above 1.50 m tall were measured for total height and diameter at breast height (dbh -1.30 m), considering (a) trees above 1.50 m tall and dbh ≤ 5 cm and (b) trees with dbh ≥ 5 cm. Also, the total height of the three largest trees per plot was measured to obtain the average upper canopy height. To avoid plot shape distortion, we applied slope correction to the plots when the slope was ≥5%. In total, 749 (3.14 m²) plots were established. For vegetation and environmental variables analyses, 36 plots were discarded as juvenile trees were absent or light measurements could not be taken. The final sample size (n) for the analysis of juveniles was 653 plots (see Table S1), TDF (n = 261), MMDF (n = 170), LDF (n =175), and LSF (n = 47) and 105 large plots (Table S2) for adult trees, TDF (n = 32), MMDF (n = 31), LDF (n = 30), and LSF (n =12).

Sampling effort

We used a systematic design with a random selection of the first plot to ensure adequate coverage of the delimited area. Systematic sampling can be more than twice as efficient as random sampling for a similar sample size (Iles 2003). These forests have a rugged topography, a dense understory with shrubs, vines, and bamboo as well as frequent rains and high humidity that limit time for field work and use of equipment. We chose a small plot size for measuring juveniles due to the large abundance of juveniles of vines, shrubs, and trees. Also, we looked at causing minimum disturbance by reducing the need to walk within the plot. On the other hand, using bigger plots could increase the potential to miss or step on individuals. Yet, the efficient selection of sampling size for evaluating plant communities remains a debatable issue (Stohlgren 2007, Sgarbi et al. 2020). The sampling intensity for each forest was between 0.17% and 0.22% for juvenile trees, whereas for the adult trees, it was between 1.61% and 6.22%. Overall, the sample fraction for juvenile trees in each forest was well-balanced (Tables S1 and S2).

To verify adequate sample coverage, we determined the sampling effort using the species accumulation curve (SAC), considered more appropriate to estimate sampling effort than the traditional species-area curve, as our sample had non-contiguous plots (Dengler 2008). Chazdon et al. 2023 define ‘sample coverage’ as the total relative abundances of the observed species, or equivalently, the proportion of the total number of individuals in an assemblage that belong to the species represented in the sample. The software iNEXT (Interpolation-Extrapolation) from Chao et al. (2014) was used to estimate the sampling effort in terms of the Hill numbers according to Hsieh et al. (2016): q ₀ (species richness), q ₁ (Shannon index or evenness, representing the typical species), and q ₂ (Simpson index of diversity representing the dominant species). We generated the SAC for each q _i number with 95% confidence intervals and the sampling effort for observed, interpolated, and extrapolated number of samples (until doubling the reference sample size).

Environmental variables

To estimate the light environment and canopy structure variables, we took hemispherical photographs in the centre of each circular plot with a Nikon FC-E8 hemispherical ‘Fisheye’ lens (optical field = 183°) attached to a 10-megapixel Nikon COOLPIX P5000 digital camera mounted on a tripod and levelled to keep the equipment completely horizontal. The photos were taken at 1.5 m aboveground under overcast conditions following the protocol of Zhang et al. (2005) for determining correctly the degree of light exposition (control of lens aperture and shutter speed), thus optimizing the estimates of light variables. Extensive field testing at sites with open sky but uniform cloud cover allowed us to find the appropriate lens aperture (f5.3 to f5.4) and shutter speeds (1/125 to 1/250) with ISO 200 (Quevedo-Rojas et al. 2015). We pre-processed the images to correct the lens distortion (183° to 180°) according to Frazer et al. (2001) based on a third-order polynomial:

$$Y = 6.638X - 0.0025X^2 - 2.401E - 0.5X^3\;\;\;\;\;0^\circ \le X \le 90^\circ$$

where Y is the radial position of a projected point measured in pixels from the optical centre of a full-resolution digital image (1600 × 1200 pixels) and X is the angular distance. We used SIDELOOK v.1.1 (Nobis & Hunziker 2005) to automatically obtain the optimum thresholds corresponding to open sky and vegetation. We processed the corrected images with Gap Light Analyzer (GLA) v. 2.0 (Frazer et al. 1999) to process and analyse the hemispherical photographs. The GLA simulates the solar radiation regime and vegetation characteristics. The software uses the latitude, longitude, growing season, elevation, and slope to estimate the parameters related to light, calculating the position of the Sun for the hours of the day and the days of the year, as well as the effect of the slope on the amount of theoretical light that would reach a surface unit at 1.5 m from the ground. Also, GLA estimates several variables related to the light environment and canopy structure. The proportions of direct and diffuse radiation in the cloudy sky were determined as indicated by Frazer et al. (1999). We estimated leaf area index (LAI), percentage of canopy openness (%CO), percentage of direct sunlight (%TDir), and diffuse sunlight (%TDif) transmitted from the canopy to the understory and variables indicating topographical position (slope and elevation). For each plot, we measured the slope (%) in the steepest direction with a clinometer CST/Berger 6-3/8, and then we transformed it to an inclination angle (degrees) concerning the horizontal plane.

Relationships between adult and juvenile tree composition

Adult species composition could influence understory composition; nevertheless, making a direct study of a cause–effect association between the presence of juvenile and adult trees at the plot level was beyond the scope of this study. To gain an insight into the possible relationships between the relative importance of adult trees of a given species and the relative importance of juveniles of that same species in the understory, we calculated the Species Importance Value Index (IVI) and the Relative Natural Regeneration (RNR). The IVI is the well-known importance value of the American Forestry School (Curtis 1959, Curtis & McIntosh 1951) which estimates the importance of each species in terms of their relative dominance (%RD), frequency (%RF), and abundance (%RA) considering trees ≥5 cm dbh. The IVI = (%RD + %RF + %RA)/3 must add up to 100 considering all tree species. On the other hand, the RNR (Finol 1971) is similar to the IVI but considers all juveniles up to 150 cm tall, where RNR = (%RA + %RF + %RSC)/3 must add up to 100; with RSC = relative size category (percentage of individuals per species in each juvenile size category). The IVI (adults) of each species by forest was compared with the %RNR (juveniles) to check whether both indicators are correlated for each species. In addition, we determined the distribution of heights by species and forest type. The shape of the height distribution indicates how abundant is the regeneration of a given species in regards to the abundance of adult trees.

Statistical analysis

For each forest, we calculated the descriptive statistics for the environmental variables and then we tested for the differences of each variable between forests through a nonparametric one-way analysis of variance based on the Kruskal–Wallis test (Kruskal & Wallis 1952) using PAST 4.03 (Hammer et al. 2001). Pairwise differences were determined with a posterior Dunn’s test to compare the medians of environmental variables between forests. Then, the p-values were corrected using the Bonferroni correction for multiple tests (p = 0.05). Further, a multivariate partial correlation test was run to analyse the correlation structure between environmental variables.

We used PC-ORD v 5.0 (McCune & Mefford 1999) to run a detrended correspondence analysis (DCA) (Hill 1979a) for comparing the floristic composition of juveniles between forest understories. In the DCA ordination plane, the closest sites (plots) reflect a greater similarity in floristic composition than sites farther away. Additionally, we checked the correlation of the environmental variables with the ordination axes to identify the dominant environmental gradients. Further, a two-way indicator species analysis (TWINSPAN) was run to corroborate the existence of groups of species assemblages (Hill 1979b). In addition, significant differences among groups were tested with the multi-response permutation procedure (MRPP) of Mielke & Berry (2007), a non-parametric method that examines differences in the assembly structure between groups defined ‘a priori’. We chose the squared Euclidean distance as the most appropriate for abundance data (McCune & Grace 2002); and then we applied indicator species analysis (ISA) to identify the species discriminating the TWINSPAN groups (Dufrêne and Legendre 1997). For each species, the ISA estimates a species indicator value (SIV) between 0 and 100, reaching the maximum when all individuals of a species occur in all sites of a group. The SIV expresses, in percentage, the capacity of a given species as an indicator of a group of sampling units chosen ‘a priori’ or as a product of a grouping. The method selects indicator species based on both high specificity and high fidelity to a specific group. This technique is robust to differences in sample sizes within groups and abundances among species. The statistical significance of the indicator value of each species was computed by a Monte Carlo randomization process with 1,000 permutations. Afterwards, we applied the weighted averaging (WA) method (Gauch 1982, ter Braak & Looman 1995) to estimate each species’ optimum regarding the environmental variables. Finally, we performed a cluster analysis to define species groups based on their environmental similarities determined from the WA using Ward’s method as a group link method. If not otherwise mentioned, PC-ORD v. 5.0. was used in all analyses.

Results

We found a total of 53 species of juvenile trees for the whole sample, but for multivariate analyses discarded 13 of them with less than 4 individuals per species to reduce the noise in the data that might distort the analysis (these species comprised only 28 individuals); therefore, for analyses, we used 40 (3,858 plants) species (Table 1). The four forests differed in structural characteristics (Table 2). Increasing canopy height is observed from the low to the tall forests, and increasing stand density is observed from sparse to dense forests. Also, the distribution of trees by height categories showed a pyramidal shape for all forests. The distribution of most species was also pyramidal suggesting an adequate number of individuals to replace the upper categories (Figure S2).

Table 1. List of identified species of juvenile trees in the cloud forest sampled area. San Eusebio University Forest, Mérida, Venezuela

Table 2. Average and standard deviations of selected forest structure variables (trees > 5 cm dbh) by forest type. TDF (tall dense forest); MMDF (medium-tall medium-dense forest); LDF (low tall dense forest); and LSF (low tall sparse forest)

Defining what is a rare species within a sample is usually a researcher’s decision (McCune and Grace 2002, Poos and Jackson 2012). Usually, species occurring at a single site are eliminated (Legendre & Legendre 1998, Poos & Jackson 2012). Other researchers suggest removing rare species that occur in <5% (McGarigal et al. 2000) or <10% (Marchant 1990, McCune & Grace 2002) of sites (see Poos & Jackson 2012).

For the juvenile trees (height = 30–150 cm), the rarefaction–extrapolation curves indicate that the sampling effort is adequate as the curves for the Hill numbers tend to reach a plateau as the number of sampled individuals increases (Figure S3). Only for q ₀ (species richness), the graph for the LSF indicates a rather low sampling effort, but for q ₁ and q ₂ emphasizing typical and dominant species, there is an ample plateau. Also, the sampling effort efficiency was between 0.975 and 1 for each forest and diversity index as shown by the coverage-based rarefaction and extrapolation for the juveniles (Figure S4) and adult trees (Figure S5), respectively.

To analyse the relationships between adult and juvenile tree composition, we calculated the IVI for each forest (Table S3) which did not show a clear dominance of any particular species in any forest, with a maximum of two species showing an IVI > 10% except the LSF with five species > 10%. There were no dominant species across sites, with none representing more than 19% of the total RNR for a given forest. Several species showed a positive relationship between IVI (adult importance) and RNR (juveniles), but others showed an inverse relationship, i.e., species with high IVI had very low RNR. For example, in the TDF, Myrcia fallax had the highest IVI (15.9) but a relatively low RNR (6.8); on the other hand, Aegiphila ternifolia with a low IVI (1.0), showed the largest share of RNR (12.2). Also, Retrophyllum rospigliosii, a species that represents a large fraction of the canopy cover due to the large crowns of emergent individuals, occupied the second place for both IVI and %RNR.

In the MMDF, M. fallax had the largest IVI (13.6) and a good portion of the RNR (9.2); however, Eschweilera monosperma, second in IVI (9.2) showed a very low RNR (2.6); whereas, Myrcia acuminata in third place with an IVI of 7.8 has the largest share in RNR (15.5). Also, M. fallax and M. acuminata are within the three most important species in terms of IVI for three of the forests and have a good share of the RNR in these forests. On the other hand, in the LSF, four species (Hyeronima moritziana, Miconia resimoides, Alchornea grandiflora, and Clusia minor) account for more than 50% of the IVI; however, together they represent only 15% of the RNR; meanwhile, less important species have a higher RNR, for example; Zanthoxylum quinduense (16.2), Eugenia tamaensis (10.4), and Prunus moritziana (10.4). So, some species that dominate the canopy in the adult stage (IVI) show low regeneration; whereas others with a low IVI have a high RNR which could lead them to dominate the canopy in the future. Although we could observe many seedlings (<5 cm tall) of several species below the crown of adults, this did not imply a high abundance of 30–150 cm juvenile trees. Most seedlings usually die after the first year of germination.

Forest types and environmental variables

There were significant differences (p < 0.01) for the median of environmental variables between forests as indicated by the Kruskal–Wallis test (i.e., for a given variable, at least one forest differed from the others). The posteriori Dunn test (Bonferroni corrected) showed the statistical significance (p < 0.05) of pairwise comparisons between forests. For example, MMDF and LDF differed significantly only in slope and elevation. The TDF and the LSF differed significantly in all variables (Figure 2). Also, the TDF differed from the other forests, showing the lowest values for variables related to transmitted light (%TDir, %TDif) and %CO, but significantly higher values for LAI. Regarding topographical variables, all forests differed significantly in slope and elevation. The multivariate partial correlation test (Table 3) showed high negative to positive (−0.80 to 0.82) correlations between light/structure variables and low correlation between topography variables (0.2). Also, correlations between light/structure and topographical variables were low (−0.09 to 0.07).

Figure 2. Descriptive statistics and Kruskal–Wallis analysis of variance test results of environmental variables: canopy openness (%CO), leaf area index (LAI), transmitted direct light (%TDir), transmitted diffuse light (%TDif), slope (Slope%), and elevation (masl) by forest: tall dense forest (TDF), medium-tall medium-dense forest (MMDF), low dense forest (LDF), and low sparse forest (LSF). Violin shapes show sample plot distribution. The same letters (a, b, c, d) for forest represent no significant differences (p = 0.05) for a Dunn post-hoc test with Bonferroni, correction for multiple tests. SEUF, Mérida, Venezuela.

Table 3. Multivariate partial correlation analysis between environmental variables. Correlation between (a) light/canopy structure variables (yellow), (b) between topography variables (orange), and (c) between a and b (green)

The DCA showed gradient lengths ≥ 4 standard deviations, indicating the appropriateness of using this method (Lepš & Šmilauer 2003). The first and second ordination axes explained a variance of 26.1% (18.8% and 7.3%, respectively). Canopy structure, transmitted light, and topographical variables were well correlated with the first two axes. In the ordination graph, the TDF plots were displayed mainly to the left of axis 1, clearly separated from the LSF plots located to the right of this axis; however, the LDF and MMDF plots were mainly at the centre of the ordination plot, showing a high overlap (Figure 3). Thus, the ordination matched only partially at the understory level with the Rangel’s forest classification. Since there was no clear separation of the understory among forest types with the DCA, we explored whether the sites could be separated with a TWINSPAN and tested for statistical differences with MRPP. TWINSPAN generated four groups (G): G1 (n = 151), G2 (n = 139), G3 (n = 291), and G4 (n =72). The value for the statistic A was 0.16 (p < 0.0001), indicating that at least one group had a different floristic composition. Notice, that for ecological communities, the A values are usually below 0.1, with A = 0.3 considered high (McCune & Grace 2002). Further, to detect which groups were different, we tested for pairwise differences between groups. Significant differences (p < 0.0001) were found for all comparisons: G1 and G2 (A = 0.09), G1 and G3 (A = 0.13), G1 and G4 (A = 0.11), G2 and G3 (A = 0.09), G2 and G4 (A = 0.14), and G3 and G4 (A = 0.09). These understory groups do not match the site composition observed in the forest types. Rather, Forests 1 and 4 match Groups G1 and G4. Whereas, Groups G2 and G3 do not match any of the forest.

Figure 3. First two axes of the detrended correspondence analysis (DCA) show the distribution of plots by Forest types. Diamonds correspond to the tall dense forest (TDF), triangles correspond to the medium-tall medium-dense forest (MMDF), circles are the low dense forest (LDF) and squares are the low sparsely-dense forest (LSF). SEUF, Mérida, Venezuela.

Then, we ran the DCA again with the same specifications, but using the groups created with TWINSPAN (understory composition) as a categorical variable. The DCA graph (Figure 4) shows the ordination axes, sites (points), environmental variables (arrows), and species optima (crosses). The correlation among environmental variables along the ordination axes indicated a significant positive correlation (p < 0.001) with axis 1 for %CO, %TDir, %TDif, and slope; while LAI and elevation had a significant negative correlation, suggesting that Axis 1 is related to an increase of both, %TDir and %TDif because of a larger %CO. There was a positive correlation between LAI and elevation along the first axis and a negative correlation with %CO, %TDir, and %TDif, suggesting increasing LAI with elevation and; therefore, lower %CO, TDir, and %TDif. On the other hand, Axis 2 showed significant correlations, although with a lower percentage of explained variance, with all variables except LAI. Correlations were negative with %CO, TDir, and %TDif; whereas they were positive with topographical variables (Table 4).

Figure 4. DCA ordination showing plots, species, and environmental variables. Site categories correspond to the TWINSPAN groups: G1 (blue diamonds), G2 (red circles), G3 (green triangles), and G4 (pink squares). Groups delimited by ellipses. Crosses show the optimal sites for juveniles of species in the understory. The arrows indicate the direction and magnitude of the environmental variables. For species abbreviations, see Table 1. San Eusebio University Forest, Mérida, Venezuela.

Table 4. Pearson correlation coefficients (r) and p-values of environmental variables with the two first axes of the detrended correspondence analysis (DCA; n = 653). Percentage of canopy openness (%CO), leaf area index (LAI), percentage of transmitted direct light (%TDir), and transmitted diffuse light (%TDif) through the canopy reaching the understory. San Eusebio University Forest, Mérida, Venezuela

In the DCA (Figure 4), Group G1 (TWINSPAN) is located towards the left end of Axis 1, meaning that the plots from this group represent environments where higher LAI and elevation predominate. In general, these are sites with the lowest light availability. Next, G2 is located close to G1, towards the left of Axis 1, indicating high values of LAI and elevation although lower than G1; whereas the second axis correlates to steeper slopes.

On the other hand, G3 is towards the center of the ordination plot showing intermediate values of LAI, %CO, %TDir, and %TDif. Nonetheless, a considerable proportion of plots in this group are placed at the higher end of Axis 2, indicating steeper slopes. Instead, G4 is located between the centre and the right-end of Axis 1 at an angle of 45° with both axes indicating larger %CO, %TDir, and %TDif.

For Group G1, the ISA (SIV Table S4) included nine species, where A. ternifolia (SIV = 62.6%; p = 0.001), Myrcianthes karsteniana (SIV = 27.6%; p = 0.001), and Casearia tachirensis (SIV = 21.7%; p = 0.001) had the largest indicator values. In G2, the ISA included four species, but only M. fallax (SIV = 64.9%; p = 0.001) was relevant. Ten species were included in G3 with M. acuminata having the highest indicator value (SIV = 42.7%; p = 0.001). Also, Ladenbergia undata (SIV = 15.2%; p = 0.001) stands out because of its position at the steepest slopes. The Group G4 included 17 species, being the most important Z. quinduense (SIV = 40.1%; p = 0.001), E. tamaensis (SIV = 39.6%; p = 0.001), P. moritziana (SIV = 29.5%; p = 0.001), Clusia multiflora (SIV = 25.0%; p = 0.001), Weinmmania jahnii (SIV = 24.2%; p = 0.001), and A. grandiflora (SIV = 22.0%; p = 0.001). The species Z. quinduense, at the right and down in the ordination plane, has its optimum in places with the highest light availability and the lowest LAI and elevation. In the same group, C. multiflora and A. grandiflora appear towards the lower right end of Axis 2, suggesting that their optimum values are found in sites with lower slopes and elevation.

The WA method estimated each species’ optimum per environmental factor (Table S5). From the WA results, the Cluster Analysis produced four clusters (C) of species. We obtained the average value for each environmental variable which determines the optimal environmental conditions for each species cluster (Table 5). Thus, cluster C-I included eight species whose estimated environmental optima are characterized by the highest %CO, %TDir, and %TDif, and the lowest LAI, slope, and elevation. All species of C-I were indicator species in G4 (TWINSPAN). These findings agree with those found in the DCA: all the species in C-I are located at the right end of Axis 1 and the lower end of Axis 2, where most LSF plots are located. In this forest, we found the greatest penetration of light, the sparser canopy, and the lowest elevation and slope, the latter being very characteristic of flood-prone sites. Opposed to C-I, C-III showed the lowest values for light, i.e., %TDir, %TDif, and %CO; but the largest LAI, which suggests darker conditions in the understory. Further, it showed the highest elevation and the steepest slopes. On the other hand, C-III included species with high SIV associated with G1 and G2 (TWINSPAN), specifically M. fallax, A. ternifolia, M. karsteniana, and Miconia meridensis. These species are located at the left end of the ordination Axis 1, especially A. ternifolia and M. karsteniana associated with sites of higher LAI, elevation, and moderate slopes. Finally, C-II and C-IV have intermediate values of %CO, %TDir, %TDif, LAI, and elevation, although C-II showed steeper slopes than the other clusters. Environmental conditions in C-II are optimal for eight species, being E. tamaensis, P. moritziana, and W. jahnii indicator species for G4 (TWINSPAN); whereas Nectandra rigida, although not identified as an indicator species, was associated with the steepest slopes as shown by their location at the upper end of DCA Axis 2.

Table 5. Environmental averages of individual optima for the species clusters obtained by cluster analysis (Ward method). Percentage of canopy openness (%CO), leaf area index (LAI), percentage of transmitted direct light (%TDir), and transmitted diffuse light (%TDif) through the canopy reaching the understory. For species abbreviations, see Table 1. San Eusebio University Forest, Mérida, Venezuela

Discussion

We hypothesized that the understory of the forest types would differ in their juvenile tree floristic composition due to differences in canopy structure/floristic composition, light penetration, and topographic variables and anticipated the existence of species or groups of indicator species in the understory of each forest, for example, by observing a predominance of shade-tolerant species in the understory in the more shaded forests.

The results partially confirmed our initial hypothesis: only two forests, TDF and LSF could be separated in terms of their floristic composition and environmental conditions; whereas the MMDF and the LDF were more similar and were not clearly separated. The DCA separated the TDF from LSF; but MMDF and LDF showed a considerable overlap of sites in the ordination space (Figure 3), suggesting few differences in the juvenile floristic composition between these latter forests. So, there were no clear-cut distinctions that could be detected from the ‘a priori’ forest classification. Nonetheless, TWINSPAN grouped the sites in terms of their internal homogeneity in species composition and species heterogeneity among sites of different groups. We found four TWINSPAN groups with statistically significant differences in species composition; two of them matching the TDF and LSF (G1 and G4, respectively); whereas G2 and G3 were a mix of the MMDF and LDF sites, but differing significantly in their floristic composition.

Overall, we found low light levels reaching all understories. Low light transmission was associated with a low %CO, agreeing with studies in tropical cloud forests. Thus, Sylvester & Avalos (2013) found a %CO between 8.53 and 17.05 in a cloud forest in Costa Rica; whereas, Acevedo et al. (2003) and Quevedo-Rojas et al. (2015) found ranges between 3.4 and 9.5 and between 0.5% and 12.8%, respectively, for the Andean cloud forests of Venezuela. Also, DeCarvalho & Oliveira-Filho (2001) reported a %CO between 5% and 18% for a cloud forest in south-eastern Brazil. These patterns suggest that low %CO is the rule in these forests, except for scattered gaps produced by tree falls where the maximum aperture in %CO does not reach 30%. In most tropical and temperate forests, the frequency distribution of light environments is highly biased towards microsites with <5%CO which is probably a consequence of the complex and multi-layered vertical structure of these forests (Lusk et al. 2006). This effect is further enhanced by the geometric effects of slope on radiation income and tree crown architecture (Montgomery & Chazdon 2002).

The average LAI varied between 2.9 and 4.9; within the range observed in other cloud forests. For the Andean cloud forests of Venezuela, Schwarzkopf et al. (2011) reported a mean LAI < 3.0 for three cloud forests situated at different average elevations. Likewise, Acevedo et al. (2003) reported an average LAI of 2.3 at 2,300 masl Quevedo et al. (2016) reported in a dense cloud forest an average LAI of 5.4 and 4.2 in gaps and undisturbed areas respectively. In the cloud forest of other regions, decreasing LAI with elevation was reported. Leuschner et al. (2007) found a decrease in LAI from 5.1 at 1,050 masl to 2.9 at 3,060 masl on an elevation transect in Ecuador. Consideration should be given to the effect of branches, trunks, and epiphytes on the values of LAI in TMCFs. Moser et al. (2007) found that the accumulated area of elements other than leaves increased considerably with altitude. An increase in epiphytes with elevation could influence the observed increase of LAI with elevation in our study; despite the relatively narrow range of elevation (2,320–2,460 m).

Light transmission was closely correlated with %CO and opposed to LAI along Axis 1. Both, %TDir (5.7–20.63) and %TDif (4.8–15.3) showed a relatively wide range; however, for these forests, both variables were <10% typical of closed canopies, except for the LSF in which both variables were >15% more typical of forests with sparser canopies. For example, Sylvester and Avalos (2013) found relatively high values of TDir and TDif for three cloud forests in Costa Rica (16–21% for both variables).

Although the second ordination axis explained less total variance, it suggested a significant linear correlation with topography, mainly with the slope and at a lower degree with elevation. Several studies confirmed the importance of the slope in the distribution and floristic diversity in the understory since it directly affects the microclimate, soil properties, water availability, and other resources under the tree canopy (Dar & Parthasarathy 2022, Huo et al. 2014). Slope variation not only influences understory light availability by determining the spatial disposition among tree crowns and the seasonal differences in sunlight penetration, but also affects water infiltration rates. For example, on higher slopes, lower water availability in drier sites may result in a higher probability of desiccation in high-light microsites, especially for drought-intolerant species (Brenes-Arguedas et al. 2011). Previous studies in the area (Hettsch 1976, Quiroz-Sandoval 2010, Valcarcel 1982) found a hydro-sequence of soil conditions associated with slope (e.g., soil water potential, hydromorphism, anoxia, organic matter decomposition rates, and nutrient availability) associated with changes in forest structure (e.g., basal area, tree density, and canopy height). Thus, the slope can be an important factor in separating groups of juvenile tree species based on their shade tolerance coupled with their tolerance to water conditions in the soil.

The ISA and the cluster analysis indicated the presence of indicator species groups and their optimal environmental conditions which provide insights into species shade tolerance. For example, Clusters III and IV included species whose juveniles were reported as shade tolerant in ecophysiological experiments: A. ternifolia and M. karsteniana did not tolerate sudden increments in radiation, suffering chronic photoinhibition; whereas C. tachirensis, Beilschmiedia sulcata, and R. rospigliosii showed larger photosynthetic plasticity responding favorably to light increments (García-Núñez et al. 1995, Quevedo-Rojas et al. 2018). On the other hand, Tetrorchidium rubrivenium and M. meridensis found in Cluster I, dominating areas with relatively high light availability (%CO), showed different degrees of acclimation to sudden increases in light levels like those occurring when ‘gaps’ are formed by tree falls. These species are resistant to photoinhibition (Quevedo-Rojas et al. 2018). For Clusters I and II, representative species were A. grandiflora, C. minor, C. multiflora, and M. resimoides identified as gap colonizing species and water deficit avoiders (Ataroff & García-Núñez 2013, García-Núñez et al. 1995, Rada et al. 2009). Planted species trials monitored up to seven years old showed that species such as T. rubrivenium, and M. meridensis behave as shade-intolerant species with fast growth and rapid crown recession; instead, Hieronyma moritziana and Billia colombiana are partially shade-tolerant needing a certain amount of shade to grow adequately. Finally, Myrciantes karsteniana, M. fallax, and A. ternifolia are shade-tolerant as expressed by slow growth and persistence of leaves under conditions of high shade below the fast-growing species (Quevedo-Rojas & Jerez-Rico 2021).

The results showed that despite the complex environment in the understory of the cloud forest, some patterns of juvenile tree composition are associated with transmitted light and topography. Several studies identified tree species according to their shade tolerance analysing their distribution along the light gradient, and their ecophysiological responses through experiments with seedlings and planted species trials. For example, Quevedo-Rojas et al. (2015) found that along a gradient of %CO most juvenile tree species exhibited the highest abundance below 6.6%, being scarce above 12.8%. At higher %CO (e.g., gaps), the survival and growth of trees were impaired by the presence of vines and bamboo. Most species showed large plasticity along this narrow range of %CO, but the distribution of individual species was, in general, non-random, with some species growing in low light and others in high light. However, the overall abundance of juveniles for the whole forest was randomly distributed (Quevedo-Rojas et al. 2015). Nonetheless, the incorporation of other gradients such as water availability, could increase the predictive power to explain these patterns. For example, Brenes-Arguedas et al. (2011) suggest that light availability interacts with other variables such as precipitation, so species with low tolerance to water stress have reduced survival and growth in dry microsites with high light. It could be argued that species composition in the understory is strongly affected by the proximity of conspecific adults; however, this is not necessarily so. For example, MacDougall & Kellman (1992), in a tropical forest of Belize, found that the floristic composition of seedlings was spatially heterogeneous, and differences in light intensity, rather than proximity of conspecific adults had the greatest influence on seedling distribution at the species level. Although median light levels were very low (∼2% of full light), there was enough light variability in patches to support species with different regeneration strategies, facilitating species coexistence. In the SEUF, Quiroz-Sandoval (2010) did not find a relationship between the closeness of conspecific adults and juveniles’ abundance. Jerez et al. (2011) observed, in several Venezuelan forests, that the regeneration of shade-intolerant species is impaired by low light, so most juveniles were of shade-tolerant species. The low proportion of large gaps in the cloud forests and their encroachment by a thick layer of lianas and bamboo might reduce the survival of juveniles of shade-intolerant species. Cloud forest species, however, could adapt to large variations in environmental conditions and may show adaptations as they go through the various growth stages by changing their degree of shade tolerance (Svenning 2000).

We could identify groups of juvenile tree species and their relation with environmental gradients could help in species selection to implement adequate regeneration techniques in the variable environments existing in total or partially degraded forests. Additional studies, like ecophysiological experiments, and integrated research in light, soils, and adult tree composition should be incorporated directly to get a more accurate picture of the determinants of understory composition. Also, active/passive restoration trials across environmental gradients will help to explain the processes of competition, facilitation, and coexistence processes that determine cloud forest dynamics. Our findings can make a significant contribution to the improvement of the ongoing Forest Landscape Restoration Programs in TMCFs in the Andes.