Predicting impacts of global climatic change on genetic and phylogeographical diversity of a Neotropical treefrog

Future climate changes may affect species distribution and their genetic diversity, hampering species adaptation to a new climate or tracking the suitable conditions. Amphibians have high sensitivity to environmental degradation and changes in temperature and humidity. Thus, the expected climatic changes by the end‐of‐century (EOC 2100) may cause local or complete extinction of some species. Here, we address the effects of climate change on genetic and phylogeographical diversity, together with the geographical distribution of the South American treefrog Scinax squalirostris Lutz, 1925. Furthermore, we assess how protected areas will conserve its genetic variation.

Global climate is changing drastically, leading to an increase of 0.6°C of the world average temperature (Jones et al., 2001) in the last century, due to the anthropogenic emission of greenhouse gases (Solomon, 2007;IPCC, 2014IPCC, , 2018. It is expected an increase of ~2°C to 4°C by mid-century (Brown & Caldeira, 2017;Fischer et al., 2018;New et al., 2011). Global warming may force species to change their current distribution range to track new suitable areas (Chen et al., 2011). However, some species may not be able to follow their optimal climatic niche as they might have limited dispersal capacity or cannot disperse across natural or anthropic barriers (Urban et al., 2012), driving them to either adapt to the new conditions or to become extinct. Species with low dispersal capacity may have higher extinction rates and population declines than vagile species (Duan et al., 2016).
Climate change can impact species distribution by reducing the climatic suitability across their range, diminishing their fitness, impacting the intraspecific genetic diversity, increasing genetic stochasticity in small populations and, therefore, fixation of poorly adapted haplotypes (Pauls et al., 2013;Urban et al., 2013).
Genetic diversity represents species' evolutionary potential to adapt to environmental changes (Urban et al., 2012(Urban et al., , 2013. High levels of genetic diversity across evolutionary independent intraspecific lineages may favour evolutionary responses to climatic changes, but low genetic diversity may reduce their evolutionary potential, making them prone to extinction (Spielman et al., 2004;Frankham, 2005;Rizvanovic et al., 2019). Thus, identifying populations and lineages that might be more vulnerable to climate change is a concern of biology and conservation genetics (D'Amen et al., 2013).
Ecological Niche Modelling (ENM) is widely used to predict changes in the geographical distribution of species under climate change (Collevatti et al., 2015;Lima et al., 2014Lima et al., , 2017Thomas et al., 2004). This framework uses the current distribution data of species to model their fundamental environmental niches and forecast potential distributions under various climatic scenarios (Peterson et al., 2011), enabling the understanding of spatial dynamics and the impact of climatic changes on species distribution in the future (Mendoza-González et al., 2013;Simon et al., 2013;Vasconcelos et al., 2018;Zhang et al., 2012). Future distributions are forecasted by projecting the current suitable climatic conditions of species onto the future climatic scenarios (Lima et al., 2017;Midgley et al., 2002). These projections assume that species conserve their current environmental niche and they will be able to track future suitable climatic conditions (Elith & Leathwick, 2009).
Coupling ENM with simulations of genetic parameters within and among populations can improve our understanding of climate change effects on genetic diversity and population genetic structure (Lima et al., 2017). Understanding the effects of climate changes on genetic diversity can help to predict how populations will cope with rapid future changes (Corn, 2005;Duan et al., 2016), and to protect the viability and genetic diversity of species in the future (Mendelson et al., 2006;Schwartz et al., 2007).
Amphibians are one of the most threatened groups by climate change, they are amongst the most vulnerable vertebrates, declining even faster than birds and mammals (Foden et al., 2013;Stuart et al., 2004). Approximately, 50% of all amphibian populations are declining or endangered (Alroy, 2015;Duan et al., 2016;Gibbons et al., 2000;Houlahan et al., 2000;Stuart et al., 2004). Because amphibians rely heavily on environmental conditions (Duellman & Trueb, 1994), climate change is expected to cause the extinction of around 40% of the amphibian species, with greater impact on endemic groups (Foden et al., 2013;Gibbons et al., 2000;Thomas et al., 2004), making them the most threatened group of animals (Mendelson et al., 2006;Stuart et al., 2004). Amphibians have low vagility, high sensitivity to changes in environmental temperature and humidity, mainly due to their permeable skin and life cycle, which includes stages in both water and land for most of them (Duellman & Trueb, 1994;Storfer et al., 2009;Stuart et al., 2004).
Changes in temperature and precipitation can lead to changes in geographical distribution and the abundance of amphibian populations (Chen et al., 2011). Furthermore, climate changes can affect the hydroperiod: the time period that a temporary pond retains water available for amphibian populations to mate and metamorphose (Carey & Alexander, 2003;Rowe & Dunson, 1995). The South American treefrog Scinax squalirostris Lutz, 1925 is a model group to study the impact of climate change on South American herpetofauna because of its wide distribution throughout the central-west, southeast and south of South America (Figure 1; Frost, 2021). It inhabits diverse habitats such as open formations, forests, grasslands and rushes. It breeds in small and temporary ponds or cattle ponds (Neves et al., 2019;Vaz-Silva et al., 2020).
Here, we address the effects of future climate change on the geographical distribution, genetic and phylogeographical diversity of Scinax squalirostris. We use Ecological Niche Modelling coupled with genetic simulations, conditioned to climatic and topographic variables, to model the dynamics of genetic clusters through time. In addition, we perform a clustering analysis of the current genetic and phylogeographical diversity to analyse their conservation within the current Protection Areas (PAs) scheme.

