Effects of fire disturbance on species and functional compositions vary with tree sizes in a tropical dry forest

Study site

The research was conducted in a seasonal dry evergreen tropical forest, the 50-ha Huai Kha Khaeng Forest Dynamics Plot (HKK FDP) from north-western Thailand (15°40′N, 90°10′E). The HKK FDP was established in 1990 to 1991. All trees with diameter at breast height (DBH, 1.3 m above ground) ≥1 cm were tagged with a unique number and identified to species (Condit, 1998; Anderson-Teixeira et al., 2015). Trees were measured and censused at a 5-year interval from 1994 to 2014. The plot contains about 300 tree species, most of which are evergreen. Dominant species include Croton roxburghii (Euphorbiaceae), Dimocarpus longan (Sapindaceae), and Polyalthia viridis (Annonaceae). The elevation of the plot ranges from 549 m to 638 m above sea level (Bunyavejchewin et al., 1998). The mean annual temperature is around 23.5 °C (Bunyavejchewin et al., 2009). Mean annual rainfall is 1,500 mm (Baker & Bunyavejchewin, 2017), with a 6-month dry season beginning in November (Bunyavejchewin et al., 2003). The major disturbances are low-intensity surface fires. There were three fire events recorded between 1990 and 2020 (Baker, Bunyavejchewin & Robinson, 2008; Trouvé, Bunyavejchewin & Baker, 2020), and fire occurred approximately once in a decade. The 1992 fire burned the south-east corner of the 50-ha plot (2 years before the first full census in 1994), the 1998 fire burned through the entire area of the plot (Bunyavejchewin et al., 2009), and the 2005 fire burned only some parts mostly around the border and very slightly in the middle of the plot (Trouvé, Bunyavejchewin & Baker, 2020). The HKK FDP had no records of human disturbances (Bunyavejchewin et al., 1998).

Functional trait measurement

We collected functional trait samples of the 49 most abundant species (Table S1) from June to September 2016. We collected leaf samples from trees across four DBH size classes: (1) ≤5 cm, (2) >5 cm to ≤10 cm, (3) >10 cm to ≤20 cm, and (4) >20 cm. We sampled at least three individual trees in each diameter class of each species. We collected three leaves and ten leaflets from the species with simple and compound leaves, respectively, except for of C. roxburghii. For C. roxburghii, only two individual trees from the smallest size class were sampled. All leaf samples were collected from trees inside the 50-ha plot. We measured leaf chlorophyll content (Chlo, SPAD) from each individual leaf by using a hand-held ‘SPAD-502 Chl meter’ (Minolta Camera Co., Osaka, Japan). We measured leaf area (LA, cm²), leaf dry matter content (LDMC, mg g⁻¹), leaf thickness (LT, mm), and specific leaf area (SLA, cm² g⁻¹) of each leaf. We calculated the mean values of each species and each size class from these individual leaves. For leaf chemical traits, we took mixed-sample leaves as represented for size class mean values to measure leaf total carbon concentration (LCC, g kg⁻¹), leaf total nitrogen concentration (LNC, g kg⁻¹), leaf total phosphorus concentration (LPC, g kg⁻¹), and leaf total potassium concentration (LKC, g kg⁻¹). All leaf functional traits were measured according to Pérez-Harguindeguy et al. (2013). We collected five wood samples from individuals outside the plot for each tree species and measured wood density following the ForestGEO protocol (www.forestgeo.si.edu/protocols). We calculated the mean wood densities of each species. Note that due to the limitation of collecting wood samples across size classes, we assigned the wood densities of all four size classes of each species the same values as takenfrom the mean wood density of that species.

Quantifying temporal turnover in species and functional composition

We compared species and functional composition of 1,250 20 m × 20 m quadrats between different censuses. We used the 1994 census as a baseline and calculated the temporal turnover of the following time intervals: 1994–1999, 1994–2004, 1994–2009, and 1994–2014. The 20 m × 20 m is an appropriate scale to investigate forest turnover for tropical forests as recommended by Swenson et al. (2006) and Swenson et al. (2012). We quantified species turnover for the community from the 49 most abundant species and for each diameter class using an abundance-based metric, the Bray-Curtis dissimilarity (‘vegdist’ function in the package ‘vegan’ in R). To calculate temporal turnover in functional traits, we performed the log10-transformation of the functional trait values to achieve normality prior to the analysis. We then standardized all functional traits by subtracting the mean and divided by the standard deviation. We calculated the pairwise dissimilarity of the functional traits using the Euclidean distance-based measurement (‘dist’ function in the package ‘stats’) of each 20 m × 20 m quadrat to represent the temporal turnover in functional composition in each quadrat. We did the same analysis for temporal functional turnover in each of the four tree diameter classes (1–5, 5–10, 10–20, >20 cm). The formula of pairwise dissimilarity calculation is described in Swenson et al. (2012). This gives:

