Introduction

Drivers of species distributions and their predictions have been a long-standing search in ecology, with approaches varying from deterministic to neutral (i.e. stochastic) and almost everything in between (e.g. near-neutral, continuum or emergent-neutral1,2). Most models are based on prior assumptions of processes that drive community dynamics. The Maximum Entropy Formalism (hereafter called MEF) makes no such, potentially unjustified, a-priori assumptions in generating predictions of species abundance distributions, as such it is a useful construct to infer processes driving community dynamics given the constraints imposed by prior knowledge (e.g. functional traits or summed regional abundances)3. Quantifying the relative importance of these distinct constraints can thus provide additional answers to understand the complexity of community dynamics (see Supporting Materials SM: boxes S1–S3). This is especially so because, although many different tests are available that link variation in taxon abundances to (1) trait variation, (2) taxon turnover between habitats or environments and (3) the distance decay of similarities between samples, none quantify the importance of these relative to each other. The MEF as applied here, however, is capable of and designed to do exactly this by decomposing variation to separate information explained by each of these aspects in a four-step model (Fig. 1 and Box S2). Its application to an unprecedented large tree inventory database on genus level taxonomy consisting of > 2,000 1-ha plots distributed over Amazonia4 and a genus trait database of 13 key functional traits representing global axes of plant strategies5 allows us to advance the study of Amazonian tree community dynamics from a new cross-disciplinary perspective.

Figure 1
figure 1

Schematic depiction of the MEF procedure. Left panel shows a genus abundances per site and a functional trait matrix per genus, bottom half outlines calculations. Middle and right panel show different scenarios of neutral and deterministic dynamics under infinite or limited migration. Figure was custom made using Adobe Illustrator (Adobe Inc., 2019. Adobe Illustrator).

Full size image

Results

Principles from information theory6,7,8 can be used in an ecological setting to predict the most likely abundance state for each taxon while simultaneously maximizing entropy based on constraints. Maximization of entropy allows quantifying the information yield for each constraint and therefor identifies which constraints reduce entropy the most. Here we specifically use Shipley’s mathematical framework (CATS) for the MEF calculations, similar to earlier studies9,10,11.

Predictive power of the four-step model

Using a uniform prior and CWM values (Community Weighted Means) as constraints accounted for 23% on average of total deviance between observed and predicted relative abundances (measured by R2KL values, see Box S2 Eq. 5). Filtered by forest type this was 34% for podzol forests, várzea 25%, igapó 23%, swamp forests 34%, 21% and 24% for Guyana Shield and Pebas terra firme respectively and 20% for Brazilian Shield terra firme forests (see Table S1 for detailed decomposition). Using observed metacommunity relative abundances as prior regardless of functional traits accounted for on average 56% for the combined dataset with for all forest types between 51 and 55%, except for the Guyana Shield terra firme with 62%. The hybrid model (including both traits as constraints and the metacommunity prior) performed slightly better for the combined dataset (average 60%) with a minimum of 57% for swamp and várzea forests and a maximum of 66% for Guyana Shield terra firme forests. To compensate for spurious relationships between regional abundances and local trait constraints, regardless of selection, explanatory power was regarded relative to model bias yielding the pure trait and metacommunity effects (Box S3, Fig. 2 and Table S1). This lowered the proportion of information accounted for and yielded average pure metacommunity effects of 40% for the overall dataset ranging between 26 and 45% for each forest type separately with pure trait effects explaining only 5% of information for the combined dataset on average with for each forest type between 3 and 9%. Although the latter was lowered substantially, the explanatory power did appear to be strongly dependent on forest type. The online supplementary material provides additional results relating to the predictive power of each model as well as the spatial gradient of the pure trait and metacommunity effects (Figs. S2–S3).

Figure 2
figure 2

Visual representation of pure trait, pure metacommunity, hybrid model and the remaining unexplained information for each separate forest type. Abbreviations indicate different types: igapó (IG), podzol (PZ), swamp (SW), Brazilian shield terra firme (TFBS), Guiana Shield terra firme (TFGS), Pebas terra firme (TFPB) and várzea (VA). Boxplots show median value of pure effects over all samples, with lower and upper hinges corresponding to 25th and 75th percentiles. Whiskers extends from hinge to largest or smallest value no further than 1.5 * IQR from hinge. Points beyond this range are plotted individually and only positive values were plotted.

Full size image

Direction and strength of selection of trait-based constraints