K E Y W O R D S
Bayesian clustering, conservation, Ecological Niche Modellings, global warming, Hylidae, phylogenetic diversity 2 | ME THODS

| Population sampling and genetic data
We sampled 219 individuals in 26 localities across the distribution of Scinax squalirostris (Brazil, Uruguay and Paraguay), with a sampling effort ranging from 1 to 10 individuals per locality ( Figure 1; Table   S1 in Appendix S1). From all individuals, we collected tissues (muscle and liver) to extract DNA using the Dneasy Blood & Tissue Kit (Qiagen ® , Chatsworth, CA).
To obtain genetic data, we sequenced two fragments of the mitochondrial DNA (hereafter mtDNA): the 12S ribosomal subunit (primers 12Sa-12Sb; Reeder, 1995) and a fragment of the cytochrome B gene (cytB, primers MVZ15L, Moritz et al., 1992 andH15149, Kocher et al., 1989; see Appendix S2 for GenBank accession numbers). We also sequenced the nuclear (hereafter nDNA) RAG-1 gene (413 bp, Heinicke et al., 2007; see Appendix S2 for GenBank accession numbers). These gene sequences are often used to access genetic diversity and phylogeography of amphibians (Barrow et al., 2020;Fusinato et  Further, we used CLUSTAL OMEGA (Sievers et al., 2011) to obtain the multiple sequence alignments. Coding sequences were tested for saturation by plotting transitions and transversions with TN93 distance (Tamura & Nei, 1993) using DAMBE software (Xia, 2013). We excluded from the final alignment the third codon position of cytB, which showed high levels of saturation ( Figure S1 in Appendix S4).

| Genetic and phylogeographical diversity
We used concatenated sequences (mtDNA and nDNA = 1,002 bp) to obtain an estimation of genetic diversity parameters using ARLEQUIN 3.11 (Excoffier et al., 2005). We estimated nucleotide (π), haplotype (h) diversities (Nei, 1987), and the number of haplotypes (k) for each population and overall populations.
Hence, even high values of genetic diversity may represent low lineage or phylogeographical diversity among populations (e.g. the diversity of lineages at different locations). Therefore, we applied a metric based on Phylogenetic Diversity (Faith, 1992), calculating the minimum spanning path of a set of haplotypes in a coalescent tree, defined hereafter as Phylogeographical Diversity (PGeoD). PGeoD was calculated using R package picante (Kembel et al., 2010).
To infer the coalescent tree, we ran a Bayesian coalescent analysis implemented in BEAST     precipitation of the coldest quarter (BIO19). All climatic layers were used in a spatial resolution of 0.5º cell size.

| Ecological niche modelling
The current and future potential distributions of the treefrog were inferred using 10 presence-only or presence-absence algorithms (Table S2 in Appendix S1). We randomly sampled pseudoabsences across the Neotropics, excluding cells with presence, maintaining the same number of absences as presence data. We fitted ENM models using fivefold cross-validation, repeated 20 times, and partitioning 70% of the data in training set (calibration) and 30% in testing set (evaluation). ENM was run with sdm R package (Naimi & Araujo, 2016). Model accuracy was assessed using the area under the receiver operator characteristic (AUC) (Allouche et al., 2006) and the true skill statistic (TSS), which take into account both omission and commission errors (Allouche et al., 2006). We used a weighted averaging of TSS statistics to ensemble the models (Araújo & New, 2007; Table S3 in Appendix S1).
The combination of all ENMs and AOGCMs resulted in 60 independent predictive maps (10 ENMs × 6 AOGCMs) for each period of time (present, 4.5 and 8.5 RCPs). We applied a hierarchical ANOVA using the predicted suitability of all models (10 ENMs × 6 AOGCMs × 3 time periods) as a response variable to identify and map the uncertainties due to the modelling components. Hence, the ENMs and AOGCMs components were nested into the time component, but crossed by a two-way factorial design within each period of time (see Terribile et al., 2012). To quantify the range sizes and range shifts among the three time periods, we transformed the suitability maps into binary maps (presence/absence) using a threshold of 0.5. Then, we calculated range size as the number of presence cells, for each one of the 60 predictive maps. Range shift was calculated as the difference in range sizes between the present-day and 4.5 or 8.5 RCP predictive maps.

| Simulations of genetic and phylogeographical diversity under climatic changes
To understand how climate change may affect genetic and phylogeographical diversities, we simulated the extinction of populations under several extinction thresholds for the two future climatic scenarios (4.5 and 8.5 RCP). We assumed that only populations occurring in areas with suitability higher than an extinction threshold will persist and contribute to the gene pool for the next generations (see Collevatti et al., 2011). As we do not know the real extinction threshold to predict population extinction, we applied a sensitivity analyses changing the extinction thresholds from 0.3 to 0.9, based on the minimum and maximum suitability range of the current populations (Table 1), to evaluate how much our conclusions are dependent on extinction thresholds values.
For each combination of extinction thresholds and climate change scenarios (4.5 RCP and 8.5 RPC), we calculated genetic diversity and PGeoD for the remaining populations. Genetic diversity was calculated using the software ARLEQUIN 3.11 (Excoffier et al., 2005). PGeoD was calculated by dropping the extinct haplotypes from a sample of 1,000 coalescent trees and re-calculating the minimum, maximum and median PGeoD of those coalescent trees in each scenario of population extinction.

| Forecasting the dynamics of genetic clusters
We performed spatially explicit simulations to predict the dynamics of genetic clusters conditioned to climatic changes using the software POPs (Jay et al., 2015). This software implements a Bayesian clustering based on genetic, geographic and environmental variables.
It assigns individuals to genetic clusters based on genetic ancestry and models the effect of climatic and landscape variables in such cluster assignments, thus predicting changes in the genetic structure in response to environmental changes (Jay et al., 2012). For the simulations, we selected two environmental variables for current and future scenarios: minimum temperature of the coldest quarter and precipitation of the wettest quarter, which explained 66% and 26%, of the variation among populations, respectively, based on the coefficient of variation. In addition, we selected topographic variables that may affect the dispersal and occurrence of anurans, such as slope and altitude (Table 1). We performed simulations for mitochondrial and nuclear haplotypes separated and excluded rare haplotypes due to convergence analyses (Table S4 in Appendix S1).
In order to define the rare haplotypes, we made an abundance rank curve and selected those with 85% of the relative frequency, thus excluding haplotypes present in less than 15% of the individuals.
We simulated genetic clusters under current climatic conditions combining genetic, climatic and topographic variables. The analysis was performed using Markov Chain Monte Carlo (MCMC) implemented in POPS (Jay et al., 2015). The MCMC runtime was set for 50,000 sweeps and the burn-in period for 5,000 sweeps using models with admixture. To define the number of clusters (K), we run the simulations four times for each K varying between 2 and 23 for both mitochondrial and nuclear data. The subset of runs minimizing the deviation information criterion (DIC, Spiegelhalter et al., 2002) and the lowest DIC values were selected. Finally, we predicted the genetic clusters for both future scenarios (4.5 and 8.5 RCP) conditioned on future climatic and topographic conditions. We measured the shifts in ancestry between contemporary and predicted ances-   Table S1 in Appendix S1.
by their number of individuals to take into account variation in sampling effort. Additionally, we excluded populations with less than 5 individuals to avoid underestimation of genetic richness due to low sampling effort. Because most of sampled populations were outside PAs, we applied a kriging interpolation of the genetic richness and PGeoD onto a raster of the current species distribution predicted by the ENM, with a resolution of 0.5º ( Figure S4 in Appendix S4).
For analysis, we used a spherical semivariogram implemented in the R package gstat (Gräler et al., 2016). We also tested the correlation between genetic richness and PGeoD considering spatial autocorrelation (Dutilleul et al., 1993), implemented in the R package SpatialPack (Vallejos et al., 2018). We calculated the median correlation of PGeoD and genetic richness for each coalescent tree from a distribution of 1,000 trees and their 95% quantile interval.
We obtained PAs polygons from Protected Planet website (www. prote ctedp lanet.net). For each PA, we extracted the interpolated genetic richness that overlapped with their polygons. Then, we performed PAs clustering using k-means method (Hartigan & Wong, 1979), implemented in k-means function of R package stats (R Core Team, 2020).
The clusters were optimized by minimizing the sum of squares within groups. The number of groups was optimized by a grid search ranging from 2 to 15 clusters. The same procedure was repeated with PGeoD. To evaluate PAs effectiveness to conserve genetic richness, we compared the mean genetic richness within and outside PAs by generating 10,000 random samplings with replacement of the genetic richness from cells outside PAs with the same number of cells within PAs. Subsequently, we calculated a t statistic and calculated the probability of observing the difference in mean genetic richness within and outside PAs with 95% confidence level.
To evaluate the effectiveness of PAs to conserve PGeoD, we cal- negative values indicate that less evolutionary information is protected in PAs than expected by chance. We ran Mean Phylogenetic Diversity analysis implemented in the R package picante (Kembel et al., 2010). All analyses in R were run in version 3.6.3 (R Core Team, 2020).  Figure S2 in Appendix S4).

