Natural history predicts patterns of thermal vulnerability in amphibians from the Atlantic Rainforest of Brazil

Abstract In the Brazilian Atlantic Rainforest (AF), amphibians (625 species) face habitat degradation leading to stressful thermal conditions that constrain animal activity (e.g., foraging and reproduction). Data on thermal ecology for these species are still scarce. We tested the hypothesis that environmental occupation affects the thermal tolerance of amphibian species more than their phylogenetic relationships. We evaluated patterns of thermal tolerance of 47 amphibian species by assessing critical thermal maxima and warming tolerances, relating these variables with ecological covariates (e.g., adult macro‐ and microhabitat and site of larval development). We used mean and maximum environmental temperature, ecological covariates, and morphological measurements in the phylogenetic generalized least squares model selection to evaluate which traits better predict thermal tolerance. We did not recover phylogenetic signal under a Brownian model; our results point to a strong association between critical thermal maxima and habitat and development site. Forest species were less tolerant to warm temperatures than open area or generalist species. Species with larvae that develop in lentic environment were more tolerant than those in lotic ones. Thus, species inhabiting forest microclimates are more vulnerable to the synergistic effect of habitat loss and climate change. We use radar charts as a quick evaluation tool for thermal risk diagnoses using aspects of natural history as axes.

| 16463 CARILO FILHO et AL. international effort to estimate the extinction risk of species) may be facing declines and/or loss of populations (Ceballos et al., 2017).
The Brazilian Atlantic Forest is one of Earth's biodiversity hot spots, having been reduced to only 12.4% of its original area (SOS Atlantic Forest, 2017). It harbors more than 50% (ca. 650 spp.) of amphibian species recognized for Brazil (Myers et al., 2000;Rossa-Feres et al., 2017;Segalla et al., 2019). In the south of Bahia state, the Atlantic Rainforest is particularly important for amphibian conservation given the high number of species, several endemic (e.g., Dias et al., 2014;Mira-Mendes et al., 2018). The overall diversity of amphibians in this region remains understudied given the numerous recent additions to knowledge of regional diversity by descriptions of new species (Carnaval et al., 2009;Dias et al., 2017Dias et al., , 2020Orrico et al., 2018). Therefore, southern Bahia can be a considerate priority area for amphibian conservation taking into account the high diversity of species, lineages, distinct life history guilds, and current pressures as climate change and fragmentation (Campos et al., 2017(Campos et al., , 2020Loyola et al., 2008;Vasconcelos et al., 2014).

Both on a local (southern Bahia) and on a biome (Brazilian Atlantic
Forest) scale, habitat fragmentation and climate may act in synergism producing disconnection between the habitat from different amphibian ontogenetic stages (Becker et al., 2007). Therefore, amphibian conservation efforts can benefit from considering life history variation in habitat use, because microhabitat, activity windows, and sensibility may vary along with the ontogenetic development in amphibians (Enriquez-Urzelai et al., 2019).
Thermal exposure is among the main threats to amphibian conservation in the Atlantic Forest (Hof et al., 2011;Silvano & Segalla, 2005). Therefore, it is paramount to investigate ecophysiological aspects linked to the thermal niche (e.g., CT Max ) to understand the consequences of climate change for these animals (Li et al., 2013). In a warming Earth, ecophysiological studies provide data to assess climatic vulnerability of species, mainly in ectotherm animals such as amphibians . The critical thermal maximum (CT Max ) is the thermal point where the animal loses its locomotor ability to respond to a stressful thermal environment (Cowles & Bogert, 1944;Taylor et al., 2020). The use of CT Max combined with the environmental thermal data allows to infer the species level of sensitivity (e.g., warming tolerance) to future global warming scenarios (Deutsch et al., 2008;Duarte et al., 2012;Gutiérrez-Pesquera et al., 2016;Simon et al., 2015). Additionally, studies have recovered phylogenetic signals in thermal niche dimensions (e.g., CT Max and CT Min ) for amphibians (Gutiérrez-Pesquera et al., 2016;Hof et al., 2010). The presence of phylogenetic signal has been tested in tadpoles of anuran species from the Atlantic Rainforest of southern Bahia, but nothing is known about adults (Gutiérrez-Pesquera et al., 2016). Since the larval and adult stages are exposed to different environmental conditions, it is paramount to understand the physiological dimensions of both development stages and elucidate any phylogenetic signatures involved, and the implications for the species fitness (i.e., population growth rate).
In Brazil, research on the impacts of global warming on amphibians is still scarce (Winter et al., 2016). Herein, we evaluate patterns of upper thermal tolerance by accessing the CT Max in adults of 47 amphibian species from the Atlantic Forest of southern Bahia. We expect that macro-and microhabitat play a role in buffering against the exposition to high temperatures. Furthermore, we investigate patterns of thermal tolerance among species groups organized by macro-and microhabitats of adults, and site of larval development. We hypothesize that forest-associated species have lower heat tolerance than generalists and those who inhabit open areas. We also expect that arboreal species will reveal lower CT Max values if compared with terrestrial, fossorial, and cryptozoic species. We used a phylogenetic model to identify whether ecological, allometric, and microclimatic covariates can predict thermal tolerance of species. Therefore, we use model selection to evaluate whether the CT Max of the studied species (obtained from experiments with adults) results from an evolutionary process predicted by Brownian motion, as already predicted for CT Max of larval forms of species in the study region (Gutiérrez-Pesquera et al., 2016). Last, we use radar charts as a rapid assessment tool for thermal risk diagnoses, using CTMax as a ranking factor to categorize ecological groups according to their thermal sensitivity.

