Amphibian spatial distribution modeling, what are we missing? Comparing adult-only and tadpole models of habitat suitability

Introduction

Species distribution models (SDMs) provide key information on the ecological requirements and potential geographic ranges of species (Peterson et al., 2011). In the simplest version of SDMs, species presence-only data are correlated with environmental predictors (Elith et al., 2011; Rojas-Soto et al., 2024). Due to their inherent association between distribution and both ecological and anthropogenic factors (Guisan et al., 2013; Eyre et al., 2022), SDMs have been central to the assessment of the conservation status of hundreds of species (Syfert et al., 2014). However, the geographic distribution is often a challenging feature to model, due to factors such as biological traits (e.g., low detection of cryptic or elusive species), data limitations (e.g., few or spatially-biased records), computing power, or even the lack of funding to conduct surveys across a species’ full range (Breiner et al., 2015; Simmonds et al., 2020; Lissovsky et al., 2021). Among these challenges, species with complex life cycles –characterized by distinct morphological stages with divergent ecological requirements (Moran, 1994; Minelli and Fusco, 2010)– are particularly difficult to model. Different life stages often occupy environments that differ markedly in structure or abiotic conditions. In this sense such species may exhibit different traits and ecological needs throughout ontogeny (Taboada et al., 2013). Yet, these ontogenetic dynamics are rarely incorporated into SDMs (but see Zurell, 2017; Stevens et al., 2023), despite their relevance for understanding variable physiological tolerances and environmental limits (Saenger and Holmes, 1992).

Anurans, comprising over 7,000 species (AmphibiaWeb) of which 40.7% are of conservation concern (Luedtke et al., 2023), represent a prominent case for this topic. Most anurans exhibit biphasic life cycles in which free-living, fully aquatic, larvae differ substantially from adults in physiological tolerances, ecological preferences, and habitat use (Turriago et al., 2015; Gutiérrez-Pesquera et al., 2016). This biological feature challenges the most commonly used approaches for distribution modelling, although SDMs have been widely used for amphibian conservation planning (Luedtke et al., 2023; Pottier et al., 2025a). Modelling a biphasic life cycle is not straight forward; to model the aquatic stage one can use water bodies boundaries, as is often used to model other aquatic species (using algorithms like AquaMaps and the Sea Around US Project; Robinson et al., 2017; Klaassen et al., 2025). On the other hand, ecoregions, altitude, temperature or other physiographical traits can be used to define the accessibility area of the terrestrial stage (Elith and Leathwick, 2009; Mammola et al., 2021). Furthermore, in climate change-related SDMs, the ability to disperse from one site to another must also be explicitly considered (Soberón and Peterson, 2005); but that ability differs strongly between life stages (Massol and Débarre, 2015; Ahmed, 2025). In this context, to correctly model the distributions of anurans with biphasic life cycles, one needs to explicitly consider how their complex life cycles constrain or widen their ecological niches along their ontogeny (Pineda and Lobo, 2009; Trumbo et al., 2012).

Physiological constraints of anurans also change in different life stages. For example, in some Hylid frogs, adults may spend long periods of time in the tree canopy far from water bodies and only reach rivers or streams during the reproductive season; whereas tadpoles rely completely on the formation and prevalence of water bodies (Lertzman-Lepofsky et al., 2020; Brooks and Kindsvater, 2022; Pareja-Mejía et al., 2023). Besides, differences in physiological optima between tadpoles and adults are widely recognized, with adults being usually more resistant to temperature changes than tadpoles (Ruthsatz et al., 2022), since water bodies have higher thermal inertia than terrestrial microhabitats (Ashcroft and Gollan, 2013). Adult mobility allows them to actively reach preferred temperatures by perching at different sites (Tracy et al., 2013), whereas tadpoles may change their activity patterns. Tadpoles can also move to shaded or open areas for thermoregulation purposes, but the buffering effect of water may make such behavioral control of body temperature less effective (Brattstrom, 1962; Ortega-Chinchilla et al., 2023). Moreover, larval development is highly sensitive to water temperature and hydroperiod, with warmer conditions reducing larval duration (Álvarez and Nicieza, 2002; Rivers-Moore and Karssing, 2014). Thus, incorporating physiological constraints and the dynamics of ephemeral water bodies is fundamental for accurately modeling the distributions of free-living anuran larvae.