| Ecological niche modelling: forecasting
Hierarchical ANOVA showed that ENM algorithms were the main factor of variation (a median of 87%) in the predicted suitability maps (Table S5 in Appendix S1), but areas with higher suitability uncertainty were on cells outside S. squalirostris' distribution ( Figure   S3 in Appendix S4). Variation of suitability within S. squalirostris' F I G U R E 3 Geographical distribution of nucleotide and haplotype genetic diversity, number of haplotypes and Phylogeographical diversity for Scinax squalirostris populations. (a) Haplotype diversity (h); (b) Nucleotide diversity (π); (c) Number of haplotypes (k); (d) Phylogeographical diversity (PGeoD). For details on population codes and localities see Table S1 in Appendix S1 [Colour figure can be viewed at wileyonlinelibrary.com] distribution was explained by ENM algorithms and changes through time ( Figure S3b,c in Appendix S4), representing the effects of climatic changes on S. squalirostris' distribution.

| Simulations of genetic and phylogeographical diversity under climatic changes
Populations of S. squalirostris showed high haplotype and nucleotide diversities (Table 1; Figure 3a-c). Highest diversity was found across South Brazil, matching the areas of higher suitability for S. squalirostris (Figures 1 and 3). Populations showed low phylogeographical diversity ranging from 8% to 23% of total PGeoD (Table 1; Figure 3d).
Our simulations of extinction thresholds (Figure 4) showed that ~50% of the genetic diversity (e.g. number of haplotypes and PGeoD) will be potentially lost with a threshold higher than 0.6 for both climate change scenarios (Figures 4a,d, 5 and 6). The decrease in genetic diversity was steeper under 8.5RCP scenario. Haplotype and nucleotide diversity had similar pattern of diversity loss, increasing diversities at higher extinction thresholds, until haplotype and nucleotide diversities collapsed due to the extinction of most populations (Figure 4b,c).