Each trait showed significant differences in lambda when compared between forest types (Fig. S1, see methods for a definition of lambda). Scatterplots of CWM trait values versus lambda show that, in general, higher lambda values correspond with higher CWM trait values (Fig. S7). A number of functional traits associated with low nutrient conditions (e.g. ectomycorrhiza) and life history strategies suited for protection against herbivores (e.g. latex, resin and high leaf C content) were clearly positively associated with abundance in nutrient poor environments (podzols), indicated by the positive lambda values. In contrast, having fleshy fruits and high leaf N and P content were clearly negatively associated with abundance on these soils. Nodulation was also negatively associated with abundance on poor soils. The ability to accumulate aluminium was positively associated with abundance on those soils commonly associated with higher aluminium content such as igapó (strong positive effects). In contrast, it was strongly negatively associated with abundance on other soils, with negative lambda values for podzol and várzea forests. Traits such as SLA or winged fruits also showed patterns dependent on forest type, although less pronounced.

Effect of regional metacommunity prior

There was a remarkable similar mean 22% decrease of the information explained purely by the metacommunity prior for each forest type (Fig. S3). For the separate forest types, although the initial pure metacommunity effect varied, the decline appeared remarkably similar with a mean 25% decrease in pure metacommunity effect for podzol, 23% for várzea, 26% for igapó and 27% for swamp forests with terra firme forests having a smaller decline of approximately 21%, averaged over the three subregions. It should be noted there is an obvious risk that when sampling size is increased, this also includes more environmental heterogeneity as samples are coming from a variety of localities potentially leading to changing composition. If this were the case, however, the regional prior (qi from Fig. 1 and Box S2) would also change, as taxa might be abundant in some places but rare or absent in others. As the metacommunity effect is the explained information that remains relative to any trait effects (i.e. information unique to the neutral prior) and the pure trait effects are the explained information remaining after correcting for pure metacommunity effects (Box S3) this effect should then be accompanied by an increase in pure trait effect for each sample. This was not observed, not even within the different forest types. Instead, the trait effect gradually went up and then remained constant (Fig. S4).

Discussion

The MEF emerges from a well-founded theoretical and empirical body of ecology and evolutionary biology, regarding natural selection, migration and population dynamics. From an ecological point of view, it can be used to quantify the relative association between directional or stabilizing selection for functional traits versus the importance of relative regional abundance regardless of these traits by imposing these as constraints. Our results show that pure trait effects, on average, explained only 5% of the information when all forest types were taken together whereas the pure metacommunity effect, however, explained eight times more taken all forest types together (40%). Greater trait dissimilarity was positively associated with higher pure trait effects, indicating trait-based selection, although the assumed influence of dispersal regardless of these traits appeared to confer more information explaining tree genus composition of the Amazon rainforest. The strength and direction of selection indicated clear directional selective pressure for life history strategies of either growth or protection, depending on forest type (see supplementary online material S-A for a more detailed exploration of ecological interpretation). Including community weighted variance as reflective of potential stabilizing selection did not provide additional information. Although this could be interpreted as indicative of weak or absent stabilizing selection, it is more likely to be an artefact of many genera not being shared among localities due to the sheer geographical scale resulting in a strong mismatch between observed, predicted and uniform relative abundances resulting in a model bias higher than information yielded by including these constraints (see also box S2).

Despite showing clear patterns in environmental selection and dispersal effects, there was a large proportion of information left unexplained (44% on average). Potentially, local demographic stochasticity could weaken any link between functional traits measured and regional abundances of genera. This would, however, mean that almost half of the information contained in relative abundances are the result of random population dynamics and are not structurally governed. Alternatively, this could be due to functional traits reflective of processes not taken into account in this study, such as traits reflective of interactions between trophic levels (e.g. traits influencing specific plant-pathogen interactions). Another and at least equally likely hypothesis for (local) unexplained information is that when scaling up, the ratio of genus richness to total abundance decreases rapidly initially but stabilizes again as relatively non-overlapping habitats are included in the regional abundance distributions and more genera are included again due to the different habitats. This would result in a change of the regional abundance distribution (i.e. the prior) to which each local community is compared, resulting in higher local unexplained information. Further study into these aspects could provide additional insight, though the data necessary for these scales is lacking for Amazonian trees.

Metacommunity importance

