Tropical tree growth sensitivity to climate is driven by species intrinsic growth rate and leaf traits

A better understanding of how climate affects growth in tree species is essential for improved predictions of forest dynamics under climate change. Long-term climate averages (mean climate) and short-term deviations from these averages (anomalies) both influence tree growth, but the rarity of long-term data integrating climatic gradients with tree censuses has so far limited our understanding of their respective role, especially in tropical systems. Here, we combined 49 years of growth data for 509 tree species across 23 tropical rainforest plots along a climatic gradient to examine how tree growth responds to both climate means and anomalies, and how species functional traits mediate these tree growth responses to climate. We showed that short-term, anomalous increases in atmospheric evaporative demand and solar radiation consistently reduced tree growth. Drier forests and fast-growing species were more sensitive to water stress anomalies. In addition, species traits related to water use and photosynthesis partly explained differences in growth sensitivity to both long-term and short-term climate variations. Our study demonstrates that both climate means and anomalies shape tree growth in tropical forests, and that species traits can be leveraged to understand these demographic responses to climate change, offering a promising way forward to forecast tropical forest dynamics under different climate trajectories.

in the long-term variation range around them, that is, an important element to detect anomaly 157 effects on tree growth across different mean climates. 158 To build the climate covariates for the tree growth models, the monthly 30-year mean and 159 standardised anomaly variables were averaged over the months between consecutive censuses (two 160 to five years). For MCWD, the maximum over the growth periods between two censuses was used 161 instead of the weighted mean. The eight resulting interannual averaged variables were used as 162 predictors to model tree growth (see Data analysis). Correlations among these variables, stand 163 structure and elevation are presented in Table S3a and the Supplementary Methods S1. 164 Stand structure 165 As stand structure can vary between plots, we include its effect on tree growth through total plot 166 basal area. Plot basal area (m²/ha) was calculated at each census, with expectations that increasing  (Table 1; Table S1 and S2 for plot and species details). Species were 172 chosen to sample those that made up 80% of the standing biomass. Trait data collection and 9 measurement are detailed in Supplementary Methods S1. We measured leaf, wood and maximum 174 size traits that relate to light, water and nutrient use (Table 1, see Table S3b for pairwise trait  175 correlations, and Fig. S1 for trait distribution along the elevation gradient). Traits were measured 176 on three individuals per species, and included photosynthesis and stomatal conductance at a 177 reference CO2 concentration of 400 µmol mol -1 and irradiance of 1500 µmol photons m -2 s -1 (Asat 178 and gsat), dark respiration (Rd) at the same CO2 concentration, the CO2-saturated photosynthesis 179 and stomatal conductance (Amax and gmax), measured at 1200 µmol mol -1 CO2. The one-point 180 method (De Kauwe et al. 2016) was used to estimate maximum carboxylation rate (Vcmax) for 181 each individual from net photosynthesis measured at 400 µmol mol -1 CO2, and maximum light-182 driven electron flux (Jmax) from net photosynthesis measured at 1200 µmol mol-1 CO2 183 (Bloomfield et al. 2018). We also measured leaf stable carbon isotope ratio (δ 13 C), nutrient 184 concentration, and leaf area, leaf mass per area (LMA), leaf thickness, wood density (from 185 branches, after bark removal). All traits were averaged at the species level for tree growth analyses. 186