| Forecasting the dynamics in genetic clusters
We found 4 mitochondrial genetic clusters (Figure 7a; K = 4, DIC = 1,040.66; Tables S6-S8 in Appendix S1) and 9 nuclear genetic clusters (Figure 7e; K = 9, DIC = 765.39; Tables S9-S11 in Appendix S1), both conditioned to environmental and topographic variables (e.g. precipitation, temperature, altitude and slope). The correlation between estimated and predicted admixture coefficients of geographical and environmental covariates was 0.99 for both mitochondrial and nuclear data, indicating that predictions of environmental variables were accurate.
Overall, spatially explicit simulations showed loss of genetic variability over time, due to environmental climate changes (Figure 7).
Three out of four climatic clusters were potentially lost for mitochondrial DNA (clusters 2 and 3) for both scenarios of climate changes ( Figure 7b,c), leading to a genetic homogenization across S. squalirostris' geographical range (Figure 7c). Nuclear data showed a geographical displacement in clusters (Figure 7e-g) and loss of clusters for both scenarios, remaining only one cluster (cluster 1; Figure 7f,g). In addition, the simulations evinced a homogenization in the spatial distribution of the genetic variability in Central, Northeast and East Brazil for 8.5 RCP (Figure 7g). Genetic cluster 1, for both mitochondrial and nuclear sequences, remained constant for both scenarios. Ancestry turnover analysis showed congruent response to climate changes and spatial genetic structure, with high spatial shifts in genetic clusters for both mitochondrial and nuclear sequences (Figure 7d,h).