Among the variety of SDMs approaches, physiology-based (mechanistic) models offer a promising avenue because they incorporate explicit environmental tolerances (Kearney and Porter, 2009; Newman et al., 2022). For ectothermic species, such as amphibians, the most studied physiological axis is temperature (Gerick et al., 2014; Cunningham et al., 2016), as it modulates the organisms’ pace of life, accelerating or slowing physiological processes. However, despite their utility, physiological datasets remain scarce, particularly for tadpoles, due to the logistical and experimental challenges associated with quantifying critical thermal limits and other physiological data (Pintanel et al., 2022; Pottier et al., 2025b). When available, these data can substantially improve SDMs by providing ecologically realistic constraints for both larval and adult distributions.

It is not surprising that most correlative SDMs of anuran species to date have been constructed with adult presence-only data. Nevertheless, the survival of recruits is fundamental for any species to persist in nature. Given the narrower geographical and physiological constraints of tadpoles compared to adults, we expect smaller suitable areas than those projected for adult individuals. To test this hypothesis, here we implemented SDMs for both life stages of ten Hylid frog species that occur in Mesoamerica under both present and future scenarios. Our goal is to improve our current understanding of the threats these species may face in the future whether such threats differ among life stages.

Methods

Study site and species

Mesoamerica is a biogeographic region located in the northern Neotropics, extending from northern Mexico to the Panama–Colombia border (Holland, 2012; Imbach et al., 2013). This region represents a convergence zone between Nearctic and Neotropical lineages, which, in a geological and evolutionary context (Cione et al., 2015; Bacon et al., 2015), has favored a high amount of endemism across multiple taxa, including anurans (Myers et al., 2000; Mittermeier et al., 2011; de Albuquerque et al., 2015; Noss et al., 2015). Although tropical species are often regarded as particularly sensitive to climate change, species inhabiting Mesoamerica remain underrepresented in global databases (Alves-Ferreira et al., 2025), and the region is similarly underrepresented in the literature on amphibian declines, despite its high biodiversity (Crawford-Ash et al, 2025). Approximately 157 Hylid species occur in Mesoamerica, of which a little more than half (88) possess sufficient presence data to be modeled. Most species lay eggs in ephemeral ponds, although some breed in permanent water bodies (Duellman and Trueb, 1994; Wells, 2007).

Given that one of our goals was to incorporate physiological limits into SDMs, we selected Mesoamerican Hylid species with available critical thermal maximum (CTmax) data in the AmphiTherm database (Pottier et al., 2025b). After that, we filtered species by occurrence availability, retaining only those with >30 unique presence points (Wisz et al, 2008), since low sample sizes can lead to unstable, unreliable, and highly overfitted SDMs (Elith and Leathwick, 2009; Owens et al., 2013). Additionally, from those species we discarded three which are microendemic and their presence data are highly spatially correlated, preventing us from adequately modeling their distribution (Boria et al., 2014). Finally, we incorporated two additional species for which we experimentally obtained thermal tolerance data (Dryophytes arenicolor and Smilisca baudinii). Our final dataset included 10 species: Boana rosenbergi, Dendropsophus ebraccatus, Dryophytes arenicolor, Rheohyla miotympanum, Scinax elaeochroa, S. ruber, Smilisca baudinii, S. fodiens, S. phaeota, and Tlalocohyla smithii, all of which reproduce in temporal waterbodies (floodplains mainly, but also in pastureland where puddles can form; Höbel, 2000).

Adult distribution modeling

We downloaded presence data for each species from the Global Biodiversity Information Facility (GBIF) and performed an exhaustive data-cleaning workflow due to the known uncertainties and limitations associated with GBIF records (Yesson et al., 2007; Nori et al., 2023). Data downloads were conducted using the “rgbif” R package (Chamberlain et al., 2025), and individual Digital Object Identifiers (DOIs) are provided in Supplementary Material (Table S1). Our cleaning procedure included: (i) removing records without geographical coordinates, (ii) excluding occurrences lacking associated bioclimatic information or falling outside the temporal range 1970–2024, and (iii) eliminating duplicate records within the same climatic grid cell (Perez-Navarro et al., 2021; Prieto-Torres, 2024).

