Possible extinction of Eversmannia subspinosa in Eastern Alborz by 2060 due to climate change: A MaxEnt study

Climate change has altered ecosystems by affecting the suitability of environments for different species. Species distribution models provide insights regarding these changes, which could be valuable in conservation planning. In this study, we model the current and projected future distribution of Eversmannia subspinosa, a shrub from the Fabaceae family, across the Eastern Alborz area. This plant is endangered in Iran and is only found in the study area.We used the Maximum Entropy (MaxEnt) method and trained the model with the least correlated and most important bioclimatic and topographic variables. CMCC-ESM2 and HadGEM3 climate models, in conjunction with IPCC Representative Concentration Pathways RCP2.6, RCP4.5, and RCP8.5, were used to predict the distribution of the species in 2040 and 2060. Our results suggest that the current habitats of E. subspinosa are mainly in southern parts of Eastern Alborz. Elevation (54.9%) and annual mean temperature (24.5%) were the main contributors to the model. CMCC-ESM2 model predicts signi�cant habitat loss up to 2040 and almost complete disappearance of high probable habitats (0.5 ≤ probability of presence) up to 2060. It also suggests a distribution shift toward higher elevations in Northern and North-Eastern areas of the study area. The model trained by HadGEM3 predicts extinction of E. subspinosa across the study area up to 2040. Filling the gaps between protected areas and national parks and limiting road construction, which blocks its migration to higher elevations in response to global warming, could delay the extinction of this rare species in Iran.


Introduction
Climate change has been identi ed as a signi cant driver of shifts in species distribution.The ongoing changes in climate have led to alterations in the habitats and geographic ranges of numerous species (Bonebrake et al., 2018; A. Dhyani et al., 2021;Gervais et al., 2021;Parmesan et al., 2015;Warren et al., 2018).Recent studies have further substantiated these ndings.While marine species are shifting poleward in response to ocean warming (Sunday et al., 2012), terrestrial species are undergoing range shifts at varying rates, depending on their speci c thermal tolerances and habitat requirements (Pecl et al., 2017).Many plant species are moving uphill and towards the poles in response to rising temperatures.This uphill and poleward migration is a survival strategy, as plants seek cooler temperatures and more suitable conditions for growth (Chen et al., 2011;Couet et al., 2022;Lenoir et al., 2015;Zu et al., 2021).Moreover, these shifts in species distribution have signi cant implications for biodiversity and ecosystem functioning.As species move to new areas, they can alter the composition and structure of communities, potentially leading to changes in ecosystem processes and services (Bellard et al. 2012;Sharma et al. 2021).
Species distribution modeling (SDM) is a powerful tool used in ecology to predict the geographic distribution of species based on environmental conditions at locations where they are known to occur (Elith et al., 2009).SDMs are based on the fundamental ecological principle that organisms are adapted to speci c environmental conditions and thus their geographic distributions can be predicted by identifying these conditions (Anderson et al., 2011).The models use presence-only or presence-absence data of species along with environmental predictor variables (such as temperature, precipitation and altitude) to predict suitable habitats for the species (Elith et al., 2009).These models can also forecast potential range shifts of species under different climate change scenarios, providing valuable information for conservation planning (Elith et al., 2006).
In this study, we focused on modeling the current and projecting the future distribution of Eversmannia subspinosa (Fisch.Ex DC.) B. Fedtsch, (Fabaceae Hedysareae), an endemic species distributed in Iran, Afghanistan, Central-Asia, North-East China, South-East Russia, Caucasia and Europe (Amirabadizadeh et al., 2023;Mironov et al., 2000).The populations of E. subspinosa in Iran, primarily observed along roads, are limited to Northern area of Semnan province, thriving in seasonal waterways and areas with shortterm rainwater accumulation.Due to the very low occurrence of these populations and the threats they face from human factors such as road construction, we classify the species E. subspinosa as endangered in Iran, based on the International Union for Conservation of Nature classi cations (Jalili et al., 1999;Mehrabian, 2013).Therefore, it is necessary to understand the changes in E. subspinosa distribution as it could affect the biodiversity and ecosystem functioning of Eastern Alborz which is one of the biologically richest areas of South-East Asia (Noroozi, 2020).Our primary objective in this SDM was to employ the MaxEnt method to map the present distribution of E. subspinosa and forecast the future distributions of this species up to 2060.To the best of our knowledge, this is the rst SDM study conducted on a species belonging to the Eversmannia genus.

