Climate change threatens micro-endemic amphibians of an important South American high-altitude center of endemism

. Species inhabiting broad altitudinal gradients are particularly exposed to the effects of global climate change (GCC). Those species reaching mountain tops are the most negatively affected. Here, using ecological niche models we estimated the climate change exposure of endemic amphibians of the most important extra-Andean mountain system of Argentina: the Sierras Pampeanas Centrales. Our results pinpoint that micro-endemic amphibians of this mountain system are heavily exposed to the effects of GCC, with important constraints of suitable climatic conditions for the six analyzed species. Among the most important ﬁndings, our models predict a high probability of a total disappearance of suitable climatic conditions for two of the species, currently restricted to mountain tops. This high exposure, in synergy with their very restricted ranges, and other important human induced threats (as ﬁsh invasion and emergent diseases), pose a serious threat to these endemic species, which can enter into the “extinction pathway” in a near future if no concrete conservation actions are taken. Our ﬁndings provide additional evidence of the great negative impact of GCC in high-altitude centers of endemism. center,

Species that inhabit steep altitudinal gradients are particularly exposed to the effect of global climate change. In fact, significant species declines and local extinctions in mountain systems have been empirically demonstrated (Pounds et al., 2006;Lenoir et al., 2008;Raxworthy et al., 2008;Rowe et al., 2015;Yousefi et al., 2015), with dramatic future forecasts as a consequence of the phenomenon (Chen et al., 2011;Engler et al., 2011;Forero-Medina, Joppa and Pimm, 2011;Enriquez-Urzelai et al., 2019). Given their physical characteristics, the environmental gradients of mountain systems are especially heterogenous. Along these environmental gradients, likely, the expected temperature increases could lead to a significant upward displacement of the appropriate climatic conditions for certain organisms (Enriquez-Urzelai et al., 2019). Logically, species inhabiting highest altitudes, and especially those reaching mountain tops, are the most negatively affected by the pattern (Nori et al., 2016). This upward displacement of the climatic niche necessarily produces an important reduction of the geographically explicit climatic niche of the species (Broennimann et al., 2006;Blaustein et al., 2010). If these organisms are not able to adapt (physiologically or ecologically) to the changes, they will become locally extinct (Bellard et al., 2012;Moritz and Agudo, 2013;Pacifici et al., 2015).
This picture becomes more worrisome if we consider that, given their particular geological histories, mountain systems are generally centers of endemism and in most cases harbor a greater diversity than neighboring regions (Jingyun, Zehao and Haiting, 2004;Körner, 2004;Brehm et al., 2008;Kozak and Wiens, 2010;Hoorn et al., 2013;Graham et al., 2014). In Argentina, due to its extension and degree of isolation, the orographic system of the Sierras Pampeanas represents the most important extra-Andean mountain formations. Within this great formation, the Sierras Pampeanas Centrales (SPC) (also known as Sierras Pampeanas de Córdoba and San Luis or Sierras Pampeanas Orientales) present many endemic plants (Cabido et al., 1998), fungi (Robledo et al., 2006;Robledo and Renison, 2010) and vertebrates (Cei, 1980;Nores, 1995). The highest endemicity levels of the system are associated with the grasslands and highland forests that develop in the summits and extensive sectors of the SPC Martínez et al., 2016). The historical processes that determined the current composition of the SPC biota were also key for the system to host a rich amphibian fauna, whose distribution is strongly linked to the altitudinal gradient of the system .
This system is inhabited by 24 amphibian species, eight of which are strictly endemic from 600 m.a.s.l. to the mountain top (i.e. 2780 m.a.s.l.). Of these endemic amphibians, four inhabit a larger altitudinal gradient: Melanophryniscus stelzneri, Melanophryniscus estebani, Odontophrynus cordobae, and Boana cordobae wich reaches the mountain top. The remaining species (Odontophrynus achalensis, Pleurodema kriegi, Pleurodema cordobae and Rhinella achalensis) are restricted to the highest altitudes. In addition, these species are the most exposed to other important extinction drivers such as Batrachochytrium dendrobatidis or the invasion of alien fishes (Lescano, Longo and Robledo, 2013). Given this scenario, our main objective was to assess and quantify the exposure to climate change (i.e. loss or gains of suitable areas due to variations in climatic conditions) to these endemic amphibian species in the near future, by using ecological niche modelling techniques, in order to generate additional relevant information to guide accurate conservation decisions that need to be taken in the nearest future.