| Conservation of genetic and phylogeographical diversity
Genetic richness and PGeoD showed positive and weak correlation (median ρ = 0.25, 95% quantile interval = 0.04-0.46). Interpolated genetic richness within PAs was classified into 5 groups (Figure 2; Figure S5 in Appendix S4; Table S12 in Appendix S1), as well as the interpolated PGeoD (Figure 2b,c; Table S12). PAs clusters were classified into 3 size classes: smaller than 0.2 km 2 , 0.2-0.4 km 2 and larger than 0.4 km 2 . PAs smaller than 0.2 km 2 were split into 3 diversity clusters. Genetic richness and PGeoD were different in smaller PAs  Table S1 in Appendix S1 [Colour figure can be viewed at wileyonlinelibrary.com] F I G U R E 6 Shifts in geographical distribution of phylogeographical diversity across the simulations under different extinction thresholds of 8.5 scenario RCP for the 26 populations of Scinax squalirostris. For details on population codes and localities, see Table S1 in Appendix S1 [Colour figure can be viewed at wileyonlinelibrary.com] F I G U R E 7 Spatial distribution of genetic clusters conditioned to environmental variables simulated for Scinax squalirostris, for mitochondrial (a, b and c) and nuclear data (e, f and g), for present-day (a and e), 4.5 RCP scenario (b and f), 8.5 RCP scenario (c and g). Each colour set represents a cluster in the figure legend. Shifts in ancestry between contemporary and predicted ancestry coefficients for both future scenarios based on the simulation of genetic clusters of all individuals of Scinax squalirostris. (d) Shifts in ancestry for the four genetic clusters based on mitochondrial data. (h) Shifts in ancestry for the nine genetic clusters based on nuclear data [Colour figure can be viewed at wileyonlinelibrary.com] clusters, but did not differ in larger PAs. Genetic richness within PAs cells was not different from cells outside PAs (t = −1.08, p = .27). On the other hand, PAs are conserving fewer PGeoD than expected by a random draw of lineages (ses.MPD = −9.7, p < 1 × 10 -4 ).