$$D_{\mathrm{P}\mathrm{W}}=f_A\sum _{i=1}^{S_A}{f}_i\overline{{\delta }_{ib}}+f_B\sum _{j=1}^{S_B}{f}_j\overline{{\delta }_{ja}}$$where S_A and S_B are the total number of tree species at time A and B in the community, respectively; f_i is the relative abundance of species i at time A, and f_j is the relative abundance of species j at time B in the community. The mean pairwise functional distance between species i at time A and all species at time B in the community is represented by $\overline{{\delta }_{ib}}$ and the mean pairwise functional distance between species j at time B and all species at time A in the community is represented by $\overline{{\delta }_{ja}}$. $f_A$ is the total number of trees at time A divided by the total number of trees at time A and B, and $f_B$ is the total number of trees at time B divided by the total number of trees at time A and B.

Null models for temporal turnover in functional trait composition

We used a null model approach to test if the observed pairwise functional dissimilarity value was significantly different from random. We constructed null distributions of temporal turnover in functional trait combinations for each quadrat using 9,999 randomizations. In each randomization, species names were randomly drawn from all studied species (49 species) for each quadrat. We calculated the pairwise functional dissimilarity of each quadrat using the randomized data. We then compared the observed value with the null distribution and estimated the quantile score (rank) of the observed pairwise dissimilarity value. The quantile scores of observed functional pairwise dissimilarity were significantly lower or higher than expected when the quantile scores of the observed values fell within the bottom 2.5^th or the top 97.5^th percent quantiles, respectively. This also implies that the temporal turnovers in functional composition were more similar or dissimilar than expected by chance.

Temporal species turnover

The Bray-Curtis (BC) dissimilarities tended to increase with the time intervals (1994–1999, 1994–2004, 1994–2009, and 1994–2014). This indicates that the species composition of most quadrats changed over time. On average the BC dissimilarities between 1994 and 1999 were the lowest in comparison to other intervals (Fig. 1 but see additional information in Fig. S1). The distribution of this interval had the highest peak compared with other time intervals.

The temporal turnover in species composition varied across diameter classes (Fig. 2 and Fig. S2). The BC dissimilarities tended to increase with the time interval between censuses for all diameter classes, and the shift was most evident in trees at the smallest diameter class (DBH ≤ 5 cm). In other words, the differences in BC dissimilarities tended to decrease with tree sizes. Furthermore, the average values were the lowest in 1994–1999 in all diameter classes. Trees in the largest size class (DBH > 20 cm) had the lowest compositional turnover between 1994 and 1999.

Temporal functional turnover

The temporal functional turnover did not increase with time intervals (Fig. 3 and Fig. S3). Overall, the quantile scores of pairwise functional dissimilarity were not significantly different from the null models in most quadrats. Within 1,250 quadrats, 77, 49, 50, and 59 quadrats exhibited functional turnover significantly different from random in 1994–1999, 1994–2004, 1994–2009, and 1994–2014, respectively, (Fig. 3). In addition, these quadrats are primarily distributed around the middle of the 50-ha study plot (Fig. 4).

While the overall temporal turnover in functional composition remained unchanged, the functional turnover decreased with tree size (Fig. 5 and Fig. S4). The number of quadrats that had functional turnover significantly lower than expected by random increased with tree size. In most quadrats, pairwise dissimilarities of saplings and treelets (DBH ≤ 5 cm and 5 cm < DBH ≤ 10) were not significantly different from null models. The quadrats with pairwise dissimilarities significantly higher than null models were distributed differently within the HKK FDP among size classes (Fig. 6). That is, in the smallest size class the quadrats with significantly higher functional turnover distributed around the southwest corner of the plot, whereas in the largest size class the quadrats with significantly higher functional turnover primarily occurred at the high elevations (Fig. 6).

Effect of fire on turnover in species composition