Although the initial explanatory power of the metacommunity prior differed between forest types, the decay pattern was very similar. As the effects of either traits or the metacommunity are measured in the goodness-of-fit predictions on local relative abundances, this implies that at small spatial scales the surrounding regional abundances provide better estimators than functional traits, while at larger spatial scales this shifts to the traits. The ecological translation would be that on small spatial scales, local communities share similar environmental conditions leaving dispersal and drift acting predominantly in changing community composition, at least for genus level taxonomy. As the potential regional pool is increased, more and more environmental heterogeneity and non-overlapping regions are likely to be introduced. The more gradual decline of terra firme forests can then arguably be attributed to these forests having the largest relative surface area of Amazonia (even for the separate subregions), potentially giving these forests an almost continuous metacommunity without gaps, resulting in a more gradual transition from metacommunity to trait relative importance. The fact that metacommunity effects do not change anymore after certain distances would indicate the effect of dispersal potentially occurs over very large distances. It should be noted that as these calculations are done at community and genus level, they do not measure single dispersal events but rather the effect of dispersal on community composition much deeper in time. In other words, this effect suggests more than a dispersal event every now and then. Instead, it argues for prolonged mixing of forests on large geographical and temporal scales, supported by recent findings demonstrating a lack of geographical phylogenetic structure of lineages for Amazonian tree genera12.

Conclusion

Using an unprecedented scale of data and applying the Maximum Entropy Formalism from information theory we show that constraints formed by regional relative abundances of genera explain eight times more of local relative abundances then constraints based on directional selection for specific functional traits, although the latter does show clear signals of environmental dependency. There is, however, still much to be explored due to the large unexplained effects and analyses on finer taxonomic (i.e. species level) and environmental (e.g. microhabitat) scales could resolve these issues. The relatively large effects of the regional pool of genera over great distances does suggest an important role for long term dispersal and mixing of Amazonian trees, especially for the Amazonian interior.

Methods

Empirical data