Using the filtered presence data, we delineated species-specific calibration areas (“M”; Soberón and Peterson, 2005) representing the geographic space that each species has likely accessed over its history. This region reflects the limits imposed by dispersal capacity and geographic and ecological barriers (Barve et al., 2011). To define these areas, we intersected species’ records with polygon-based shapefiles of Terrestrial Ecoregions (Dinerstein et al., 2017) and Neotropical Biogeographic Provinces (Morrone et al., 2022), excluding any grid cells outside these two zones. We chose this approach following the performance assessment of Rojas-Soto et al. (2024).

To characterize the species’ ecological niche and distribution patterns, we used environmental data from WorldClim 2.1 at ~10 km² cell size resolution (Fick and Hijmans, 2017). We excluded BIO8 (Mean Temperature of Wettest Quarter), BIO9 (Mean Temperature of Driest Quarter), BIO18 (Precipitation of Warmest Quarter), and BIO19 (​​Precipitation of Coldest Quarter) due to known artifacts that generate spatial discontinuities (Escobar et al., 2014; Booth, 2022). Predictor variables were selected for each species using a Variance Inflation Factor criterion <10 to reduce multicollinearity, as implemented in the “usdm” R package (Naimi et al., 2014). Detailed information on the final set of variables per species is provided in Appendix S1. Each model was calibrated using environmental data from the entire M, to capture the full breadth of the species’ ecological niche (Chevalier et al, 2022). This prevents ‘niche truncation’, a source of bias where models calibrated on a geographic subset misrepresent the species’ environmental tolerances, leading to unreliable projections (Guisan et al, 2025). For species whose ranges extend beyond Mesoamerica, we used their full-range and then cropped the final suitability projection to our focal study area.

We implemented SDMs using MaxEnt 3.4.3 (Phillips et al, 2017) through the “kuenm” R package (Cobos et al., 2019) to evaluate multiple model configurations and select optimal parametrizations based on low omission rates (<5%) of training data (Anderson et al., 2003), the highest partial ROC values (>1.0; Peterson et al., 2008), and Akaike information criterion (AICc; Merow et al., 2014). MaxEnt estimates habitat suitability for each pixel (0 = unsuitable, 1 = highly suitable) based on the contrast between environmental values at presence records and background points (Elith et al., 2006, 2011). For each species, we generated 510 candidate models combining 17 regularization multipliers (β=0.1–10) and 15 feature-class combinations (combinations of linear, quadratic, product, and threshold). Final models were generated using 80% of the occurrence and evaluated with the remaining 20%, running 10 bootstrap replicates with 500 iterations and cloglog output format (Phillips et al., 2017). In addition, the defined environmentally suitable conditions were then projected onto future global climate scenarios (see below). We enabled the ext.type=“ext_clam” option in the kuenm package to detect areas with potentially novel environmental conditions under future climates (Elith et al., 2011; Peterson et al., 2018). For each climate scenario, we calculated the median suitability across replicates (Campbell et al., 2015) and then converted the continuous suitability maps (logistic output: 0–1) into binary presence–absence predictions using the 10th percentile of training presence threshold (Liu et al., 2013). All selected model characteristics, including sample size, GBIF DOI, and ROC statistics are presented in the Supplementary Material (Supplementary Table S1).

Tadpole distribution modeling

To model tadpole distributions, we incorporated the aquatic environments required during this life stage. Due to the lack of high-resolution information for the studied area, we developed a mechanistic layer representing pond availability and temperature across Mesoamerica, using the microclimatic model NicheMapR (Kearney and Porter, 2017). This model reconstructs local microclimates based on climatic and terrain information. Following the approach of Enriquez-Urzelai et al. (2019), we modeled potential ponds for every pixel in the study area at ~10 km² resolution, assuming a soil type composed of 90% clay to ensure water retention and using the pedotransfer() function to estimate hydraulic soil properties. Microclimatic simulations were run with the micro_global() function, using its soil moisture functionality and increasing the rain multiplier parameter (rainmult = 4) to simulate water accumulation. The model provides hourly estimates of pool depth and temperature (approximating saturated soil temperature) for an entire year.

