Family imprint reveals basin-wide patterns of Amazon forest embolism resistance

Domain: chemistry
Published: 2026-02-26
Editorial outcome: major_revision

Abstract
Amazon rainforests face intensifying water stress due to increases in vapour pressure deficit and changing hydrological regimes. Embolism resistance (Ψ₅₀) is a critical metric of tree survival under drought conditions, it is defined as a plant’s capacity to resist disruption of xylem water flow due to air bubble formation from water stress. However, measurements of Ψ₅₀ are only available for a limited number of Amazon locations and species. Conversely, data on forest taxonomic composition are abundant across Amazonia, and if Ψ₅₀ is conserved phylogenetically, these data could provide a way to scale-up drought resistance patterns. Here we evaluate Ψ₅₀ measurements across non-flooded Amazonian tree taxa and reveal a moderate phylogenetic signal, with phylogenetic conservatism evident at the family-level. Notably, Fabaceae is amongst the most embolism-resistant tree families in Amazonia. Leveraging the phylogenetic signal we use species composition and tree size data from 448 forest plots across Amazonia to produce a macroecological assessment of Amazonian vulnerability to embolism. The resulting estimate spatial pattern reveals that forests in the Brazilian and Guiana Shield regions, where Fabaceae abundance is high, show strong resistance to embolism. In contrast, tree communities in Western Amazonia appear more vulnerable to embolism, suggesting a reduced capacity to withstand future drought conditions.

Introduction
The Amazon region is home to the largest and most diverse tropical forest in the world, and plays an important role in planetary biogeochemical cycles. Recent findings have documented substantial changes in non-flooded Amazonian forests (terra-firme), including floristic and functional composition, structure, and dynamics, potentially associated with ongoing changes in climate and atmospheric composition. In recent decades, the Amazon has been subjected to a sequence of large-scale drought events (1998, 2005, 2010, 2015–16, 2023 and 2024), as well as a continued warming of 0.6°–0.7°C since 1950, exposing plants to higher water stress due to increased vapour pressure deficit (VPD). Climate model projections also suggest that the severity of drought effects on forests will continue to increase and that temperatures will likely rise to levels without historical analogues. Together, these climatic changes are expected to exacerbate water stress in Amazon rainforests.

Xylem embolism resistance is a key structural trait that determines the ability of plants to tolerate water stress. This is because water stress is associated with increasingly negative xylem water potentials, which may result in the formation of water vapour/air bubbles (emboli) in the xylem and compromise water transport to the canopy. Typically, embolism resistance is quantified as the xylem water potential at which a tree’s hydraulic conductivity declines to 50% of its maximum value (Ψ₅₀), with a more negative Ψ₅₀ implying higher embolism resistance. In Amazonia, embolism resistance has been shown to explain transpiration and canopy conductance responses to extreme drought, as well as patterns of species distributions and differential mortality patterns under imposed drought.

In recent years, advances have been made in understanding how embolism resistance varies locally and along basin-wide precipitation gradients. These studies reveal that embolism resistance varies in response to broader climatological differences but also local topographical variation in water availability (e.g., refs. 20, 30). Still, a key challenge remains in that embolism resistance is a complex and time-demanding trait to measure, thus far only measured at a small number of locations across Amazonia and therefore making inference of large-scale patterns difficult. Inventory data on tree species composition and dominance, by contrast, are now widely available across Amazonia (e.g., the RAINFOR network encompasses 600 forest plots). If embolism resistance has a phylogenetic signal, then the widespread availability of forest inventory data could enable the use of phylogenetic imputation to scale up trait observations. This approach allows missing values to be estimated based on the assumption that closely related species tend to share similar traits.

Embolism resistance has indeed been shown to be phylogenetically constrained at a global scale, indicating that a species’ resistance to embolism is shaped by its ancestral lineage. Yet, Amazonian tree species are still very poorly represented in such analyses, which span global gradients of water availability and a vast amount of evolutionary clades, potentially leading to a stronger phylogenetic signal than what can be expected for regional-scale analyses. Efforts to determine the phylogenetic conservatism of embolism resistance in Amazonian trees have been made by studies in Central Amazonia at two sites. These studies considered 28 and 16 congeneric species found on valley/plateau or flooded/non-flooded Amazonian habitats, respectively, and found no phylogenetic signal in embolism resistance. Rather, Ψ₅₀ values of species growing in flooded areas with ample water availability were significantly larger than those of congeneric species found in non-flooded areas. This suggests repeated convergent evolution of resistance to embolism based on topo-hydrological position, likely driven by the strong environmental pressure that these environments impose.

While apparently inconsistent with global studies, existing results for Amazonia are based on hydraulic traits datasets covering only a few taxa (e.g., ref. 25 considered six families and eight genera). Given that the Amazon contains more than 6000 known tree species belonging to 803 genera and more than one hundred families, a wider evaluation is essential to fully determine whether or not a phylogenetic signal in embolism- and drought-resistance exists in the Amazon tree flora.

Here we make use of the pan-Amazonian hydraulic trait database, encompassing data for 129 species, 88 genera, 36 families and 14 orders, combined with the most up-to-date Amazonian molecular-genus-level phylogeny to show the extent to which embolism resistance across non-flooded (terra-firme) Amazonian tree taxa is phylogenetically conserved. We then integrate data on hydraulic traits and floristics from lowland non-flooded forests across the entire Amazon Basin to create the most comprehensive macroecological assessment of Amazonian xylem embolism resistance to date. Our dataset includes the most important Amazonian tree families (e.g., Fabaceae, Moraceae, Lecythidaceae and Lauraceae), both in terms of number of stems and aboveground biomass and wood productivity, allowing us to evaluate the relative importance of family vs. genus-level controls on Ψ₅₀ variation and to investigate differences in Ψ₅₀ across the major Amazonian families.

Results
Evidence of phylogenetic control on embolism resistance of Amazonian trees

We found evidence of a moderately strong phylogenetic signal (PS) for embolism resistance of Amazonian tree genera (Fig. 1A, Tab. 1), as demonstrated by a value of Blomberg’s k of 0.44 (p = 0.02) for Ψ₅₀ across the entire dataset. This value is at the top-end of reported k values for tropical forest tree traits; by comparison, wood density across Amazonian taxa was found to show a phylogenetic signal with a Blomberg’s k value of 0.30. This is, however, still lower than the value expected under a Brownian motion model of evolution. The observed phylogenetic signal for embolism resistance remained similar when also accounting for environmental variation (i.e. excluding forests with long dry season length) (SI Fig. 1, Tab. 1). Nested analysis of variance on a sub-sample of our dataset (n = 10 families, 25 genera, 68 species), further revealed that family (F = 2.411, df = 9, p = 0.003) is a stronger predictor of Ψ₅₀ than genus (F = 1.501, df = 15, p = 0.171) or species (F = 3.389, df = 39, p = 0.261), implying that the phylogenetic control of hydraulic properties in Amazonian trees is not a statistical artefact. The sub-sampling was performed to ensure adequate replication within groups and to meet the test requirements.

Fabaceae have particularly high embolism resistance

Phylogenetic randomisation analysis revealed that the Myrtales and Fabales orders were particularly important in influencing the magnitude of the phylogenetic signal within our dataset, both having marked resistance to embolism (Fig. 1B). Among the families, Fabaceae and Rubiaceae stand out as particularly embolism-resistant. An analysis of variance (ANOVA) constrained to 12 widely-occurring Amazon tree families-those present in at least three sites (SI Figs. 2, 3), revealed significant differences in embolism resistance across major Amazon tree families (ANOVA and Tukey HSD post hoc: F = 2.45, df = 6, p = 0.05), with Fabaceae standing out as having, on average, particularly low values of Ψ₅₀ (i.e., more resistant, Fig. 2A). In contrast, Myristicaceae and Euphorbiaceae had lower embolism resistance (less negative Ψ₅₀), while Moraceae had intermediate resistance. Notably, Fabaceae consistently ranked among the most embolism-resistant families across the entire dataset, regardless of sampling density across sites and species (SI Fig. 4).

Fabaceae are the most abundant and ecologically dominant plant family in Amazonia, as well as being broadly distributed and dominant across South and Central American and African tropical rainforests. Given their overarching ecological importance, we then tested whether Fabaceae, on average, are more emboli-resistant than non-Fabaceae within Amazonia. We found that Fabaceae (mean Ψ₅₀: −2.6 ± std 1.0) have a significantly more resistant xylem than non-Fabaceae (mean Ψ₅₀: −1.8 ± std 0.8) across our entire dataset (Fig. 2B, W = 938, p = 0.00015), with the difference being especially marked in the Amazonian forests with an intermediate dry season length that occupy most of the Amazon basin (SI Fig. 5, W = 228, p = 0.0007). A similar pattern was found for ever-wet aseasonal forests (SI Fig. 5, W = 32, p = 0.035), although the sample size of Fabaceae in this region was much smaller than in other regions. In long dry season length forests where species are more adapted to water stress, Fabaceae and non-Fabaceae have equally resistant xylems (SI Fig. 5, W = 83, p = 0.922). Overall, Fabaceae span a broad range of Ψ₅₀ values, reflecting the environmental and taxonomic diversity of the family (SI Fig. 4). Yet, despite this wide variation, within each site across the precipitation gradient, Fabaceae tend to consistently rank among the more embolism-resistant taxa, occupying the resistant end of the hydraulic spectrum across Amazonia forest plots (SI Fig. 6). Analysis of a broader global dataset also yielded similar results, with Fabaceae showing greater xylem resistance when compared to other families, and being among the top most emboli-resistant families across tropical seasonal forests, tropical rainforests and temperate seasonal forests (SI Fig. 7).

Embolism resistance of Amazonian tree communities

Based on the phylogenetic signal observed in our database, we used inventory data from 448 permanent plots across Amazonia (from the RAINFOR forest data network, curated by ForestPlots.net) to impute Ψ₅₀ values to forest plots for which community-level Ψ₅₀ measurements were not available. We gap-filled missing data using the approach of ref. 1, assigning values from the next available taxonomic level (see “methods”). After gap-filling, we calculated the community weighted mean (CWM) of each plot to capture the central tendency of embolism resistance at each site, weighting each species’ trait value by its relative dominance (species basal area/total plot basal area). As these plots include information on species composition, abundance, and the size of each tree, we were able to estimate community-level embolism resistance across Amazonia, accounting for the relative local dominance of each taxon. Our estimates reveal a broad macroecological pattern of basal-area-weighted mean Ψ₅₀ (Fig. 3) where forests in the Brazilian and Guiana shields have tree communities with particularly high embolism resistance. Conversely, Western Amazon forests generally have communities with the least resistant xylem in the Amazon Basin. These patterns were shown to be consistent across two spatial analysis methods (interpolation and spatial-environmental clustering) and regardless of interpolation parameters used (SI Fig. 8, SI Figs. 13–15, and see “methods”).