Our results showed that the abundance-based species compositional changes tended to increase with the time intervals in HKK in the past two decades (Fig. 1 and Fig. S1). One possible reason is that fire events occurred prior to both 1994 and 1999 tree censuses, and the recovery time from fires was similar for both censuses (Table S2). In 1992, a fire burned some parts of the study site, and in 1998 a fire burned around 1,500 km² of the surrounding area and thoroughly across the 50-ha study plot (Baker, Bunyavejchewin & Robinson, 2008). There are only two years apart between the 1992 fire and the 1994 census, and one year apart between the 1998 fire and the 1999 census. As surface fires killed mostly small individuals (Baker, Bunyavejchewin & Robinson, 2008), species composition 1 or 2 years after a fire was probably composed of surviving trees and recruits of fast-growing species. The number of trees in the 50-ha plot declined from 78,923 to 74,347 individuals between 1994 and 1999, probably due to the 1998 fire which killed a large number of individuals. Since then, only one fire occurred in 2005 and the abundance in the 50-ha plot increased to more than 127,000 trees in 2014.

Even though surface fires have been observed to reduce the number of trees, the abundance of many species can quickly recover with a few years (Cochrane & Schulze, 1999; Baker, Bunyavejchewin & Robinson, 2008). This is because the reduction in tree density opens up more space for recruitment. In addition, a decrease in tree density may weaken density dependent effect due to competition or natural enemies, which may lead to higher seedling survival rate (Augspurger, 1984). Fast-growing species, such as C. roxburghii, had high mortality during fire events but can quickly recover after fire disturbances (Baker, Bunyavejchewin & Robinson, 2008). Unlike other slow-growing species, these species may benefit from the reduction in tree density caused by fires. Abundances of C. roxburghii increased greatly from 8,896 in 1994 to 33,990 in 2014. On the other hand, Baccaurea ramiflora, a thin bark species (Baker & Bunyavejchewin, 2006) that fires greatly affect its mortality (Baker, Bunyavejchewin & Robinson, 2008), has slightly decreased in abundance from 2,567 in 1994 to 2,048 in 2014.

The patterns of temporal turnover in species composition across size classes (Fig. 2 and Fig. S2) were similar to the overall pattern (Fig. 1 and Fig. S1). Our results particularly showed a strong species turnover pattern in small trees through time. This finding highlights that species composition in big trees was more stable than in small trees over time. The effect of surface fire may contribute to these patterns. Big trees are more effective in preventing the heat of fires from entering their living tissues than small trees (Gutsell & Johnson, 1996; Brando et al., 2012). While fires kill a huge number of small trees (Cochrane & Schulze, 1999; Baker, Bunyavejchewin & Robinson, 2008; Brando et al., 2012), new recruits of diverse species enter the post-fire community. According to our results, species turnover is lower in larger trees (DBH > 5 cm), perhaps due to low mortality rates and limited recruitment in large size classes. Even though a previous study in this forest has reported that fires could effectively kill large trees of some common species such as C. roxburghii, B. ramiflora, and P. viridis (Baker, Bunyavejchewin & Robinson, 2008), the mortality of these large trees was relatively low. For example, the 1998 fire burned through the study plot, and the mortality of the larger C. roxburghii trees was 34 times lower than that of the smaller trees in the 1999 census. We also found that the largest-sized trees had the lowest compositional turnover between 1994 and 1999 (Fig. 2). The result may reveal the relatively minor effects of 1992 and 1998 fires on species composition in big trees. Those surviving big trees that persist in 1994 and 1999 might contribute to the composition similarities in 1994 and 1999.

Effects of fire on turnover in functional composition

Our results showed that functional turnover in many quadrats was not different from the expectations of stochastic processes (Fig. 3 and Fig. S3). The results indicate that the temporal shifts in functional composition in the whole plot are predominantly random. In other words, the changes in tree species composition in this forest are mostly random through time in terms of functional trait combinations. However, the functional turnover was significantly lower than null expectations in some quadrats (Figs. 3 and 4). This implies that the observed temporal turnover patterns in these quadrats were a result of deterministic processes, such as environmental filtering, that favor species with optimal functional traits in that environment (Kraft, Valencia & Ackerly, 2008; Cornwell & Ackerly, 2009; Jung et al., 2010). In our case, high canopies of dipterocarps in the flat areas around the middle of the plot create a low-light condition that is suitable for shade-tolerant species. These species usually have dense wood, providing these trees with resistance from physical damages under the shade (van Gelder, Poorter & Sterck, 2006). After these shade-tolerant trees die, new recruits under the same dipterocarp canopies will be constrained by the light condition, resulting in another set of shade-tolerant species with dense woods. Consequently, functional composition in this area will remain similar through time.