Literature on the microenvironments used by tadpoles has revealed that their preferred depth is species-specific (Melo et al., 2018; Tavares-Junior et al., 2020). As such, according to our fieldwork experience, we determined a minimum depth of 10 cm that can be generalized for our ten study species. To determine the ideal ponds within the study area, we considered a temporal threshold of pond availability that would guarantee the developmental time leading to metamorphosis. Since such data are unavailable for most species, we used experimental developmental times of three representative Hylids: D. arenicolor (three populations: -104.08 N, 20.34 W; -97.56 N, 18.32 W; and -99.25 N, 19.77 W), S. baudinii (one population: -98.90 N, 18.02 W) and E. xera (one population: -97.55, 18.32). These species represent distinct range patterns and ecological strategies within Hylidae: D. arenicolor exhibits broad environmental tolerance across its range from Utah to Southern Mexico (Hernández-Herrera et al., 2024); S. baudinii is widely distributed throughout Mesoamerica (Duellman, 1970); and E. xera is a microendemic species, serving as a proxy for rare and narrowly distributed hylids (Canseco-Márquez and Gutiérrez-Mayén, 2010).

To assess time to metamorphosis, we collected 150 tadpoles in the developmental Gosner stage 25 (Gosner, 1960) of each of the three mentioned species and locations. After two days of acclimatation, we raised tadpoles in plastic containers (27.5×18×13.5 cm) with constant water level of 4 L and 10 individuals per container, under a 12:12 photoperiod and under three temperature regimes: 20 °C, 26 °C and 30 °C. We fed tadpoles with a 4:1 mixture of rabbit chow and fish flakes ad libitum and changed water weekly to ensure water quality (Hernández-Herrera et al., 2024). We registered the number of days required for each individual to reach metamorphosis.

To account for the potential intra-specific variation of developmental time, while also excluding outliers, we pooled the data from all organisms, divided the data into quartiles and selected the first and third quartiles as the limits of each species’ threshold (17–41 days). We next used the minimum value of this threshold to filter the simulated pond data. We selected the pixels in which a simulated pond remained filled for this amount of time, from March to November (excluding the winter period, when activity drops to a minimum). Finally, we created a binary pond raster where 1 represents ponds adequate to sustain larval development and 0 represents unsuitable conditions, which we used to constrain adult SDMs, thereby incorporating aquatic requirements into the tadpole distribution models.

Future climate scenarios and species distribution projections

We utilized the GCM_CompareR tool (Fajardo et al., 2020) to select adequate global climate models for our climate-change analyses. This tool compares projections from all available global climate models and identifies which models project conditions closest to the ensemble-mean and, how much each model deviates from present-day climates. For this assessment, we used climate data from the Coupled Model Intercomparison Project 6 (CMIP6; Stoerk et al., 2018) available at the Worldclim web portal. Based on these results, we selected two global climate models (GFDL-ESM4 and IPSL-CM6A-LR) for SDMs in future years (“2041-2060” vs. “2061-2080”). These two climate models have shown improved performance in representing precipitation, zonal-mean atmospheric fields, equatorial ocean subsurface structure, and El Niño–Southern Oscillation dynamics in the Americas (Cook et al., 2020; Zelinka et al., 2020). It is important to note that all projections were performed using the intermediate Shared Socio-economic Pathway scenario SSP3-7.0, which assumes high greenhouse gas emissions and limited climate-change mitigation (Riahi et al., 2017). Although some studies suggest that extreme scenarios (SSP5-8.5) may be more likely (Schwalm et al., 2020), we adopted a more conservative pathway consistent with recent IPCC assessments (IPCC, 2022).

We forecasted the potential impacts of global climate change for adults on each species by comparing present and future geographic projections for each species under the aforementioned scenarios. Such comparisons were classified into three categories: (1) stable suitable areas, where cells remained suitable under both current and future climates; (2) loss of suitable sites, where cells suitable at present became unsuitable in the future; and (3) gains of suitable areas, where initially unsuitable cells became suitable under future conditions. The areas of stability, gain, and loss were calculated under two dispersal assumptions: “contiguous dispersion” (in which species can colonize new suitable areas within their “M” region but cannot cross biogeographic barriers) vs. “non-dispersion” (in which species may retain or lose currently occupied areas but cannot expand into newly suitable sites, representing the most unfavorable scenario; Peterson et al., 2002). To simplify the interpretation only the results from contiguous dispersion are presented here, but the non-dispersion can be seen on the Supplementary Material.

