Climate change threatens native potential agroforestry plant species in Brazil

Climate change is one of the main drivers of species extinction in the twentyfirst-century. Here, we (1) quantify potential changes in species' bioclimatic area of habitat (BAH) of 135 native potential agroforestry species from the Brazilian flora, using two different climate change scenarios (SSP2-4.5 and SSP5-8.5) and dispersal scenarios, where species have no ability to disperse and reach new areas (non-dispersal) and where species can migrate within the estimated BAH (full dispersal) for 2041–2060 and 2061–2080. We then (2) assess the preliminary conservation status of each species based on IUCN criteria. Current and future potential habitats for species were predicted using MaxEnt, a machine-learning algorithm used to estimate species' probability distribution. Future climate is predicted to trigger a mean decline in BAH between 38.5–56.3% under the non-dispersal scenario and between 22.3–41.9% under the full dispersal scenario for 135 native potential agroforestry species. Additionally, we found that only 4.3% of the studied species could be threatened under the IUCN Red List criteria B1 and B2. However, when considering the predicted quantitative habitat loss due to climate change (A3c criterion) the percentages increased between 68.8–84.4% under the non-dispersal scenario and between 40.7–64.4% under the full dispersal scenario. To lessen such threats, we argue that encouraging the use of these species in rural and peri-urban agroecosystems are promising, complementary strategies for their long-term conservation.

Hundreds of species are unquestionably promising for future human welfare 1 , yet a large number of these species may potentially be threatened by impacts of climate change in the coming decades. Climate change is one of the important drivers affecting species survival, causing global biodiversity loss [2][3][4][5] . Climate change affects species in different ways, such as altering the suitability of current habitat of species, resulting in accelerated extinction rates 2,6,7 . The United Nations Intergovernmental Panel on Climate Change (IPCC) estimates that if Earth's average temperature rises between 2 and 3 °C, about 20 to 30% of all terrestrial biodiversity will be at high risk of extinction by the end of the century 8 . In the last century, land and ocean temperature showed a warming of approximately 1.0 °C 3 , which may increase another 1.4 °C to 5.0 °C by 2100, if we do not reduce greenhouse gas emissions 9 . The latest IPCC special report reinforces the importance of keeping the temperature increase below 1.5 °C, in order to keep negative effects on natural resources, ecosystem functioning, food security and biodiversity to a minimum 10 . Considering the potential climate change scenarios with additional temperature increases, both widespread species and narrow-ranged endemic species will likely suffer irreparable consequences with regard to their distribution range and abundance 5 .
Brazil is the world's most biodiversity-rich country, with 33,161 known species of vascular plants 11 , and harboring some of the largest remnants of tropical old-growth forests 12 . Despite the large number of native plant species, many of which with major untapped socioeconomic potential, the Brazilian agricultural industry exploits only a few, and largely, exotic crops 13 . Agroforestry species are those that have the function of simultaneously benefiting people's livelihoods and the ecological systems while showing great potential for multispecies intercropping 14 . These species are often characterized by their multiple uses, different harvest seasons and potential for market adoption [15][16][17][18]  www.nature.com/scientificreports/ pasture production, silviculture, orchards, bioenergy, green manuring, as well as in integrated, biodiverse, multifunctional agroforestry 13,19 . Besides enhancing biodiversity and promoting the socio-economic development of local communities [20][21][22][23] , agroforestry systems can play a pivotal role in mitigating the effects of climate change: they sequester more atmospheric carbon than conventional farming [24][25][26] . Although agroforestry practices can ameliorate the impacts of climate change in Brazil, these agroecological systems are also vulnerable 27 . Considering the rapidly increasing human demand for plant products, native plant species from megadiverse countries undoubtedly represent a reservoir of genetic diversity, providing beneficial alleles for crop improvement and higher adaptive potential to face global changes 13,28 . Changes in land use may not be the main driver impacting these species as they are widely distributed among the neotropics and are easily found along streets and city squares across Brazil, some almost ruderal, regenerating in open areas of cities 13,29 . The future impact of climate change on species distributions should be taken into account for setting conservation priorities, as well as for promoting species conservation through their sustainable use 30 . Spatial and temporal changes of species' suitable habitat can be predicted with ecological niche models (ENMs), the most widely used tool to assess species vulnerability to changing climatic conditions [36][37][38] . Besides that, modeled habitats based on climatic variables allow us to consider the impacts of climate change on the species' area of habitat, which is the habitat available to a particular species within its range 39 . Here, to consider those impacts, we modeled the species' bioclimatic area of habitat (BAH). Species dispersal is pivotal to the survival of species in the face of rapid climate change 40 . Thus, to better understand species responses, this central process that determines the potential spread of a population needs to be addressed in conservation assessments 7 . Although several studies have sought to better understand the impact of climate change on the distribution of plant species with narrow-ranged distribution or threatened with extinction, we note that no study has yet focused on species of agroforestry interest in Brazil, which generally have widespread distribution. These species are promising for conservation-by-use, an approach used by people communities for millennia in different ecosystems in Brazil 41,42 .
Here, we apply an ENM approach to (1) quantify potential changes in BAH of 135 native potential agroforestry species from the Brazilian flora using two climate change scenarios (SSP2-4.5 and SSP5-8.5) and two dispersal (non-dispersal and full dispersal) scenarios for 2041-2060 and 2061-2080. We then (2) assess the preliminary conservation status of each species using IUCN Red List of Threatened Species criteria 43 .