The ATDN4,13,14 consists of over 2000 tree inventory plots distributed over the Amazon basin and the Guiana Shield, collectively referred to as Amazonia (a map of all current plots can be found at https://atdn.myspecies.info/). Only those plots with trees ≥ 10 cm diameter at breast height were used, leaving 2011 plots with a mean of 558 individuals per plot identified to at least genus level. Most plots used are 1 ha in size (1414) with 492 being smaller (minimum size of 0.1 ha) and 105 larger (maximum size of 80 ha). Genera have been standardized to the W3 Tropicos database15 using the Taxonomic Name Resolution Service (TNRS, see16). After filtering based on above criteria and solving nomenclature issues, 1,121,935 individuals belonging to over 828 genera remained. Plots were distributed over seven abiotically different forest types: Podzol forests (PZ), Igapó (IG, black water flood forests), Várzea (VA white water flood forests), Swamp (SW) and Terra firme forests (TF) with subregions BS (Brazilian Shield), GS (Guyana Shield) and PB (Pebas) (see also17 for details regarding these forest types). Trait data were extracted from several sources. Wood density was mostly derived from18. Traits related to leaf characteristics mostly came from four large datasets19,20,21,22,23,24, including additional data from other sources25,26,27 as well as unpublished data (J. Lloyd, A.A. de Oliveira, L. Poorter, M. van de Sande & Mazzei, M. van de Sande & L. Poorter). Data on seed mass came from28,29,30 as well as different flora’s and tree guides. As this particular trait can vary over several orders of magnitude, this was included on a log-scale29,31. Ectomycorrhizal aspects were derived from literature32, the same applies to nodulation33,34. Traits involved in aluminum accumulation were based on35,36 and references therein. For binary traits (yes/no), a genus was considered having a certain trait only when > 50% of the genus was positive for that specific trait.

Functional traits and trait imputation

Constraints were formed by Community Weighted genus Means (CWM) of functional traits (Table 1), related to key ecological life history aspects. According to principles of natural selection, CWM values will likely be biased towards favourable trait values for that particular environment in the case of directional selection, as taxa with these traits will be more abundant due to environmental selection. Previous studies included community weighted variance (CWV) as well as indicative of potential stabilizing selection11,37. In our case, however, including CWV as constraints resulted in a model bias that was consistently higher than information including trait or metacommunity aspects, CWV was therefore not included as constraints in the final analysis. As for many traits it has been shown earlier that the interspecific variability was larger than the intraspecific variability, this allowed the use of data from different sources to at least calculate a mean species trait value. Genus trait values were subsequently computed as genus-level means of species values if known within the genus and considered constant for each genus. Genus level of taxonomy was used as the available trait database had the most information on this taxonomic level (see Table 1). Unknown values for traits were estimated by Multiple Imputation with Chained Equations (MICE, see15) by delta adjustment, subtracting a fixed amount (delta), with sensitivity of this adjustment to the imputations of the observed versus imputed data analysed using density plots (Fig. S8) and a linear regression model. This procedure was done using the mice package38, available on the R repository, under predictive mean matching (pmm setting, 50 iterations). Results showed imputations were stable and showed near identical patterns with each imputation scenario (see Figs. S5–S6 and Table S2). After imputation, all trait values were transformed to Community Weighted Means (CWM) of each trait (J) for each plot (K) (\({\overline{T} }_{JK}\)) as \({\overline{T} }_{JK}= {\sum }_{i=1}^{S}{t}_{ij}{ra}_{ik}\) with ra the relative abundance of the ith genus in the kth plot following earlier uses37.

Table 1 Overview of used functional traits.
Full size table

MEF procedure predictions and ecological inference

Figure 1 provides a schematic procedure overview, box S1 provides an overview of important terms and Boxes S2–S3 further mathematical details. Initially, a maximally uninformative prior is specified, where qi (Box S1 Eq. 1) equals 1/S, indicative of each species having equal abundances, and trait constraints are randomly permuted multiple times (n = 50) among genera to test whether inclusion of specified constraints significantly changes derived probability distributions (see also39). Subsequently, the same prior is used but now observed trait CWM values belonging to specific genera are used as constraints (following earlier applications using simulated communities11). Third, observed regional abundances are used as prior with permutated trait constraints and finally both observed regional abundances and observed trait CWM are used as prior and constraints. Maxent2, an updated version of the maxent function currently in the FD library of R provided the computational platform. Proportions of uncertainty explained by each model are given by the Kullback–Leibler divergence R2KL, a generalization of the classic R2 goodness of fit. In contrast with standard linear regression models having squared goodness-of-fits measurements, the R2KL is much more related to the concept of relative entropy, quantifying the information lost when one distribution is compared to another by means of quantifying the statistical distance between two distributions40. Pure trait, pure metacommunity, joint metacommunity-trait and unexplained effects are calculated as proportions of total biologically relevant information (Box S1 and Box S2). Data was rarefied to smallest sample size (swamp forests; 28) and calculations bootstrapped 25 times. Results indicated no significant change compared to using all data, hence the total dataset was used for all analyses.

Strength and direction of selection

Predictions of genus relative abundances are computed as a function of traits reflected in the CWM values and a series of constants (λjk: the Lagrange Multipliers). Each multiplier quantifies the association between a unit of change for a particular trait j and a proportional change in predicted relative abundance pik (the ith genus in the kth community) considering all other traits are constant, formally described as: \({\raise0.7ex\hbox{${ \partial p_{ik} }$} \!\mathord{\left/ {\vphantom {{ \partial p_{ik} } {\partial t_{ij} }}}\right.\kern-0pt} \!\lower0.7ex\hbox{${\partial t_{ij} }$}} = \lambda_{jk} p_{ik} \left( {1 - p_{ik} } \right)\) (see appendix 1 from41). Positive values indicate larger trait values associated with higher abundances (positive selection), negative values indicate the opposite (negative selection) with changes proportional to lambda. Values approximating zero indicate no association between specific traits and relative abundances of species. Decomposing λjk and comparing by means of a One-Way Analysis of Variance for each trait separately between forest types allows studying both the strength and direction of selection in different habitats. Note that this is done for the same constraint between forest types, as lambda values for each constraint do not scale linearly between different constraints.

Estimation of metacommunity size

Iteratively increasing the regional species pool considered as prior in concentric circles of a fixed radius of 50 km allows estimating the spatial effect of metacommunity size. Due to computational limits, the number of permutations for the MEF calculations (see above) was reduced to two, shuffling the combinations of genera and traits at least once. Comparison of results from the analyses using all plots indicated small effects of a smaller perturbations (average of 5% difference for metacommunity effect between 5 and 50 permutations). The relationship between pure metacommunity effect and radius of metacommunity size was fitted using a smoothing loess regression (function loess and predict; R-package stats with span set at 0.1). Fits subsequently were used to predict values of metacommunity effect based on geographical distance to visualize general patterns for each forest type. Exponential decay of pure metacommunity effect was described using a self-start asymptotic regression function (SSasymp) of the form y(t)∼yf + (y0 − yf)e−exp(log(α))t (nls from stats42. A list of all packages used in R in addition to those preloaded can be found in the supplementary online material (SA2).