To project tadpole future distribution, we downloaded current and future monthly averages of maximum air temperatures and precipitation from the WorldClim database (using the same future scenarios as before). We calculated the monthly differences between future and current conditions and averaged them to obtain a single change factor per pixel. These values were incorporated as offsets in the warm parameter of the micro_global() function in NicheMapR, allowing microclimatic simulations to account for future environmental changes. We reran the microclimatic model for each scenario to generate future pond-availability rasters (using the same temporal threshold of metamorphosis time), which were then used to trim the corresponding future binary SDMs, thereby incorporating larval aquatic requirements into the projections. As with adults, we compared the geographic projections between present and future tadpole distributions to identify areas of stability/gains and areas predicted to experience local extirpation.

Integrating CTmax to species distribution models

Once we obtained the binary maps of current and future distribution for each species, we detected areas with temperatures exceeding the critical thermal maximum (CTmax) using empirical values ​​from the AmphiTherm database (Pottier et al., 2025b). To do this, we retrieved maximum temperature layers from WorldClim at the same resolution and using the same global climate models applied in the SDMs. We aligned these layers with the final presence maps in each scenario to detect cells in which predicted temperatures exceeding species-specific CTmax. We then calculated the proportion of each species’ estimated distribution exposed to temperatures above Ctmax. The species-specific CTmax values are shown in Table 1.

Table 1. Projected changes in suitable habitat area and exposure to critical thermal maxima (CTmax) for Mesoamerican amphibian species under future climate scenarios and dispersion assumptions.

Habitat areas are given in km². Future projections are based on two Global Climate Models (GFDL-ESM4 and IPSL-CM6A-LR) for the periods 2041–2060 and 2061–2080. Two scenarios are shown: Dispersion, where species can colonize new climatically suitable habitats, and Non-dispersion, where species are confined to their current suitable habitats.

Physiological limits for tadpoles were also incorporated using CTmax (Kearney and Porter, 2009). Considering that these limits may vary across ontogenetic stages in amphibians (Ruthsatz et al., 2022), tadpole CTmax values may differ from those of adults (Sheridan and Bickford, 2011; Duarte et al., 2012). Direct measurements for most species were unavailable; however, given the relatively low evolutionary lability of CTmax within lineages (Araújo et al., 2013; Muñoz and Bodensteiner, 2019), we opted to use published tadpole CTmax values from 30 neotropical hylid species (Pintanel et al., 2022) as a surrogate, yielding a mean threshold of 40.65 °C. Using this value, each pixel of the simulated water-body layer (see above) was evaluated to determine whether CTmax was reached at any time during the larval developmental window. It is relevant to consider that such temperatures represent the upper limit for a limited time; therefore, if the temperature of the pixel reached or surpassed the limit, it was considered as inhabitable for tadpoles, restricting their distribution (Calosi et al., 2008; Sinclair et al., 2016).

Results

SDMs under current and future climate scenarios

Adult-based models showed wide differences in their projected current distributions, with S. eleachroa presenting the smallest estimated extent (4,360 km²) and S. baudinii the largest (117,770 km²). These differences remained similar when the present distribution of tadpoles were modeled (ranging from 170 to 6,330 km²). All species showed smaller predicted distributions for tadpoles than for adults (only 4% of the adult distribution on average), which is expected as they are constrained by water bodies availability by the layer constructed with NicheMapR. These results showed an average stage specific difference of 96% in the distribution of tadpoles when using only adult data for modelling.