Study area
The study was conducted in Semnan province, Iran, spanning an area of approximately 96,816 square kilometers, located within the geographical coordinates of approximately 36.7°N (northernmost point), 34.6° N (southernmost point), 52.3° E (westernmost point) and 56.4° E (easternmost point) (Fig. 1).This region encompasses the southern slopes, hills and plains of the Eastern Alborz, a mountain range in northern part of Iran (Noroozi et al., 2008).
The climate of the study area is predominantly shaped by the hot and dry air currents originating from the desert plain (Ganjaeian, 2023).However, several other factors contribute to the region's climate, including its distance from the sea, the orientation and extent of the mountains, the altitude and prevailing wind patterns.
The precipitation is notably scarce, predominantly falling as rain during the colder seasons.The annual average precipitation is approximately 145 mm.The humidity levels increase from west to east across the area and from south to north.Species E. subspinosa is a type of short spinose shrub with a height of 10 to 40 centimeters (Fig. 2).It has upright, thin divergent branches, with almost opposite thorns up to 19 millimeters long.The leaves are imparipinnate, 1.5-5 cm long and 7-15-foliolate.The owers are arranged in many-owered, oblong racemes measuring up to 15 cm long.The dark brown calyx is covered with trichomes, measured 3.5-5.5 mm long with lanceolate teeth and an acute apex.The corolla is purple.The standard is oblong-obovate, measuring 13.5 to 15.5 mm long, with a retuse apex.The wings are dark brown and lanceolate, measuring around 5 mm long, with elongate oblong auricles at the base.The legume is brown, elongate, measuring approximately 4 cm, curved and tuberculate-granular.The seeds are 3 mm long, ovate, brown, smooth, shiny, more or less at.The owering season is from late April to mid-June (Amirabadizadeh et al., 2023).
In Iran, populations of E. subspinosa are found in seasonal waterways and puddles, mostly along the main roads in the northern part of Semnan province (Fig. 1).They grow in slightly alkaline (pH = 7.0-8.0)loamy soils with low to moderate salinity (EC = 2-8 dS/m).