We also tested whether distance to water table depth (obtained from ref. 46, a key factor driving local water availability, would be related to the pan-Amazonian variation of embolism resistance (SI Fig. 9C). In contrast to what has been found across Amazonia and Cerrado biomes, our results indicate that variation in water table depth did not contribute to explain patterns of Ψ₅₀ distribution at the basin-wide level in Amazonia (statistical results are shown in the figure caption). However, the role of distance to water table on tree drought resistance manifests mainly at a local scale and may not be accurately represented at the levels of generalisation and uncertainty associated with coarse water table depth estimations. These limitations could potentially mask the combined effects of water table depth and regional precipitation gradients that have been shown in previous studies.

Methods
Amazonian tree embolism resistance (Ψ₅₀) dataset

In this study, we used species and genus mean values of Ψ₅₀ from 129 species, 88 genera, belonging to 36 families and 14 orders, obtained from the pan-Amazonian hydraulic traits dataset, which integrates newly collected samples from Western and Southern Amazonia together with published data from Central–Eastern Amazonia. The dataset is deposited as a ForestPlots.net data package, which can be accessed using the ref. 37. In this dataset, embolism resistance was characterised by constructing xylem vulnerability curves for 129 species among 11 sites encompassing effectively the entire Amazon climatological gradient, incorporating forests with: long dry season length, intermediate dry season length, and ever-wet aseasonal forests. Xylem vulnerability curves were constructed for each species at each site by fitting a curve based on pooled data from individuals of the same species, quantifying xylem embolism formation as a function of branch dehydration. To do so, ref. 21 used the pneumatic method, which consists of measuring the air discharge from terminal branch ends to assess embolism formation, and the bench dehydration technique to induce water stress in a given branch.

Phylogenetic signal

We used the most up-to-date Amazon tree genus-level phylogeny to explore the existence of phylogenetic signal (PS) in Ψ₅₀ across different Amazonian genera. To quantify the strength of PS, we used Blomberg’s k values, which provide a measure of the strength of phylogenetic signal based on comparing the observed variance for a given trait against the variance that would be expected under a Brownian motion (BM) model of trait evolution. While ‘k’ values close to zero imply evolutionary independence (random trait distribution among the branches of the phylogenetic tree), values close to one denote phylogenetic non-independence, indicating high trait similarity among closely related clades. The significance (p-value) of k was estimated through a randomisation exercise in which the tips of the phylogenetic tree were randomised 1000 times and the resulting distribution of k values was compared to the observed value of k. We defined k values to be significant if they fell outside the 2.5% − 97.5% percentile range of the simulated distribution. Blomberg’s k test has been shown to be able to detect PS on trees with at least 20 observations and is thus appropriate for the size of our dataset. We further tested the sensitivity of the observed phylogenetic signal (PS) by repeating the analyses on a restricted subset of the dataset. As the spatial datasets spanned dry-adapted long dry season length forests as well as core Amazon forests, we re-ran the analyses excluding these long dry season forests to verify that observed PS patterns were not driven simply by differential sampling of genera across different climate regimes. Data were transformed to meet the normality assumptions of the test (Table 1). Due to the skewed distribution, a cube root transformation was applied when a square root transformation was not sufficient to normalise the data. This choice was guided by visual inspection of histograms and the Shapiro-Wilk normality statistical test. To test which specific taxonomic groups may strongly contribute to PS, we employed the randomisation approach of ref. 43. In brief, the approach consists of first estimating the ancestral value for each node in the phylogeny using ancestral state reconstruction and then randomising the tips of the phylogeny 1000 times to generate a random distribution of ancestral nodes, which are compared to the observed reconstructed node value.

Complementarily, to better understand the variation in Ψ₅₀ explained by different taxonomic levels (species, genus, family), we performed a nested analysis of variance (ANOVA) [family/genus/species]. To ensure the minimum replication within groups and meet the requirements of the test, this analysis was limited to genera for which we had data from more than one species and families with at least two genera or genera/families with the same species/genera sampled in more than 1 site. Data were square root transformed to meet the normal distribution criteria of the test. The reduced dataset in total accounted for 10 families, 25 genera, and 68 species.

Embolism resistance variation across families

To evaluate embolism resistance variation across widely occurring families in our dataset, we conducted one-way ANOVA followed by Tukey’s Honestly Significant Difference (HSD) tests. In this analysis, we selected only families that occurred in at least 3 sites along a wide mean MCWD (Maximum Climatological Water Deficit) gradient to strike a balance between including a broad representation of Amazonian tree diversity while maintaining statistical robustness in family-level comparisons across the precipitation gradient (Fig. 2A, SI Fig. 2).

Given the importance of Fabaceae both in this dataset (SI Figs. 3, 4) and across Amazonian forests more generally we also investigated whether there were differences in embolism resistance between Fabaceae and non-Fabaceae by performing a Wilcoxon rank sum test with continuity correction (Fig. 2B). We did this across the entire dataset but also for specific forest types (long dry season length, intermediate dry season length, and ever-wet aseasonal forests), following ref. 21 (SI Fig. 5). We used the global dataset from ref. 45 to investigate the distribution of Fabaceae embolism resistance in other systems. We focused our analyses on plant families represented by at least two individuals per biome, excluding the temperate rainforests due to limited Fabaceae representation. Winteraceae was excluded due to its lack of vessels (SI Fig. 7).

Geographic patterns of Amazon embolism resistance

Based on the phylogenetic signal observed in our database, we computed the community-weighted mean Ψ₅₀ for a further 448 Pan-Amazonian plots from the RAINFOR inventory network. Individual tree data for each plot were obtained from the ForestPlots.net database. We selected single censuses with the best available level of taxonomic identification at the species level, for each structurally mature lowland tropical forest plot within the Amazonia sensu latissimo definition, therefore excluding dry forests, flooded forests, and plots with elevation > 1000 m above sea level, as well as those affected by direct human disturbance. As the results of the nested ANOVA indicated a dominant family control on embolism resistance, we restricted the mapping exercise to plots that shared at least 60% of their dicotyledonous arboreal family composition, in basal area terms, with our Ψ₅₀ database.

We gap-filled missing data following the methods of ref. 1, whereby we used values for the next available taxonomic level. For example, for species for which we had not measured Ψ₅₀, we used genus-level means when available and family-level means when genus-level data were not available. When family-level data were not available, we used the plot mean value. After this procedure, we then calculated the community weighted mean (CWM) to describe the central tendency of tree community trait values at each site. The CWM₅₀ was calculated by weighting the trait value of each species by its relative dominance (species basal area divided by plot total basal area) in each study plot. To test the gap-filling procedure, we replaced the original species values by family-mean values (simulating gap-filling) for the plots on which Ψ₅₀ were directly sampled (SI Fig. 11B linear model: p < 0.0001; R2 = 0.78) and for the surrounding plots with similar species composition (SI Fig. 11A linear model: p < 0.0001; R2 = 0.72) and compared the resulting CWM Ψ₅₀ values.

Once the plot dataset was gap-filled, we used two combined procedures to estimate macroecological/geographic patterns of Ψ₅₀ across Amazonia. First, we performed a geographically constrained clustering of the gap-filled data points, using the ClustGeo package of R programming language version 4.4.2. This method implements a Ward-like hierarchical clustering algorithm that includes spatial and environmental constraints. The algorithm takes two dissimilarity matrices, D0 representing environmental distances (i.e., ‘feature space’) and D1 representing geographic differences (i.e., ‘constraint space’), and a mixing parameter alpha varying between 0 and 1. The criterion minimised at each stage is a convex combination of the homogeneity criterion calculated from D0 and the homogeneity criterion calculated from D1, with the relative weight of each given by alpha. The value of alpha is optimised by graphical analysis of the relative proportion of explained inertia associated with D0 and D1, so that the ideal alpha value increases spatial contiguity without deteriorating the quality of the solution based on the feature space.

We evaluated the effect of alpha values ranging from 0 to 1 at 0.1 intervals, for a total number of clusters ranging from 2 to 8, using plot-level Ψ₅₀, Mean Cumulative Water Deficit (MCWD) and annual and monthly Water Table Depth (WTD) as our feature space. Regardless of the number of clusters tested, the optimal alpha value was consistently 0.3. We then visually compared clustering results from 2 to 8 clusters, and decided as ideal the largest number of clusters that did not produce clusters with excessive spatial overlap between clusters, nor clusters with too few data points. Our final choice of clustering parameters was then alpha = 0.3 and k = 5 clusters, with the resulting clusters and corresponding average Ψ₅₀ values for each clusters shown on SI Fig. 8.

We then used the interpolation method of Inverse Distance Weighting (IDW) to better visualise this basin-wide geographic pattern as a continuous surface, using the gstat package running on R version 4.4.2. IDW results are mainly contingent on two parameters, the inverse distance power coefficient (idp), which determines the strength of the decay of the weighting with respect to distance, and the maximum number of neighbour points to be considered when calculating a new interpolated value (nmax). To optimise the choice of these parameters, we used a Leave One-Out Cross Validation (LOOCV) strategy where several parameter combinations were tested by interpolating N−1 observations in the dataset, and then calculating the error between the interpolated value and the actual value of the left-out observation. This process is repeated N times to generate a set of N error values for each parameter combination. The parameter combination that minimised the Root Mean Squared Error (RMSE) was selected as optimal. After a three-step parameter search procedure, we converged on the best IDW parametrisation of idp = 0.3 and nmax = 32, and used these parameters to generate interpolated surfaces for all five folds of the cross validation (see below), averaged these five results to produce a final interpolation map over a raster grid with a cell size of 10 × 10 km, using the WGS-84 grid and the Equal Earth Americas projection (EPSG 8858). To reduce linear artefacts resulting from the uneven spatial distribution of the original plots, we smoothed each resulting interpolation using a 9 × 9 mean convolution filter prior to averaging.

To estimate accuracy, we performed a spatially-constrained k-fold cross-validation using the five spatial clusters as data folds to quantify fine-grained spatial uncertainty. For each fold, interpolation was performed for N−1 folds (training data), while the remaining fold (testing data) was used to calculate the Root Mean Squared Error and the cross-validation based Variance Explained (VEcv) uncertainty metrics. Overall, quantitative agreement between field values and interpolated values had an RMSE of 0.24 MPa, with the RMSE value for each fold varying from 0.19 MPa to 0.32 MPa (SI Fig. 13). We can see on the figure that most of the mismatch between observed and interpolated values comes from the smoothing inherent to the IDW interpolation method. For this reason, the calculated Variance Explained (VEcv) was only 6.7%. We further tested the adequacy of IDW as a visualisation tool by assessing the overall agreement of the IDW interpolated trends with the average values of the spatial clusters generated. First, we calculated the mean Ψ₅₀ value for each spatial cluster using the source point data (point averages). We then delineated the minimum convex hull enclosing all data points in each cluster, and used these polygons to extract all Ψ₅₀ pixel values from the interpolated raster, which were averaged to provide a second mean Ψ₅₀ value (interpolated averages). We then compared the agreement between field-based and interpolation-based averages for all five clusters, showing that both had a very good agreement (SI Fig. 14).

Finally, to test the robustness of our generalisation, we performed Multivariate Environmental Similarity Surface (MESS) analysis, using the dismo package, version 1.3.16. This analysis tests how well each location in the prediction (interlation) space matches the environmental conditions of the sourced data (field points). Values between 0 and 100 indicate that the location is within the min-max range of all environmental variables, with 50 indicating that all environmental variables are right in the centre of their min-max ranges. Environmental suitability is calculated individually for each environmental variable at each location, and the final MESS value for the locations is given as the worst similarity score across all variables. MESS values below zero indicate that at least one of the environmental variables at that location exceeds the environmental envelope of the sourced data, and therefore, spatial predictions should be taken with caution. We calculated MESS statistics for our interpolated surface using the Minimum Cumulative Water Deficit and Water Table Depth raster layers as environmental variables, and extracted their corresponding values for each of the 448 field plots as reference values for MESS calculation. The analysis indicated that over 93% of the study area (Amazonia latissimo sensu) had MESS values within the range of 0 to 100, with areas below zero corresponding to isolated pixels, the highlands Andean region and a few regions in the lowland western Amazonia (SI Fig. 15). When plotting the data, we have added for context a terrestrial forest mask, produced by combining the Amazon wetlands mask, Brazilian land cover types from MapBiomas, and the latest deforestation layer from the Brazilian PRODES monitoring programme. Main rivers were obtained from the World Bank Data Catalogue.

Taken as a whole, the multiple approaches to accuracy analysis emphasise that our interpolation and clustering approaches do succeed at showing the broad, basin-wide macroecological patterns of vulnerability to embolism that are already present in the estimated plot-level values. Importantly, however, neither of these procedures produces pixel-level accurate spatial predictions of Ψ₅₀ across Amazonia - and thus we strongly discourage attempts to use the IDW results as a ‘data layer’ for further analysis.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This is an interesting study that tests the degree of phylogenetic signal in xylem vulnerability to cavitation among a sample fo Amazonian tree species. The results suggest a significnat phylogentic signal with Fabaceae exhibiting more negative values as a group.\nMajor issues- The study is a little too Amazon-focussed. Certainly this is an important system and might be used to model rainforest evolutionary patterns, but these data need to be presented in the context of other systems. Many important questions raised by this study were not considered. Such as:\n1. How do these P50 values compare with vlaues from other clades and communities?\n2. P50 and “drought sensitivity” are assumed to be the same thing. If this is the case then a more detailed consideration of the combined effect of phylogeny and climate on the distribution of P50 might be expected.\n3. The range of P50 overall seems very small compared with some other studies of communities (uncited work by Blackman, Brodribb Powers and others). This makes the community of the Amazon trees seem rather convergent in general, but his compnent is missed by failing to consider oither systems.\n4. Are Fabaceae “tougher” than other families in other systems?\n5. How does the P50 of Amazonian Fabaceae compare with ither data?\nThese are a few obvious and important questions that highlight the problem with being too focused on Amazonia. I suggest a significant broadening of perspecitve would greatly benefit the impact of the work.", "reviewer_2": "Review Report for Nature Communications\nJ. V. Tavares, et al.\n“Family imprint reveals basin-wide patterns of Amazon forest embolism resistance”\nDegradation of Amazonian tropical forests due to rapid climate change and land use pressure is an essential concern of the Earth system. Whether and how tree species diversity plays the role of tropical forest resilience remain challenging question.\nBy this paper, authors paid attention to the branch xylem resistance to embolism among angiosperm tree taxa across Amazonian forests, using the recent database of tree species’ Ψ50 values compiled by Tavares et al. (2023) in relation to phylogeny and biogeographic distributions based on the region-wide database of tree inventory plots. Authors highlighted that genera of Fabaceae in particular exhibit high xylem embolism resistance in terms of low Ψ50 values than genera of other families, and low-Ψ50 species tended to be distributed in drier forest types and climatic subregions. Together with the previous papers (Oliveira et al. 2019; Tavares et al. 2023) indicating that Ψ50 as a reasonable indicator of drought resistance variation across tree species and forest types on the Amazonian basin, this paper offers interesting information for understanding the role of phylogenetic diversity and history in forest functioning and predicting the tree community change in Amazonian forests under anthropogenic pressure.\nHere I would pose some points to be considered for easier understanding of wide-ranged readers of this paper.\nGeneral points:\n(1) On non-Fabaceous taxa\nThe Ψ50 dataset indicates a wide variation across taxa. Some families and genera are characterized by high or low Ψ50. On line 283-284, you demonstrate the two orders, Fabales and Myrtales, “significantly” and “markedly” distinct from others based on the phylogenetic linkage analysis. Fig. 1A also indicates genera of other than these two orders show high (and low) resistance as well. Figures 1 and 2 also demonstrate there is large cross-taxa overlap in Ψ50. Under the title of “family imprint”, it is better not (only) mentioning two characteristic APG orders, but indicate typical families in Ψ50 distribution, e.g. Rubiaceae by genus Coussarea. ‘Opposite’ high Ψ50 taxa such as Myristicaceae would be worth mentioned, because they would be more sensitive to ongoing rapid climate change.\n(2) On other traits of Fabaceae\nAs is briefly discussed (line 386-387), Fabaceae species are characterized by symbiotic nitrogen fixation and compound leaves. Authors discuss that small leaflet size reflecting high nitrogen content would underlie high drought resistance of this family. Meantime, it has been suggested that leaves with high nitrogen content tend to have low leaf mass per area (LMA), short leaf lifespan, thus quick turnover of leaves (e.g., Reich, 2014, J. Ecol.). Such leaf trait would contribute to drought resistance of legume tree species. In relation to quick leaf turnover, drought deciduousness is also key to drought resistance (e.g., Oliveila et al. 2021). These would not be necessarily linked to, but act complementary to xylem embolism resistance. It looks authors have collected dataset of LMA of those examined taxa (Tavares et al. 2023), I wonder it worth to compare LMA variation between Fabaceae and non-Fabaceae taxa.\nSpecific points:\nline 152, line 170, ...: The term “terrestrial forest” (for Terra firme) sounds confusing. Better rephrase, e.g. “non-flooded forest”.\nline 204-209: Only abundance or composition? May abundance change caused by tree deaths during severe drought events?\nTable 1: It looks you made transformations for “absolute” Ψ50 values, because water potential is always negative. Explain clearly. Also provide reasoning of cube root transformation rather than square root of “absolute” Ψ50 for drier-ecotone excluded dataset (more skewed distribution of Ψ50).\nFig. 1 (and SM Dig. 1): Fix inset color bar caption “sqrt(p50)”. Either sqrt|Ψ50| or “square-root of absolute Ψ50”, and remove “-” (minus). For SM Fig. 1, “cube root of Ψ50” (not “absolute”).\nline 291: Fabaceae is “broadly dominant across most of the world’s tropical forests” —— to be rephrased. This statement is not the case for insular Southeast Asian tropical rainforests. There are some important species, but we do not observe any family-level dominance or abundance there, in contrast to continental American and African rainforests. I wonder your finding of high drought resistance of Fabaceae may provide a partial explanation of Fabaceae abundance in continental tropical climate experiencing drought periods.\nFig. 3B: In inset “P50” is “Ψ50”.\nline 391: “… and tropical forests more generally” better rephrase or remove as above.\nline 436: Provide the elevation range of sites? (later it appears in line 502).\nline 467-468: Define how to classify “dry ecotonal forests” from others. It appears unclear “ecotonal” in Results. Every time note “dry” or “drier”. Definition is also to be repeated in Results (e.g. line 268; Table 1).\nline 469: Is “climate regimes” biome types? You can also incorporate climate measures (e.g. MCWD groups) instead of forest-biome types.\nAll figures are informative and collectively summarize well the paper.\nSM Fig. 1: It would be better to provide comparable panel to Fig. 1B.\nSM Figs 2 and 3 are interesting. Are some of these forest types “drier ecotonal forests”? If so, indicate them in figure legends.", "reviewer_3": "In this paper, the authors explore phylogenetic and spatial patterns in embolism resistance among tropical tree species.\nThey find a moderate phylogenetic signal, with Fabaceae displaying an especially high embolism risk. The phylogenetically conserved signal allows them to produce a map of embolism risk across the Amazon basin, showing that Fabaceae-dominated systems likewise exhibit low risk.\nThis is a nicely written manuscript with a clear set of analyses and goal. I think it is timely given the focus on tree traits, and the relative difficulty measuring and quantifying traits related to plant embolism risk and moisture regulation more broadly.\nHowever, from a methodological perspective, I think some of the decisions could be better justified, and some of the inferences are not fully supported by the current analyses.\nFirst, the authors restricted their main analysis (e.g., Fig 2) to the widely occurring tree families. But doing so essentially limits this work to being only relevant to the 12 widely-occurring families. For example, what do the results look like if you use all families? Is Fabaceae just low risk among widely occurring families, or is it low even when the others are considered? If part of the goal of this work is to look across a wide number of taxa (as suggested by the authors, L 226-236), then this decision to restrict the data seems counterproductive, and the requirement of \"six sites\" used to determine widely occurring species seems arbitrary. It's also unclear if this subsetting is used throughout the results (e.g., L 294, L 378) or just for specific analyses. The choice to exclude families seems particularly relevant when discussing the relative difference of Fabaceae across forest types (L 294-300). These results suggest that Fabaceae is simply a uniquely widespread dry-adapted species. By excluding these dry-adapted specialists (who I presume are not widespread), then of course Fabaceae might look to be very dry-adapted among the generalists. A more careful treatment of this issue in the results (e.g., stratifying the ANOVA by forest type) could help shed light on this.\nThe authors also note that they exclude specific plots based on forest type and composition (L 503). But yet the resulting map is projected across the entire Amazon basin, without regard to forest type or composition. These decisions seem at odds with one another. Are the excluded plots regionally clustered or located on specific environmental conditions? If regions are excluded, then these should be masked from the extrapolation (or, more rigorously, something like Multivariate Environmental Similarity Surfaces should be used to prevent environmental extrapolation).\nIn line with this, the location of the forest plots is very highly clustered (Fig. 3), and there's no attempt to quantify the relative extent of spatial extrapolation (or environmental extrapolation). In such a setting, using inverse distance to interpolate the space seems suspect, and at a minimum, I would expect to see how predictive accuracy decreases as a function of decay distance. There's no reason to expect that regions 500 km from observed data points can be interpolated accurately, and the resulting map (Fig 3a) supports this, as it essentially just outlines the strong clusters with no ability to truly interpolate the space between.\nTo obtain a realistic goodness of fit for the mapping results, I'd like see how these models perform using more rigorous spatially buffered leave-one-out CV or similar, using a kernel distance at which spatial autocorrelation drops to near zero (see Roberts et al. 2017, doi:10.1111/ecog.02881; Ploton et al 2020, doi: 10.1038/s41467-020-18321-y). I would also like to see the full raw cross-validation predictions for all LOO points, rather than the R2 of means (as in Fig. S5, which surely overestimates R2). I'd also be interested in seeing this cross-validation accuracy split up between the non-gap-filled points (L 508) and the gap-filled points, as I suspect the model performs worse if only worse due to the gap filling. I was also surprised not to see any spatial uncertainty, e.g., at a minimum using resampling or bootstrap approaches to quantify the role of sampling noise. As such, it's difficult to interpret the map with any confidence or make reliable inferences, apart from what is given in the raw plot data (e.g., Fig. 3B).\nLastly, given that the authors wish to explore the role of environmental drivers like water depth, is there a reason they did not use other mapping/modelling approaches that explicitly link embolism to covariates? The current approach of a simple linear model (Fig 8C) is likely confounded by a variety of underlying abiotic factors, such as precipitation and temperature regimes, soil types, etc. If the authors want to make robust inferences about water table depth, I would suggest a more structured approach to better control for confounders and interactions.\nA few minor comments:\nSupplementary Table 1 appears to be missing\nWhat dataset did the water table depth come from? Is this Mattos et al?\nThe data description could have a bit more detail. Are all of the analyses and plots species-level averages? Did the original dataset contain multiple measurements within each species? If so, what the the level of intra vs. interspecific variation in the full dataset?"}, "author_response": "Responses to referees\nManuscript “Family imprint reveals basin-wide patterns of Amazon forest embolism\nresistance”\nTavares et al.\nTracking Number NCOMMS-24-52403\nReferees' comments:\nReviewer #1 (Remarks to the Author):\nComment 1.1. This is an interesting study that tests the degree of phylogenetic signal in xylem\nvulnerability to cavitation among a sample of Amazonian tree species. The results suggest a\nsignificant phylogenetic signal with Fabaceae exhibiting more negative values as a group.\nMajor issues- The study is a little too Amazon-focussed. Certainly this is an important system\nand might be used to model rainforest evolutionary patterns, but these data need to be\npresented in the context of other systems. Many important questions raised by this study were\nnot considered. Such as:\nResponse: We thank the reviewer for these comments on the novelty and significance of our\nmanuscript. Following the referee’s suggestion, we broadened the overall discussion of our\npaper, by adding information about Fabaceae in other systems (Choat et al., 2012) - please\nsee responses to comments 1.4 and 1.5.\nReference: Choat, B., Jansen, S., Brodribb, T.J., Cochard, H., Delzon, S., Bhaskar, R., Bucci,\nS.J., Feild, T.S., Gleason, S.M., Hacke, U.G. and Jacobsen, A.L., 2012. Global convergence\nin the vulnerability of forests to drought. Nature, 491(7426), pp.752-755.\nComment 1.2. 1. How do these P50 values compare with values from other clades and\ncommunities?\nResponse: We believe this topic is answered in the responses to comments 1.4 and 1.5, on\nwhich we added global Ψ₅₀ data from Choat et al (2012).\nComment 1.3. 2. P50 and “drought sensitivity” are assumed to be the same thing. If this is the\ncase then a more detailed consideration of the combined effect of phylogeny and climate on\nthe distribution of P50 might be expected.\nResponse: We agree with the referee that an improvement in the terminology was needed in\nthe text. We therefore removed ‘drought sensitivity’ and just used the terminology ‘resistance\nto embolism’ throughout the text.\nComment 1.4. 3. The range of P50 overall seems very small compared with some other\nstudies of communities (uncited work by Blackman, Brodribb Powers and others). This makes\nthe community of the Amazon trees seem rather convergent in general, but his component is\nmissed by failing to consider other systems.\nResponse: We respectfully disagree with the referee’s assertion that the overall Ψ₅₀ range in\nour study is significantly smaller than in previous studies. When placing our dataset in a\nbroader context, we find that the range of Ψ₅₀ values reported here - derived from Tavares et\nal. (2023) - broadly overlaps with values reported for tropical forests in Choat et al. (2012),\nplease, see figure 1 below.\nIt is also important to highlight key methodological and ecological differences between\nour study and the cited works. Blackman, Brodribb, and Jordan (2012) and Powers et al.\n(2020) focused on leaf embolism resistance, which is expected to differ from branch xylem\nembolism resistance (Levionnois et al., 2020). Furthermore, their studies were conducted in\ntropical montane and dry forests, with some forests experiencing prolonged dry seasons (5–\n6 months dry season length - DSL) and receiving ~1600–1700 mm of mean annual\nprecipitation. This contrasts with our sampling design, which incorporates aseasonal,\nintermediate, and long DSL lowland Amazonian rainforests. Given these climatic differences,\nthe broader Ψ₅₀ range in their dataset is expected. For comparison, Powers et al. (2020)\nreported Ψ₅₀ values ranging from -1.55 to -6.52 MPa, while Blackman, Brodribb, and Jordan\n(2012), studying montane rainforest and dry sclerophyll forest, found values between -1.03\nand -4.31 MPa. As shown in the figure below, our dataset reaches values as low as -5.1 MPa,\ndemonstrating substantial overlap with previously reported ranges.\nReference: Levionnois, S., Ziegler, C., Jansen, S., Calvet, E., Coste, S., Stahl, C., Salmon,\nC., Delzon, S., Guichard, C. and Heuret, P., 2020. Vulnerability and hydraulic segmentations\nat the stem–leaf transition: coordination across Neotropical trees. New Phytologist, 228(2),\npp.512-524.\nFig 1. Tree embolism resistance distribution across Amazonian terrestrial forests (Tavares et\nal., 2023) in grey compared with global vegetation (Choat et al., 2012) in red.\nComment 1.5. 4. Are Fabaceae “tougher” than other families in other systems?\nResponse: To address the referee’s question, we analyzed data from Choat et al. (2012),\nwhich evaluated xylem resistance to embolism (Ψ₅₀) across six biomes: boreal tundra,\ntemperate rainforest, temperate seasonal forest, tropical rainforest, tropical seasonal forest,\nand woodland/shrubland. We focused on plant families represented by at least two individuals\nper biome, excluding temperate rainforests due to limited Fabaceae representation.\nWinteraceae was excluded due to its lack of vessels.\nFabaceae does not rank among the most resistant families (based on median Ψ₅₀) in\nmediterranean biomes (woodland/shrubland), which aligns with expectations given the\nhydraulic strategies of trees in these ecosystems. However, in tropical seasonal forests and\ntropical rainforests and even temperate seasonal forests, Fabaceae tend to exhibit greater\nxylem resistance compared to other families. The broader dataset confirms our findings that\nFabaceae is among the top most emboli- resistant families.\nFig 2. (New SM Figure 6). Xylem embolism resistance (Ψ₅₀) across angiosperm families in\ndifferent biomes, based on the global dataset from Choat et al. (2012). Panels show data for:\n(A) woodland/shrubland, (B) temperate seasonal forest, (C) tropical seasonal forest, and (D)\ntropical rainforest. Red points indicate Fabaceae species; black points represent species of\nother families. Dashed lines represent biome-wide means.\nComment 1.6. 5. How does the P50 of Amazonian Fabaceae compare with other data?\nResponse: Please see response to comments 1.4 and 1.5\nComment 1.7. These are a few obvious and important questions that highlight the problem\nwith being too focused on Amazonia. I suggest a significant broadening of perspective would\ngreatly benefit the impact of the work.\nResponse: We thank the referee for raising this point, which we believe indeed benefits the\nimpact of the work. We attempted to incorporate this suggestion into the text (lines 310-314\nand 525-529) and by adding SM Fig. 6.\nReviewer #2 (Remarks to the Author):\nReview Report for Nature Communications J. V. Tavares, et al. “Family imprint reveals basin-\nwide patterns of Amazon forest embolism resistance”\nComment 2.1. Degradation of Amazonian tropical forests due to rapid climate change and\nland use pressure is an essential concern of the Earth system. Whether and how tree species\ndiversity plays the role of tropical forest resilience remains a challenging question. By this\npaper, authors paid attention to the branch xylem resistance to embolism among angiosperm\ntree taxa across Amazonian forests, using the recent database of tree species’ Ψ50 values\ncompiled by Tavares et al. (2023) in relation to phylogeny and biogeographic distributions\nbased on the region-wide database of tree inventory plots. Authors highlighted that genera of\nFabaceae in particular exhibit high xylem embolism resistance in terms of low Ψ50 values\nthan genera of other families, and low-Ψ50 species tended to be distributed in drier forest\ntypes and climatic subregions. Together with the previous papers (Oliveira et al. 2019; Tavares\net al. 2023) indicating that Ψ50 as a reasonable indicator of drought resistance variation\nacross tree species and forest types on the Amazonian basin, this paper offers interesting\ninformation for understanding the role of phylogenetic diversity and history in forest functioning\nand predicting the tree community change in Amazonian forests under anthropogenic\npressure. Here I would pose some points to be considered for easier understanding of wide-\nrange readers of this paper.\nResponse: We thank the referee for the comments and seek to address the referee’s concerns\nin subsequent comments.\nComment 2.2. General points: (1) On non-Fabaceous taxa the Ψ50 dataset indicates a wide\nvariation across taxa. Some families and genera are characterized by high or low Ψ50. On\nline 283-284, you demonstrate the two orders, Fabales and Myrtales, “significantly” and\n“markedly” distinct from others based on the phylogenetic linkage analysis. Fig. 1A also\nindicates genera of other than these two orders show high (and low) resistance as well. Figures 1 and 2 also demonstrate there is large cross-taxa overlap in Ψ50. Under the title of\n“family imprint”, it is better not (only) mentioning two characteristic APG orders, but indicate\ntypical families in Ψ50 distribution, e.g. Rubiaceae by genus Coussarea. ‘Opposite’ high Ψ50\ntaxa such as Myristicaceae would be worth mentioned, because they would be more sensitive\nto ongoing rapid climate change.\nResponse: Comment accepted. We added the referee's suggestion by mentioning more\nfamilies explicitly on lines 288-297.\nComment 2.3. (2) On other traits of Fabaceaes is briefly discussed (line 386-387), Fabaceae\nspecies are characterized by symbiotic nitrogen fixation and compound leaves. Authors\ndiscuss that small leaflet size reflecting high nitrogen content would underlie high drought\nresistance of this family. Meantime, it has been suggested that leaves with high nitrogen\ncontent tend to have low leaf mass per area (LMA), short leaf lifespan, thus quick turnover of\nleaves (e.g., Reich, 2014, J. Ecol.). Such leaf trait would contribute to drought resistance of\nlegume tree species. In relation to quick leaf turnover, drought deciduousness is also key to\ndrought resistance (e.g., Oliveila et al. 2021). These would not be necessarily linked to, but\nact complementary to xylem embolism resistance. It looks authors have collected dataset of\nLMA of those examined taxa (Tavares et al. 2023), I wonder it worth to compare LMA variation\nbetween Fabaceae and non-Fabaceae taxa.\nResponse: Following the referee’s suggestion, we compared LMA variation between\nFabaceae and non-Fabaceae species. Our results show no significant differences in LMA\nvalues between these groups (please see Figure 3 below, W = 1229, p-value = 0.17). Within\nFabaceae, we further examined whether specific subgroups (Caesalpinioideae, Dialioideae,\nand Papilionoideae), leaf architecture (bipinnate, pinnate, unifoliate leaves), or nitrogen-\nfixation ability influenced the high embolism resistance observed in this family. However, our\nanalysis indicates that embolism resistance is a general trait across Fabaceae, with no\nevidence that these factors drive variation in Ψ₅₀ values. While these factors may contribute\nto Fabaceae’s overall success, our results indicate that embolism resistance may be another\nlayer contributing to the success of this Family and is independent of leaf traits such as LMA\nor nitrogen-fixation ability. We added this discussion into the text (lines 401-414) and added\nthe figure below into the Supplementary Material (SM Fig. 9)\nFig 3. (New SM Figure 9). Fabaceae embolism resistance in variation across leaf traits,\nphylogenetic subgroups and nitrogen fixation ability. (A) Comparison of leaf mass per area\n(LMA) between Fabaceae (red) and non-Fabaceae (black) species across our entire dataset.\n(B–D) Variation in xylem embolism resistance (Ψ₅₀) within Fabaceae subgroups: (B) by\nphylogenetic subfamily (Caesalpinioideae, Dialioideae, Papilionoideae), (C) by leaf\narchitecture (bipinnate, pinnate, unifoliate), and (D) by nitrogen-fixation ability (yes/no; based\non Gei et al., 2018). Red diamonds indicate group means; black dashed lines indicate the\noverall Fabaceae mean Ψ₅₀; blue dashed lines represent the overall dataset mean.\nReference: Gei, M. et al. (2018) ‘Legume abundance along successional and rainfall gradients\nin Neotropical forests’, Nature Ecology and Evolution, 2(7), pp. 1104–1111.\nSpecific points:\nComment 2.4. line 152, line 170, ...: The term “terrestrial forest” (for Terra firme) sounds\nconfusing. Better rephrase, e.g. “non-flooded forest”.\nResponse: We incorporate the referee’s suggestion by replacing “terrestrial forest” for “non-\nflooded forest” throughout the text.\nComment 2.5. line 204-209: Only abundance or composition? May abundance change caused\nby tree deaths during severe drought events?\nResponse: We understand that the referee has assumed that we are talking about 'abundance'\nin an ecological sense, instead, we used ‘abundance’ in a generic sense. To avoid confusion,\nwe replaced 'abundance of forest inventory measurements' with 'widespread availability'\nof forest inventory data.\nComment 2.6. Table 1: It looks you made transformations for “absolute” Ψ50 values, because\nwater potential is always negative. Explain clearly. Also provide reasoning of cube root\ntransformation rather than square root of “absolute” Ψ50 for drier-ecotone excluded dataset\n(more skewed distribution of Ψ50).\nResponse: We have clarified in the text (lines 279-281 and 496–499) that transformations\nwere performed to the absolute values of Ψ₅₀ to satisfy the test requirement of normality.\nRegarding the choice of cube root transformation rather than square root, when we excluded\nthe drier ecotonal forests (long DSL) the distribution of Ψ₅₀ values became more strongly right-\nskewed and a square root transformation was insufficient to achieve normality for statistical\ntesting. We therefore applied a cube root transformation that better normalized the distribution.\nThis choice was guided by visual inspection of histograms and Q-Q plots, as well as the\nShapiro-Wilk normality statistical test.\nComment 2.7. Fig. 1 (and SM Dig. 1): Fix inset color bar caption “sqrt(p50)”. Either sqrt|Ψ50|\nor “square-root of absolute Ψ50”, and remove “-” (minus). For SM Fig. 1, “cube root of Ψ50”\n(not “absolute”).\nResponse: We thank the referee for spotting this mistake, which is now fixed.\nComment 2.8. line 291: Fabaceae is “broadly dominant across most of the world’s tropical\nforests” —— to be rephrased. This statement is not the case for insular Southeast Asian\ntropical rainforests. There are some important species, but we do not observe any family-level\ndominance or abundance there, in contrast to continental American and African rainforests. I\nwonder if your finding of high drought resistance of Fabaceae may provide a partial\nexplanation of Fabaceae abundance in continental tropical climates experiencing drought\nperiods.\nResponse: We rephrased the sentence as “Fabaceae are the most abundant and ecologically\ndominant plant family in Amazonia (Ter Steege et al., 2013; Fauset et al., 2015, Cardoso et\nal., 2017), as well as being broadly distributed and dominant across South and Central\nAmerican and African tropical rainforests (Gentry, 1988)”. (lines 300-301).\nComment 2.9. Fig. 3B: In inset “P50” is “Ψ50”.\nResponse: We thank the referee for spotting this mistake, which is now fixed.\nComment 2.10. line 391: “… and tropical forests more generally” better rephrase or remove\nas above.\nResponse: Referee’s suggestion is now incorporated.\nComment 2.11. line 436: Provide the elevation range of sites? (later it appears in line 502).\nResponse: We think that the elevation range of the sites is better place in line 502 (first\nsubmitted version) now numbered as line 538 (current revised version), because it\nincorporates all the 448 forest plots used in the macroecological assessment, instead of\nmentioning this information for only the 11 sites from Tavares et al., 2023, to which line 436\n(first submitted version) is referring to.\nComment 2.12. line 467-468: Define how to classify “dry ecotonal forests” from others. It\nappears unclear “ecotonal” in Results. Every time note “dry” or “drier”. Definition is also to be\nrepeated in Results (e.g. line 268; Table 1).\nResponse: Referee’s suggestion is now incorporated. We replaced “ecotonal” from the results\nand kept the same terminology used in Tavares et al., (2023) through the text: long dry season\nlength, intermediate dry season length, and ever-wet aseasonal forests.\nComment 2.13. line 469: Is “climate regimes” biome types? You can also incorporate climate\nmeasures (e.g. MCWD groups) instead of forest-biome types.\nResponse: Referee’s suggestion is now incorporated, please see response to the comment\nabove.\nComment 2.14. All figures are informative and collectively summarize the paper well.\nResponse: We thank the referee for the comment.\nComment 2.15. SM Fig. 1: It would be better to provide a comparable panel to Fig. 1B.\nResponse: Referee’s suggestion is now incorporated.\nComment 2.16. SM Figs 2 and 3 are interesting. Are some of these forest types “drier ecotonal\nforests”? If so, indicate them in figure legends.\nResponse: Referee’s suggestion is now incorporated. We now indicate in the legend which\nare the long DSL, intermediate DSL and ever-wet aseasonal forests.\nReviewer #3 (Remarks to the Author):\nComment 3.1. In this paper, the authors explore phylogenetic and spatial patterns in embolism\nresistance among tropical tree species. They find a moderate phylogenetic signal, with\nFabaceae displaying an especially high embolism risk. The phylogenetically conserved signal\nallows them to produce a map of embolism risk across the Amazon basin, showing that\nFabaceae-dominated systems likewise exhibit low risk.\nThis is a nicely written manuscript with a clear set of analyses and goal. I think it is\ntimely given the focus on tree traits, and the relative difficulty measuring and quantifying traits\nrelated to plant embolism risk and moisture regulation more broadly. However, from a\nmethodological perspective, I think some of the decisions could be better justified, and some\nof the inferences are not fully supported by the current analyses.\nResponse: We thank the referee for these comments on the novelty and significance of our\nmanuscript and seek to address the referee’s concerns in subsequent comments.\nComment 3.2. First, the authors restricted their main analysis (e.g., Fig 2) to the widely\noccurring tree families. But doing so essentially limits this work to being only relevant to the\n12 widely-occurring families. For example, what do the results look like if you use all families?\nIs Fabaceae just low risk among widely occurring families, or is it low even when the others\nare considered?\nResponse: We restricted the main analysis in Fig. 2A, specifically subsetting families that\noccurred in 3 or more sites in our dataset, totalling 12 families. However, we agree with the\nreviewer that it is also valuable to examine the full distribution of all sampled families. We have\nnow included this extended comparison in SM Fig. 4 (also shown below), which presents the\nΨ₅₀ values across all families in the dataset, regardless of sampling density across site and\nspecies. This allows readers to assess the relative position of Fabaceae within the broader\nphylogenetic context.\nFinally, we would like to emphasize that our Fabaceae vs. non-Fabaceae comparisons (e.g.,\nin Fig. 2B) used data from all non-Fabaceae taxa in the dataset—not only those from the 12\nmost widely occurring families. The restriction to 12 families applies only to the left panel of\nFig. 2A, which was intended to highlight family-level differences across the full Amazonian\nprecipitation gradient. We made the subsetting rationale and when it is applied more clear to\nthe readers in lines (289-290, 301-307 and 513-529)\nFig. 4 (New SM Figure 4). Embolism resistance variation across all sampled families in the\ndataset, regardless of sampling density across site and species. Dashed horizontal lines show\nmean trait value across families. Boxplots show the 25th percentile, median and 75th\npercentile. Vertical bars show the interquartile range ×1.5 and data points beyond these bars\nare potential outliers. The dashed vertical line shows the mean 𝛹 across the pan-Amazonian\n50\ndataset.\nComment 3.3 If part of the goal of this work is to look across a wide number of taxa (as\nsuggested by the authors, L 226-236), then this decision to restrict the data seems\ncounterproductive, and the requirement of \"six sites\" used to determine widely occurring\nspecies seems arbitrary.\nResponse: We thank the reviewer for raising this point and we revised the text accordingly\n(lines 289-290, 301-307 and 513-529). Through the referee’s question we realised that our\ntext contained a mistake and also it was not clear enough to the reader. Our decision to restrict\nthe main family-level comparison (Fig. 2A) to widely occurring families was made to ensure\nthat differences in embolism resistance reflected intrinsic taxonomic or functional variation\nrather than environmental filtering linked to narrow geographic or climatic distributions. For\nthis, we selected species that occurred in at least three sites (not six - this was a mistake in\nthe text), please see SM Fig. 2. We have chosen this to strike a balance between including a\nbroad representation of Amazonian tree diversity while maintaining statistical robustness in\nfamily-level comparisons across the precipitation gradient.\nComment 3.4. It's also unclear if this subsetting is used throughout the results (e.g., L 294, L\n378) or just for specific analyses. The choice to exclude families seems particularly relevant\nwhen discussing the relative difference of Fabaceae across forest types (L 294-300). These\nresults suggest that Fabaceae is simply a uniquely widespread dry-adapted species. By\nexcluding these dry-adapted specialists (who I presume are not widespread), then of course\nFabaceae might look to be very dry-adapted among the generalists. A more careful treatment\nof this issue in the results (e.g., stratifying the ANOVA by forest type) could help shed light on\nthis.\nResponse: We appreciate the referee’s thoughtful comment and acknowledge that the\ndistinction between analyses based on the full dataset versus the subset of widely occurring\nfamilies was not always made clear in the original text. To improve clarity, we have now revised\nthe Results section to explicitly state, at each relevant point (e.g., Lines 288-297 and figure\ncaptions), whether the full dataset or a subset was used, along with the rationale behind each\nchoice. Specifically, for the analysis presented in Lines 294–300 (submitted version), we used\nthe entire dataset, not only the subset of widely occurring families. The subset of 12 families\nwas used exclusively for the analysis shown in Fig. 2A (as discussed in our response to\nComment 3.2), to minimize the confounding influence of families with narrow environmental\ndistributions on broad taxonomic comparisons.\nRegarding the reviewer’s point on Fabaceae appearing to be uniquely dry-adapted, we\nwould like to clarify that Fabaceae species in our dataset are broadly distributed across the\nfull Amazon precipitation gradient (SM Fig 3). As seen in that figure, Fabaceae is not restricted\nto dry or highly seasonal forests but occurs widely, including in ever-wet regions. To further\naddress this concern, we conducted ANOVA tests within each forest type, comparing Ψ₅₀\nvalues between Fabaceae and non-Fabaceae species (SM Fig 5). These analyses revealed\nno significant differences in Ψ₅₀ within forests experiencing long dry seasons, supporting the\nidea that Fabaceae’s high embolism resistance is not merely an artifact of excluding dry-\nadapted specialists with narrow ranges. These additional analyses and their interpretations\nare now presented in the Supplementary Material and referenced accordingly in the main text\n(Lines 308–314).\nComment 3.5. The authors also note that they exclude specific plots based on forest type and\ncomposition (L 503). But yet the resulting map is projected across the entire Amazon basin,\nwithout regard to forest type or composition. These decisions seem at odds with one another.\nAre the excluded plots regionally clustered or located on specific environmental conditions? If\nregions are excluded, then these should be masked from the extrapolation (or, more\nrigorously, something like Multivariate Environmental Similarity Surfaces should be used to\nprevent environmental extrapolation).\nResponse: We appreciate the referee’s thoughtful comment. We would like to clarify that the\ninterpolation used to produce Figure 3 was employed as a visualization device only, to\nemphasize the broader, large scale pattern of Ψ₅₀ in the basin. That is why we don’t even\ninclude numeric values on the plot legend for the interpolated values. We fully agree that if the\ngoal was to produce a high-resolution, spatially explicit model or ‘layer’ of embolism resistance\nacross the Amazon, a more rigorous approach would be required, particularly to avoid\nextrapolation into unsampled environmental space. In that case, incorporating tools such as\nMultivariate Environmental Similarity Surfaces (MESS) or other environmental masking\napproaches would be essential. But while generating such a data layer would be an interesting\nand relevant research problem (which has been partially addressed by the recent publication\nby Chen, S.et al. Nature 631, 111–117 (2024). https://doi.org/10.1038/s41586-024-07568-w\nusing satellite and environmental data to predict hydraulic safety margins), we feel it is beyond\nthe scope of this paper, especially since, as recognized by the reviewer, the spatial distribution\nof the points is far from ideal to support quantitative modelling.\nMasking forest type and composition is also a good suggestion. In response to it, we\nhave improved our treatment of the interpolation by applying a mask to restrict the visualization\nto terra-firme forests (Revised Fig 3, which is also shown below). Specifically, we now exclude\nflooded forests, white-sand forests, and deforested areas, which differ substantially in edaphic\nand ecological conditions and were not included in our core analysis. These areas were\nmasked using data from MapBiomas Collection 8 and the wetlands dataset by Hess et al.\n(2003). We believe these adjustments help ensure that the visualization aligns better with the\ndata used in the analysis.\nFig 5. (Revised Figure 3). Estimated basin-wide spatial variation of Amazonian vulnerability\nto embolism. A. Macroecological patterns of community weighted mean value of 𝛹 , created\n50\nusing Inverse Weighted Distance interpolation of 448 upland moist forest plots distributed\nacross Amazonia sensu latissimo (see methods). Our analyses exclude dry forests, flooded\nforests, and plots with elevation >1000 m above sea level, as well as those affected by direct\nhuman disturbance. B. Estimated community weighted mean values for all 448 non-flooded\nforest plots across the basin. Each dot represents a forest plot and the colour shows its\ncommunity-weighted mean 𝛹 Due to the spatial scale of the map, plot locations in panel B\n50.\nare displaced (jittered) to remove overlap and improve visualization; exact plot locations are\nshown in panel. On both panels, we masked out areas of Amazonia with very different\nenvironments, notably flooded forests, white sand forests and deforested areas, to emphasise\nour focus on upland (terra-firme) moist forests.\nComment 3.6. In line with this, the location of the forest plots is very highly clustered (Fig. 3),\nand there's no attempt to quantify the relative extent of spatial extrapolation (or environmental\nextrapolation). In such a setting, using inverse distance to interpolate the space seems\nsuspect, and at a minimum, I would expect to see how predictive accuracy decreases as a\nfunction of decay distance. There's no reason to expect that regions 500 km from observed\ndata points can be interpolated accurately, and the resulting map (Fig 3a) supports this, as it\nessentially just outlines the strong clusters with no ability to truly interpolate the space\nbetween.\nResponse: The reviewer is correct in stating that the absolute accuracy of IDW interpolation\nwould not be high in this case. However, we emphasize that our use of IDW was solely as a\nvisualization tool to illustrate broader basin-wide patterns and highlight spatial clustering,\nrather than for precise estimations at unsampled locations. This approach has been previously\napplied in ecological and remote sensing studies, such as Mitchard et al. (2014), who used\nIDW to generate a biomass map from forestplots.net clusters. Their study demonstrated the\nutility of IDW in capturing broad-scale spatial trends from irregularly distributed field data while\nacknowledging its limitations for precise spatial predictions. Similarly, our application of IDW\nserves to outline the major clusters in space rather than to provide absolute values with high\nlocal accuracy. See our response below for a further exploration of the dataset to\naccommodate this.\nReference: Mitchard, E.T., Feldpausch, T.R., Brienen, R.J., Lopez-Gonzalez, G.,\nMonteagudo, A., Baker, T.R., Lewis, S.L., Lloyd, J., Quesada, C.A., Gloor, M. and Ter Steege,\nH., 2014. Markedly divergent estimates of Amazon forest carbon density from ground plots\nand satellites. Global ecology and biogeography, 23(8), pp.935-946.\nComment 3.7. To obtain a realistic goodness of fit for the mapping results, I'd like see how\nthese models perform using more rigorous spatially buffered leave-one-out CV or similar,\nusing a kernel distance at which spatial autocorrelation drops to near zero (see Roberts et al.\n2017, doi:10.1111/ecog.02881; Ploton et al 2020, doi: 10.1038/s41467-020-18321-y). I would\nalso like to see the full raw cross-validation predictions for all LOO points, rather than the R2\nof means (as in Fig. S5, which surely overestimates R2). I'd also be interested in seeing this\ncross-validation accuracy split up between the non-gap-filled points (L 508) and the gap-filled\npoints, as I suspect the model performs worse if only worse due to the gap filling. I was also\nsurprised not to see any spatial uncertainty, e.g., at a minimum using resampling or bootstrap\napproaches to quantify the role of sampling noise. As such, it's difficult to interpret the map\nwith any confidence or make reliable inferences, apart from what is given in the raw plot data\n(e.g., Fig. 3B).\nResponse: To accurately estimate spatial autocorrelation (i.e. a semivariogram) we would\nneed a uniform distribution of the distances among points, from near to very far. With our\nclustered point distribution, no acceptable semivariance model could be fit, a problem also\nencountered by Mitchard et al. (2014), who used IDW to generate a biomass map from\nforestplots.net clusters. What we opted to do instead was to apply a spatially-constrained\nenvironmental clustering algorithm to the dataset (using the ClustGeo R package,\nhttps://cran.r-project.org/web/packages/ClustGeo/index.html). Using this algorithm, we first\nclustered the data points into spatially coherent groups of similar Ψ₅₀, MCWD and Water Table\nDepth values. We then used a k-fold validation strategy excluding each cluster in turn to\ngenerate 5 separate IDW interpolations, quantifying the errors using the remaining cluster\ndata. This ensures that the validation data is not spatially correlated with the training data at\neach fold. Our estimated RMSE averaged over the 5 folds using this strategy was 0.24 (s.d.\n0.05), and we included the range of RMSE values in the main text on lines 596-599. To further\nsmooth out the final visualization and better incorporate the error variability, we produced a\nfinal interpolation map by averaging the five interpolated maps resulting from each fold. We\nalso added SM Fig 7 (presented below) into the manuscript to show the results of the\nclustering, to further emphasize the focus on the broader geographic pattern. Finally, we\nmodified the text to make it more clear to users that the resulting ‘map’ should not be\nconsidered an accurate representation of absolute Ψ₅₀ values over the Amazon, but as a\nvisualization of a broad geographic pattern. In fact, we removed all mentions of the word map\nand refer to it instead as a macroecological/geographic pattern.\nFig. 6 (New SM Figure 7). Top: results of the geographically constrained clustering of plots\nbased on the combination of 𝛹 , MCWD (Mean Cumulative Water Deficit) and WTD (Water\nAmazonia sensu latissimo (see methods). Areas with very different edaphic environments -\nflooded forests, white sand forests and deforested areas – are masked from the background.\nBottom: Distribution of 𝛹 values within each spatial cluster.\nlike water depth, is there a reason they did not use other mapping/modelling approaches that\nexplicitly link embolism to covariates? The current approach of a simple linear model (Fig 8C)\nis likely confounded by a variety of underlying abiotic factors, such as precipitation and\ntemperature regimes, soil types, etc. If the authors want to make robust inferences about water\ntable depth, I would suggest a more structured approach to better control for confounders and\ninteractions.\nResponse: We agree with the reviewer, but as answered above for comment 3.5 , our intent\nwas not to quantitatively and accurately predict Ψ₅₀ values spatially across the basin nor to\ngenerate a ‘data layer’, but rather to demonstrate the existence of a broad\nmacroecological/geographic pattern. This level of complex modelling is worth investigating in\nthe future (and in fact has been already approached, as we mentioned, in terms of hydraulic\nsafety margin by Chen et al. (2024)), but this would comprise an entirely different study and\nmanuscript.\nNevertheless, we would like to emphasize the fundamental role of our expanded community-\nlevel dataset (achieved through phylogenetic imputation) in, for example, supporting remote\nsensing efforts to assess forest resistance and resilience from space. Remotely sensed\nvegetation indices have been widely used to investigate resilience loss in forests (e.g. Forzieri\net al., 2022; Tai et al., 2022), but ground-based validation of these approaches remains\nscarce. Our study provides crucial empirical insights into vegetation resistance to embolism,\nbased on an extensive dataset spanning 448 forest sites across Amazonia. This large-scale\ndataset can be harnessed to validate and improve future remote-sensing based models\nquantifying forest resilience, and to enhance the interpretation of remote sensing indicators,\nthus improving our ability to detect early warning signals of forest functioning decline.\nReference:\nChen, S., Stark, S.C., Nobre, A.D., Cuartas, L.A., de Jesus Amore, D., Restrepo-Coupe, N.,\nSmith, M.N., Chitra-Tarak, R., Ko, H., Nelson, B.W. and Saleska, S.R., 2024. Amazon forest\nbiogeography predicts resilience and vulnerability to drought. Nature, 631(8019), pp.111-117.\nTai, X., Trugman, A.T. and Anderegg, W.R., 2023. Linking remotely sensed ecosystem\nresilience with forest mortality across the continental United States. Global Change Biology,\n29(4), pp.1096-1105.\nForzieri, G., Dakos, V., McDowell, N.G., Ramdane, A. and Cescatti, A., 2022. Emerging\nsignals of declining forest resilience under climate change. Nature, 608(7923), pp.534-539.\nA few minor comments:\nComment 3.9. Supplementary Table 1 appears to be missing\nResponse: We thank the referee for spotting this mistake. Supplementary Table 1 refers to an\nanalysis that was part of an earlier version of the manuscript but has since been removed\nduring revisions. We have now updated the supplementary materials to reflect the current\nversion of the manuscript and removed the reference to this table to avoid confusion.\nComment 3.10. What dataset did the water table depth come from? Is this Mattos et al?\nResponse: We thank the referee for spotting this mistake. Information about water table depth\nwas obtained from a raster file available at Fan et al., (2013). We added this information in the\nlegend of SM Fig 8 and methods.\nReference: Fan, Y., Li, H. and Miguez-Macho, G., 2013. Global patterns of groundwater table\ndepth. Science, 339(6122), pp.940-943.\nhttp://thredds-gfnl.usc.es/thredds/catalog/GLOBALWTDFTP/catalog.html\nComment 3.11. The data description could have a bit more detail. Are all of the analyses and\nplots species-level averages? Did the original dataset contain multiple measurements within\neach species? If so, what the the level of intra vs. interspecific variation in the full dataset?\nResponse: We have now added a more detailed description of the dataset and analytical\nstructure to the Methods section (please, see lines 467–470). The Ψ₅₀ data available from\nTavares et al., (2023) consists of one fitted vulnerability curve per species per plot, based on\npooled data from individuals of the same species at each site. Because of this structure, our\ndataset does not include multiple individual-level measurements within species at the same\nsite, which limits our ability to explicitly quantify intra- versus interspecific variation in embolism\nresistance."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The MS first claims to show strong phylogenetic signal in xylem vulnerability, which is then used to impute P50 in forest plots where species were not measured, then interpolated across the landscape to create a “resistance to embolism?” map. This is an interesting method that creates a nice overview of possible landscape function, but the accuracy of this approach is highly questionable in my opinion.\nA key finding is stated as the demonstration that Fabaceae are more embolism resistant than other families in Amazonia. Combined with the finding that “family (F=2.74, p=0.02) is a stronger predictor of (cid:1)50 than genus (F=1.88, p=0.07) or species (F=1.1, p=0.4)”. This leads the reader to assume that much of the imputation is driven by family-level means. The problem with this approach is that the variation in P50 in Fabaceae spanned almost the entire range of all other taxa. When combined with the fact that the Fabaceae mean P50 was not significantly different from other families, one begins to see that the (potentially adaptive) variation that is likely to be captured by the predictions is extremely low.\nOf the 3 outputs of this paper, none seem very robust.\n1. That phylogeny has an influence on P50 is already known\n2. That Fabaceae is a particularly tough does not seem particularly well supported considering the enormous diversity of\nP50 in the family (presumably correlated with the large diversity of the family).\n3. That coarse scale mapping on P50 provides a useful insight into regional-scale “vulnerability” seems unlikely or at least\ndemanding of a far more robust ground-truthing.\nThus, I feel that a much more cautious approach to the interpretation of these data is warranted. Having said this, I do\nbelieve there is promise in such an approach, so perhaps the manuscript simply requires a more circumspect rather than\ncategorical viewpoint.", "reviewer_2": "Dear Authors,\nI found authors revised the paper throughout in response to reviewers' comments and questions inclusive of my own.\nAdded information in review responses added further information. Revised paper presents clear and new insights from the\nregion-wide phylogenetic and geographic analysis focusing on stem embolism tolerance are worth to be published.\nReviewer 2", "reviewer_3": "I appreciate the clarification the authors have provided regarding the subsetting of data and the comparison of Fabaceaes\nrelative to other families. The additional analysis and clarification of the text address my key concerns in this area, providing\nconvincing evidence that Fabaciea is not simply a uniquely widespread drought-tolerant specialist. I likewise appreciate the\nexclusion of specific forest types in the spatial projection, as well as the clustered cross-validation for overall model fit.\nThough, for this last point, what is the R2 of observed vs. predicted rather than the RMSE? Depending on how the data are\nnormalised, an RMSE of 0.24 could equate to an R2 as low as RMSE^2 = 0.05 (see Li 2017 for a discussion on this point\nhttps://doi.org/10.1371/journal.pone.0183250).\nWith this in mind, there is still inadequate assessment of the uncertainty in the spatial projections. The authors state that they\nare unable to do traditional spatial uncertainty/extrapolation analyses because \"the spatial distribution of the points is far\nfrom ideal to support quantitative modelling,\" and that, \"with our clustered point distribution, no acceptable semivariance\nmodel could be fit.\" But this is exactly the issue -- if the data do not allow for appropriate statistical analysis, then how is the\nanalysis justified? The authors state that, \"our intent was not to quantitatively and accurately predict Ψ₅₀ values spatially\nacross the basin...\" and that determining the extent of spatial extrapolation is, \"beyond the scope of this paper...\" Yet they are\nvisually depicting these spatial patterns and referencing these patterns in the title, so one would hope that doing so\n\"accurately\" is an aim. These calculations (the MESS statistics, spatial uncertainty via boostrapping, etc) are relatively easy\nto do with so few points, and would provide important context for understanding which regions are being estimated with\nhigher/lower uncertainty.\nThe justification the authors provide for their approach is Mitchard et al. (2014). But in the decade since this paper was\npublished, there has been growing consensus about the appropriate analysis of spatial patterns and what is necessary\nwhen producing spatial projections. It's certainly possible that the spatial projections are reasonable, but without estimates of\nuncertainty/extrapolation, it's equally likely that the patterns are heavily biased. The uncertainty/extrapolation analysis\ndoesn't need to be incredibly laborious, but even some indication of which regions have higher/lower uncertainty would\nhelp. If the data do not allow for proper error assessment, or if the spatial projections have a small area of applicability or low\ngoodness of fit, then I'm not sure that a spatial analysis is feasible given the current data."}, "author_response": "Responses to referees\nManuscript “Family imprint reveals basin-wide patterns of Amazon forest embolism\nresistance”\nTavares et al.\nTracking Number NCOMMS-24-52403A\nReferees' comments:\nReviewer #1 (Remarks to the Author):\nComment 1.1. The MS first claims to show strong phylogenetic signal in xylem vulnerability,\nwhich is then used to impute P50 in forest plots where species were not measured, then\ninterpolated across the landscape to create a “resistance to embolism?” map. This is an\ninteresting method that creates a nice overview of possible landscape function, but the\naccuracy of this approach is highly questionable in my opinion.\nResponse: We thank the referee for highlighting the importance of the study and we have tried\nto clarify the referee’s concerns over the uncertainty of our approach in the next comments.\nComment 1.2. A key finding is stated as the demonstration that Fabaceae are more embolism\nresistant than other families in Amazonia. Combined with the finding that “family (F=2.74,\np=0.02) is a stronger predictor of 𝛹50 than genus (F=1.88, p=0.07) or species (F=1.1,\np=0.4)”. This leads the reader to assume that much of the imputation is driven by family-level\nmeans. The problem with this approach is that the variation in P50 in Fabaceae spanned\nalmost the entire range of all other taxa. When combined with the fact that the Fabaceae\nmean P50 was not significantly different from other families, one begins to see that the\n(potentially adaptive) variation that is likely to be captured by the predictions is extremely\nlow.\nResponse: We agree with the referee that Fabaceae exhibits a wide range of Ψ₅₀ values in our\ndataset. However, it is important to note that Fabaceae species are broadly distributed across\nthe entire Amazon precipitation gradient (SM Fig. 3), including ever-wet and intermediate and\nlong DSL regions, and consistently rank among the more embolism-resistant taxa within each\nsite (Fig. 1; now also added to the SM Fig. 6). This pattern indicates that, despite local\nenvironmental variation, Fabaceae tend to occupy the resistant end of the hydraulic spectrum\nacross Amazonia forest plots (lines 310-313).\nBecause of this within-family variation, mean Ψ₅₀ values do not differ significantly\nwhen comparing Fabaceae against other families with wide occurrences (those present in ≥3\nsites; Manuscript Fig. 2A; see lines 289–297 and 518–522). Nevertheless, when Fabaceae are\ncompared to all non-Fabaceae taxa combined, they emerge as markedly more embolism\nresistant (Manuscript Fig. 2B). This underscores that while the family spans a diversity of\nhydraulic strategies, Fabaceae as a whole tends to consistently be represented among the most\nresistant species in Amazonian forests.\nFig. 1 (New SM Figure 6) | Embolism resistance variation across sites (Tavares et al., 2023),\nwhich are are ordered left to right as drier (NVX) to wetter (SUC) environments: long dry\nseason length - DSL (NVX, KEN1, KEN2); intermediate DSL (FEC, TAP, CAX, TAM, MAN)\nand ever-wet aseasonal forests (ALP1, ALP2, SUC). Boxplots show the 25th percentile,\nmedian and 75th percentile. Vertical bars show the interquartile range ×1.5 and data points\nbeyond these bars are potential outliers. The dashed horizontal line shows the mean 𝛹 across\n50\nthe pan-Amazonian dataset. Red points represent Fabaceae species.\nComment 1.3. Of the 3 outputs of this paper, none seem very robust.\n1. That phylogeny has an influence on P50 is already known\n2. That Fabaceae is particularly tough does not seem particularly well supported considering\nthe enormous diversity of P50 in the family (presumably correlated with the large diversity of\nthe family).\n3. That coarse scale mapping on P50 provides a useful insight into regional-scale\n“vulnerability” seems unlikely or at least demanding of a far more robust ground-truthing.\nResponse: We appreciate the referee’s comments and the opportunity to clarify the major\ncontributions of our study. While previous global studies have demonstrated phylogenetic\nconservatism in Ψ₅₀ (e.g: Laughlin et al. 2023; Sanchez-Martinez et al. 2023), these analyses\nincluded extremely limited representation of Amazonian taxa (please see figure below).\nPublished Data Figures Redacted\nFig. 2 | Geographic distribution of plots used to evaluate in global assessments of embolism\nresistance phylogenetic conservatism. Images from (Laughlin et al. 2023; Sanchez-Martinez et\nal. 2023).\nOur work fills this gap by showing that embolism resistance is also phylogenetically\nconserved within Amazonia, across 129 species spanning the full precipitation gradient.\nPhylogenetic signal often weaken at regional scales compared to global analyses (Swenson et\nal. 2006; Krasnov, Poulin, and Mouillot 2011; Graham et al. 2018; Losos 2008), so\ndemonstrating this pattern within a single biome adds important evidence for the evolutionary\nstructuring of hydraulic traits in the largest and most diverse tropical forest of the world.\nWe agree that Fabaceae spans a broad range of Ψ₅₀ values, reflecting the environmental\nand taxonomic diversity of the family. However, a key finding is that Fabaceae species\nconsistently occupy the resistant end of the spectrum within sites across the entire precipitation\ngradient (see response to Comment 1.2 and Fig. 1). Thus, despite the within-family variation,\nFabaceae are disproportionately represented among the most embolism-resistant species within\nAmazonian forest plots, making their hydraulic strategy ecologically significant.\nWe agree that the IDW-based map should not be interpreted as a predictive “data layer”,\nbut rather as a visual guide to understand large-scale biogeographical patterns. Please see the\ndetailed responses given to referee 3 below.\nReferences:\nGraham, Catherine H., David Storch, Antonin Machac, and Nick Isaac. 2018. “Phylogenetic\nScale in Ecology and Evolution.” Global Ecology and Biogeography: A Journal of\nMacroecology 27 (2): 175–87.\nKrasnov, Boris R., Robert Poulin, and David Mouillot. 2011. “Scale-Dependence of\nPhylogenetic Signal in Ecological Traits of Ectoparasites.” Ecography 34 (1): 114–22.\nLaughlin, Daniel C., Andrew Siefert, Jesse R. Fleri, Shersingh Joseph Tumber-d, William M.\nHammond, Francesco Maria Sabatini, Gabriella Damasceno, et al. 2023. “Rooting Depth\nand Xylem Vulnerability Are Independent Woody Plant Traits Jointly Selected by\nAridity , Seasonality , and Water Table Depth.” https://doi.org/10.1111/nph.19276.\nLosos, Jonathan B. 2008. “Phylogenetic Niche Conservatism, Phylogenetic Signal and the\nRelationship between Phylogenetic Relatedness and Ecological Similarity among\nSpecies.” Ecology Letters 11 (10): 995–1003.\nSanchez-Martinez, Pablo, Maurizio Mencuccini, Raúl García-Valdés, William M. Hammond,\nJosep M. Serra-Diaz, Wen Yong Guo, Ricardo A. Segovia, et al. 2023. “Increased\nHydraulic Risk in Assemblages of Woody Plant Species Predicts Spatial Patterns of\nDrought-Induced Mortality.” Nature Ecology & Evolution.\nhttps://doi.org/10.1038/s41559-023-02180-z.\nSwenson, Nathan G., Brian J. Enquist, Jason Pither, Jill Thompson, and Jess K. Zimmerman.\n2006. “The Problem and Promise of Scale Dependency in Community Phylogenetics.”\nEcology 87 (10): 2418–24.\nComment 1.4. Thus, I feel that a much more cautious approach to the interpretation of these\ndata is warranted. Having said this, I do believe there is promise in such an approach, so\nperhaps the manuscript simply requires a more circumspect rather than categorical viewpoint.\nResponse: We thank the referee for this thoughtful comment. We agree that a cautious\ninterpretation is necessary given the inherent uncertainties of trait extrapolation and spatial\ninterpolation, and we have revised the manuscript to adopt a more circumspect tone throughout\nin the Discussion (lines, 426-428, 431-433, 452-459). Specifically, we now emphasise that our\nfindings highlight patterns and hypotheses emerging from the largest Amazon-wide hydraulic\ntrait dataset to date, rather than categorical statements about the absolute values of embolism\nresistance across the basin. At the same time, we believe that our work demonstrates the\npromise of combining phylogenetic imputation with community inventory data to explore\nhydraulic strategies at macroecological scales. By grounding the analysis in 129 species across\n11 sites spanning the full Amazonian precipitation gradient, and explicitly discussing the\nlimitations of our approach (lines 581-582, 600-612), we aim to provide a robust empirical\nframework that can guide future data collection and modelling efforts.\nReviewer #2 (Remarks to the Author):\nComment 2.1. Dear Authors,\nI found authors revised the paper throughout in response to reviewers' comments and\nquestions inclusive of my own.\nAdded information in review responses added further information. Revised paper presents\nclear and new insights from the region-wide phylogenetic and geographic analysis focusing\non stem embolism tolerance are worth to be published. Reviewer 2\nResponse: We thank the referee for their positive assessment and encouraging feedback. We\nare pleased that the revisions and additional information clarified our approach. We greatly\nappreciate the constructive comments provided throughout the review process, which have\nsignificantly strengthened the manuscript.\nReviewer #3 (Remarks to the Author):\nComment 3.1. I appreciate the clarification the authors have provided regarding the\nsubsetting of data and the comparison of Fabaceaes relative to other families. The additional\nanalysis and clarification of the text address my key concerns in this area, providing\nconvincing evidence that Fabacea is not simply a uniquely widespread drought-tolerant\nspecialist. I likewise appreciate the exclusion of specific forest types in the spatial projection,\nas well as the clustered cross-validation for overall model fit. Though, for this last point, what\nis the R2 of observed vs. predicted rather than the RMSE? Depending on how the data are\nnormalised, an RMSE of 0.24 could equate to an R2 as low as RMSE^2 = 0.05 (see Li 2017\nfor a discussion on this point https://doi.org/10.1371/journal.pone.0183250).\nResponse: Thank you for your comments. We do not understand however why the reviewer is\nasking to see the R2 values, when the discussion in the suggested reference clearly states that\nr and r2 are not appropriate measures of prediction accuracy when comparing observed vs\npredicted values. As the authors state in the suggested reference: “The r and r2 do not measure\nthe accuracy and are incorrect accuracy measures. The existing error measures suffer\nlimitations. VEcv and E1 are recommended for assessing the accuracy. The applications of\nthese accuracy measures would encourage accuracy-improved predictive models to be\ndeveloped to generate predictions for evidence-informed decision-making”. We have thus\ncalculated VEcv (Variance Explained) values as per the recommendation of the suggested\npaper, and present these results below.\nComment 3.2. With this in mind, there is still inadequate assessment of the uncertainty in the\nspatial projections. The authors state that they are unable to do traditional spatial\nuncertainty/extrapolation analyses because \"the spatial distribution of the points is far from\nideal to support quantitative modelling,\" and that, \"with our clustered point distribution, no\nacceptable semivariance model could be fit.\" But this is exactly the issue -- if the data do not\nallow for appropriate statistical analysis, then how is the analysis justified?\nComment 3.3. The authors state that, \"our intent was not to quantitatively and accurately\npredict Ψ₅₀ values spatially across the basin...\" and that determining the extent of spatial\nextrapolation is, \"beyond the scope of this paper...\" Yet they are visually depicting these\nspatial patterns and referencing these patterns in the title, so one would hope that doing so\n\"accurately\" is an aim.\nResponse to 3.2 and 3.3: We apologise if our statement was not clear before. Our data is\nclustered spatially, which prevents the use of the specific statistical method of kriging, as\nsuggested before. Other methods of quantitative modelling (e.g. regression, random forests or\nother types of continuous predictive modelling based on other explanatory variables will be\nfurther limited by the mismatch between the spatial resolution of environmental layers (often\n1km pixels or more) and the ecological scale that would be relevant for such modelling. Tree\nphysiology is strongly determined by microclimate and microtopography, usually varying on\na scale of a few meters to a few hundreds of meters (https://doi.org/10.1088/1748-\n9326/ad0064). This local scale variability is not well captured by existing environmental layers,\nwhich justifies why such modelling approaches are yet to be used successfully for predicting\nΨ₅₀.\nHowever, we are not attempting to generate accurate spatial predictions of Ψ₅₀ for every\nmapped location/pixel (i.e. a Ψ₅₀ ‘data layer’). We are using mathematical (not statistical)\nspatial interpolation methods to further visualise the broad biogeographical patterns in the\ndistribution of Ψ₅₀ across the Amazon basin, already shown by the data. This is why we have\nopted to use Inverse Distance Weighted interpolation, a simpler and recognised mathematical\nmethod of spatial interpolation that makes no statistical assumptions, unlike kriging. IDW aids\nin visualisation of the general spatial patterns as we intend. And for this very reason, we have\nrefrained from showing actual Ψ₅₀ values on the legend of Figure 3A, and only provide actual\nvalues on Figure 3B (which shows actual plot-level values, not interpolated values.\nThere is a large difference between predicting Ψ₅₀ values accurately at each location\n(i.e. creating a ‘data layer’) and using spatial interpolation methods to visualise broad spatial\npatterns and regions with higher and lower Ψ₅₀ values. For example, Bergmann’s rule is a very\nwell known and recognised biogeographical spatial pattern where endotherms tend to have\nlarger body mass towards the poles and less body mass towards the tropics. If we interpolated\npoint data on body mass across the globe using IDW we would be able to visualise this pattern\n- yet no one would expect Bergmann’s rule to accurately predict body size at a given location\nfor any taxa (https://doi.org/10.1111/gcb.16860). In the same vein, all we are affirming in our\nmanuscript is that some regions in the Amazon tend to have larger Ψ₅₀ values than others - and\nIDW helps us better visualise these patterns than plotting actual values or the clustering results\nalone - but we cannot with our data predict actual Ψ₅₀ values with enough accuracy to\nconfidently present them as a ‘spatial data layer’. Still, we strongly believe evidencing these\nbroad patterns provides a valuable contribution to understanding how biogeography plays a\nrole in the vulnerability of Amazonian forests to more frequent and stronger droughts. We have\nfurther emphasised the nature of the analysis as a visualisation only on lines 581-582 and 610-\n612.\nWith that in mind, we have performed additional uncertainty assessment as requested\nby the reviewer and emphasised by the associate editor, but we would ask you to please keep\nthe above in mind as we reply to the comments. It is expected that validation scores based on\npoint-based observed vs ‘predicted’ (not really a prediction, as stated above) scores will be low\n- especially as we spatially smooth the IDW results for better visualization. therefore reducing\nthe range of Ψ₅₀ in the interpolated map. However, we also provide further evidence that the\nvisual pattern of Ψ₅₀ distribution is correct (see answer below).\nComment 3.4. These calculations (the MESS statistics, spatial uncertainty via boostrapping,\netc) are relatively easy to do with so few points, and would provide important context for\nunderstanding which regions are being estimated with higher/lower uncertainty.\nResponse: In the previous review we assessed uncertainty via k-fold cross-validation using\nenvironmental-spatial clustering to define the folds for training and test samples. We used this\nstrategy instead of bootstrapping as it controls for the spatial autocorrelation between training\nand test data - at each fold the validation data is from a different geographical region than the\ntraining data, ensuring independence between training and testing. For this review, we included\nthe calculation of Variance Explained (VEcv) as per the reference suggested by the reviewer.\nWe have also calculated the MESS statistic as requested, using the MCWD and Water Table\nDepth layers as environmental variables and extracting their corresponding values for each of\nthe 448 plots as the reference values for MESS calculation. We have modified the text on lines\n600 to 638 to further describe these analysis and the validation results, and included three\nadditional figures as Supplementary Material to further support it.\nOverall, the the IDW interpolation reproduces the patterns already revealed by\nclustering, with more negative Ψ₅₀ values for the northern and southeastern Amazon. This can\nbe seen on Supplementary Material Figure 7, reproduced below:\nFig. 3 (SM Figure 8) | Top: results of the geographically constrained clustering of plots based\non the combination of 𝛹 , MCWD (Mean Cumulative Water Deficit) and WTD (Water Table\nAmazonia sensu latissimo (see methods). Areas with very different edaphic environments -\nflooded forests, white sand forests and deforested areas – are masked from the background.\nBottom: Distribution of 𝛹 values within each spatial cluster..\nagreement with observed data (for the reasons explained in the response to comments 3.2 and\n3.3), as seen on the new Supplementary Material Figure 12, below:\nFig. 4 (New SM Figure 13) | Observed versus interpolated (predicted) 𝛹 values based on\n50\nfive-fold cross-validation, showing how the smoothing and averaging that results from the IDW\ninterpolation method used reduced the range of 𝛹 values to between -2 MPa and -2.4 MPa,\n50\nversus -1.3 MPa to -3.8 MPa in the observed data, and explaining the low quantitative\nagreement at the local (point/pixel) level.\nRoot Mean Square Error for the k-fold crossvalidation was of 0.24 MPa, with the value for\neach fold varying from 0.19 MPa to 0.32 MPa. We can see from the plot above that most of\nthe mismatch between observed and ‘predicted’ values comes from the smoothing of\n‘predicted’ values that is part of the IDW interpolation. While observed data ranged from -1.3\nMPa to -3.8 MPa, the ‘predicted’ data only ranged from -2 MPa to -2.4 MPa. For this reason,\nthe calculated Variance Explained (VEcv) was only 6.7%.\nWhile this value would be considered very low for a quantitative predictive model, we\nemphasize again that our use of IDW is just to provide a better and more continuous visual\nrepresentation of the broad scale Ψ₅₀ patterns already revealed by clustering. To demonstrate\nthat, we have first calculated the mean Ψ₅₀ value for each spatial cluster using the source point\ndata (point averages). We then delineated the minimum convex hull enclosing all data points\nin each cluster, and used these polygons to extract all Ψ₅₀ values from the final interpolated\nraster, which were averaged to provide an second mean Ψ₅₀ value calculated from the\ninterpolated values (interpolated averages). If the IDW interpolation is indeed reproducing the\ngeneral spatial pattern already shown by the clustering, we would expect these two averages to\nagree well. And this is what we in fact observe (new Figure SM 13, reproduced below).\nFig. 5 (New SM Figure 14) | Agreement between average 𝛹 values calculated from the field-\n50\nbased estimates, for each spatial cluster (SM Figure 7), and average 𝛹 values calculated from\n50\ninterpolated data, for the area determined by the minimum convex hull enclosing each spatial\ncluster.\nThe MESS analysis did not reveal any strong environmental disagreement between the\nreference locations and the full interpolation domain (Amazonia latissimo sensu). Other than\nisolated pixels, only the highland Andean region and a few regions in lowland western\nAmazonia were outside of the range of the reference data (i.e. MESS values below 0 or above\n100, new Figure SM 14, reproduced below). Overall, 93% of the entire area had MESS values\nwithin the 0 to 100 range, further indicating environmental conditions were mostly within the\nreference range. MESS values ranged from -846.5 to 95.7, with a median of 13.6 and an\naverage of 15.6.\nFig. 5 (New SM Figure 15) | Multivariate Environmental Similarity Surfaces (MESS) analysis\nof the interpolated surface, using Maximum Cumulative Water Depth and Water Table Depth\nlayers as environmental variables, and their corresponding values for each of the 448 plots as\nthe reference values for MESS calculation.\nComment 3.5. The justification the authors provide for their approach is Mitchard et al.\n(2014). But in the decade since this paper was published, there has been growing consensus\nabout the appropriate analysis of spatial patterns and what is necessary when producing\nspatial projections. It's certainly possible that the spatial projections are reasonable, but\nwithout estimates of uncertainty/extrapolation, it's equally likely that the patterns are heavily\nbiased. The uncertainty/extrapolation analysis doesn't need to be incredibly laborious, but\neven some indication of which regions have higher/lower uncertainty would help. If the data\ndo not allow for proper error assessment, or if the spatial projections have a small area of\napplicability or low goodness of fit, then I'm not sure that a spatial analysis is feasible given\nthe current data.\nResponse: Please see our responses above clarifying the nature of our analysis and the\nexpanded analysis of uncertainty / error assessment as requested."}, {"round": 3, "reviewer_comments": {"reviewer_1": "Authors have responded to my previous comments with some justifications in their rebuttal, but without making any substantive changes to the text. My suggestion that there should be more consideration of potential weaknesses in the proposed arguments, has been ignored in the text. This is OK I suppose, but it does seem unfortunate to make bold statements about climate vulnerability, when it is possible that P50 patterns may more strongly reflect phylogenetic patterns than tolerance of drought, Consider the case if Fabaceae tend to be more shallow rooted than other families, their more negative P50 will not correspond to a greater capacity to withstand drought, but rather a compensation for their shallow roots. Furthermore, if the presence if Fabaceae exerts a major influence on community “vulnerability”, then this may not reflect a sensitivity of the community to rainfall shortage in the way that is being sold here. This is a theoretical scenario, but it is based on classical observations in the field, and by the general patterns observed in Choat et al 2012 (Nature).", "reviewer_3": "I appreciate the substantial effort the authors have invested in the additional spatial modelling and the clarity of analysis. I have no critiques with the overall question or analysis, nor do I have an issue with showing a spatial map per se, but I continue to remain skeptical about the robustness and usefulness of the current mapping approach.\nAs stated previously, my main concern is extrapolation beyond data support. Large parts of the central basin lack observations; thus the maps necessarily predict into unsampled regions. Despite the authors' reframing this as simply a \"visualization\", this is misleading---it is still an extrapolation into regions where there are no data---and needs to be evaluated as such. Ultimately, we need a way to distinguish reliable spatial projections from poor ones, or at least provide the readers with the necessary information for them to judge on their own.\nRegarding predictive performance, I fully agree with the authors that spatial projections need not (and rarely do) deliver precise point estimates, but they must be approximately unbiased (i.e., track the 1:1 line) even if they have high point-to-point variance. That is, the discussion invoking Bergmann’s rule appears to conflate inference about the mean trend with the width of prediction intervals; wide prediction intervals (low predictive accuracy at the point level) do not justify biased trends on aggregate. The new cross-validation (SM Fig. 13) suggests negligible predictive accuracy with substantial bias. And while SM Fig. 14 indicates good within-fit performance on aggregate, to justify extrapolating trends beyond clustered regions, this likewise needs to be conducted out-of-fit. However, even if the authors were to aggregate the clusters in SM Fig. 13 and plot obs vs. pred for each of the 5 folds, it appears that there would be no real signal, suggesting that even on aggregate, the Kriging approach does not yield approximately unbiased results.\nThe MESS analysis is a useful direction, but the variable choice needs clarification. Much of the text emphasizes Mean Cumulative Water Deficit and Water Table Depth, whereas the MESS uses “Maximum Cumulative Water Depth.” If this is not a typo, please justify the choice and the exclusion of other plausible drivers (e.g., temperature, precipitation/BioClim, topography). Because MESS depends on the included variables, restricting to two requires a rationale; alternatively, show that these two suffice by demonstrating relatively unbiased spatial CV performance with them.\nGiven the current CV results, I remain hesitant about including Fig. 3 as a basin-wide product. It remains very possible that it is either misleading or essentially random outside of the sampled areas. If the editor prefers to treat the figure as a qualitative visualization rather than a quantitatively accurate projection, I’m happy to defer to that framing. But in that case, I think more caveats would be required, and perhaps adding a much coarser categorical classification (high/med/low), which presumably can be predicted with greater accuracy."}, "author_response": "Responses to referees\nManuscript “Family imprint reveals basin-wide patterns of\nAmazon forest embolism resistance”\nTavares et al.\nTracking Number NCOMMS-24-52403B\nReviewers Comments\nReviewer #1 (Remarks to the Author):\nComment 1.1. Authors have responded to my previous comments with some justifications in\ntheir rebuttal, but without making any substantive changes to the text. My suggestion that\nthere should be more consideration of potential weaknesses in the proposed arguments, has\nbeen ignored in the text. This is OK I suppose, but it does seem unfortunate to make bold\nstatements about climate vulnerability, when it is possible that P50 patterns may more\nstrongly reflect phylogenetic patterns than tolerance of drought, Consider the case if\nFabaceae tend to be more shallow rooted than other families, their more negative P50 will\nnot correspond to a greater capacity to withstand drought, but rather a compensation for\ntheir shallow roots. Furthermore, if the presence if Fabaceae exerts a major influence on\ncommunity “vulnerability”, then this may not reflect a sensitivity of the community to rainfall\nshortage in the way that is being sold here. This is a theoretical scenario, but it is based on\nclassical observations in the field, and by the general patterns observed in Choat et al 2012\n(Nature).\nResponse: We thank the referee for this thoughtful comment and for highlighting the\nimportance of carefully distinguishing between different mechanisms underlying drought\nresponse. We agree that a more cautious interpretation is warranted, and we have revised\nthe manuscript to adopt a less categorical tone throughout the Results and Discussion.\nRegarding the referee's theoretical scenario: “if Fabaceae tend to be more shallow\nrooted than other families, their more negative Ψ will not correspond to a greater capacity\n50\nto withstand drought, but rather a compensation for their shallow roots” - we respectfully\nsuggest that this scenario represents a drought avoidance vs. tolerance distinction rather\nthan a contradiction of our interpretation (Santiago et al., 2018; Laughlin et al., 2023).\nDeep-rooted species rely on drought avoidance, while shallow-rooted species likely require\nmore embolism-resistant xylem to tolerate low water potentials (Brum et al., 2018; Oliveira et\nal., 2021). These strategies are not mutually exclusive and do not undermine the ecological\nmeaning of lower Ψ values, instead, they highlight that hydraulic traits must be interpreted\n50\nwithin a broader whole-plant framework. That said, we fully agree with the referee that Ψ\n50\nalone cannot capture the full complexity of drought resistance (Choat et al., 2018). In\nresponse, we have now strengthened the manuscript by adding explicit caveats about the\nlimitations of interpreting Ψ as a direct proxy for drought vulnerability and clarifying that our\n50\nconclusions concern basin-wide patterns of embolism resistance, not absolute drought\nresistance across all possible mechanisms.\nWe agree that community-level patterns of Ψ may partially reflect phylogenetic\n50\ncomposition. However, this does not preclude ecological interpretation, instead, it\nemphasises the central role of evolutionary history in shaping hydraulic strategies. To make\nthis clear, we now explicitly state that our spatial patterns represent embolism resistance\ndistributions, not direct predictions of drought mortality.\nFinally, we appreciate the referee’s recognition that our approach holds promise, as\nnoted in their previous round of comments. In the revised manuscript, we now explicitly\nemphasise that our results identify emergent macroecological patterns and hypotheses,\nrather than providing definitive predictions of climate-change vulnerability. We agree that\nadopting a more circumspect tone strengthens the paper. By grounding our analysis in an\nAmazon-wide empirical dataset and openly acknowledging its limitations, we aim to offer a\nrobust foundation for future data collection and modelling efforts, as well as for integrative\nresearch linking ecophysiology, demography, and forest dynamics. Please, the changes in\nthe revised manuscript are marked as yellow.\nReferences:\nSantiago, Louis S., et al. \"Coordination and trade‐offs among hydraulic safety, efficiency and\ndrought avoidance traits in Amazonian rainforest canopy tree species.\" New Phytologist\n218.3 (2018): 1015-1024.\nLaughlin, Daniel C., et al. \"Rooting depth and xylem vulnerability are independent woody\nplant traits jointly selected by aridity, seasonality, and water table depth.\" New Phytologist\n240.5 (2023): 1774-1787.\nBrum, Mauro, et al. \"Hydrological niche segregation defines forest structure and drought\ntolerance strategies in a seasonal Amazon forest.\" Journal of Ecology 107.1 (2019):\n318-333.\nOliveira, Rafael S., et al. \"Linking plant hydraulics and the fast–slow continuum to\nunderstand resilience to drought in tropical ecosystems.\" New Phytologist 230.3 (2021):\n904-923.\nChoat, Brendan, et al. \"Triggers of tree mortality under drought.\" Nature 558.7711 (2018):\n531-539.\nReviewer #3 (Remarks to the Author):\nComment 2.1. I appreciate the substantial effort the authors have invested in the additional\nspatial modelling and the clarity of analysis. I have no critiques with the overall question or\nanalysis, nor do I have an issue with showing a spatial map per se, but I continue to remain\nskeptical about the robustness and usefulness of the current mapping approach.\nResponse: We thank the referee and acknowledge their personal position. We feel that after\naddressing their concerns and including the additional analyses requested, the manuscript\nnow offers enough additional information and support to allow readers to reach their own\nconclusions (as requested below) about our general spatial patterns. Moreover, we added\nadditional caveats regarding the analysis and figure in three different sections of the\nmanuscript, in response to editor requests - please see the attached Author Checklist and\nthe changes in the revised manuscript, which are marked as yellow.\nComment 2.2. As stated previously, my main concern is extrapolation beyond data support.\nLarge parts of the central basin lack observations; thus the maps necessarily predict into\nunsampled regions. Despite the authors' reframing this as simply a \"visualization\", this is\nmisleading---it is still an extrapolation into regions where there are no data---and needs to be\nevaluated as such. Ultimately, we need a way to distinguish reliable spatial projections from\npoor ones, or at least provide the readers with the necessary information for them to judge\non their own.\nResponse: We understand the referee position but still disagree fundamentally with the\nperception that there is no distinction between visualising/spatializing general patterns and\nproviding accurate and fully quantitative spatial predictions. We agree that our work\nextrapolates into regions where there are no data - this is the fundamental aspect of the\npaper. But we believe at this point we have made it clear that while we are confident we are\ncapturing broad biogeographic/macroecological patterns, there is still much work to be done\nin further confirming and ‘unpacking’ these patterns. We have refrained since the first\nversion of the paper to make any quantitative claims or support any quantitative inference\nfrom our spatialization, exactly because we agree with the referee that it could be misleading\n- and it is in fact why although we report the code used to generate the data in the figure (as\nrequired by the journal), we are not making the actual results available as a ‘ data layer’. It\ndoes not have the required level of quantitative accuracy we would deem necessary to\nrelease it as such. But the visualization and discussion of the broad biogeographic patterns\nis still valid, and very relevant and important to drive forward the discussion about the\nvulnerability and the future of the Amazon biome. As stated by the referee, we are confident\nwe now provide all the necessary information for readers to judge on their own.\nComment 2.3. Regarding predictive performance, I fully agree with the authors that spatial\nprojections need not (and rarely do) deliver precise point estimates, but they must be\napproximately unbiased (i.e., track the 1:1 line) even if they have high point-to-point\nvariance. That is, the discussion invoking Bergmann’s rule appears to conflate inference\nabout the mean trend with the width of prediction intervals; wide prediction intervals (low\npredictive accuracy at the point level) do not justify biased trends on aggregate. The new\ncross-validation (SM Fig. 13) suggests negligible predictive accuracy with substantial bias.\nAnd while SM Fig. 14 indicates good within-fit performance on aggregate, to justify\nextrapolating trends beyond clustered regions, this likewise needs to be conducted out-of-fit.\nHowever, even if the authors were to aggregate the clusters in SM Fig. 13 and plot obs vs.\npred for each of the 5 folds, it appears that there would be no real signal, suggesting that\neven on aggregate, the Kriging approach does not yield approximately unbiased results.\nResponse: First, we must reiterate we did not use a Kriging approach in this work, exactly\nbecause it requires a strong spatial autocorrelation model that is not present in our dataset\ndue to geographic sampling bias. Our first version of this visualisation included only the\ninterpolated data points based on distance. Since then, we introduced the clustering analysis\nto our method - addressing the referee’s suggestion that the environment should be better\nrepresented. We also performed the MESS validation as requested to prove that although\nwe are extrapolating to regions with no data, these regions are still within the environmental\nenvelope covered under the point data. We have also already clarified (using Fig SM13) that\nmost of the bias results from the inherent smoothing of the IDW method, which we further\nsmooth to avoid unnatural sharp transitions in values. This bias would have been an issue if\nwe were making any quantitative predictions from the data, but all we are showing is the\nvery broad variations in the underlying averages (i.e. inferences about the mean as stated\nabove). We disagree that there is no real signal, and we believe figure SM 14 clearly shows\nthis - plotting the average observations per cluster versus the average interpolated values\nper cluster show a very close agreement. We thus hope that the multiple additions to the\nmethod we have made during the review process are now sufficient for readers to make their\nown conclusions.\nComment 2.4. The MESS analysis is a useful direction, but the variable choice needs\nclarification. Much of the text emphasizes Mean Cumulative Water Deficit and Water Table\nDepth, whereas the MESS uses “Maximum Cumulative Water Depth.” If this is not a typo,\nplease justify the choice and the exclusion of other plausible drivers (e.g., temperature,\nprecipitation/BioClim, topography). Because MESS depends on the included variables,\nrestricting to two requires a rationale; alternatively, show that these two suffice by\ndemonstrating relatively unbiased spatial CV performance with them.\nResponse: First, we apologise for the typo. It should have said Maximum Cumulative Water\nDeficit on the Figure 15 SM legend and this has now been corrected.\nWe have chosen to use MCWD and WTD as our predictor variables because these\nare two of the variables with the most direct causal relation with plant water stress and\ncauses of embolism. They encompass several basal climatic factors (temperature, rainfall,\nevapotranspiration, VPD), and we have shown previously (Tavares et al. 2023 Nature\n617:11–117) that MCWD has a strong correlation to Ψ values.\nprevious round of review responses, detailed modelling of the environmental controls on Ψ\n50\nwould be a worthy scientific problem but is outside the scope and purpose of this paper.\nAgain, our aim is to better understand spatial patterns in Ψ , and not to quantitatively\n50\nestimate Ψ on a per-pixel basis across the Basin. We therefore believe we are using MESS\n50\nas intended, to demonstrate that as far as WTD and MCWD are concerned, the plot\nlocations from where data was interpolated (not statistically predicted) do cover the entire\nexpected range of these variables in the Amazon. We feel at this point there is not much\nroom left for further rebuttal, and that the referee and ourselves have very different\nperceptions and expectations - we could continue this debate for a long time. As stated by\nthe referees themselves, the most important issue is that our analysis is transparent and well\nsupported enough for its merit to be judged individually by readers, and we are confident we\nhave reached this stage now after all the requested reviews.\nComment 2.5. Given the current CV results, I remain hesitant about including Fig. 3 as a\nbasin-wide product. It remains very possible that it is either misleading or essentially random\noutside of the sampled areas. If the editor prefers to treat the figure as a qualitative\nvisualization rather than a quantitatively accurate projection, I’m happy to defer to that\nframing. But in that case, I think more caveats would be required, and perhaps adding a\nmuch coarser categorical classification (high/med/low), which presumably can be predicted\nwith greater accuracy.\nResponse: we appreciate the referee’s stance, but we disagree that it is misleading or\nessentially random outside of the sample areas. As the Ψ values approximated for the\n50\nunsampled areas are entirely driven by a distance-decay function that averages across the\nexisting samples - so by design correlated with the existing observations, it is not\nmathematically possible for these to be random. Moreover, we also have further confirmed\nthe patterns using two independent methods: interpolation and clustering. We reiterate that\nwe are not positioning Figure 3 as a data product of any kind, but rather presenting it as a\nvisual aid to better show the patterns already present in the point data (which can be\nassessed by looking at Figure 3B, albeit not so immediately clear and readable - hence the\nimproved visualization of Figure 3A). We have answered all previous and current comments\nof the referee and have on their request modified our method and included several metrics of\nmodel performance and generality. We believe at this point there is sufficient information in\nthe document and supplementary material for readers to come to their own conclusions\nabout our results, and to revisit/rebut them with new data to further advance our knowledge\non the vulnerability of Amazonian forests to ongoing climatic changes."}]

Decision letter