Figures 1–3 show a graphical representation of the distributional changes for the IPSL-CM6A-LR Global Climate Model in the 2061–2080 period. Although all other maps are shown in Supplementary Figures S1–S9, we present and discuss the numerical results of all scenarios. Under future climate scenarios, two species (S. ruber and S. fodiens) consistently showed contractions under all scenarios, although these remained below 4%. Four species showed mixed results, with an increase under some Global Climate Models or temporality and a reduction with others (B. rosenbergi, S. eleachroa, S. phaeota, and T. smithii Figures 1–3). Finally, four species showed consistent increases in their distribution, regardless of Global Climate Model or temporality, which can be as little as 2.29% (D. arenicolor, IPSL-CM6A-LR 2041-2060), and as much as 130.02% (D. ebraccatus, IPSL-CM6A-LR 2061-2080; Table 1). In contrast, tadpole-based SDMs showed that, under climate change, all species experienced distributional reductions, with losses ranging from 5.00% (S. ruber, IPSL-CM6A-LR 2041-2060) to 62.03% (S. fodiens, GFD-ESM4 2061-2080; Table 2). No species showed expansion in larval habitat suitability under any future scenario.

Table 2. Projected changes in suitable habitat area for Mesoamerican amphibian species (tadpoles) under future climate scenarios, considering pond availability and temperature.

Habitat areas are given in km². Future projections are based on two Global Climate Models (GFDL-ESM4 and IPSL-CM6A-LR) for the periods 2041–2060 and 2061–2080.

Mechanistical distributions under CTmax constraints

Overall, when applying CTmax as a physiological constraint to adult distributions, four species (D. arenicolor, S. baudinii, S. fodiens, and T. smithii) lost portions of their present suitable range. These reductions were minor for some species (as low as 0.83%) but reached nearly one-third of the projected distribution for S. fodiens (27.97%). Besides, under future climate scenarios, mechanistic constraints produce heterogeneous responses. Four species (B. rosenbergi, S. ruber, S. eleachroa, and S. phaeota) remained unaffected across all global climate models and periods, as none of the projected units reached or surpassed CTmax. In contrast, six species exhibited additional range losses under CTmax constraints, ranging from 0.09% (for R. myotympanum under GFD-ESM4 2041-2060) to 77.44% (for Smilisca fodiens under IPSL-CM6A-LR 2041-2060). For species such as S. baudinii and T. smithii, these mechanistic reductions altered the patterns observed in correlative models, reversing previously projected expansions or mixed trends and resulting instead in net distributional contractions (Table 1). For tadpoles, applying CTmax had no effect on the projected distributions of S. ruber (consistent with adult-stage results) or R. myotympanum. For all remaining species, CTmax further reduced larval distributions, with losses reaching up to 51.9% in S. fodiens (Table 2).

Discussion

Understanding how species with complex life cycles respond to environmental change is essential for accurate biodiversity forecasting (Possingham and Roughgarden, 1990). Yet most SDMs still rely solely on adult records and ignore stage-specific ecological constraints and physiological limits. Our goal was to test how incorporating information on the biphasic life cycle of Mesoamerican hylid frogs —particularly, the more restrictive environmental and physiological requirements of tadpoles— alters current and future distribution estimates under climate change. By combining correlative SDMs with mechanistic constraints linked to pond availability and thermal limits, we found that most examined species may experience reductions in suitable areas, and that these areas may shrink further when projected temperatures exceed species-specific CTmax values. Together, these findings support our central hypothesis that failing to consider larval stages and their physiological thresholds may lead to an incomplete characterization of species’ potential distributions. Nevertheless, it is necessary to gather further data on tadpole distribution and improve water bodies layers to assess if the observed pattern is not an artifact of the data limitations.

Populations facing rapid climatic change may either adapt, track their niche through migration, or become extirpated (Bell and Collins, 2008; Pyron et al., 2015). Although informative as worst-case estimates, non-dispersal scenarios may be ecologically unrealistic for most anurans, as many species could colonize new suitable areas if landscape connectivity is maintained (Kindlmann et al, 2025). However, effective dispersal is uncertain. Climate change velocity is increasing globally, and if environmental change outpaces the intrinsically low dispersal capacity of amphibians, range tracking becomes unlikely (Smith and Green, 2005). Their limited home ranges, low vagility, and dependence on humid microhabitats further restrict movement and increase sensitivity to barriers and habitat discontinuities. Thus, even if dispersal is theoretically possible, geographic barriers, fragmentation, biotic interactions, and physiological limits often hinder effective movement (Thierry et al, 2024). In this context, evaluating both dispersal and non-dispersal assumptions is essential to capture species’ vulnerability under climate change.