Species and records
We compiled and organized a database containing occurrence records for the entire distribution of all the endemic amphibians from the SPC (Boana cordobae, Melanophryniscus stelzneri, Melanophryniscus estebani, Odontophrynus cordobae, Odontophrynus achalensis, Rhinella achalensis, Pleurodema kriegi and Pleurodema cordobae). Sibling and cryptic species of Melanophryniscus and Pleurodema were considered as a single taxon, hereafter named as Melanophryniscus spp. and Pleurodema spp. respectively (see details in Lescano et al., 2015). Records were obtained from relevant literature (Filipello and Crespo, 1994;Di Tada, 1994;Rosset et al., 2009;Baraquet et al., 2012;Lescano et al., 2015; and cites within), specific material from herpetological collections (Fundación Miguel Lillo and Museo Argentino de Ciencias Naturales) and field sampling ( fig. 1, supplementary table S1). Data were georeferenced using Google Earth software and mapped using ArcGis 10.3. Our final database contained 156 individual records (43 records for B. cordobae, 32 records for Melanophryniscus spp., 17 records for O. achalensis, 22 records for O. cordobae, 16 records for Pleurodema spp. and 18 records for R. achalensis). The complete database is available in supplementary table S1.

Climatic variables
We ran Pearson's correlation between 19 bioclimatic variables at a spatial resolution of 30 arc seconds retrieved from WorldClim, and selected five of them based on their collinearity (r < 0.75) and biological relevance for species survival. The selected variables were: Mean Diurnal Range, Isothermality, Temperature Seasonality, Mean Temperature of Wettest Quarter, Precipitation Seasonality, and Precipitation of Wettest Quarter. All selected variables were masked with the extent of the Sierras Pampeanas Centrales of Cordoba and San Luis, Argentina, considering 400 meters above sea level to the maximum altitude of the mountain system. We considered three different climatic scenarios for these variables: (a) current climatic conditions, and (b) two different Global Circulation models: IPSL-CM54-LR and INMCM4 (hereafter IPSL and INM respectively), for the year 2050. For these future hypothetic scenarios, we considered a Representative Concentration Pathway of 4.5 at intermediate gas concentration. All the bioclimatic variables were downloaded from WorldClim 1.4 (www.worldclim.org; Hijmans et al., 2005).

Ecological niche models (ENMs)
ENMs were built to estimate climatically suitable places for the selected species in SPC. These ENMs provide spatially explicit habitat suitability maps, and can therefore be used to inform management and conservation decisions (Guisan and Thuiller, 2005). There are currently many alternative ENM algorithms, and each one can have a differential performance and accuracy under different conditions (Diniz-Filho et al., 2009). A proposed solution to take this variability into account is to combine different algorithms into ensembles in order to find areas of consensus (as suggested by Araújo and New, 2007), and an optimal model complexity (Merow et al., 2014). This approach has successfully applied in many ecological studies to develop optimal habitat suitability maps for vertebrates (Coetzee et al., 2009;Grenouillet et al., 2011;Scales et al., 2016;Nori et al., 2017;Austin et al., 2019). Here, ensembles were generated by using three different modelling algorithms: Generalized Linear Models (GLM), Artificial Neural Networks (ANN), and Random Forest (RF). All the modelling analysis were generated in the R package (R Core Team, 2014) biomod2 (Thuiller et al., 2014).
We calibrated models using a subset of 70% of the records of each species and used the remaining 30% to test the models. Additionally, we generated a random set of pseudo-absences for model calibration and evaluation. We defined the number of pseudo-absences as twice the number of presences. Pseudo-absences were generated within the limits of the study area and outside a 20 km buffer zone of the presences, in order to reduce the possibility of including pseudo-absences in places climatically identical to the presence pixels (Barve et al., 2011).
After calibration, each model was projected for the three climatic scenarios. We selected the maximum suitability with an omission rate of 0.05 of total of presence as threshold value. For model evaluation, we used the remaining 30% of the data in order to calculate two widely used metrics: the True Skill Statistics (TSS, Liu et al., 2011) and the Area Under the "Receiver Operating Characteristics" Curve (hereafter AUC / ROC, Fielding and Bell, 1997). These metrics allows us to compare the relative performance of each model. We repeated the analysis three times, thus providing a threefold internal crossvalidation of models. True Skill Statistics values vary from −1 to +1, with values equal or lower than zero implying statistical fits that are no better than those generated at random. The AUC / ROC values vary from 0 to 1, with values higher that 0.5 being better than random. Final ensembles only included those projections with TSS higher than 0.6 and AUC / ROC higher than 0.8. Ensembles were calculated as an average of all these projections, weighted by their AUC / ROC scores.
In order to assess a species' exposure to global climate change, we overlapped final ensembles at current and future climatic conditions and quantified spatial change in the suitable climatic surface. Additionally, for each species, we quantified the number of pixels at different classes of altitude in order to assess the displacements of suitable surfaces in the altitudinal gradient.

