Compositional response of Amazon forests to climate change

Abstract Most of the planet's diversity is concentrated in the tropics, which includes many regions undergoing rapid climate change. Yet, while climate‐induced biodiversity changes are widely documented elsewhere, few studies have addressed this issue for lowland tropical ecosystems. Here we investigate whether the floristic and functional composition of intact lowland Amazonian forests have been changing by evaluating records from 106 long‐term inventory plots spanning 30 years. We analyse three traits that have been hypothesized to respond to different environmental drivers (increase in moisture stress and atmospheric CO 2 concentrations): maximum tree size, biogeographic water‐deficit affiliation and wood density. Tree communities have become increasingly dominated by large‐statured taxa, but to date there has been no detectable change in mean wood density or water deficit affiliation at the community level, despite most forest plots having experienced an intensification of the dry season. However, among newly recruited trees, dry‐affiliated genera have become more abundant, while the mortality of wet‐affiliated genera has increased in those plots where the dry season has intensified most. Thus, a slow shift to a more dry‐affiliated Amazonia is underway, with changes in compositional dynamics (recruits and mortality) consistent with climate‐change drivers, but yet to significantly impact whole‐community composition. The Amazon observational record suggests that the increase in atmospheric CO 2 is driving a shift within tree communities to large‐statured species and that climate changes to date will impact forest composition, but long generation times of tropical trees mean that biodiversity change is lagging behind climate change.