Eight of the ten analyzed species showed projected range expansions under at least one global climate model or period, emphasizing the strong influence of dispersal assumptions on SDM outcomes (Peterson et al., 2011). However, incorporating CTmax altered these patterns for two species. Species unaffected by CTmax filtering were mostly concentrated in southern Mesoamerica, reflecting regional climatic stability and their relatively small ranges (Sandel et al., 2011), which reduce exposure to thermally unsuitable conditions. In contrast, widespread species exhibited considerable losses —up to 20% of their modeled range— once physiological thresholds were applied, and thermally sensitive taxa such as S. fodiens exhibited declines exceeding 70%. These results align with critiques that correlative SDMs often overlook physiological mechanisms that ultimately constrain species persistence, particularly thermal tolerance.

In widespread species, suitability is frequently overestimated because models include areas where individuals cannot physiologically tolerate extreme temperatures (Sunday et al., 2014). For those species with a restricted distribution, this overestimation is less severe because presence data are less spatially dispersed, allowing models to perform more conservatively. Integrating physiological information therefore improves the ecological realism of forecasts (Evans et al., 2015). Using CTmax as a mechanistic constraint is especially valuable because this trait is less evolutionarily labile than other physiological parameters that may shift over short time scales via acclimation or genetic accommodation (Muñoz et al., 2014; Muñoz and Bodensteiner, 2019). CTmax is also tightly linked to enzymatic and metabolic performance, and upper thermal limits typically produce steep performance declines beyond the tipping point they represent (Huey and Kingsolver, 1989). Although CTmin may also influence distributions, upper limits generally impose stronger constraints under warming climates (Araújo et al., 2013; Bennett et al., 2021).

Comparisons between adult- and tadpole-based SDMs reflect the contrasting ecological constraints of each life stage and highlight marked differences in data availability. Adult models captured the broader mobility of post-metamorphic individuals, which can access terrestrial habitats far from breeding sites and may temporarily tolerate suboptimal but non-lethal conditions while dispersing (Smith and Green, 2005; Sinsch, 2014; Fonte et al., 2019). In the case of tadpoles, distributions were tightly restricted by the availability, persistence, and thermal properties of water bodies, which cannot be circumvented until metamorphosis (McDiarmid and Altig, 1999). However, we like to point out the limitations associated with the lack of reliable layers to model ephemeral water bodies and the low number of tadpole records. Differences in adult and tadpole SDMs emphasize the relevance of explicitly modeling ontogenetic stages in species with complex life cycles, as each stage is exposed to distinct environmental pressures and habitat limitations (Petitgas et al., 2013; Vera-Candioti et al., 2023). Importantly, reproduction cannot occur outside suitable aquatic habitats, regardless of adult mobility, making larval constraints a key determinant of species’ climatic vulnerability. This is particularly relevant since amphibians are especially susceptible to land-use change and habitat fragmentation—two factors that reduce connectivity and eliminate breeding sites—integrating larval ecology into SDMs, is critical for accurate conservation assessments (Teixeira et al., 2025).

Data availability limitation is one of the main challenges for robust SDM predictions, current and future (Simmonds et al., 2020; Lissovsky et al., 2021). In the case of amphibians with complex life cycles, the first major limitation is the lack of reliable layers of ephemeral water bodies inhabited by tadpoles. Mechanistic models like the one we used to generate such layers are a good approximation; however, detailed modelling also has its own limitations, which we discuss later. Regarding presence records, most of them are from adults (Nori et al., 2025; Vera-Candioti et al., 2023). However, our results show that the omission of complete life-cycle data may lead to an overestimation of species distribution in both widespread and restricted-range species, although this needs to be confirmed when more reliable data is available. These findings highlight two urgent needs: tadpole records in different localities and full descriptions of tadpoles for accurate recordings (several species do not have such descriptions; Vera-Candioti et al., 2023). We argue that tadpoles, as a life cycle stage, keep similarities with the “rare species modeling paradox”, which states that rare species are usually the most at need for ecological modelling to better understand their distributions, but the hardest ones to model due to limited information (Lomba et al., 2010). Tadpoles have fragmented and temporally constrained habitat, which emphasizes the need for ecological modeling, considering the risk amphibians are facing currently (Luedtke et al., 2023). Tadpole data shortfalls will not be easily overcome in the near future; therefore, it is useful to use the information of as many species as possible, as well as their relationships, as an indicator for what other species may face in present and in the future (Mondanaro et al., 2023; Smith et al., 2019).