| Specimen collection and study area
We captured specimens from Atlantic Forest fragments in south-

Research Highlights
• The thermal vulnerability among the Atlantic Rainforest amphibians is correlated to ecological pattern of habitat occupancy.
• Forest species that reproduce in streams and whose adults are arboreal are the most susceptible to climate change.
Road, Ilhéus municipality), and Camamú Bay Islands (Maraú municipality). We performed visual encounter surveys through habitat and collected adult individuals of different amphibian species (orders Anura, Gymnophiona) by hand.

| Critical thermal maximum and warming tolerance
We subjected individuals to Hutchison's dynamic method (Lutterschmidt & Hutchison, 1997), through a constant and gradual increase in temperature until reaching the Critical thermal maximum (CT Max ) (see Taylor et al., 2020 for a review on thermal tolerance methods). We acclimated individuals for 72 hr at a room temperature of 25°C prior to assessing CT Max , following a 12-hr photoperiod of light and dark regimes. We conducted the experiments in a specimen-size experimental chamber with a 5-mm high dechlorinated water layer (controlling the risk of desiccation), covered with a net to prevent individuals from escaping. We placed the experimental chambers inside a water bath to create a homogeneous heating system for experimental trials.
We started the bath's temperature at 25°C and gradually increased it by 0.25°C min −1 (Gutiérrez-Pesquera et al., 2016). We used the "belly-up" condition (Brattstrom, 1968;Taylor et al., 2020) to diagnose CT Max for specimens by the lack of response to periodically employed stimuli (5 touches per min with a glass rod). Once the specimens were flipped over, exposing their venters, but remained immobile (not returning to normal position in 10 s), we considered the temperature of the chamber to have reached the CT Max and the gradual warming was stopped. Subsequently, the chamber bottom temperature was measured using a contact thermometer (Miller & Weber, Inc.; 0.1°C accuracy). The specimens were then placed in a 23-25°C water container for cooling. We monitored all individuals for 24 hr, to ensure the lethal temperature was not reached during experiments, and only data from those who survived this period were included in the analyses. The CT Max for each species was defined as the average of its sampled individuals. For each specimen, we measured snout-vent length (SVL) and head width (HW) with a 1-mm precision caliper after the experimental protocol (including the observation time) in order to minimize the influence of handling in their survivorship. We also weigh (W) each individual just before experimental procedures using a 0.01-g precision scale.
We obtained microclimate data (average temperature-T Mean ; maximum temperature-T Max ) by installing temperature data loggers (HOBO Pendant Temp/Light, UA-002-64) in an area of regenerating Atlantic Forest remnants (secondary forest) during the three-month (November-February) period with highest annual temperatures in the studied region. We installed a data logger in each microhabitat category used by amphibians (see section below). The warming tolerance (WT) was calculated from the difference between the CT Max of each species and the maximum temperature of the microhabitat they use (sensu Deutsch et al., 2008), and we projected the future vulnerability by using the most recent International Panel on Climate Change (IPCC) indices for more and less acute warming scenarios (RCP 8.5-4.8°C; RCP 4.5-2.6°C, respectively) available on