Results
All algorithms showed good performance. We found a mean TSS = 0.78, sd = 0.06; mean AUC / ROC = 0.92, sd = 0.03 for Boana cordobae, and nine projections were included in the final ensemble. We found a mean TSS = 0.84, sd = 0.09; mean AUC / ROC = 0.94, sd = 0.05 for Melanophryniscus spp., with all projections being included in the final ensemble. Odontophrynus achalensis obtained a mean TSS = 0.95, sd = 0.13; mean AUC / ROC = 0.98, sd = 0.07, with all projections included in the final ensemble for the species. Odontophrynus cordobae obtained a mean TSS = 0.64, sd = 0.20; mean AUC / ROC = 0.82, sd = 0.11, with six projections included in the final ensemble. We found a mean TSS = 0.93, sd = 0.11; mean AUC / ROC = 0.97, sd = 0.04 for Pleurodema spp., with all projections included in the final ensemble. Finally, we found a mean value of TSS of 0.99, sd = 0.01; mean AUC / ROC = 0.96, sd = 0.02 for Rhinella achalensis, and again, all nine projections were included in For the future forecasts, we found a clear general pattern characterized by an upward displacement of suitable climatic conditions for the species with a consequent reduction of the surface occupied for the climatic niche, especially for those species reaching the mountaintops ( fig. 3A, B). The species with the highest exposure were O. achalensis and Pleurodema spp., which showed a projected total loss of climate niche in the geographic space for both future scenarios. These were followed by R. achalensis with a total loss for the IPSL scenario and a reduction of 72. 2% fig. 3A, B).