A deeper understanding of adult movement patterns and physiological limits is also needed to determine whether presence records represent residents or transient dispersers (Smith and Green, 2005; Sinsch, 2014), as conventional SDMs treat all occurrences equally. Anurans can move up to 15 km during mating migrations or dispersal events (Semlitsch, 2008), but that is an extreme scenario, and most species movements are constrained between 0.1 to 1 km (Pittman et al., 2014), whereas other studies have modeled that anurans move between 11 to 13 km (Smith and Green, 2005). Therefore, a good selection of the accessible area and data cleaning is essential, as differences in this feature could lead to differential results, mainly overestimating potential distribution areas under dispersion scenarios (Barve et al., 2011; Soberón, 2007). Yet, most correlative approaches still treat species as static entities and ignore population-level plasticity, acclimation, and demographic processes such as migration or biotic interactions (Benito-Garzón et al., 2019; Zurell et al., 2023). Instead, these models typically rely on simplified assumptions of either “full dispersal” or “no dispersal” abilities, as shown above (Peterson et al., 2002; Thuiller et al., 2019). These omissions can bias predictions, especially for species with free-living larvae, where developmental timing strongly influences survival and the ability to track changing climates (Newman, 1992; Buckley and Kingsolver, 2012; Svenning et al., 2014; Swab et al., 2015). Incorporating such traits is essential for improving projections under global change and for more accurately evaluating species’ conservation status (Pacifici et al., 2015; Urban et al., 2016).

The use of large-scale microclimatic modeling presents computational challenges, as simulations at fine temporal and spatial resolutions require substantial processing time. Also, microclimatic traits usually affect at micro, site, and even local scales (<10 m to 10 km, according with Pearson and Dawson, 2003), which may bias the results at big scales, like ecoregions. Nonetheless, microclimatic approaches are invaluable in regions like Mesoamerica, where detailed hydrological layers and high-resolution maps of temporary water bodies are lacking (Kearney and Porter, 2009, 2017). Because ephemeral ponds, floodplains, and temporary wetlands are essential breeding habitats for many hylids, their explicit modelling could be considered a necessary component of SDMs. Recognizing breeding sites as fundamental to the reproductive process —and accounting for their implications for species survival when modeling distributions (Cobos and Alonso Bosch, 2018)– is of the utmost importance, particularly given the widespread use of SDMs for conservation planning and protected-area prioritization (Vera-Candioti et al., 2023; Nori et al., 2025).

For some anuran and urodele groups, adult-based SDMs may still provide reliable approximations of species distributions, since environmental conditions experienced across life-cycle stages remain similar. This is the case for fully aquatic taxa (Pippa, Xenopus, Telmatobius, newts and axolotl) and direct developers (Brachycephaloidea and many Plethodontids species). In these groups, larval or juvenile stages occupy the same habitats as adults, and niche divergence across ontogeny is minimal. However, for other species with free living larvae, even those laying eggs in permanent water bodies (including many Ranids, Buffonids and other Hylids), adult-only models may overestimate geographic and environmental suitability. Therefore, species distribution and their projections must be taken cautiously, especially when SDMs are used to guide conservation prioritization or restoration planning. Reconsidering the physiological constraints of different life-cycle stages can be a necessary step to continue improving our understanding of species distribution, and a powerful tool for future conservation planning efforts in species with complex life cycles.

Generative AI statement

The author(s) declared that generative AI was used in the creation of this manuscript. The authors state that ChatGPT was used to assist with English translation and grammar improvement. All scientific content, interpretations, and conclusions are entirely the responsibility of the authors.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/famrs.2026.1755900/full#supplementary-material

Supplementary material Supplementary material at the source

Footnotes

Edited by: Luisa Maria Diele Viegas, Federal University of Bahia (UFBA), Brazil

Reviewed by: Jill Awkerman, United States Environmental Protection Agency (EPA), United States

Dillon Monroe, Southern Utah University, United States