187
We addressed our questions through three sets of Bayesian multilevel models (M1 to M3; details 188 in Supplementary Methods S1). 189 M1: Tree growth response to climate means and anomalies, and species differences in their 190 sensitivities to climate 191 In M1, we used 12,853 individuals from all 509 species to test the effects of climate on tree growth, 192 and to investigate tradeoffs among species between intrinsic growth rate and growth sensitivity to 193 climate covariates. We built a two-level hierarchical Bayesian model of AGR, where the hierarchy 194 included an upper level of response (hereafter grand coefficients or effects, affecting AGR across 195 species) above a lower, species-level response. The higher level modelled AGR responses to 196 covariates via hyperparameters (i.e. statistical distributions from which species-level intercepts and 197 slope coefficients arose), while the lower level captured species-specific growth sensitivities to 198 model covariates, and species-level intercepts (hereafter intrinsic AGR) captured unexplained 199 growth variation across individuals, growth periods, and plots. 200 More specifically, we modelled individual log(AGR) as a species-specific function of (i) initial 201 tree size (approximated by log(DBH) at the beginning of a growth period), (ii) the local 30-year 202 mean of a climate variable, (iii) the anomalies of the same climate variable averaged over the 203 studied growth period, and (iv) stand structure (approximated by plot basal area at the beginning 204 of a growth period), using varying slopes (also known as random slopes) and a covariance matrix 205 to estimate correlations among species-specific AGR sensitivities to the covariates, as:   Supplementary Methods S1). Parameter α1 is the species-level departure from the grand intercept 256 (α0) for an increase of one standard deviation in the log(trait Tj) value of species j (direct effect of 257 trait on AGR), while β2-4, 1 are the departures from the grand slope of the corresponding model 258 covariates for an increase of one standard deviation in the log(trait Tj) value of species j (trait 259 mediation of AGR response to climate and stand structure; see Supplementary Methods S1 for 260 ecological interpretations of trait coefficient signs). We did not include the role of species traits in 261 AGR response to tree size because some traits can change through tree ontogeny (Fortunel et al. 262 13 2020) and our trait data does not encompass species tree size ranges. M2 models were run 263 separately for each of the four climate variables and for each of the 15 functional traits to manage 264 model complexity (representing a total of 60 M2 models). 265 In both M1 and M2 models, we standardised the response variable log(AGR) and all covariates -266 but climate anomaliesto mean zero and unit standard deviation, to allow relative importance 267 comparisons between covariates through slope coefficients (Schielzeth 2010), and to ease plausible 268 weakly-informative prior assignment to the parameters (McElreath 2020) (see Supplementary 269 Methods S1). We did not standardise averaged monthly anomalies to maintain their interpretability 270 as deviations from long-term means in terms of plot-specific units of standard deviation (see eq. 2; 271 i.e. mean anomaly covariate slope coefficients are not directly comparable to other covariate mean 272 slopes). 273 M3: Plot-level tree growth response to climate anomalies and interaction with mean climate 274 M3 models evaluate plot-level growth response to climate anomalies, and whether it varies 275 depending on local mean climates (e.g. whether plot-level AGR sensitivity to VPD anomalies is 276 higher in drier sites). We focused on the tree growth at the plot level, and modelled the expected 277 log(AGR) as a linear function of mean climate and climate anomalies. We used a similar Bayesian 278 hierarchical model as described for M2, where plot-specific average AGR depended on climate 279 anomalies, whose effect on AGR itself depended on the plot mean climate, as: 280 where αk is the average growth rate in plot k, and β1k characterises the growth response of plot k to 288 standardised climate anomalies for time interval t. α0 is the mean intercept value (i.e. mean absolute 289 growth rate) across plots, and α1 is the departure from the grand mean for one unit increase in mean 290 climate (see detailed equations and priors in Supplementary Methods S1). β1,0 is the grand slope of 291 climate anomalies, and β1,1 is the departure from this grand mean for a one unit increase in mean 292 climate (mediation of the effect of anomalies on growth by the plot mean climate). Parameters γj, 293 δt, λi are varying intercepts for species, census periods, and individual stems, respectively. 294 We run M3 models only for two climate variables (VPD and SRAD), as we found they were the 295 most important climate variables for tree growth in M1 and M2 models (see Results). 296 Standardisation of variables was carried out as for M1. 297 Trends in climate over time 298 To explore the implications of the effects of climate anomalies on tree growth, we built a separate 299 set of hierarchical Bayesian models to test for linear temporal trends in mean annual climate 300 variables between 1971 and 2019. We used varying year slopes per plots to allow plot-specific 301 trends (model details in Supplementary Methods S1). We also run the models for the  Table S4). There was no general temporal trend for MCWD or SRAD 307 (Fig. 1c). 308