Discussion
The results of this study evidence a clear pattern of differential exposure to global climate change for endemic amphibians of the Sierras Pampeanas Centrales. In general, the conservation status of the species is consistent with the potential distribution predicted by our ENMs (fig. 2)  with the current distribution in background (red). On the foot of each map it is quantified the exposure to GCC for each future scenario relative to the current. The exposure (Exp) to GCC for each future scenario relative to the current situation is calculated as the percentage of N°pixels of future scenario / N°pixels of current distribution. B) Histograms showing the altitudinal range with suitable climatic space for each species at each respective scenario. has a very restricted range, and categorized as NT, was the exception. These findings reinforce a previously described pattern (e.g. Lenoir et al., 2008;Chen et al., 2011;Engler et al., 2011;Comte et al., 2013;Enriquez-Urzelai et al., 2019). The species with wider distributions in SPC are currently categorized under non-threatened categories [Boana cordobae (LC), Melanophryniscus spp. (LC) and Odontophrynus cordobae (LC)].
As expected, our models predict an upward migration of suitable climatic conditions for the studied species ( fig. 3A, B). Logically, those species not reaching mountaintops (Melanophryniscus spp., O. cordobae) or inhabiting the entire altitudinal gradient (B. cordobae) were proportionally the least affected (although most of them experienced reductions of suitable climate surface) ( fig. 3A, B). The most exposed species are micro-endemic and heavily threatened by other human activities, such as introduced fishes for recreational activities, diseases or land cover changes (Lescano, Longo and Robledo, 2013;. In fact, our results predict a total loss of suitable climatic conditions for these micro-endemic species (O. achalensis and Pleurodema spp.) by 2050 ( fig. 3C).
Global Climate Change will generate strong negative effects in amphibian species and populations of different parts of the world, such as Near East Africa (Kafash et al., 2018), North America (Wright et al., 2016), and Southeast Asia (Bickford et al., 2010). The same strong pattern of reduction of adequate areas was also found in other mountain ranges systems, including the Australian Wet Tropics bioregion (Williams, Bolitho and Fox, 2003) or Europe (Lenoir et al., 2008;Engler et al., 2011, Araujo Thuiller andPearson, 2011). A strong negative impact of GCC is also expected for highland amphibian species in China (Duan et al., 2016) and Ethiopia (Gower et al., 2013). The negative response of amphibians to GCC in highland areas is particularly important because these environments harbor high diversity of endemic amphibians all across the world (Fu et al., 2006;Bergl, Oates and Fotso, 2007). Additionally, given their isolation and particular geological history, mountains systems are particularly important in terms of conservation (Roelants, Jiang and Bossuyt, 2004).
The projected exposure to GCC estimated in this study ( fig. 3B) must not be considered as a future geographical reduction of the distributional ranges of species (see Dawson et al., 2011), but as a hypothetic distribution. Our results should only be interpreted as displacements (and eventually, reductions) of the suitable climatic conditions for the species in these hypothetical future climatic scenarios for 2050. These models could be improved by incorporating other kind of information, such as the physiological or ecological capacity of species to tolerate and adapt to climatic changes. Nevertheless, these ENMs are very helpful to generate hypotheses on the potential response to the species to GCC, being a useful tool to guide conservation efforts.
Based on our findings (among others; Lescano, 2018), R. achalensis, one of the three most threatened micro-endemic species of this study (and one of the most exposed to GCC; figs 2 and 3C), was recently recategorized from NT to EN (IUCN, 2018). Of the remaining two micro-endemic species, one of the most exposed species (O. achalensis) is currently within a threatened category of IUCN (VU), while Pleurodema spp. remains categorized as NT, despite its very small distributional range (approximately 2000 km 2 ; fig. 1, supplementary table S3), in synergy with other important factors which threaten the species (Lescano, Longo and Robledo, 2013;Lescano et al., 2015). According to criterion B1 of the IUCN a species is categorized as threatened, if besides of the restricted extent of occurrence, at least two additional necessary conditions are met (being one of them an inferred or projected reduction in the extent of occurrence; see IUCN, 2015). According to our results, for Pleurodema spp., at least one additional condition is met (i.e. extent of occurrence > 20 000 km 2 with a drastic projected reduction). Thus, Pleurodema spp. is a possible candidate to be categorized as Vulnerable. It would be very important to carry out additional studies of the species related to assess the other two conditions (related to fragmentation and population fluctuations), in order to determine the appropriate conservation status for this species.
Our findings are not only evidence of the great potential impact of the GCC on microendemic amphibians that inhabit mountain ranges, but also provide additional evidence of the large negative impact of GCC on highland endemism centers (Dirnböck, Essl and Rabitsch, 2011). While this general problem seems to be understood by decision makers, it is clear that there is a delay or even neglect in the conservation efforts. Conservation measures are needed to mitigate the effect of GCC on biodiversity. For this particular system, some concrete actions could be important, as: (i) protect and monitor headwaters in where these endemic species breed, (ii) control the invasion and avoid new introductions of salmonids, (iii) promote the environmental education activities to make aware stockholders (e.g. local population, decision makers, tourists) to the worrying picture which these species are undergoing. All these palliative measures, could help to reduce the pressure on endemic amphibian species of SPC allowing them to cope with future climate changes.