Furthermore, we found functional composition resulting from species with highly functional similarities was more stable in 1994–1999. After the 1992 and 1998 fire disturbances, the recolonizing species possessed similar functional traits to those before the fire. Consequently, the functional turnover in 1994–1999 was lower than in other years (Fig. 3). Another reason could be that, based on species turnover results, we found species compositions were similar between 1994 and 1999. Species turnover during this time was possibly mainly determined by trees that survived both the 1992 and 1998 fires.

As species composition varied through time, the functional composition remained unchanged. One possible reason for this apparent decoupling between species and functional turnover is that species that died during this study period were replaced by species with similar traits through some deterministic processes. Our findings are consistent with previous research in the 50-ha Barro Colorado Island Forest Dynamics Plot (BCI FDP) in Panama. Functional composition in BCI FDP was relatively stable through time, despite there being shifts in species composition (Swenson et al., 2012).

We further observed functional trait distributions of each trait over time in the community, which could help us better explain the temporal trends in functional composition. Only a few functional traits shifted through time. WD slightly decreased, while LPC increased with the time interval (Fig. 7). Accordingly, species that grow fast commonly have high values of LPC, LNC (Chai et al., 2016), SLA (Poorter et al., 2008; Prado Júnior et al., 2015), and low WD (Enquist et al., 1999; Muller-Landau, 2004). If fires kill many trees and create suitable habitats for recruits of fast-growing species, one would expect functional composition after fire disturbances to be contributed mainly by fast-growing species. Consistent with this expectation, we found that the dominant and fast-growing species C. roxburghii, which has low WD and high LPC, greatly increased in abundance over time. This implies that a prominent change in the abundance of this species may have caused the increase in LPC and decrease in WD. Other functional traits varied very little through time, suggesting that this forest is generally composed of species with similar functional traits over time. Thus, the compositional changes in species in this forest are mainly caused by shifts in the abundance of a dominant fast-growing species and species with similar functional traits over time.

Temporal functional turnover varied with tree sizes and exhibited similar patterns in species turnover. We found that big trees (DBH > 10 cm) mostly showed significantly lower functional turnover than small trees (Fig. 5). The functional turnover in many quadrats of big trees was significantly lower than expected by null models, indicating that the shifts in species composition of big trees were due to replacement by species with similar functional traits. In other words, most large trees that survived fires through time determined the functional composition after burning. The temporal functional turnover of small trees appeared to be more random than big trees. The composition of small trees may be equally contributed by recruits of fast-growing species and surviving trees. The areas which contained most quadrats with random patterns in functional turnover were distributed throughout the plot. Remarkably, high temporal functional turnover was observed in some quadrats which differed according to areas and varied among diameter size classes (Fig. 6).

Higher temporal functional turnovers than expected by null models mainly occurred in many quadrats around the southwest corner of HKK FDP for trees of the smallest size. Functional turnovers in these areas are possibly determined by the increase in C. roxburghii abundance. On the other hand, quadrats with significantly high functional turnover in big trees were distributed mainly at high elevations (Fig. 6). In fact, these areas are relatively dry and have a high percentage of rock with low soil moisture compared to lower elevations. There are few big trees in comparison to low elevations, and only few species, such as Lagerstroemia tomentosa, live specifically in these areas (Bunyavejchewin et al., 2009). Big trees of L. tomentosa disappeared in some quadrats while other species such as C. roxburghii and Mallotus philippensis established in these areas through time. Consequently, temporal functional turnover in these areas was perhaps driven by species with different functional traits.

Although surface fires could kill some large trees for some species (Baker, Bunyavejchewin & Robinson, 2008), other factors, such as drought, may also be important. Drought not only induces fire occurrences (Cleary & Genner, 2004; Cochrane, 2009) but also increases the mortality risk of large trees (Bennett et al., 2015). During the study period (1991–2014), several ENSO events were recorded (1991–1992 (Harger, 1995), 1993–1994, 1997–1998 (Wooster, Perry & Zoumas, 2012)). Droughts associated with these ENSO events may also contribute to the death of the big trees in HKK.

Although temporal turnover in functional composition changed very little, functional turnover varied greatly throughout the entire 50-ha plot (Figs. 4 and 6). This variation across quadrats might be affected by the spatial heterogeneity of fire intensity (Trouvé, Bunyavejchewin & Baker, 2020). Factors such as light conditions, soil and topographic variables at local scales might also contribute to the spatial variation in the functional turnover patterns. For instance, the soil moisture in the quadrats at low elevations was relatively high, influencing the fire occurrence and intensity (Trouvé, Bunyavejchewin & Baker, 2020), thus contributing to the temporal forest turnover in HKK.