Results
Model performance. We evaluated the model performance through a null-model for significance testing of presence-only ENMs and retained 135 significant ENMs, corresponding to 97.1% of all species. Overall, final models showed high accuracy, indicated by AUC values ranging from 0.850 ± 0.139 to 0.985 ± 0.058, demonstrating a clear ability to distinguish suitable from unsuitable habitats. We detected no spurious correlations through inspection of species-response-curves. Impacts on species BAH. Under the non-dispersal scenario, the average decline in BAH was predicted to be between 38.5% (SSP2-4.5) and 43.5% SSP5-8.5 by 2041-2060 and between 43.4% (SSP2-4.5) and 56.3% (SSP5-8.5) by 2061-2080. For the full dispersal scenario, however, the average decline of BAH was predicted to be between 22.3% (SSP2-4.5) and 29.7% (SSP5-8.5) by 2041-2060 and between 27.4% (SSP2-4.5) and 41.9% (SSP5-8.5) by 2061-2080 (Supplementary Table S1). Although the majority of species predicted BAH losses over different scenarios, some are predicted to experience BAH gains ( Fig. 1., Table 1, Supplementary Table S1). We noticed that some species were predicted to lose their entire BAH by 2041-2060 and 2061-2080, such as the medicinal species Cunila microcephala in the non-dispersal scenario and the forage species Ornithopus micranthus, in all scenarios, except for the SSP5-8.5 (2061SSP5-8.5 ( -2080 in the non-dispersal scenario and SSP2-4.5 and SSP5- 8.5 (2061-2080) in the full dispersal scenario. Species with the greatest increase in BAH were the forage species Indigofera sabulicola with an increase of 388%, the ornamental species Epidendrum fulgens (263%), and the forage species Echinochloa polystachya (259%) in the SSP5-8.5 for 2061-2080 considering the full dispersal scenario (Table 1, Supplementary Table S1).

IUCN Red list preliminary assessment.
Assessing the geographic range (B1a + B2a criteria), we found that only 4.3% of the native species was qualified for a Threat category (Supplementary Table S2). However, when considering the predicted quantitative habitat loss due to climate change (A3c criterion) the percentages increased. We observed that 68  www.nature.com/scientificreports/ Rapid range changes of species may, in turn, impact material, non-material and regulating contributions of nature to people 2,48,57,58 . The global productivity of farms may be negatively affected by climate change throughout most of the tropics 59 . Smallholders and traditional Brazilian communities such as the caiçaras, quilombolas and indigenous peoples make use of these species in agroforestry practices 13 , a land use management system that increases carbon sequestration in biomass and soils 24,60 . These systems improve people's livelihoods by simultaneously providing income, food security, fuel, medicine, forage, and/or other goods and services. Moreover, increased tree cover due to agroforestry may help to mitigate climate change 23,57,61 . Here we present evidence of how currently suitable areas for cultivation of these species may become unsuitable in the future. This may lead to a severe decline in people's livelihoods and regional food security 23,62 . Moreover, under climate shifts species survival is likely to be threatened by biogeographic barriers such as an agricultural or otherwise ecologically degraded landscape matrix which would prevent species migration to climatically more suitable areas 63 .
The outcomes of climate change predicted by biogeographical and ecological studies have been neglected and have barely been integrated into conservation planning 30,64 . Prioritization of conservation efforts is often based on a species' extinction risk 65 . Determining whether a taxon is threatened with extinction depends on biological indicators, such as rapid population decline, and qualifying species in a threat category may assist in decision making 43    www.nature.com/scientificreports/ for threatened species in Brazil, agroforestry systems can act as an alternative to overcome part of the conflicts between conventional agricultural production and the conservation of natural resources, as can be seen in the new forest code (Law 12.651/2012), which provides explicit provisions for sustainable agroforestry. Additionally, a specific Law (12.854/2013) promotes forest recovery activities and the implementation of agroforestry systems 34 . It is well known that deforestation is a primary driver leading species to extinction 67 . However, climate change is expected to overtake this driver in a few decades 44 . The majority of the species prioritized in the "Plants for the future initiative", found mainly in the Atlantic Forest have great potential to be conserved through sustainable practices, particularly by smallholders and traditional communities 13 . The Atlantic Forest remains endangered as a result of continued deforestation and the future of this forest relies on well-structured conservation plans based on reliable information 55,68,69 . The Atlantic Forest cover has been reduced to less than 20% of its original size 70,71 , distributed mainly in small and disturbed fragments of less than 50 hectares 56,72 . Thus, most of the endemic species in the Atlantic Forest biome could already be qualified as critically endangered according to IUCN criteria 43 . On the other hand, despite our predicted catastrophic scenarios for native Atlantic Forest species, we observe that species such as Bracatinga (M. scabrella), Brazilian peppertree (Schinus terebinthifolia)  www.nature.com/scientificreports/ and Peruvian groundcherry (Physalis peruviana), are distributed across the neotropics and easily found along streets and city squares all over Brazil. Many of these are pioneers, some almost ruderal, regenerating easily in open city areas 13,29 . Many of these forest species have endured over 500 years of deforestation, and still remain abundant in the Atlantic Forest even after losing approximately 80% of their natural habitat 29,56,71 . Hence, we note the need for studies that address species response to global changes to better understand the resilience potential of these species. Protecting people's livelihoods in a rapidly changing climate may be one of the great challenges of the twentyfirst century. Although it is not shared by all of the scientific community, as discussed by Loreau 73 , species with economic value seem to have advantages for conservation over those with non-economic value as can be seen in long-term human activities such as protection, transport and planting of useful species and removal of non-useful species by local communities 31,42,74 . Indeed, socioeconomic underutilization of plant resources may in some cases even jeopardize socioecological synergies of tropical forest resilience 58 . All species analyzed here have a great potential to be conserved through a conservation-by-use approach, because of their different uses that do not necessarily jeopardize reproduction and persistence. They play an important role among local communities 41,48,50,75,76 and farmers' livelihoods 31,42 . Searching for evidence of conservation among species with economic and cultural values, Reis et al. 42 noticed that the species Ilex paraguariensis, A. angustifolia and B. antiacantha were intentionally favored through protection, transplantation and selection by farmers. Furthermore, Donazzolo et al. 32 noted that management of Acca sellowiana populations, retained high level of genetic diversity and tended to increase the species genetic variability. We argue that adopting measures, such as the establishment of new agroforestry systems to increase carbon sequestration, the selection of varieties capable of withstanding new climates and the improvement of habitat connectivity to facilitate species migration/dispersal should be a strategy for short-term and long-term conservation and people's livelihoods.
ENMs are widely used to forecast the distribution of species across geographic space and time. Building meaningful models to estimate the future distribution of species for an uncertain future requires very specific decisions and interpretations with extreme caution 51,77,78 . Several uncertainties and complexities are related to our study. Modeling a large number of species can make the species-specific selection of predictors methodologically and practically complex 77,79 . To mitigate this, we selected the most suitable environmental predictors for different plant growth forms. Model performance evaluation is a key step for ENM studies and probably the most problematic one owing to its complexity 79 . The random cross-validation approach is the most common practice, adopted by modelers to evaluate model performance, where datasets are split into k folds, using one part to test the model and the remaining (k-1 folds) to calibrate the model 80,81 . To reduce the over-optimistic nature of cross-validation, we applied a null-model for significance testing of presence-only ENMs 82 . Binarization of continuous probabilities output is commonly employed by modelers to quantify species range changes and build species richness over time 83 . Nevertheless, Santini et al. 51 recently concluded that this practice reduces the predictive probability of models. Although we binarized ENMs outputs to quantify the climate change impacts, we applied a threshold highly indicated for conservation purposes for showing high performance in the identification of suitable areas and commonly used 50,[84][85][86] . We assumed that the species are at equilibrium with the environment 87 and occurrence records were sampled randomly 88 . Furthermore, we included no biotic interactions 40,89 , adaptations and evolution 90 into our modeling approach. For dispersal 37,40 , we considered two scenarios (non-dispersal and full dispersal), but we limited the full dispersal scenario to species BAH, since we understand that plant dispersal rates over 0.1 km/year might not occur for vascular and non-vascular plants 5,64,91 or over 100 km for Brazilian tree species as result of climate change 44,92 . Although humans fundamentally affect dispersal and alter landscapes by transporting individuals 93-95 , we did not include human-mediated dispersal data in our models due to the lack of information related to human migrations as well as for each specific species. Another limitation of our study is to restrict our analysis to the estimated BAH of species, which may mask some macroecological patterns, yet adopting this conservative approach allows us to observe more concise species responses and diminish model overfitting 36,44,96 .
In summary, we showed that future climate will likely trigger a decline in BAH between 38.5-56.3% under the non-dispersal scenario and between 22.3-41.9% under the full dispersal scenario of several native potential agroforestry species from the Brazilian flora. Additionally, we found that only 4.3% of the studied species could be threatened under the IUCN criteria B1 and B2. However, when considering the IUCN criterion A3, 68.8-84.4% (non-dispersal) and 40.7-64.4% (full dispersal) of our species of interest could be qualified as threatened. Although accessing genetic material with quality for native species might be difficult and the scenarios used here estimate considerable losses for 2041-2060 and 2061-2080, we argue that actions such as the promotion of these species in agroecosystems are promising alternatives to increase their population sizes. We urge that public policies involving farmers and local communities be adopted, as practices and management systems implemented by them have proven to maintain landscapes with productive forest fragments, and consequently favors species and forest conservation. Lastly, we highly recommend the development of scientific research towards biotechnological applications to select promising genotypes for a changing global climate. www.nature.com/scientificreports/ Convention on Biological Diversity (CBD) and International Treaty on Plant Genetic Resources for Food and Agriculture, particularly with regard to promoting the sustainable use of biodiversity components 13 , these species provide food security for local communities and have commercial value in national and foreign markets.

Methods
In spite of the fact that not all species have already been found in current agricultural systems or been managed by farmers, they all have one or multiple uses and can be combined in mixed cropping. Some species such as A. angustifolia 31 , A. sellowiana 32 , E. edulis 33 , Ilex paraguariensis 34 , and Mimosa scabrella 35 have already traditionally been combined with other agricultural crops in managed landscapes in southern Brazil. All taxonomic authorities and species common names are listed in Supplementary Table S3.
Species occurrence data. Occurrence data for the 139 species evaluated here was downloaded from the Global Biodiversity Information Facility (https:// doi. org/ 10. 15468/ dl. vjezvb) 97 . We collected a total of 28,860 unique records. The sample size for species ranged from 12 (Ornithopus micrantus) to 5464 (Casearia sylvestris). We standardized botanical names using the R package 'flora' 98 , which uses the nomenclature accepted by the Brazilian Flora 2020 project (http:// flora dobra sil. jbrj. gov. br/). To avoid modeling truncated niches, we extracted all records from an extent, defined by latitudes 60°S-15°N and longitudes 90°-30°W 99 . We checked the geographical consistencies of all records using the cleaning pipeline proposed by Gomes et al. 36 . Firstly, we removed all occurrences outside the Neotropics. Then, we removed all records with missing latitude and longitude, using the function 'cleancoordinates' from the R package 'CoordinateCleaner' 100 . Finally, we estimated the kernel density for each species in order to remove spatial outliers using the density function from the R package 'stats' 101 . As geographic sampling biases are common among biological collections 102,103 , which can lead to overrepresentation of environmental conditions, we spatially filtered the species occurrence data over a distance of 20 km using the R package 'spThin' in order to diminish spatial autocorrelation 104 . We did not model species for which there were less than ten records, as models fit with few data may not be reliable 105 Table S4). These predictors are critical in determining the distribution limits of a wide range of plant growth forms and are highly related to plant physiological responses 109,110 . We then (2) checked for multicollinearity by examining the correlation structure of the predictor variables through the variance inflation factor (VIF) for epiphyte, fern, graminoid, herb, hydrophyte, lithophyte, shrub, tree and vine species 111 . This measure evaluates how much the variance of an estimated regression coefficient increases if their predictors are correlated 112 . We kept only predictors with VIF values below 5 113 . The VIFs were checked using the function 'vifstep' in the R package 'usdm' 114 . The retained predictors are shown in Supplementary Table S5.
Modeling approach. We used bioclimatic habitat suitability to assess the potential impacts of climate change on species' BAH and inform IUCN Red List assessments 115 . To delimit BAH, we incorporated a 100 km buffer around species extent of occurrence (EOO), as we understand that plant dispersal rates over 0.1 km/ year might not occur for vascular and non-vascular plants 5,64,91 , and adopting a conservative approach reduces model overfitting 36,44,96 . EOOs were quantified by drawing a minimum convex polygon (MCP) around known species records as recommend by IUCN 43 . Current and future potential habitats for species were predicted using MaxEnt v.3.4.1 k, a machine-learning algorithm used to estimate species' probability distribution 116 . MaxEnt is the presence-based method widely used for having high performance when compared to other available algorithms 36,105,[117][118][119] . ENMs were fitted using the following parameters in the MaxEnt: bootstrap method with 100 replicates, 500 maximum iterations, 10,000 points of background, and Cloglog output format. We only kept linear and quadratic features to avoid overfitting of the models and as recommended by Merow et al. because of the absence of a biological justification with the variables used 120,121 . Furthermore, we inspected speciesresponse-curves to avoid spurious calibrations, following the evaluation strip method proposed by Elith et al. 122 . This method investigates the effect of one variable at a time, keeping the others constant at their mean values 122 .
To assess robustness and alert policy-makers for the uncertainties typically associated with these methods, each ENM was tested against a bias corrected null-model as proposed by Raes and ter Steege 82 . The AUC values of the ENMs built with n occurrence records were tested against the upper AUC values of the lower quantile of 95% of the AUC values obtained from 100 × n points drawn and predicted randomly. Only significant ENMs were projected to future climatic conditions. Future projections. The climate projections were carried out according to the Sixth Assessment Report (AR6) of the IPCC, using two Shared Socioeconomic Pathways (SSPs) as reference (SSP2-4.5 and SSP5-8.5). SSPs are projections of future climates, based on different socioeconomic assumptions such as population, technological, and economic growth. The SSP2-4.5 (2041-2060) is a stabilization scenario, assuming global temperature increases ranging from 2.1 to 4.3 °C and mean warming of 3.0 °C. The SSP5-8.5 (2061-2080) represents the worst-case scenario, assuming absence of climate change policies, with global temperatures continue to rise throughout the twenty-first century, with estimates ranging from 3.8 to 7.4 °C and mean warming of 5.0 °C 9 . We averaged eight different global climate models: BCC-CSM2-MR, CNRM-CM6-1, CNRM-ESM2-1, CanESM5, IPSL-CM6A-LR, MIROC-ES2L, MIROC6 and MRI-ESM2-0 to take into account the uncertainties related to future climate conditions 81 . The fitted consensus ENMs were projected to these two datasets to obtain predicted future maps of habitat suitability for each species. www.nature.com/scientificreports/ To map changes in future ranges of species, we converted the continuous habitat suitability into binaries using the maximum training sensitivity plus specificity threshold 85,123 . This threshold is indicated for conservation purposes for its high performance in the identification of suitable areas 50,[84][85][86] . To assess whether species will face a decline or expansion in BAH under future climate conditions, we quantified the difference between the relative number of pixels occupied in current and future BAHs. We assumed that species have no dispersal capacity and full dispersal capacity for 2041-2060 and 2061-2080 timeframes. In the first scenario, species do not have the ability to disperse and reach new areas (pixels) in future climate scenarios. In the second, species can migrate within the estimated BAH of each species. Here, to ensure transparency and reproducibility for reporting ENMs, we adhered to the ODMAP (Overview, Data, Model, Assessment, Prediction) protocol v1.0. (Supplementary  Table S6), as proposed by Zurell et al. 78 . All analyses were conducted within R environment version R 4.1.1 101 .
IUCN Red list preliminary assessment. The IUCN Red List assessments provide important information related to species status, trends and threats for the establishment of conservation planning and improvement of decision-making 43,124,125 . Criterion B is linked to the geographic range and has two sub-criteria (B1 and B2), which are based on the extent of occurrence (EOO) and B2 on the area of occupancy (AOO), respectively 43 . Further, three other conditions (a, b, and c) describe aspects of the biology and potential decline of the taxon in response to threats 43 . At least one sub-criterion and two conditions must be met to qualify a given species as threatened 43 . Following the guidelines for using the IUCN Red List Categories and Criteria version 14, we calculated the geographic range (B1a + B2a criteria) using the R package 'ConR' 126 . Additionally, we evaluated predicted quantitative habitat loss due to climate change by assessing the decline in habitat quality (A3c criterion), suspected to be met in the future, to qualify whether a particular species would be in a threat category 43 . The classification of threat includes: Vulnerable (VU), Endangered (EN) and Critically Endangered (CR). To qualify a species as Vulnerable, qualitative habitat loss must be ≥ 30%, Endangered ≥ 50% and Critically Endangered ≥ 80%. These categories are related to the risk of extinction of species in the wild 43 .