Figure S2 -Relationship between genus-level traits and dominance.
In red on the lower panels graphs show relationship including 679 genera from Extended Amazonia dataset. Top graphs, in blue, show the relationships for the 57 genera with more than 500 individuals across the plots analysed here.
Abundance and basal area are presented in the natural logarithm scale. Water deficit affiliation values were calculated in Esquivel-Muelbert et al. (2017), wood density was obtained from the Wood Density Database (Chave et al., 2009, Zanne et al., 2010, and potential size from Fauset et al. (2015).
Basal Area (log scale)

Appendix S3 Trends in climatological water deficit across Amazonian plots
Maximum cumulative water deficit (MCWD) is a metric calculated from precipitation and evapotranspiration data (Aragao et al., 2007, Chave et al., 2014. There are two major sources of information of precipitation for the Amazon Basin: ground-based data from weather stations and satellite-based data from the Tropical Rainfall Measuring Mission (TRMM -Huffman et al., 2007).
Data from weather stations are interpolated at 0.5 o resolution by the climatological research unit (CRU - Harris et al., 2014) are available from 1901 to 2014. Ground-based data are available for a longer period of time, however the number of weather stations within the Amazon can compromise the quality of these data. Thereby, we verify the plot-level yearly values of MCWD calculated using only CRU data against MCWD using precipitation data from TRMM.
For evapotranspiration we used CRU data, which is calculated using the Penman-Monteith equation (Allen et al., 1994). To understand the whole of changes in precipitation vs. temporal variation in evapotranspiration on the observed changes in MCWD we calculate MCWD trends maintaining evapotranspiration constant at 100 mm. For these examinations of the components of MCWD we used the plots from the Extended Amazonia dataset.
Both datasets show differences on how trends in MCWD are distributed across space, especially for South East Amazon ( Figure S3.2), which is not included in the main analyses, i.e. Core Amazonia.
Overall CRU-based MCWD showed similar results from TRMM-based values and the major difference in the direction of trends for some plots seems to be related to the time window analysed: when the same time window is analysed (1998-2013) both datasets show very similar trends for the basin: TRMM -1.2 (95% CI -2; -1); CRU -0.8 (-1.1; -0.5). Thereby, we used CRU-data in the analyses, as it spans over the same time window as the vegetation dataset.
Values of trends in MCWD for constant evapotranspiration (-0.7 mm y -1 , 95% CI = -0.9, -0.5; results for Core Amazonia) did differ marginally significantly from trends when allowing  Plot area and monitoring period are expected to affect the plot-level trends as forest stands monitored over shorter periods or smaller area are more likely to be affected by stochastic phenomena, such as tree fall. Thus, monitoring plots with greater sample effort should represent better the trends in functional and floristic composition. When analysing a combination of inventory plots with different sample sizes and monitoring periods, smaller plots and those monitored over a shorter length of time should bring more variation to the overall trend increasing the error to estimate ratio. In these cases, greater weights are given to plots with greater sample effort (Brienen et al., 2015.
In order to decide which are the best weights to be implemented in our analyses we evaluate the effect of monitoring period and plot area on the deviation from the mean estimate for each model, i.e. the absolute difference between the change in functional composition in each plot and the Pan-Amazonian mean slope for each model. The effect of sample effort was assessed through the slope of the relationship between deviation from mean estimate and area or monitoring period (Figure S4.1a and b) and compared to the effect when weights are applied ( Figure S4.1 c-f). To verify whether the weighting procedure is appropriate we tested for the relationship between weights and residuals vs.
When including the weights in the models their fit improved overall for our main analyses, the bootstrap mean of slopes ( Figure S4.1). The relationship between residuals and area or monitoring period flatten when using the squared root of these parameters as weights for each plot. Following this exploratory analysis the squared root of area * squared root of monitoring period seems to be an adequate weighting procedure for this type of analyses ( Figure S4.2). The results differ for the LMM analyses where weights seem not to be necessary ( Figure S4.3). Even light weights, such as the cubic root of area and monitoring period seem to overweight these parameters ( Figure S4.3 b and c). This indicates that the random effect including in the model seem to be enough to account for the variation among plots in terms of area and monitoring period. Thus, no weights were applied for LMM analyses.

-Mean linear slopes in individual-based functional composition across the 165 plots
in the Extended Amazon dataset. For each trait, the bootstrap mean annual changes in community weighted mean (CWM) weighted by the squared root of plot size x monitoring period. In brackets: 95% confidence intervals. CWM was calculated for: water deficit affiliation (WDA), potential size (PS) and wood density (WD). The analyses were repeated for recruits, dead trees and the difference between recruits and dead trees (net fluxes). In bold significant trends, i.e. where CIs do not overlap zero.

Community
Potential size (cm y -1 ) Water Deficit affiliation (mm y -1 ) Wood Density (g cm -3 y -1 ) Whole community -0.001 (  Intercept, slope, and percentage change per year (in brackets), of plot-level water deficit affiliation, potential size and wood density between 1985 and 2015. Trends were calculated by fitting a linear mixed model (LMM) to a time-series of census level community-weighted mean for each trait. The LMM consider plot slope and intercept as random effects. The analyses were repeated for recruits, losses and the difference between recruits and losses (net fluxes). Values in bold and followed by + represent slopes that significantly differ from zero, considering respectively α = 0.05 and α = 0.1.    Table S7.1 -Mantel correlation between the trends in functional composition and the spatial distance (latitude and longitude) between our plots and groups of plots (cluster, three letter codes in figure S1). Values in bold represent Mantel correlations that significantly differ from zero, considering respectively α = 0.05 rejecting the null hypothesis of absence of spatial structure in our data. Note that the correlation between the spatial distance between the plots and the similarity between the trends is often low and is not significant when analysed are performed at the cluster level.  Table S7.

-Mean linear slopes in individual-based functional composition across the 38 clusters
in the Core Amazon subset. For each trait, the bootstrap mean annual changes in community weighted mean (CWM) weighted by the squared root of cluster area x mean monitoring period across plots of a cluster. In brackets: 95% confidence intervals. CWM was calculated for: water deficit affiliation (WDA), potential size (PS) and wood density (WD). The analyses were repeated for recruits, dead trees and the difference between recruits and dead trees (net fluxes). In bold significant trends, i.e. where CIs do not overlap zero. Note that the do not differ from the result at the plot-level analyses (Table 1 and Table 1 and 2 in the main text.

Appendix S10 -Trends of Amazonian plant functional types
We investigate the trend in abundance for the four Amazonian functional types defined by Fyllas et al. (2012). This approach is an independent test for our hypotheses as it integrates a series of foliar traits not used in our core analyses. Fyllas et al. (2012) classified 276 Amazonian tree species within four functional types based on foliar and structural data: small pioneer, small statured nonpioneer, tall pioneer and tall non-pioneers, see Fyllas et al. (2012) for details on each functional group and how they are defined. We applied this classification to the 243 species described by Fyllas et al.
(2012) that were also found in our data and test for the trend in abundance within each of these functional types. Then we calculated the trend in abundance of each functional type following the same procedure described in the main text for individual taxa. Trends were calculated by fitting a linear mixed model (LMM) to a time-series of census-level community-weighted mean for each trait. The models consider plot intercept and slope as a random effects. In bold slopes that significantly differ from zero, considering α = 0.05.

Appendix S11 -Trends for domesticated taxa
It has been hypothesized that large parts of Amazonia are influenced by 'legacy' effects of indigenous forest management, with marked impacts still visible in contemporary forest composition including the relative prevalence of domesticated species . If so, then in the absence of such management we might expect some ecological relaxation, with community change toward reduced dominance of those species which were favoured by indigenous forest management. Thus, if most plots are still recovering from previous humanuse, we would expected domesticated taxa to be declining and for any declines to be greater than those of non-domesticated taxa. To test for this prediction, we investigated whether there has been a change in the proportion of domesticated species sensu Levis et al. (2017) within the plots.
We use the same approach used in the main text when analysing changes in functional composition (see Trends in functional composition section for details) but here testing for temporal trends in the proportion of domesticated species within the inventory plots. Firstly, we calculated the bootstrapped mean and 95% CI of linear slopes of the percentage of domesticated species as a function of time across all plots. Secondly, we used linear mixed effect model (LMM) where the percentage of domesticated species were a function of time using function lmer from the R package lme4 (Bates et al., 2014). Finally, we analysed each domesticated species individually based on their trends in abundance (details in the section Trends in floristic composition section in the main text).
Overall, the proportion of individuals of domesticated species did not change over time.
This was consistent regardless the analytical method: the LMM showed 0.7% intercept, slope = 5 x 10 -4 % y -1 ; P-value = 0.35; bootstrap indicated a non-significant trend of 4x10 -4 % y -1 (CI: -5 x 10 -4 ; 1 x 10 -3 % y -1 ). When analysed individually we also failed to detect a decrease in domesticated taxa. Levis et al reported 85 woody species to be domesticated by pre-Columbian peoples, of which 63 are in our dataset. The results from the bootstrap analyses indicate 27 species increasing in abundance (only 5 significantly) (Table S11.1). The same analyses show 36 species decreasing in abundance, only 4 significantly (Table S11.2). These results were consistent with the outputs form the LMM analyses.  Levis et al. (2017) that increased significantly in abundance across the Amazon Basin. Columns represent: species identity, position in the ranking of dominance across Amazonia, number of plots in which the species was found within Core Amazonia, bootstrapped mean and 95% CI (in brackets) of trends in relative abundance, intercept, slope, and percentage change per year (in brackets) from linear mixed model (LMM) fitted to a time-series of census level relative abundance of each species. The LMM considers plot slope and intercept as random effects. Values in bold and followed by + represent slopes that significantly differ from zero, considering respectively α = 0.05 and α = 0.1. Sample size was too small to precisely generate LMM for Platonia insignis.