| D ISCUSS I ON
Our findings suggest that climate change can potentially affect S.
squalirostris, leading to an expansion of the geographical range towards southern South America, but a retraction in the northern part of its range. Although the range expansion may seem to be a positive response, this outcome will cause losses of both genetic and phylo- Brazil and Central Argentina (Frost, 2021), show a reduction of 31%, 43% and 52%, respectively, of suitable areas by 2050 (Vasconcelos & Nascimento, 2016) suggesting that a reduction in suitable areas in the future is a common outcome for amphibians in response to climate change, most likely due to their sensitivity to environmental conditions (Lopez-Alcaide & Macip-Ríos, 2011;Schivo et al., 2019;Velasco et al., 2021;Vilela et al., 2018;Zank et al., 2014). However, it may have hidden genetic consequences.
Our results pointed out the potential increase in the nucleotide and haplotype diversity as populations with low genetic diversity become extinct, remaining populations from the southern distribution, which have high genetic diversity. However, the number of haplotypes and PGeoD will decrease with the extinction of populations, evidencing a potential genetic filtering at the species level, leading to the persistence of few lineages and haplotypes, and homogenization of genetic ancestry remaining only one genetic cluster for mitochondrial and nuclear DNA. Thus, although some populations will have high genetic diversity, S. squalirostris will lose haplotypes, diversity of evolutionary lineages and genetic differentiation among populations, which may hinder the species response to climatic changes and increase its risk of extinction (Foden et al., 2013;Wright et al., 2008).
Temperature regimes are expected to change dramatically under the 4.5 and 8.5 RCP scenarios resulting in the homogenization of genetic clusters throughout the S. squalirostris' range. S. squalirostris is adapted to cold and dry climates. Hence the temperature increases and changes in precipitation regimes predicted for the EOC might affect the geographical distribution and the spatial patterns of genetic diversity distribution. S. squalirostris reproduces and its tadpoles develop in several types of waterbodies, mainly in temporary ponds (Vaz-Silva et al., 2020). Temperature increase may modify the hydroperiod, increasing evaporation rates thus decreasing the time for reproduction and development of the tadpole (Mathews, 2010). Amphibian tadpoles may show plasticity in relation to changes in hydroperiod, by accelerating metamorphosis rates, which may cause the adult size and life span to shrink (Liess et al., 2013), leading to population decline (McMenamin et al., 2008). Furthermore, ectothermic species tend to have limited physiological ability to adjust to climate changes (Sheldon, 2019), due to the low plasticity in thermal tolerance (Gunderson & Stillman, 2015).
However, some species may circumvent this physiological constrain by changing their activity period (Sheldon, 2019).
Despite the fact that genetic diversity is a surrogate of species response to environmental changes (Urban et al., 2013), fitnessrelated traits are expressed by a polygenic regulatory network (Franks & Hoffmann, 2012), which needs a mapping at the genomic level to properly assign the loss of genetic variability to a reduction in climate change adaptability. Also, interpolated macro-scale variables do not capture micro-climatic and microhabitat conditions that may provide suitable local environments acting as buffers against extreme climate changes (Jucker et al., 2020;Scheffers et al., 2014).
The success of microhabitat buffering may vary according to the species' life history. Free-living larval stage species, such as Scinax squalirostris, may be less vulnerable to climate change than directly developed amphibians (Scheffers et al., 2013). Nonetheless, microhabitat need to be protected against land use conversion to be effective. Microhabitat fragments will not support suitable populations for the long term as evidenced by macroecological and macroevolutionary studies (Rolland et al., 2018;Souza et al., 2019).
Scinax squalirostris is distributed in a region with a high number of PAs, which may protect its microhabitats. However, the current network scheme of PAs is not protecting its phylogeographical diversity.
As a consequence, S. squalirostris may lose the genetic information accumulated along with its evolution in response to climate changes (Hoffmann & Willi, 2008). Furthermore, there is no guarantee that S.
squalirostris will track the suitable climatic conditions, due to the high speed of climate changes (Arenas et al., 2012;Chen et al., 2011), the low dispersal capacity (Blaustein et al., 1994;Hillman et al., 2014), or failing to disperse through natural or anthropogenic barriers (Urban et al., 2012). The forecasted expansion of the species distribution is towards a region with high urbanization (Lopez et al., 2015), intensive crop farming (Laufer & Gobel, 2017;Lopez et al., 2015;Moreira & Maltchik, 2015) and cattle grazing (Lopez et al., 2015;Medan et al., 2011). Such anthropogenic disturbances may compromise population persistence in the future.  Magalhães et al., 2017;Quan et al., 2017). Therefore, establishing more PAs in S. squalirostris southern distribution is highly important to maintain its genetic and phylogeographical diversity in the future.
Given the low density of PAs in the South of Brazil, the expansion and implementation of new Priority Conservation Areas will benefit other amphibian species (either widespread or restricted) that will also be negatively affected by climate change (see also Schivo et al., 2019;Vasconcelos & Nascimento, 2016;Zank et al., 2014).
The reduction in suitable areas and the decline of amphibian species within PAs due to climate change are of major concern in conservation planning (Loyola et al., 2008(Loyola et al., , 2014Vasconcelos et al., 2018). Scinax squalirostris' populations from Central-West and Southeast Brazil have distinct morphological and acoustic characteristics compared with populations from South Brazil (Faria et al., 2013;Giaretta et al., 2020;Pombal et al., 2011). Our findings show that populations from the Central-West and Southeast Brazil have lower genetic diversity and will have lower suitability in the EOC. Therefore, lineages from these areas will be highly threatened and have a higher probability of going extinct. The loss of cryptic lineages may also affect ongoing diversification processes and hence future biodiversity (Bálint et al., 2011). Therefore, under more unpredictable and rapid climate changes, protection of the genetic and phylogeographical diversity becomes increasingly important for the survival of populations and species (Coates et al., 2018;Lourençode-Moraes et al., 2019;Moritz & Faith, 1998). It is important to be aware that there is a gap in genetic sampling from Southwest populations that could be underestimating the number of genetic clusters and phylogeographic diversity. We also acknowledge the limitations of our results because we used a low number of loci. For instance, a higher number of loci would provide more lineages, increasing the phylogeographical diversity. However, although we sequenced 1,000 bp, including both mitochondrial and nuclear DNA, we found 63 haplotypes for mtDNA and 73 haplotypes for nuclear DNA, in 219 individuals. With relatively short sequences, we obtained high polymorphism and genetic diversity. Furthermore, despite the number of loci, these regions have been successfully used to access genetic and phylogeographic diversity in amphibians (Barrow et al., 2020;Fusinato et al., 2013;Mota et al., 2020). We believe that the genetic data, which are widely used in the literature, together with the sampling effort of S. squalirostris throughout its geographical distribution (26 locations), can provide a favourable dataset to answer our questions about its occurrence and impacts on the genetic diversity due to climate change.
In conclusion, our results show that populations of S. squalirostris can, at first sight, be positively affected by climatic changes, by increasing their range size and shifting suitable areas in the EOC.
However, the retraction of suitable areas in the northern distribution boundaries may decrease its genetic and phylogeographical diversity, along with the number of genetic clusters, compromising the response of the species to the fast environmental changes.
Further studies covering a larger genome sampling will improve our understanding of S. squalirostris adaptation to climate change. The current geographic distribution of the species encompasses a high number of PAs, however, the predicted range change in the future might shift its distribution away from the currently established PAs, encompassing a lower number within its range, thus protecting less of their genetic and phylogeographic diversity. Finally, our findings point out the need of establishing PAs in South Brazil and Northeast Argentina as a strategy for long-term conservation of S. squalirostris.

ACK N OWLED G EM ENTS
We are very grateful to editor Dr Henri Thomassen and the three anonymous reviewers, whose suggestions provided major improvements to the paper. This work was supported by grants from CNPq (project no. 475333/2011-0) and CAPES/PROCAD (project no. and RGC with grants and fellowships, which we gratefully acknowledge.

CO N FLI C T O F I NTE R E S T
The authors declare no conflict of interest.

PE E R R E V I E W
The peer review history for this article is available at https://publo ns.com/publo n/10.1111/ddi.13299.

DATA AVA I L A B I L I T Y S TAT E M E N T
The GenBank accession numbers of genetic sequences used for analyses can be found in the online Supporting Information (Appendix S2). The maximum clade credibility tree is in Appendix S5.
The geographical coordinates of Scinax squalirostris occurrences are in Appendix S6.