| Groupings and statistical analyses
We analyzed how ecological features of species covariate with CT Max through a perspective of macrohabitat use (forest-Fo; generalist-Ge; and open areas-Op) and according to the site of larval development (lentic; lotic; marsupial; and terrestrial). To analyze possible trends in WT, we divided the species into groups by microhabitat (i.e., arboreal, cryptozoic, fossorial, and terrestrial). We based the functional group allocation of each individual species (see Table S2), following the ecological aspects available in Haddad et al. (2013).
Within functional groups, CT Max was calculated from the average CT Max of all individuals of the species assigned to each group.
We compared our data on CT Max of the adult's specimens graphically with those available in the literature for their respective conspecific tadpoles looking for ontogenetic tendencies on thermal tolerance We constructed radar charts to evaluate the thermal vulnerability of species given ecological covariates (e.g., macrohabitat, site of larval development, and microhabitat). We used the mean values of CT Max obtained for each category (e.g., forest, generalist, and open areas) within the ecological covariates as a categorization parameter on a scale ranging from 1 to 3, where the value 1 was given to categories with the lowest average CT Max within each covariate and the value 3 was given to those with the highest average CT Max .

| Phylogenetic comparative methods
Comparative methods that take phylogenetic information under consideration are broadly used when data on several species violate assumptions of statistical independence. Phylogenetic generalized least squares (PGLS) models use phylogenies to make estimates on expected covariance among species data (Garamszegi, 2014). We performed regression analysis by PGLS in a Brownian evolutionary motion model, using the CAPER package (Orme et al., 2013) in the statistical Program R (R Core Team, 2020) to test the effects of allometric (SVL, HL, and W), ecological (macrohabitat and microhabitat), and microclimatic (mean and maximum temperature-T Mean and T Max ) variables on CT Max . To account for phylogeny in analyses, we used the topology of Pyron and Wiens (2011), restricting the phylogeny to 23 species common with our sampling (Figure 1). We admitted the position of the nearest taxon, similar to Gutiérrez-Pesquera et al. (2016), for three taxa (Aplastodiscus gr. albosignatus; Scinax gr. ruber; and Scinax cf. x-signatus) that were not included in the phylogeny of Pyron and Wiens (2011). We used the lambda value to verify whether the covariance between the variables used follows the evolutionary pattern predicted by a Brownian motion model (where λ = 1) (Garamszegi, 2014). We ranked and selected among competing models using Akaike's information criteria (AIC). The model AIC weight value (wi) was used to evaluate the relative likelihood of model support (Burnham & Anderson, 2002). Analysis was conducted in R (R Core Team, 2020).

| Environment and CT Max
We evaluated the environmental influences of CT Max of 47 species (Table S1). CT Max differed among species' macrohabitat groups (H = 26.77, df =2, p < 0.01, n = 308) (Tables S1 and S2  Table 2. CT Max values between groups with different larval development sites were also different (H = 27.7, df = 3, p = 0.00, n = 303) (Tables S1 and S2). Species with larval development in lotic habitats had lower CT Max values than species with larvae developing in lentic habitats (p < 0.01; Table 2). Our only marsupial frog species had lower CT Max than the other groups with exception of the lotic species (Table 1).
We sampled a single species of Caecilian (Siphonops annulatus, the only species in the terrestrial group). Its CT Max was higher than that of marsupial anurans and species with lotic development and lower than that of lentic species. Our results regarding the terrestrial and marsupial group must be interpreted with caution, once in our sampling we were only able to collect individuals of only one species for each of these groups.
We rank-scaled our ecological covariates according to thermal tolerance within radar graphs using each covariate as one axis and recovered ten general patterns for our sampled species. Eight patterns had at least one covariate limiting their tolerance and increasing the vulnerability (value 1 in chart scale) (Figure 2). The most vulnerable pattern is of terrestrial species, with lotic larval development, that inhabit forests. In contrast, species that combine fossorial adult microhabitat with lentic larval development, and are habitat generalists, presented the most tolerant profile.

| Warming tolerance
The WT of the 47 studied amphibian species differed according to  (Table 3). TA B L E 1 Variation of critical thermal maximum (CT Max ) and warming tolerance (WT) by ecological covariates. Means and their respective standard deviations (mean ± SD), Sample size (N), Lower and upper ranges of the 95% confidence intervals (lower 95 and upper 95, respectively), The Kruskal-Wallis chi-squared index (H) and analyses p-value (p)

| PGLS model selection
The best-supported model (AIC = 886, Wi = 0.63) indicated that amphibian CT Max is explained by adult habitat and larval development site (Table 4). We did not observe a phylogenetic signal on CT Max within the three best models.

