Positive feedbacks and alternative stable states in forest leaf types

Abstract

The emergence of alternative stable states in forest systems has significant implications for the functioning and structure of the terrestrial biosphere, yet empirical evidence remains scarce. Here, we combine global forest biodiversity observations and simulations to test for alternative stable states in the presence of evergreen and deciduous forest types. We reveal a bimodal distribution of forest leaf types across temperate regions of the Northern Hemisphere that cannot be explained by the environment alone, suggesting signatures of alternative forest states. Moreover, we empirically demonstrate the existence of positive feedbacks in tree growth, recruitment and mortality, with trees having 4-43% higher growth rates, 14-17% higher survival rates and 4-7 times higher recruitment rates when they are surrounded by trees of their own leaf type. Simulations show that the observed positive feedbacks are necessary and sufficient to generate alternative forest states, which also lead to dependency on history (hysteresis) during ecosystem transition from evergreen to deciduous forests and vice versa. We identify hotspots of bistable forest types in evergreen-deciduous ecotones, which are likely driven by soil-related positive feedbacks. These findings are integral to predicting the distribution of forest biomes, and aid to our understanding of biodiversity, carbon turnover, and terrestrial climate feedbacks.

Fig. 1. Bimodal patterns of forest types at the continental and global scale.

A The spatial distribution of 45,276 forest inventory sites from the FIA was used for the continental analysis. The continuous color scale represents the relative abundance of evergreen trees within a plot (red = 100% deciduous; blue = 100% evergreen). B Histogram of the observed plot-level evergreen percentage across the mainland US. The black dots and error bars show the medians and 2.5–97.5% quantiles of the null model predictions driven by environmental filtering (zero adjusted Poisson distribution). C Spearman’s rank correlation coefficient between evergreen abundance and deciduous abundance in the observed data (red bar) versus the simulated results of the null model (black histogram). D The location of 815,578 forest plots from the global GFBI database. E and F are the same as panels B, and C, but for the global data. Hartigan’s dip test showed significant multi-modality (here is bimodality) in the observed values in panels B (n = 45,276, one-sided p-value < 0.001) and E (n = 815,578, one-sided p-value < 0.001). Source data are provided as a Source Data file.

Fig. 2. Observed positive feedback in con-phenological demographics in the mainland US.

A Survival probability of an individual deciduous tree (DE, red) or evergreen tree (EV, blue) within a purely evergreen or deciduous forest stand. B Individual deciduous or evergreen tree growth (stem diameter increment in cm per year) when the surrounding trees are purely evergreen or deciduous. C Recruitment rates of deciduous or evergreen trees in deciduous or evergreen dominated forest plots. All plotted data are 100 samples drawn from the 95% CI of the corresponding full model (n = 45,276), controlling for environmental conditions and stand structure. Results are presented as boxplots (medians as centre with 25th and 75th percentile as bounds). The differences between all compared pairs are highly significant (one-sided t-test p-value < 0.001, n = 100). Source data are provided as a Source Data file.

Fig. 3. Three pairs of feedback simulations versus null simulations.

Histograms in A–D represent the relative evergreen abundance within 1000 simulated forest plots after 2000 years. The color scale represents the percentage of evergreen forest within a plot (red, 100% deciduous; blue, 100% evergreen). A, Outcome of null demographic simulations with uniform initialization, from growth, recruitment, and mortality models fit without con-phenological feedback predictors. The inset in (A) shows the initial uniform distribution of relative evergreen abundance across the 1000 plots. The uniform distribution shifts to evergreen-dominated forests after 2000 years. B, Outcome of demographic simulations with uniform initialization from demographic models fit with con-phenological predictors. Most plots are dominated by either of the two leaf phenology strategies (bimodal distribution). C Outcome of null demographic simulations with bimodal initialization, from growth, recruitment, and mortality models fit without con-phenological feedback predictors. The inset in panel C shows the initial bimodal distribution of forest composition. As in panel A, the bimodal distribution shifts to evergreen-dominated forests after 2000 years. D Outcome of demographic simulations with bimodal initialization from demographic models fit with con-phenological predictors. Most plots are dominated by either of the two leaf phenology strategies. E, F Hysteresis simulations along mean annual temperature gradients, where 80% of forest plots (800 of 1000) were initially dominated by evergreen trees (EV, blue line) or deciduous trees (DE, red line). Demographic simulations were run either without feedback predictors (E) or with feedback predictors (F). Shaded regions represent the 95% CI of the mean relative abundance of evergreen trees across each set of simulations (n = 1000). Source data are provided as a Source Data file.

Fig. 4. Hemisphere-wide patterns of forest types.

A Empirical map of bimodality in forest types generated from GFBi data. Colors reflect the value of the bimodality index (BI), with red representing BIs < −0.22 (deciduous-dominated forest clusters), blue representing BIs > 0.22 (evergreen-dominated forest clusters), yellow representing BIs from −0.22 – 0 (bistable-deciduous forest clusters) and cyan representing BIs from 0–0.22 (bistable-evergreen forest clusters). B Projected map of bimodality in forest types across the Northern Hemisphere based on random forest modelling. Colors reflect the projected value of the bimodality index (BI), and share the same scale as in A. Predictions in B were made for forest regions (1) above 15 degrees northern latitude, where > 98% of the GFBi data are located, (2) whose environmental conditions well represented by our training data (> 90% interpolation, see Fig. S9A). C Relative frequencies of plot-level relative evergreen abundance in deciduous-dominated, bistable-deciduous, bistable-evergreen and evergreen-dominated regions. Source data are provided as a Source Data file.

Fig. 5. Determinant analysis showing the variable importance based on random forest models.

Seven key covariates of forest leaf phenology composition were used to train the random forest models, whereby we ran each model 100 times on 100 bootstrapping training sets and then computed the mean and standard deviation of the variable permutation importance. A, Variable permutation importance (mean ± std, n = 100, individual data points are overlayed on the bar charts) for global random forest analysis using the bimodality index as response variable. Variables along y-axis in A are ordered by their mean importance. The continuous color scale represents the variable importance from high (yellow) to low (dark blue). B, Variable permutation importance (mean ± std, n = 100) for random forest analysis of plots within deciduous-dominated forests (left panel), bimodal forests (middle panel, including both bistable-deciduous and bistable evergreen forests) and evergreen-dominated forests (right panel). Analyses in panel B all use plot-level relative evergreen abundance as response variable. Source data are provided as a Source Data file.