Analysis of model outcomes 309
All model parameter posteriors were summarised through their median and 95%-highest posterior 310 density interval (HPDI) (i.e. the narrowest posterior interval encompassing 95% of the probability 311 mass, corresponding to the coefficient values most consistent with the data; (McElreath 2020)). 312 Model covariates were considered important at two high levels of confidence, when their 313 coefficient had a posterior probability of over 95% or 90% of being either positive or negative 314 (HPDI not encompassing zero). 315 The goodness-of-fit of the models was assessed through the squared Pearson correlation between 316 the observed AGR and the AGR predicted by the fitted model (R 2 ). M1 and M2 models had high 317 explanatory power, with R² of 0.46 and 0.52 on average, respectively. M3 models, with VPD and 318 SRAD as climate variables, had an R² of 0.67 and 0.63, respectively. 319 Bayesian updating of parameters was performed via the No-U-Turn Sampler (NUTS) in Stan 320 (Carpenter et al. 2017), using three chains and 3000 steps (1500 warmings). All models mixed well 321 and converged (Rhat within < 0.01 of 1). Models were run in the R environment (Team 2020)  Results 326 Contribution of climate means and anomalies to tree growth

327
The main climate drivers affecting tree growth across species were the climate means and 328 anomalies in Tmean, SRAD and VPD (Fig. 2, Fig. S3 Fig. 2; Fig. S2; Table S5). Tree growth sensitivity to climate, stand 335 structure and tree size varied widely among species (illustration in Fig. S3). Similar results were 336 obtained from the M2 models (subset of 75 species with trait data) ( Fig. S5a-d, Table S5), though 337 we no longer detected the effects of temperature anomalies and VPD and solar radiation means in 338 this reduced dataset. 339 Coordinated tree growth sensitivities to climate means and anomalies 340 Using the fitted matrix of correlations among species-level intercepts and slopes from the M1 341 models (matrix R, see eq. 3.5) allowed testing whether fast-and slow-growing species, and species 342 growing better at opposite extremes of the range of mean climates showed different sensitivities to 343 climate anomalies. Fast-growing species (i.e. with high intrinsic AGR) were more sensitive than 344 slow-growing species to the negative effects of both VPD anomalies and plot basal area on tree 345 growth ( Fig. 3c Table S5). Drier tropical rainforests showed steeper decrease in 356 plot-level growth to positive VPD anomalies ( Fig. 4a; Table S5). Cloudier forests exhibited 357 stronger decrease in plot-level growth with positive SRAD anomalies ( Fig. 4b; Table S5). 358 Functional traits influence species intrinsic tree growth and their response to climate 359 drivers 360 Based on M2 models, species intrinsic growth increased with dark respiration rate (Rd), DBHmax, 361 leaf P content, Asat, Vcmax, leaf δ 13 C and LMA. (Fig. 5; Fig. S5e; details in Table S5). Species 362 traits also mediated the effects of climate and forest structure on tree growth, either by accentuating  Table S5). Leaf δ 13 C 365 and P content exacerbated the negative effects of positive anomalies in SRAD on tree growth, while 366 Amax, gmax, gsat and Jmax attenuated them ( Fig. 5; Fig. S6f, Table S5). The negative effects of 367 anomalies in VPD on tree growth were exacerbated in species with high leaf δ 13 C, DBHmax, leaf P, 368 and LMA, further confirming that VPD anomalies had the most negative effects on fast-growing 369 species (Fig. 3c), but also those with low gmax or leaf area ( Fig. 5; Fig. S6g). Tree growth was less 370 reduced by denser forest environments (high plot basal areas) in species with high wood density, 371 low Rd and low leaf δ 13 C ( to 2019 in our dataset, we detected a 3.8-fold stronger VPD increase rate across all plots (0.045 408 hPa / yr, 90%-HPDI: 0.019, 0.066; R² = 0.80; details in Table S4; e.g. Fig. 1b). This trend itself 409 was stronger than the 1971-2019 trend in our dataset (0.020 hPa / yr; R² = 0.84; Table S4), 410 indicating a sharper-than-previously-thought VPD increase in the past two decades. This rapid 411 increase of VPD anomalies through time combined with the generalised ensuing decrease in tree 412 growth and growth sensitivity variability to VPD among species ( Fig. S3; Table S5) suggests that 413 tropical forest composition and functions may be strongly altered by ongoing climate change, 414 especially by VPD. It is worth noting that soil water deficit also depends on evapotranspiration 415 estimates accuracy and variables unaccounted for, here (e.g. soil water retention capacity, 416 topography), so that the importance of soil-related water stresses should be interpreted with 417 caution. 418 In spite of the suppressing effects of increasing anomalies in SRAD, VPD, and Tmean, average 419 growth rates were higher in warmer and sunnier forests (i.e. higher long-term means), across 420 species (Fig. 2) and within many species (Table S5). While long-term Tmean was highly correlated 421 with elevation (r = -0.95; We showed that two aspects allowed understanding the broad range of species differences in 431 growth response to VPD anomalies: the long-term mean VPD where species grew better, and the 432 contrast between slow-and fast-growing species (Fig. 3a, c). The models including plot-specific 433 responses to climate anomalies additionally showed that forest growth sensitivity to VPD 434 anomalies was stronger in drier forests, mostly at the higher end of the VPD range (Fig. 4a) water stress and therefore closer to a threshold of further growth decrease than moist rainforests. 441 This effect may not be linear and will need to be further tested with more plots encompassing 442 diverse water-stress conditions. 443 Similarly, Sullivan et al. (2020) recently showed that warmer forests may be closer to a temperature 444 threshold beyond which woody productivity would decrease. In our study, this would translate into 445 expectations that forests and species adapted to warmer conditions would respond more negatively 446 to further temperature increases. Our results are consistent with this expectation but suggest that 447 the temperature effect manifests itself indirectly through VPD. 448 Species that grew faster in cloudier forests showed the strongest growth reduction due to positive 449 SRAD anomalies (Fig. 3b). This may reflect species differences in light-use strategies, with species 450 that grow well under low direct-sunlight conditions not benefitting from brighter conditions. This 451 was supported by the stronger negative effects of SRAD anomalies in species with lower maximum 452 photosynthetic capacity, stomatal conductance and electron transport capacity (Fig. 5) An overarching finding is that species traits can enhance our understanding of differences in species 470 growth response to the anomalies of SRAD and VPD, and to forest stand structure. Our results 471 confirmed that resource-acquisitive species overall had higher intrinsic growth rate and that their 472 growth was more sensitive to positive anomalies in SRAD and VPD. This highlights a tradeoff 473 between fast growth (via high allocation to acquisitive tissues) and sensitivity to atmospheric water 474 stress, consistent with expectations from the 'fast-slow' plant economics spectrum (Reich 2014). 475 Most physiological traits directly related to photosynthesis (Table 1)  values of these traits attenuated the tree growth reduction following increasing SRAD anomalies 480 ( Fig. 5; Fig. S7), suggesting that species investing in a more responsive and flexible photosynthetic 481 machinery may cope better with unusually-high direct exposure to sunlight. While most traits that 482 increased species intrinsic growth rate also exacerbated the negative effects of VPD anomalies on 483 tree growth, the mediation of SRAD anomalies by species traits was mostly independent of the 484 fast-slow spectrum ( Fig. 5; Fig. S5, S6). For example, while leaf P concentration, stable carbon 485 isotope ratio and the maximum photosynthetic capacity tended to increase intrinsic growth rate, 486 the two former accentuated while the latter attenuated the negative effects of SRAD anomalies on 487 tree growth (Fig. 5). 488 Stand structure as driver of tree growth variation 489 Plot basal area consistently strongly reduced tree growth across species and explained more growth 490 variation than mean climate for all four climate variables. Although plot basal area was partly 491 correlated with elevation, the 30-year average of Tmean and VPD (r = 0.63, -0.59, and -0.47, 492 respectively; Table S3a), the slope coefficient of basal area remained virtually unchanged across 493 models including Tmean, VPD, or the other less correlated covariates (and was much steeper than 494 the slopes of long-term Tmean or VPD), so that the stand structure effect detected here is unlikely 495 to indirectly reflect Tmean or VPD. Furthermore, faster growth in less dense environments across 496 forest plots suggests a release from competition for light. This is supported by the general light-497 limitation suggested by the faster growth in sunnier sites. Slower growth in denser environments 498 may also suggest an increase in competition for resources or attacks by natural enemies. 499 Neighbourhood crowding has indeed been shown to strongly reduce tree growth in tropical and 500 these studies, we found that conservative species with high wood density suffered less growth 502 reduction from increasing plot basal area, while acquisitive species with high dark respiration rate 503 and leaf δ 13 C were more sensitive to basal area (Fig. S5, S6). 504 23 505 In summary, we have shown how long-term demographic data across multiple plots encompassing 506 environmental gradients, combined with functional traits collection can yield insights into how 507 climate affects interannual variation of tree growth at different temporal scales, and give important 508 clues into which species and forests may be particularly vulnerable to climate change, and why. 509 Our findings emphasise the importance of functional traits -and notably those directly related to 510 photosynthesis and water use efficiency -to understand species differences in demographic 511 sensitivity to abiotic and biotic drivers. Future efforts to further characterise how climate and 512 neighbourhood crowding affect tree growth, survival, and population growth across environmental 513 gradients, and how these effects are mediated by species traits will help improve predictions of 514 forest response and future ecosystem functions to climate change under different trajectories. and processed data will be archived in a Dryad repository, whose DOI will be added at the end of 754 the article, and will be available from the corresponding author upon request.  Figure S1: Trait turnover along the elevation gradient. 759 Figure S2: Effects of long-term climate average, short-term anomalies, tree size and stand 760 structure on tree growth rate. 761 Figure S3: Illustration of the variability of tree growth responses to the climatic drivers 762 among species. 763 Figure S4 Coordination among the species-level growth responses to stand structure and 764 intrinsic growth rate. 765 Figure S5: Effects of tree size, mean climate, climate anomalies, stand structure, and species 766 functional traits on intrinsic growth rate. 767 Figure S6: Mediation of climatic and stand structure effects on tree growth by species 768 functional traits. 769 Figure S7:   Table S5). 801    Table S5). 841   Although most traits see a significant increase or decrease in their values, the whole trait range 882 remains well represented across the whole gradient. In addition, it is worth noting that tree 883 growth decreased with elevation ( Fig. 2; positive effect of historical mean Tmean on growth, 884 which is strongly negatively correlated to elevation, r = -0.95, see Table S3a), while several 885 traits that increased with elevation (Fig. S1) also had a positive effect on intrinsic growth rate 886 ( Fig. 5 and Fig. S5e) (e.g. Vcmax, LMA). This indicates that the effects of traits on growth are 887 actual trait effects rather than indirect elevation effects. This is confirmed by the fact that the 888 trait effects on intrinsic growth rate (from the model including Tmean as climate variable) (   Fig. S4 shows that tree growth is more reduced by high plot basal area among fast-growing 927 species than among slow-growing species, suggesting a trade-off between a fast average growth 928 rate but higher sensitivity to competition for light or other resources (or natural enemies), and 929 less sensitivity to stand structure but a slower average growth rate (see Discussion). 930 represented here, for clarity, but are in Fig. S6. Circles, thick and thin intervals are median, 90%-and 95%-HPDI of coefficient 941 posterior probability distributions. Red and blue circles indicate negative and positive effects on tree growth, respectively, for 942 the covariates with clear effects (95%-HPDI not encompassing zero); white circles indicate coefficients whose 95%-HPDI 943 include zero. e: Trait effects on intrinsic growth rate (i.e. on the intercept; α1 coefficient) from the model including VPD (Fig.   944 S5c). The direct traits effects on growth rate from other models (a, b, d) were similar and not shown here for clarity.
The models run on the 75 species for which trait data were measured yielded similar results 946 than the models run on all 509 species and without trait effects, regarding the climatic and stand 947 structure covariate effects, with some differences. The effect of tree size, negative overall based 948 on all 509 species (Fig. 2), became positive. This indicates that while tree growth rate still 949 decreased with tree size in some species (Table S5) Tables   978  The supplementary Table S1-S7 are in a separate Excel document. 979

Supplementary Methods S1
980 Study sites and demographic data 981 Individual tree annual absolute growth rates were calculated for 12,853 trees in 23 permanent 982 forest plots of tropical rainforest located in northern Queensland, Australia, between 12°44' S 983 to 21°15' S and 143°15' E to 148°33' E, and encompassing an elevation gradient between 15 984 and 1200 m a.s.l. (Fig. 1a). Twenty of these plots (0.5-ha, 100 × 50 m) were established between 985 1971 and 1980 to provide long-term ecological and demographic data (Bradford et al. 2014), 986 while three plots were established more recently along the same elevation gradient (Table S1). to 3563 mm), solar radiation (17.8 to 19.4 Mj m -2 day -1 ) and vapour pressure deficit (6.5 to 11.8 993 hPa) (Table S1). At plot establishment, all trees with stems ≥ 10 cm diameter at breast height 994 (DBH) were mapped, identified to species level and measured for diameter. The 20 long-term 995 plots were re-measured every two years for ten years, and then at three-to four-year intervals, 996 with diameter, recruits and deaths recorded, summing up to 10 to 16 censuses per plot. The 997 remaining four plots were established more recently, between 2001 and 2012, and were 998 resampled one to three times (Table S1). 999 All available censuses were used to calculate individual annualised absolute growth rate (AGR) 1000 based on DBH at times 1 and 2 (dates; t1 and t2), as: 1001  The matrix of fitted correlation coefficients among species-level intercepts and slopes (αj, β1j, 1193 β2j, β3j and β4j) allows evaluating correlations among species intrinsic growth rate (intercepts 1194 αj) and species AGR sensitivity to model covariates (β1-4j). For instance, a model with a negative 1195 ραj, β3j parameter and a negative β3,0 slope would indicate that species with higher intrinsic growth rate (αj) tend to have higher sensitivity (i.e. more negative slopes) to climate anomalies. 1197 Using covariance matrix to pull information across species-level intercepts and slopes through 1198 the multinormal distribution improves the accuracy of posterior likelihood estimates both across 1199 and within species (hierarchical levels 1 and 2, respectively) while limiting risks of overfitting 1200 through adaptive regularising priors, or partial pooling (e.g. McElreath 2020). In both M1 and M2 models, we standardised the response variable log(AGRi, j, k, t) and all 1250 covariatesbut climate anomaliesto mean zero and unit standard deviation, to allow relative 1251 importance comparisons between covariates through slope coefficients (Schielzeth 2010), and 1252 to ease the assignment of plausible priors to the parameters (McElreath 2020) (eqs. 4.7-4.9). 1253 We did not standardise the averaged monthly anomalies to maintain their interpretability as 1254 deviations from long-term means in terms of plot-specific units of standard deviation (see eq. 1255 2; i.e. mean anomaly covariate slope coefficients are not directly comparable to other covariate 1256 mean slopes). Individual trait measurements were averaged per species and log-transformed 1257 prior to standardisation, thus implying that parameter βc, j corresponds to βc, 0 at the mean trait 1258 value of the dataset. 1259 M3: Plot-level tree growth response to climate anomalies and interaction with mean climate 1260 1261 M3 models evaluate plot-level growth response to climate anomalies, and whether it varies 1262 depending on local mean climates (e.g. whether plot-level AGR sensitivity to VPD anomalies 1263 is higher in drier sites). We focused on the tree growth at the plot level, and modelled the 1264 expected log(AGR) as a linear function of mean climate and climate anomalies. We used a 1265 where αk is the average growth rate in plot k, and β1k characterises the growth response of plot 1283 k to standardised climate anomalies for time interval t. α0 is the mean intercept value (i.e. mean 1284 absolute growth rate) across plots, and α1 is the departure from the grand mean for one unit 1285 increase in mean climate. β1,0 is the grand slope of climate anomalies, and β1,1 is the departure 1286 from this grand mean for a one unit increase in mean climate (mediation of the effect of 1287 anomalies on growth by the plot mean climate; i.e., crosslevel interaction between the plot-level 1288 climate anomaly effect and the population-level mean climate effect). Parameters γj, δt, λi are 1289 varying intercepts for species, census periods, and individual stems, respectively.
We run M3 models only for two climate variables (VPD and SRAD), as we found they were 1291 the most important climate variables for tree growth in M1 and M2 models (see Results). 1292 Standardisation of variables was carried out as for M1. 1293 Trends in climate over time 1294 To explore the implications of the effects of climate anomalies on tree growth, we built a 1295 separate set of hierarchical Bayesian models to test for linear temporal trends in mean annual 1296 climate variables between 1971 and 2019. We used varying year slopes per plots to allow plot-1297 specific trends (model details in Supplementary Methods S1). We also run the models for the  Table S4). There was no general temporal 1303 trend for MCWD or SRAD (Fig. 1c). 1304

561
# There were no major differences here, but a slightly better fit was obtained using the second option 562 # in the tree growth models, so that mcwd_uet-based variables were used.