| D ISCUSS I ON
Literature suggests that available macro-and microhabitats play an important role in thermoregulation and ecophysiology of forestassociated tropical amphibians (Nowakowski et al., 2018;Scheffers et al., 2013Scheffers et al., , 2014. Our results corroborate this notion because values of adult amphibian CT Max of forest species were significantly lower than those from generalist or open environment species. Canopy protection and the availability of thermal refuges seem to be crucial to avoid overheating of forest species due to their low heat tolerance (Scheffers et al., 2014). Thus, the species of amphibians that inhabit forests will depend on even less macro-and microhabitats available in a global context marked by habitat loss (Stuart et al., 2004) and the exposure of thermal refuges to rising temperatures (Ficetola et al., 2015).
Among the three best models, there was no association between the CT Max and the T Mean and T Max of the microhabitat. Literature suggests that the CT Max of aquatic ectothermic organisms (e.g., crustaceans, fish, and aquatic insects) is more influenced by average and maximum environment temperatures because terrestrial ectothermic organisms are able to use habitat heterogeneity to exploit different microenvironments in order to maintain their body temperatures independent of thermal averages from the air (Bogert, 1949;Gunderson & Stillman, 2015;Kearney et al., 2009;Stevenson, 1985;Sunday et al., 2011).
In a study conducted in Costa Rica, Frishkoff et al. (2015) found an association between the occurrences of lentic larval development species (i.e., puddles) with deforested areas, while those with lotic larval development (i.e., streams) or direct development seemed to prefer forests. In our study, species that reproduce in ponds tolerated higher temperatures than lotic species, suggesting that lentic species that can eventually use lotic water bodies for reproduction (e.g., some species from genus Boana, Aplastodiscus, and Rhinella) are better competitors than lotic breeders in warmer environments (Haddad et al., 2013). Due to the low sample size of marsupial and direct development species, we were unable to characterize variation in CT Max with confidence for these groups. However, species that do not depend on water bodies for their development depend on the humidity of the forest to prevent desiccation of eggs (Frishkoff et al., 2015;Scheffers et al., 2013).
The Atlantic Forest faces a historical crisis of deforestation and fragmentation (Moura et al., 2018;SOS Mata Atlântica, 2017;Wanger et al., 2020). As degradation usually implies increased edge effect on forest matrices (Kapos, 1989), and consequently a higher incidence of light, forest specialist taxa are exposed to microclimatic modifications of their thermal refuges (Nowakowski et al., 2018;Tuff et al., 2016). Once exposed to anomalous and stressful thermal conditions, forest specialist species amphibians in southern Bahia face F I G U R E 2 Radar charts showing ten patterns found for vulnerability according to ecological variables of macrohabitat, microhabitat, and larval development site. Each chart (a and b) is composed by five representative species for the general ecological patterns studied in the current study. Ecological categories were classified (1, 2, or 3) by CT Max (ºC) values (1-less tolerant species; 3-more tolerant species) according to the results of ecological covariates in the group analyses. Charts with smaller areas denote greater thermal vulnerability than those with larger areas F I G U R E 3 Effect of microhabitat on the warming tolerance of amphibian species from the Atlantic Forest of southern Bahia. Vertical bars denote 95% confidence intervals around the mean (point) of each functional group the risk of isolation and population disturbance from forest fragmentation, which associated with possible thermal stress could lead to local extinction events (Becker et al., 2007;Tuff et al., 2016). In a region with high amphibian richness and endemism rates (Carnaval et al., 2009;Dias et al., 2014;Mira-Mendes et al., 2018;Vasconcelos et al., 2014), such a scenario raises an alarm concerning the extinction risk of endemic species.
Amphibians are particularly sensitive to habitat loss and fragmentation (Becker et al., 2007), and regenerating environments (secondary forests) usually present lower species richness when compared to pristine environments (Thompson & Donelly, 2018 In the present study, the most tolerant species were those with cryptozoic or fossorial microhabitat. Nevertheless, other categories (arboreal and terrestrial) maintained their WT within tolerable physiological limits, according to current thermal signatures. We found 2-12 species that will be at risk of extinction given likely climate change scenarios in the Atlantic Forest considering the projections of temperature increase between 2.6 and 4.8°C (RCP 4.5 and RCP 8.5, respectively-IPCC, 2014). Terrestrial and arboreal species were more sensitive (lower WT) than cryptozoic and fossorial ones.
Given that microhabitat warming can also be a consequence of forest cover loss, cryptozoic and fossorial species are also endangered TA B L E 3 Species most sensitive to temperature increase in view of more and less pessimistic scenarios proposed by the IPCC (2014) (RCP8.5 and RCP4.5, respectively). Values in the table refer to the maximum critical temperature (CT Max ), species microhabitat, maximum critical temperature of microhabitats (T Max ), current heating tolerance (WT), and those based on warming projections of IPCC (2014) (WT RCP4.5 and WT RCP8.5 , respectively) by the increase in temperature and risk of desiccation (Kapos, 1989;Nowakowski et al., 2018;Tuff et al., 2016).
Even though our results point to habitat as the variable that best explains CT Max , the influence of the environment is not restricted to the upper tolerance limits. Diversity within families of anuran amphibians is mainly explained by species microhabitat, and to a lesser extent by the thermal niche they occupy (Moen & Wiens, 2017).
However, several studies have recovered phylogenetic signals in niche thermal dimensions for various groups, including amphibians (e.g., Gutiérrez-Pesquera et al., 2016;Hof et al., 2010;Olalla-Tárraga et al., 2011). Our results agree with Moen and Wiens (2017)  Given that our analysis did not recover any phylogenetic signal, to what extent does the adaptive potential of the analyzed Amphibia groups influence survival in a warming landscape? According to Moritz et al. (2012), although no significant variation was found in CT Max between peripheral and central lineages of the same species, there are differences among upper thermal limits between edge species and those within forest. Hypothetically, populations at the periphery of fragments may present genotypes that give them some resilience to future warming scenarios (see Moritz et al., 2012). However, phenotypic plasticity seems to have little ability to buffer harmful effects in a progressive context of rising global temperatures (Bellard et al., 2012). Terrestrial ectotherms have low acclimation capacity when it comes to temperature rises, and therefore are less tolerant to temperature changes than aquatic ectotherms (Gunderson & Stillman, 2015). We observed that differences between the adult CT Max (this study) and tadpoles of 14 species (Gutiérrez-Pesquera et al., 2016) confirm the previously presented pattern-tadpoles (aquatic environment) tend to present higher maximal CT Max than adults (environment specific) (Figure 4).
For tadpoles, Tejedo et al. (2012) point out that the ability to adapt to temperature increase is linked to the thermal environments experienced by lineages. Thus, generalist species (tolerating more variable thermal environments) would be better adapted to cope with temperature rises than thermally specialized taxa (e.g., forest and open-air species) . Although our data for adult amphibians also indicate such a pattern, it is still necessary to investigate whether larval forms exhibit similar behavior for thermal tolerance as a function of aquatic microhabitat. Although preliminary, we note that the tolerance patterns between larvae and adults signaled here point to a thermal separation between the adult and larval life stages (Becker et al., 2007).
By combining the data available from Gutiérrez-Pesquera et al. (2016) for species with low plasticity in tolerance limits (Gunderson & Stillman, 2015;Tejedo et al., 2012), progressive reduction of forest patches (SOS Mata Atlântica), and the projections of temperature increase (IPCC, 2014), we point to threat scenarios close to those projected by previous studies (Ceballos et al., 2017).
The consequences for such scenarios pass through population reductions (Becker et al., 2007) that can lead to loss of genetic diversity for populations (Ceballos et al., 2017;Moritz et al., 2012) and species (Carnaval et al., 2009), as well as extinction events (Alroy, 2015;Ceballos et al., 2017).
Overall, our results reinforce the general understanding that global warming is certainly a threat to biodiversity but produces different pressures according to the natural history of each species. The results obtained from the model selection confirmed the same trend. Therefore, the use of graphic components based on ecological variables, as well as the application of rigorously established scales, proved to be an efficient risk/vulnerability assessment strategy for the species studied here. The same approach may be used to estimate the susceptibility of other species, or even groups, to anthropogenic impacts.

F I G U R E 4 Distribution of CT
The restrictions imposed by the thermal niches must be taken into account for conservation actions (Araújo et al., 2013;Damasceno et al., 2014). Thermal limit data are only available for about one tenth of Atlantic Forest amphibians-ca. 65 species. (Gutiérrez-Pesquera et al., 2016;Simon et al., 2015;Tejedo et al., 2012;present study).
In addition to CT Max , other ecophysiological variables (e.g., volun-

tary limits) could improve vulnerability projections for Atlantic
Forest species if implemented in mechanistic niche models (Taylor et al., 2020). Our research shows that natural history plays a key role in thermal tolerance and thermal vulnerability in amphibians of the Brazilian Atlantic Forest and may be a proxy for thermal niche.
Considering the low latitudinal variation of CT Max and conservationism in the upper thermal limits (Araújo et al., 2013), our results can be used for extrapolations within this biome.

CO N FLI C T O F I NTE R E S T
None declared.