Data Collection and Preprocessing
Presence Data: To initiate the study, eld surveys were conducted in the study area between the years 2007 and 2021 to record occurrences of E.subspinosa (Fig. 1).Geographical coordinates of these occurrences were recorded using a GPS device.
Environmental Layers: Bioclimatic data including 19 bioclimatic variables, primarily utilized for ecological purposes (Merkenschlager et al., 2023) and elevation data (Table 1), at a 30-second resolution were obtained from WorldClim2.1 (https://www.worldclim.org/)(Fick et al., 2017).Additionally, "land use" and "soil quality" data, which are considered static in SDM studies (Stanton et al., 2012), were sourced from Harmonized World Soil Database v 1.2 (https://www.fao.org/)(Fischer et al., 2008).QGIS software v 3.30.3(www.qgis.org)was employed to extract study area-speci c data from all the above-mentioned data in projection system EPSG:4326, to serve as environmental layers.Static variables (e.g.topographic data) including elevation, slope and aspect derived from the elevation data that was obtained from WorldClim2.1 using raster analysis tool in QGIS software in projection system EPSG:4326, further enhancing the environmental dataset (Table 1) (Stanton et al., 2012).
Table 1 All the environmental variables initially considered in the modeling process and the selected variables considered for the nal modeling.The model was trained using 16 out of the total 19 presence records, with the remaining 3 records used for validation.In addition, we utilized 10,000 background points, randomly selected by the model across the study area.The replication type was bootstrap and we generated ten replicates of the model with random selection of the presence data.This approach was taken to reduce uncertainty and create an ensemble model.

Environmental variable
The Area Under the Curve (AUC) was used to assess accuracy of the model as it is a widely used metric for measuring the accuracy of models in MaxEnt modeling (Morales et al., 2020).The AUC value ranges from 0 to 1, where an AUC of 1 indicates a perfect model and an AUC of 0.

Results and discussion
Based on the occurrence data, E. subspinosa is found to favor at lands as the value over occurrence is 0 degrees (Table 2).This is consistent with the known ecology of the species, which is found in seasonal waterways and areas where rainwater accumulates for a short period of time.
Our ndings also indicate that E. subspinosa prefers the annual mean temperature ranges from 12.78°C to 15.49°C and it thrives in areas where the precipitation during the coldest quarter of the year is between 61-78 mm.Precipitation seasonality, expressed as the coe cient of variation, is deemed suitable for this species when it lies between 53.93% and 70.68% (Table 2).
The species distribution modeling was performed using the MaxEnt approach, which has been widely recognized as one of the most accurate and reliable techniques for presence-only studies such as this one (Valavi et al., 2022).We also chose MaxEnt due to its popularity in distribution modeling of rare and endemic plants (Qazi et al., 2022).
The ensemble model achieved an AUC value of 0.981 which suggests high accuracy of the model (Fig. 3) (Anderson et al., 2011;Araújo et al., 2005;Yackulic et al., 2013).The maximum training sensitivity plus speci city categorical threshold was determined to be 0.256.Table 3 and Fig. 4 present the areas of potential habitats (p ≥ threshold), less probable habitats (threshold ≤ p < 0.5) and highly probable habitats (0.5 ≤ p ≤ 1).The present time model suggests that this species could be found mainly in the North and Northern-East areas of the study area in elevations spanning from 1117 m to 2154 m.The probability of presence increases from Southern desert lands with lower elevations, lower precipitation and higher temperatures toward Alborz mountains in the North but it starts to decrease again when the elevation gets higher (Fig. 2 and Fig. 4) and at lands become scarce.According to this model, the area of potential habitats for E. subspinosa could be 8070.496km 2 (Table 3 and Fig. 4).

Projections based on the CMCC-ESM2 climate model, in conjunction with all three IPCC Representative
Concentration Pathways used in this study, indicate a signi cant decrease in potential habitats for E. subspinosa by 2040 compared to the present time.These projections indicate that by 2040, less probable habitats are not expected to change signi cantly.However, highly probable habitats are predicted to decrease to 88.71% and 59.50%.According to the SDM predicted using CMCC-ESM2 and RCP8.5, even less probable habitats would decrease to 83.66% by 2040 (Table 3 and Fig. 4).
Consistent with the optimism of RCP2.6, the moderation of RCP4.5 and the pessimism of RCP8.5, the potential habitats predicted by the CMCC-ESM2-RCP4.5 model would cover fewer area than those predicted by the CMCC-ESM2-RCP2.6model, but more area compared to the potential habitats predicted by the CMCC-ESM2-RCP8.5 model (Table 3, Fig. 4).In addition, high probable habitats would become more rare in lower elevations, where high potential habitats would be replaced by low probable habitats and non-suitable areas.Conversely, higher elevations would become more suitable for E. subspinosa (Fig. 2 and Fig. 4), which could be interpreted as the impact of global warming.However, probably as the at lands suited for the species are less common in higher elevations of the Eastern-Alborz mountain range, the potential habitats would cover fewer area in 2040 in comparison to the present time.
Projections from the CMCC-ESM2 climate model, in conjunction with the RCP2.6,RCP4.5 and RCP8.5, suggest that by 2060, the potential habitat areas could decrease to 0.84%, 8.55% and 12.12% of the current potential habitat area, respectively.In addition, models based on CMCC-ESM2 predict that high probable habitats would decrease to less than 1% and less probable habitats would reduce to under 22.82% by 2060, if any of these pathways come to pass (Table 3 and Fig. 4).These projections show that almost all the high probable habitats would disappear and only some parts of the less probable habitats would remain in higher elevation in the north of the study area.
Interestingly, based on the CMCC-ESM2-RCP4.5 model, the area of potential habitats of E. subspinosa in 2060 would be signi cantly more than that projected by CMCC-ESM2-RCP2.6 model and slightly less than that projected by CMCC-ESM2-RCP8.5 (Table 3 and Fig. 4), which seems to be inconsistent with expectancy of the RCPs.Since the species would migerate to higher eleveations, one reason for that inconsistency could be the presence of suitable at lands (plateaus) on highest elevations of Eastern-Alborz mountain range, which would fall into suitable temperature range of E. subspinosa if harsher pathways come to pass and temperatures gets higher.
The projections based on the HadGEM2 climate model indicate a more severe potential impact, suggesting the complete extinction of this legume species in the study area by 2040, if either RCP2.6,RCP4.5 or RCP8.5 come true (Table 3).So, these projections leave no room for further discussion.
Without a direct comparison of the two models over the study area, it is challenging to de nitively determine which model is more suitable for this speci c study and which projection is more accurate.To the best of our knowledge, our study is the rst Species Distribution Modeling (SDM) study conducted in Iran using CMCC-ESM2 climate model.shown varied responses to climate change, with some species increasing and others decreasing in population size.This nding highlights the complex interplay between climate change and species distribution, underscoring the need for continued research and conservation efforts (Rostö, 2020).

Older versions of
The Jackknife of AUC test yielded an AUC of over 0.93 for the Bio1 environmental variable when used in isolation in the Maxent model.This was followed by elevation, which resulted in an AUC between 0.92 and 0.93 and slope, which had an AUC of less than 0.8.The omission of both Bio1 and Bio15 variables led to an AUC of less than 0.96.However, the omission of other variables did not signi cantly affect the model's AUC (Fig. 5).
These results suggests that temperature plays a signi cant role in determining habitat suitability of the species.Interestingly, despite having a lower contribution percentage (24.5%)compared to elevation (54.9%), the Bio1 variable showed a higher permutation importance (71.2%) that underscores the sensitivity of the species to temperature variations, implying that it may have speci c temperature requirements for survival and reproduction, which is in alignment of narrow range of Bio1 variable (12.78°C-15.49°C)over occurrence locations (Table 2) and Bio1 response curve of the model (Fig. 6).
On the other hand, the elevation variable, while contributing the most (54.9%) to the model, did not demonstrate as high a permutation importance as Bio1.This could suggest that while elevation plays a role in distribution of the species, its in uence is more generic and less sensitive to speci c changes in value.Elevation has the most percent contribution to the model with 54.9%, followed by Bio1 with 24.5% contribution to the model.However, Bio1 variable with a permutation importance of 71.2% is the most in uential variable when considering the effect of permuting the values of each variable across the training presence and background data (Table 2).
Many factors contribute to the way the distribution of plants shifts in response to raising tempreatures.In their comprehensive 2021 study, Steven R. Beissinger and Eric A. Riddell highlight the crucial role of physiological traits of species and the availability of suitable habitats in determining the speed and direction of species range shifts (Beissinger et al., 2021).Chardon et al. (2023) found that microhabitat variables weakly predict recruitment and seedling survival, with varied responses to environmental factors across species and life stages (Chardon et al., 2023).Adding to this complexity, a recent research has highlighted the importance of considering intraspeci c ecological variation when predicting species distributions for conservation and climate change impact assessment (Varaldo et al., 2023).

Conclusion
The present study aimed to model the current and future distribution of E. subspinosa which is an endangered species in Iran, found only in rare locations on the Southern side of the Eastern Alborz mountain range in Semnan province.Alborz mountains are in high priority for conservation efforts since they are home to many rare species as one of the areas in the South-East Asia with highest biodiversity.However, they have been under intensive human activities like road construction, recently (Noroozi, 2020).
Our presence-only SDM study suggested a signi cant loss in the area of the habitats of E. subspinosa, an endangered species in Iran, across the study area based on the CMCC-ESM2 climate model by 2060 due to global warming.In addition, the CMCC-ESM2 model suggests that this legume plant would migrate to higher elevations seeking at lands with suitable temperatures.The HadGEM3 climate model suggests an even more severe impact of climate change on the future distribution of E. subspinosa, predicting its complete extinction by 2040 under RCP2.6,RCP4.5 and RCP8.5 pathways in Semnan province, Iran.E. subspinosa is not only vulnerable to global warming across our study area, but human activities like road construction could also affect the future distribution of this endangered plant.Despite the tendency of this species to migrate to higher elevations due to rising temperatures, new roads may limit its capability to do so.While there are some protected areas and natural parks across the main habitats of E. subspinosa, the gaps between these areas could hinder the response of the species to global warming and accelerate its probable extinction over the study area (Noroozi, 2020).Therefore, lling these gaps with de ning new protected areas and national parks is highly recommended to not only delay the extinction of this species, but also to protect other rare and endemic species of the study area (Noroozi, 2020;Noroozi et al., 2023).
Finally, we concluded that temperature, elevation and slope are the most important environmental variables affecting E. subspinosa distribution.Moreover, we should not ignor the role of other factors like ecological relationships with other species, the main way of propagation and seed germination and also human activities like road constructions, to plan conservation efforts more effectively.We suggest further studies using more occurrence locations or presence-absence data over a wider study area which could span from Western China to Eastern Europe.Using citizen science data (Hutchinson et al., 2017) could be useful due to the rarity of this species.We also suggest future studies to consider reproduction traits of the plant and human activities affecting its habitats.
groups for their production and provision of model output, the Earth System Grid Federation (ESGF) for their role in data storage and access and the numerous funding bodies that back CMIP6 and ESGF.
We also wish to acknowledge the organizations that contributed to the Harmonized World Soil Database v 1.2.This includes the Food and Agriculture Organization of the United Nations (FAO), the International Institute for Applied Systems Analysis (IIASA), ISRIC-World Soil Information, the Institute of Soil Science, Chinese Academy of Sciences (ISSCAS), and the Joint Research Centre of the European Commission (JRC).Their collaborative efforts have been instrumental in our research.

Figures
Figure 2 The HadGEM have been used in several studies in Iran.Ahmadi et al. (2020) conducted a study on the potential impact of future climate on the distribution of Taxus baccata L. in Hircanian forests of Northern Iran, using six climate models including HadGEM2, HadGEM2-CC and MIRI-CGCM3.They concluded that T. baccata would lose the majority of its habitats by 2050 due to global warming(Ahmadi et al., 2020).In another study, Ghehsareh Ardestani and Heidari Ghahfarrokhi used HadGEM2-CC to project the future distribution of Salvia hydrangea in the Central Zagros mountain range, Western Iran.They found out that the area of good habitat suitability (0.5 < p < 0.75) and moderate habitat suitability (0.25 < p < 0.50) of S. hydrangea in the study area would signi cantly increase but area of low habitat suitability (0 < p < 0.25) and high habitat suitability (0.75 < p < 1) would decrease under HadGEM2-CC-RCP4.5  and HadGEM2-CC-RCP8.5 models by 2050 and 2070.Conversely, MIRI-CGCM3 climate model suggested signi cant loss of the habitats of this species (Ghehsareh Ardestani et al., 2021).Tarnian et al. (2021) used HadGEM2-AO, CCSM4 and MIROC5 to predict the future distribution of Daphne mucronata across Iran, under RCP4.5 and RCP8.5.They reached the conclusion that this species would lose 35% and 61% of its habitats in Iran and almost would disappear in elevations lower than 1500 m and 2000 m, under RCP4.5 and RCP8.5 pathways, respectively (Tarnian et al., 2021).To the best of our knowledge, our study is the rst to use HadGEM3 for SDM in Iran.There are not many SDM studies conducted on plants across Eastern-Alborz mountain range in Iran.In one study similar to ours,Fatemi et al. (2018) predicted the future distribution of Juniperus excelsaacross Central-and Eastern-Alborz mountain range by MaxEnt modeling up to 2070.They found that J. excelsa is projected to lose some of its habitats and is expected to gradually migrate to higher elevations in response to rising temperatures(Fatemi et al., 2018).Many SDM studies have projected a shift in distribution of plants toward higher elevations due to global warming in recent years.A 2022 study by Auld et al. found that climate change is causing rapid elevational shifts among Australian alpine plants.Some species are moving upslope, but continued warming could deplete suitable habitats within a century, emphasizing the need for ongoing conservation efforts(Auld et al., 2022).A recent study on alpine species in Norrbotten County, Northern Sweden, has

Table 2
Phillips et al., 2006)na et al., 2021;Stanton et al., 2012)based on CMCC-ESM2(Cherchi et al., 2019;Vichi et al., 2011)and HadGEM3(Andrews et al., 2020) climate models paired with three IPCC Representative Concentration Pathways RCP2.6, RCP4.5 and RCP8.5(Meinshausen et al., 2020)for years 2040 and 2060 sourced from WorldClim2.1.Again, QGIS software was applied to extract study area-speci c data in projection system EPSG:4326 to serve as projection layers in MaxEnt modeling.Variable Selection and Correlation Analysis: To optimize the model, the R package v 4.3.2 was used to generate a correlation matrix using the Pearson method, identifying and removing highly correlated variables (correlation ≤ 0.8)(Hanberry, 2024;Harisena et al., 2021;Stanton et al., 2012).Finally, only ve variables including elevation, bio1, bio19, bio15, and slope were considered for the modeling process as they were the most in uential environmental variables with correlation ≤ 0.8 (Table2).Environmental variables selected for the nal modeling.Values over occurrence locations, percent contribution and permutation importance are also provided.It expresses a probability distribution where each grid cell has a predicted suitability of conditions for the species (S.B.Phillips et al., 2006).Maxent distribution maps are interpreted as the predicted probability of species presence across the geographical area of interest(Elith et al., 2011).Areas with higher predicted probabilities are considered more suitable for the species based on the environmental conditions (S. B.Phillips et al., 2006).
Species distribution modeling was performed using the Maximum Entropy (MaxEnt) approach, speci cally utilizing version 3.4.4 of the MaxEnt software (S.J. Phillips et al., 2019).MaxEnt is a popular tool for modeling species distributions (Elith et al., 2011; S. B. Phillips et al., 2006).To model the habitat suitability of E. subspinosa in the study area, the selected set of environmental variables including elevation, slope and bioclimatic variables bio1, bio15 and bio19 were loaded in MaxEnt 3.4.4visual interface as environmental layers.To anticipate the potential impacts of climate change on the distribution of E. subspinosa, projection layers were loaded in MaxEnt 3.4.4visual interface.
Dhyani et al., 2021)6)del performs no better than random guessing.In MaxEnt modeling, maximum AUC of a well-tuned model is 1 -a/2, where "a" is the prevalence of the species (proportion of the study area occupied by the species) (S.B.Phillips et al., 2006).Maximum training sensitivity plus speci city was used as the threshold rule as the occurrence locations were scarce and the study conducted presence-only (A.Dhyani et al., 2021).
Phillips et al., 2006;Tesfamariam et al., 2022) based on the probabilities (p) determined by the ensemble model.These are: (a) Potential habitats where p ≥ threshold, (b) Less probable habitats where threshold ≤ p < 0.5 (c) High probable habitats where 0.5 ≤ p ≤ 1(A.Dhyani et al., 2021; S.Dhyani et al., 2020).The area of these habitats was calculated for both current and future scenarios across the study area using R package v 4.3.2.The jackknife test was applied to show the importance of each environmental variable in the modeling process(Kalle et al., 2013;Wang et al., 2024).In addition, response curves were generated to determine the effect of each environmental variable on the predicted probability of presence of E. subspinosa over the study area (Labarca-Rojas et al., 2022; S. B.Phillips et al., 2006;Tesfamariam et al., 2022).

Table 3
The predicted potential, less probable and high probable habitat area of E. subspinosa over the area of Semnan province, Iran (the study area) in the present time compared to that in 2040 and 2060 based on CMCC-ESM2 and HadGEM2 climate models under RCP2.6,RCP4.5 and RCP8.5 pathways.