Environmental, land cover and land use constraints on the distributional patterns of anurans: Leptodacylus species (Anura, Leptodactylidae) from Dry Chaco

Subtropical dry forests are among the most vulnerable biomes to land transformation at a global scale. Among them, the Dry Chaco suffers an accelerated change due to agriculture expansion and intensification. The Dry Chaco ecoregion is characterized by high levels of endemisms and species diversity, which are the result of a variety of climates and reliefs, allowing a wide variety of environments. The amphibian group exhibits a high richness in the Dry Chaco, which has been barely studied in relation to land cover changes. We used ecological niche models (ENMs) to assess the potential geographic distribution of 10 Leptodactylus species (Anura, Leptodactylidae), which are mainly distributed within the Dry Chaco. We characterized these distributions environmentally, analyzed their overlap with land cover classes, and assessed their diversity of ecoregions. Also, we evaluated how these species potential distribution is affected by the transformation of land, and quantified the proportional area of the potential distribution in protected areas. We found that temperature seasonality is the main constraint to the occurrence of the species studied, whose main habitats are savannas, grasslands and croplands. The main threats to these species are the effects of climate change over spatial patterns of seasonality, which could affect their breeding and reproduction mode; the loss of their natural habitat; the exposure to contaminants used by intensive agriculture and their underrepresentation in protected areas.


INTRODUCTION
Habitat destruction and fragmentation produced by changes in land use and land cover (LULC), and climate change are major factors influencing the global decline of populations and species (Bennett & Saunders, 2010). The humid tropics and amazon basin were the main focus of research and debate in relation to the effects of land transformation over biodiversity (Aide et al., 2013). Nevertheless, in Latin America, dry forests and savanna/shrub biomes are experiencing the second highest rate of absolute deforestation, near threatened by the International Union for Conservation of Nature (2016). The genus shows the highest diversity of reproductive modes within Leptodactylidae family, the most diverse family of Neotropical anurans (200 spp). Conservation priorities of a certain taxonomic group may also inform about the conservation requirements of other groups (Rodrigues & Brooks, 2007). Under this assumption and taking into account that the most diverse group may be representative of a wide group of taxa, the genus Leptodactylus might be a good indicator to assess the responses of sensitive species to changes in LULC in the Dry Chaco ecoregion.
Current techniques of ecological niche models (ENMs) allow relating species distribution data (species occurrence at known locations) with information about the environmental and/or spatial characteristics of the locations (abiotic factors). The environmental conditions of the localities where a species occur provide of just a partial image of the niche which can be represented in the geographic space; thus, they may be informative about the potential distribution of the species (Lobo, 2015). This association leads to confusion between the concepts of Species Distribution Models (SDM) and ENM, since it is natural to talk about SDM when inferences from the occupied area are involved. However, if we try to model potential areas, which essentially involves geographic localities with favorable conditions for the occurrence of a species (i.e. conditions contained in its existing fundamental or realized niche, but which may be present in other unoccupied regions), we must use ENM concept (Peterson & Soberón, 2012). ENM is a robust method to characterize regional species distributions (Seoane et al., 2006;Ficetola et al., 2010), offering reliable information regarding environmental constraints. These spatial analyses are the basis for assessing the effects of LULC changes over the potential distribution of species, and for proposing conservation strategies.
In this study, we estimated the potential geographic distribution of 10 Leptodactylus species, which are mainly distributed within the Dry Chaco; and we characterized them environmentally, i.e. by assessing which environment variables are the most relevant to determine the occurrence of these species. To spatially characterize the distributions of Leptodactylus species, we analyzed the proportion of the ecoregions and the different land cover types in the full range of the potential geographic distributions. We used existing maps of cultivated areas from 1976 to 2013 to evaluate how the potential distributions of these species have been affected by the expansion of the agricultural frontier. To evaluate the state of protection of Leptodactylus chacoan species, we quantified the proportion of the potential distribution of these species included within protected areas, focusing on species included in International Union for Conservation of Nature (IUCN) conservation categories.

Occurrence data
Presence points databases were constructed based on revised specimens of collections, databases of herpetological collections, bibliography and free databases. The material revised proceeded from the following herpetological collections: Argentina: Fundació n Miguel Lillo (FML), Museo Argentino de Ciencias Naturales (MACN), Museo de la Plata (MLP), private collections of Dr. María Laura Ponssa (L) and Dr. Julián Lescano (JL); Brazil: Museu Nacional de Rio de Janeiro (MNRJ), Museu de Zoologia de São Paulo (MZUSP), and Paraguay: Instituto de Investigaciones Bioló gicas del Paraguay (IIBP), Museo Nacional de Historia Natural de Paraguay (MNHNP). The database of Laboratorio de Genética Evolutiva (LGE, Instituto de Biología Subtropical, Misiones, Argentina) was also considered in the analysis. One hundred thirty five scientific articles were analyzed to extract presence points. When presence points were scarce, we completed species records with data obtained from Global Biodiversity Facility Information GBIF (accessed in March 2016). We cleaned this information to avoid mistakes from outdated taxonomic arrangements. All presence records are available in Supplemental Information 1.
We confirmed the species identifications when it was possible to avoid primary source of error in niche modeling. We identified the species recognizing external morphology features proposed in literature as descriptions of species, taxonomic and phylogenetic revisions (e.g. Heyer, 1970;Heyer, 1978;Heyer, 2014;de Sá et al., 2014). We used a binocular stereoscope when it was necessary.
We used the geographic coordinates of each specimen to georeference species occurrences, using Google Earth whenever these coordinates were missing. We estimated the uncertainties of these locations by considering the extent of the locality: each locality was defined as a circle, with a point marking the most likely position (the geographic center of the named place), and a radius representing the maximum distance from the point within which the locality is expected to occur (Wieczorek, Guo & Hijmans, 2004). The database was cleaned of records which were outside the known altitudinal range of each species and whose uncertainty was higher than 10,000 m. Sampling bias in geographic space (i.e. localities which are overrepresented in the collections) may lead to biases in the inferences about the environmental requirements of the species. To avoid this spatial bias and its consequent autocorrelation, we disaggregated occurrences of each species at 10,000 m with 'ecospat' package (Broenniman et al., 2014;Boria et al., 2014) implemented in R software. The conserved occurrence points were plotted in maps to check the consistency with the known ranges of the species. We analyzed a total of 2,306 localities, and the number of occurrences per species varied between 33 for L. laticeps and 467 for L. latrans (Table 1).

Ecological niche model and environmental characterization
To model the presence of a species, we used the theoretical Biotic-Abiotic-Mobility approach (BAM) of Soberon & Peterson (2005), which captures and links the geographic and environmental dimensions of species distributions (Peterson & Soberón, 2012). The BAM model is composed of three components: 1) "B" components are the biotic conditions, i.e. the appropriate suite of both present species (e.g. food) and certain absent species (e.g. strong competitors, diseases). In general, this component is not included in the modeling process due to the fact that it is very difficult to make accurate spatial quantification of this kind of data (Barve et al., 2011). According to this and considering that Leptodactylus species exhibit a broad food spectrum, for which considering the component B becomes much more difficult, we do not included it in the present work.
2) The "A" components represent the abiotic conditions, e.g. bioclimatic variables; and 3) The "M" component is the region of the world which has been accessible to the species via dispersal over relevant periods of time. This heuristic scheme assumes that stable populations of a species will be found only in the intersection of the B, A, and M components, (B∩A∩M) (Soberon & Peterson, 2005).
To delineate the "A" component, we selected environmental variables (i.e. climate and soil variables) matching the temporal and geographic resolution of the studied species occurrences. Climate variables refer to temperature and precipitation. Since anurans are ectotherms, they are more vulnerable to environmental changes than other tetrapods (Duellman & Trueb, 1994). Temperature is the most pervasive factor affecting the rate of biological reactions and physiological processes of amphibians (Rome, Don Stevens & John Alder, 1992); while precipitation is the most important extrinsic factor controlling the seasonality of reproduction in anurans (Duellman & Trueb, 1994). As a matter fact, in most of the species of temperate latitudes, reproductive activity depends on temperature Table 1 Compiled records: number of compiled records from different sources (herpetological collections, private anuran collections, bibliographic records and Global Biodiversity Information data base) with geographic coordinates. Unique records: number of unique locality records. Checked records: number of locality records whose specimens were checked from the total number of unique records. Records not checked: locality records whose specimens were not checked from the total number of unique records. Collections: not checked specimens of scientific collections. GBIF: locality records obtained from Global Biodiversity Information Facility, Bibliographic records: locality records obtained from literature. Disaggregate records: number of disaggregate records to 10,000 m. and precipitation (Duellman & Trueb, 1994). Soil was hypothesized to be an important predictor of the distribution of arid-adapted anurans, since it can affect water loss rates during aestivation (Schalk, Montaña & Springer, 2015). Soil pH influences hatching success in Pseudophryne bibronii (Chambers, Wilson & Williamson, 2006), thus, soil characteristics could affect the success of the survival rate of terrestrial egg-laying species, such as most Leptodactylus species. Climate information was obtained from WorldClim database (Hijmans et al., 2005), which includes 19 variables of temperature and precipitation. Soil information was obtained from SoilGrids1km database (Hengl et al., 2014, through ISRIC-WDC Soils), which includes eight global soil variables that summarize different aspects of soil at six depths. Most of these variables are not independent; consequently, Pearson correlations were used to detect and exclude highly correlated variables (r ! 0.8). For soil variables, correlation analyses were made within each depth and between depths (2.5, 22.5 and 45 cm depth). When high correlations were detected, we kept the variables that most likely affect the ecology of anurans. Eight bioclimatic variables were retained: BIO1 = Annual Mean Temperature, BIO2 = Mean Diurnal Range (Mean of monthly (max tempmin temp)), BIO4 = Temperature Seasonality (standard deviation of mean month temperature Â 100), BIO5 = Maximum Temperature of Warmest Month, BIO6 = Minimum Temperature of Coldest Month, BIO12 = Annual Precipitation, BIO13 = Precipitation of Wettest Month, BIO14 = Precipitation of Driest Month. Because the correlations of each soil variable between depths were high (r > 0.9), only the variables of the most superficial depth (2.5 cm) were kept. Six soils variables were retained: 1) Sand content, 2) Coarse Fragments, 3) Soil organic carbon, 4) pH index, 5) Bulk density, and 6) Cation-exchange capacity. The environmental databases used presented a spatial resolution of 2.5′ (approx. 5 Â 5 km per pixel at the equator). A critical step for calibration, validation and comparison of ENMs is the definition of the "M" component; i.e. the extent of the region where the model calibration will be performed. The "M" component must be customized for each species and must be related to an explicit a priori hypothesis (Barve et al., 2011). It is reasonable to take the set of biotic regions within which a species is known to occur as the hypothesis of the M region (Barve et al., 2011). Here, the M region was defined for each species by identifying the ecoregions (sensu Olson et al., 2001) with known presences and by buffering the resulting region to generate a soft edge. The buffer size was between 0.1 and 0.5 , depending on each species data. A buffer of 0.5 was applied when the points were located marginally in one ecoregion, and a smaller buffer (less than 0.5 ), when the points were located more uniformly in the ecoregions. This prevents from including areas where the species had no access. We pruned the selected areas when they included some regions where species obviously did not occur. For example, since Leptodactylus species are mainly distributed in the lowlands, we excluded higher altitudes from the "M" regions according to the altitudinal range known for each species.

Species
ENMs were calibrated using a maximum entropy method, MAXENT v3.3.3K (Phillips, Anderson & Schapire, 2006) in batch mode implemented in R-package "ENMGadgets" (Barve & Barve, 2016). Maxent fits models based on the probability distributions that show maximum entropy (i.e. closest to uniform), subject to the correlations of known presences of species with the environmental conditions across the study area (Phillips, Anderson & Schapire, 2006). Most default settings were kept but we defined five bootstrap replications and raw outputs (i.e. continuous probability values). We randomly partitioned occurrence data in two subsets: 50% of the occurrences were used to calibrate the model, and the remaining 50% were used for model evaluation. Given the observed limitations of the Receiver Operating Characteristics Curve (ROC) approach we used a partial ROC area under the curve (AUC) to evaluate the performance of the model (Lobo, Jiménez-Valverde & Real, 2008). Partial ROC allows for differential weighting of omission and commission errors and focuses on meaningful predictions for model evaluation (Peterson, Papeş & Soberón, 2008). We restricted ROC space to predictions corresponding to omission error 5%, and we then randomly sub-sampled 50% of the available evaluation data 1,000 times and estimated the ratios of simulated AUCs and null expectations of AUC (Peterson, Papeş & Soberón, 2008) using Partial-ROC software (Barve & Barve, 2016). Probabilities were determined by direct count of null replicate frequencies of AUC ratios falling below the observed value.
To identify the variables which contributed the most to the model, we used a jackknife test on the original environmental variables. This test removes one environmental variable from the full list at a time and recalculates the model to quantify the contribution of each variable to the overall model performance.

Spatial analyses
The resulting ENM can be projected in the geographic space, thereby depicting the potential distribution area of each species (Araújo & Peterson, 2012). We assumed that the species analyzed are in ecological equilibrium (i.e. their populations are not migrating) and that the potential distribution obtained corresponds to their historical distribution.
The assumption of equilibrium of a population within its ecological niche (or single environment) implies ignoring mutation pressures from the environment and migration processes (Levene, 1953). State of stable equilibrium should not be viewed as a primary property of ecological systems, but it is a feature that can emerge from extrapolation to large spatial scales (DeAngelis & Waterhouse, 1987). The estimates of the geographic range extent obtained by ENM techniques have proved to be more representative than those obtained by traditional methods (e.g. minimum convex polygon (MCP)), particularly when data are scarce or when species are rare (Marcer et al., 2013;Pena et al., 2014;Syfert et al., 2014); and also avoid potential subjective bias of experts (Fourcade et al., 2013). We used a minimum threshold (Anderson, Lew & Peterson, 2003) to convert raw model outputs (i.e. continuous probability values) to distribution estimates (i.e. a binary representation). We modified this threshold-which, by default, includes 100% of the training data-to consider the potential error in the occurrence data, by choosing the highest threshold which included (100-E)% of the training data, assuming E = 5% (Table 2). E is an estimate of the proportion of the occurrence data that is likely to include the georeferenciation error (Peterson, Papeş & Soberón, 2008). This is a conservative method that minimizes the commission rate. This threshold defines that 95% of the observed presences are predicted as such, representing, in other words, the sensitivity of the model (Allouche, Tsoar & Kadmon, 2006). We then projected the distributional estimates of the ENM in space to obtain a quantitative potential geographic distribution. The spatial reference used for the analyses was Lambert Azimuthal Equal Area (LAEA) projection at a spatial resolution of approximately 5 Â 5 km per pixel at the equator. The analyses and maps were developed with QGis software v. 2.16.2 (QGIS Development Team, 2016).

Land cover and ecoregion characterization
To spatially characterize the potential distribution of each Leptodactylus chacoan species, we quantified the proportion of their extent which overlapped with ecoregions (sensu Olson et al., 2001) and with the different categories of land cover for Latin America produced by Blanco et al. (2013). The 24 land cover categories (e.g. broadleaf forests, needle leaf forests, mixed forests, shrublands, grasslands, water bodies, urban areas and croplands) are based on the FAO/UNEP Land Cover Classification System (LCCS). We also assessed the relationship between Leptodactylus occurrences and major habitat types, by assigning land cover categories in two groups "Forest habitats" (e.g. Subtropical broadleaf deciduous and evergreen forests, Tropical broadleaf deciduous and evergreen forests) and "open habitats" (e.g. Croplands, Tropical, Subtropical and Temperate shrublands, Tropical, Subtropical and temperate grasslands).

Land transformation in the Dry Chaco
We assessed whether the potential distribution of Leptodactylus frogs matches that of chacoan cultivable areas, i.e. suitable agriculture areas which are likely to be transformed to either pastures or croplands. In order to evaluate this, we estimated the proportion of each species distribution area affected by agriculture expansion, taking into account cultivated plots of the entire Dry Chaco region from 1976 to 2013 (Vallejos et al., 2014).

Protected areas
To explore the degree of protection of Leptodactylus species, we quantified the proportional protected area within the full range of the potential geographic extent of each species, taking into account the protected areas categories I to VI assigned by the IUCN (

Ecological niche model and environmental characterization
The partial ROC tests indicated significant predictive ability of the models (p < 0.001). For all species, bioclimatic variables (temperature and precipitation) were the most important predictors. The temperature variables (temperature seasonality, the maximum temperature of warmest month and the minimum temperature of coldest month) were the most important variables for seven out of 10 species. Only for one species' potential distribution (L. podicipinus) the precipitation (precipitation of wettest month) was important. For two species' potential distribution a soil variable (PH index) was the most important (Table 2).

Spatial analysis
The potential distribution areas ranged from less than 40,000 km 2 for L. laticeps to 7.5 million km 2 for L. fuscus, with values between 1 and 3.5 million km 2 for most species (Fig. 1). Leptodactylus fuscus and L. latrans were the species with the largest potential distribution ranges (above 5 million km 2 ) (Fig. 1).

Land cover and ecoregion characterization
The level of overlaps of the potential distribution areas of Leptodactylus species ranged from 30 to 72% with "forest habitats" and from 29 to 71% with "open habitats" (Table 3). Among forest habitats categories, the highest overlaps were found in Subtropical broadleaf deciduous forests, Tropical broadleaf evergreen forests and Subtropical broadleaf evergreen forests. Among open habitat categories, the highest overlap occurred with croplands.
All the species considered in this study distributed beyond the Dry Chaco ecoregion (Fig. 2). Leptodactylus laticeps, L. bufonius, L. latinasus, L. elenae and L. chaquensis are mostly Dry chacoan species, with at least 25% of their geographic range occurring in the Dry Chaco. These species potentially inhabit less ecoregions than the other species (Supplemental Information 2). On the other hand, the remaining five species (L. latrans, L. fuscus, L. podicipinus, L. mystacinus, L. gracilis) potentially inhabit in a greater number of ecoregions (Supplemental Information 2).

Land transformation in Dry Chaco
By the end of 2013, 1.8 million km 2 of Dry Chaco (19% of the natural area in the entire ecoregion) were transformed to agriculture (Vallejos et al., 2014). The percentage of remaining natural areas inhabited by Leptodactylus chacoan species in the Dry Chaco decreased every year from 1976 to 2013 due to agricultural expansion, although the greatest transformation occurred after 1996 (Fig. 3). The largest percentage of potential chacoan area loss was 16% (L. mystacinus and L. podicipinus), while for seven species the loss of chacoan area were between 15% (L. bufonius, L. gracilis, L. latinasus and L. latrans) and 14% (L. chaquensis, L. fuscus and L. laticeps). The minimum loss of potential chacoan area was 13% (L. elenae) (Fig. 3).

DISCUSSION
We found that the main constraining variable of the potential distribution of Leptodactylus chacoan species is temperature (temperature seasonality, maximum temperature of warmest month and minimum temperature of coldest month), and that these species are found in both forest and open areas. In the Dry Chaco ecoregion, they inhabit cultivable areas which are being transformed from natural vegetation to pastures  and croplands. Leptodactylus species are underrepresented under the current system of protected areas. This combination of features suggests that Leptodactylus chacoan species are vulnerable to both climate change and habitat transformation. This situation is especially critical for L. laticeps, which besides being already threatened is distributed mainly in the Dry Chaco and is poorly represented in the system of protected areas.

Ecological niche model and environmental characterization
Our results show that the most important bioclimatic variable explaining species occurrence (in 70% of the studied species) was temperature. Other studies have also shown that climate (especially temperature and rainfall) is a determinant factor for amphibians (e.g. Parris, 2004;Urbina-Cardona, Olivares-Pérez & Reynoso, 2006). Since amphibians are ectothermic organisms, environmental temperature is probably a major constraint to their ecophysiological traits; for example, their mobility and their energy balance may be affected by temperature (Bennett, 1990). Physiological traits turn amphibians strongly sensitive to environmental conditions (Ferder & Burrgren, 1992), and consequently to climatic change (Pounds, Fogden & Campbell, 1999;Pounds et al., 2006). Projected climate changes over the 21st century include rises in surface temperature, and in the frequency and duration of heat waves and the occurrence of extreme precipitation events (IPCC, 2014). Also, seasonality will be affected by more frequent hot temperature extremes and fewer cold temperature extremes over most land areas (IPCC, 2014). Particularly, in south-eastern South America, which includes Dry Chaco ecoregion, changes in precipitation produce wetter conditions (Christensen et al., 2007). Under these scenarios, it is expected that populations change their geographic ranges to coldest areas, i.e. to higher altitudes or latitudes (Harsch et al., 2009;Thomas, 2010; Chen et al., 2011), although in the case of amphibians, this is unlikely due to their low dispersal ability. Measurements of the home range size of amphibians show a mean home range size of 40 m 2 (Wells, 2007). Moreover, the capacity for amphibian behavioral thermoregulation is limited because the cooling effect of evaporative water loss from skin counteracts heat gain by basking (Hutchinson & Dupe, 1992). Studies revealed that the effects of climate change over amphibians involve changes in the timing of breeding of some species (Blaustein et al., 2010;Carey & Alexander, 2003), instead of changes in their distributions (Corn, 2005). For example, no colonization of higher altitude areas associated with declines of Anaxyrus boreas has been found (Livo & Yeakley, 1997). In subtropical areas, the period of activity of amphibians and their population dynamics are strongly influenced by climate seasonality, e.g. reproduction behavior of most species occurs in the warm season of the year (Conte & Machado, 2005). For many species of Leptodactylus genus, the reproductive season matches with wet and/or warm season, e.g. L. latinasus, (Gallardo, 1964;Canavero et al., 2008), L. bufonius (Cei, 1980), L. gracilis, (Canavero et al., 2008;Ximenez & Tozetti, 2015), L. mystacinus, (Gallardo, 1964), L. fuscus, (Martins, 1988;Prado, Uetanabaro & Haddad, 2005;Lucas et al., 2008), L. elenae (Prado, Uetanabaro & Haddad, 2005), L. chaquensis (Gallardo, 1987;Prado, Uetanabaro & Haddad, 2005;Canavero et al., 2008) and L. podicipinus (Cei, 1980;Vizotto, 1967;Rossa-Feres & Jim, 1994); although it was reported that populations of L. podicipinus from the Brazilian pantanal exhibit a continuous reproductive cycle throughout the year (de Almeida Prado, Uetanabaro & Lopes, 2000;Prado, Uetanabaro & Haddad, 2005). The Dry Chaco presents a monsoonal climate, in which rains are concentrated in the warm season. Since the reproduction of most amphibians is associated to water, modifications in water balance in the area could have significant effects on their breeding timing and behavior. Thus, our results suggest that climate change could affect Leptodactylus chacoan species in different and unpredictable ways.

Land cover and ecoregion characterization
In general, the Leptodactylus species studied inhabit both forest and open areas in similar proportions, with even a greater proportion in open areas. This result matches that of Torres et al. (2014), who found that the response of the probability of presence of all the analyzed amphibian species in Argentinean northwestern Dry Chaco to woody biomass was negative. The importance of open habitats for Leptodactylus species suggests that these species evolved in landscapes which were more similar to open habitats than to forests, a scenario proposed for the pre-European landscape of the Chaco (Morello & Saravia Toledo, 1959a;Morello & Saravia Toledo, 1959b;Adámoli et al., 1972;Adámoli et al., 1990;Bucher & Huszar, 1999), in which the substantial contemporary woodland advance has been a result of domestic livestock overgrazing (Morello & Saravia Toledo, 1959a;Morello & Saravia Toledo, 1959b;Adámoli et al., 1990). Thus, savannas might have been an important cover type where part of the chacoan biota evolved (Torres et al., 2014).
Among the open habitat categories, "croplands" was the most important for all species, which is also consistent with observational field data (e.g. Peltzer et al., 2006;Ponssa & Barrionuevo, 2008;Attademo et al., 2014;Prado & Rossa-Feres, 2014;Guerra & Aráoz, 2015). The IUCN categorization (2016) mentions that L. bufonius, L. chaquensis, L. elenae, L. fuscus, L. latinasus, L. latrans, L. mystacinus and L. podicipinus are or appear to be well adapted to anthropogenic disturbances, while Leptodactylus gracilis is found in anthropic areas only in southern Brazil (Heyer et al., 2004). The present study also shows that L. gracilis may inhabit anthropogenic areas, since 71% of their potential geographic range is covered by open habitats, from which 41% is covered by croplands. Moreover, L. gracilis, L. latinasus and L. mystacinus have been suggested to be invading species, due to their abundance in agricultural ponds, and their ability to constitute stable populations in these environments (Sanchez et al., 2013).

Land transformation in Dry Chaco
The dynamics of land cover transformation to agriculture shows that one half of analyzed Leptodactylus species may inhabit mainly cultivable areas. The percent loss of their potential distribution area within the Dry Chaco is near of the transformation rate of natural habitats to cultivated plots from 1976 to 2013. The largest percentages of decrease were found for L. mystacinus and L. podicipinus. This is particularly worrying for L. laticeps and L. bufonius, species which are distributed mainly in the Dry Chaco, which the loss of chacoan area represent the 14% and 10% of their entire potential distribution respectively. The entire dry diagonal region is suffering the same process of natural habitat loss described for the Dry Chaco, e.g. the Cerrado biome is being transformed for soybean production (Fearnside, 2001;Warnken, 1999). The Caatinga and the Cerrado (Beuchle et al., 2015), are also experiencing a continuous net loss of natural vegetation, and an expansion of sugarcane crops in the ''dry diagonal'' (Manzatto et al., 2009). This scenario represents potential threats to species with high percent loss in the Dry Chaco (mentioned above) and to those inhabiting the "Dry Diagonal," such as L. latrans and L. podicipinus. The projected increase in annual rainfall may accelerate deforestation rates due to agriculture expansion, which would have a significant effect in the conservation of Leptodactylus chacoan species. It has been observed that deforestation in the Chaco is to some extent controlled by annual rainfall (Grau, Gasparri & Aide, 2005). This may be especially relevant for the conservation of L. laticeps. At least 20% of the potential distribution area of this species, which is already threatened, has already been transformed to agriculture in the last 40 years. Additionally to habitat loss and/or fragmentation in the Dry Chaco, the fact that Leptodactylus anurans potential geographic range currently corresponds to croplands and cultivable areas raises concern about the potential threats to which they may be exposed, e.g. contamination. It has been shown that anuran species breeding within or around agricultural areas are usually exposed to pesticides (Peltzer, Lajmanovich & Beltzer, 2003;Peltzer et al., 2006), which produce teratogenic effects, growth and development retardation, and consequently, decreasing survivorship (Hayes et al., 2006). This could be a serious risk for these well adapted to anthropogenic disturbances Leptodactylus species, and particularly for Leptodactylus latinasus, L. chaquensis, L. latrans, L. mystacinus, L. gracilis and L. fuscus, whose reproduction habits involve the exploitation of soil depressions in cultivated areas (Vasconcelos & Rossa-Feres, 2008;Silva & Rossa-Feres, 2007;Guerra & Aráoz, 2015;Peltzer et al., 2006).
Most of the anuran species undergo ontogenetic niche shifts (aquatic larvae and terrestrial adults); which is the case of Leptodactylus anurans, whose clutch deposition mode involves aquatic or terrestrial foam nests, aquatic larvae, and terrestrial adults. Thus, anthropic landscape transformation may produce a "habitat-split" (Becker et al., 2007), i.e. a spatial separation between remnants of terrestrial habitat and breeding sites (Dunning, Danielson & Pulliam, 1992). Habitat split presents a strong negative effect on anuran species with aquatic larvae, and is a determinant factor of population size, structure, and distribution, acting within a single generation (Becker et al., 2007).

Protected areas
We found that the percentages of Leptodactylus species distributions represented in protected areas are scarce (11% or less). Consequently, most of their geographical ranges (i.e. approximately 90%) may be vulnerable to human disturbance. A similar pattern was found for threatened anurans of Northeastern Brazil (Adelophryne baturitensis, Adelophryne maranguapensis, Allobates olfersioides and Agalychnis granulosa) (Campos, Brito & Solé, 2013). Added to the underrepresentation of chacoan species of Leptodactylus in protected areas, it has been found that the current assignment of protected areas to IUCN categories does not correspond to the expected gradient of naturalness (Leroux et al., 2010;Bishop et al., 2004;Chape et al., 2005;Dudley, 2008). Leroux et al. (2010) showed that the global protected areas network lacks of strictly-protected areas with low human influence.
Leptodactylus laticeps is under Near Threatened and Vulnerable conservation categories by Cortez et al. (2004) and the Argentinian assessment (Vaira et al., 2012), respectively. Our results show that L. laticeps inhabits the Dry Chaco ecoregion almost exclusively, mainly the Subtropical broadleaf deciduous forest (a typical natural land cover of Dry Chaco). The land transformation analysis showed that this species inhabits cultivable lands inside the Dry Chaco ecoregion, which, in combination with the underrepresentation of L. laticeps in protected areas and the expansion of agriculture around chacoan protected areas, increases the risk of ecological isolation (Matteucci & Camino, 2012). In Copo National Park, an Argentinean protected area of Dry Chaco, the increase in human corridors (e.g. roads) between 1976 and 1988 has been followed by an increase of the surrounding parceled land between 1988 and 2007 (Matteucci & Camino, 2012). Furthermore, habitat connectivity is severely compromised in the Argentinean Dry Chaco, where there are no ecological corridors connecting protected areas (Burkart, 2007). In Paraguay, actions to keep habitat connections have been limited to proposals such as the Biodiversity Corridor of the Dry Chaco (Comisión Mundial de Áreas Protegidas Paraguay, 2007). This proposal includes protected areas located at the northern end of the potential distribution range of L. laticeps: Teniente Agripino Enciso, Médanos del Chaco and Defensores del Chaco National Parks, the Natural Monument Cerro Timane Cabrera, Guasú Natural Reserve; and the Kaa Illa National Park in Bolivia. Also, the protected area systems in Latin America are deficient in the administration of economic resources, equipment, human resources and legal and regulatory frameworks (Brown et al., 2006;Castaño-Uribe, 2008). For instance, in Latin America and the Caribbean region: 1) Paraguay presents the lowest number of human resources per unit of protected area and the lowest rate of investment resources of the nation by protected hectare; 2) in Argentina, the staff has no specialized functions, lacking of own equipment; and 3) from 2001 to 2006, both in Argentina and in Paraguay, the amount of money allocated to protected areas decreased (Castaño-Uribe, 2008). In addition to these global concerns, it has been reported that L. laticeps species suffer of international commercial exploitation for pet trade (Schaefer & Céspedez, 2012;Cortez et al., 2004).

CONCLUSIONS
To conclude, the main constraint to the potential distribution range of Leptodactylus chacoan species is temperature. Consequently, climate change associated to modifications of seasonality patterns could affect the breeding time and reproductive mode of these anurans. Land transformation for agriculture activities in the Dry Chaco exhibits an increasing trend in natural habitats rate of loss, where Leptodactylus chacoan species inhabit. This implies three potential adverse effects on amphibians: 1) habitat loss; 2) "habitat split," produced by landscape fragmentation; and 3) exposure to contaminants used by intensive agriculture activities. The state of Leptodactylus protection is deficient, due to the scarcity of protected natural habitats in the region, and the underrepresentation of the potential distribution range of Leptodactylus species within protected areas (on average, around 5% of their potential distributions). This situation implies a strong concern about L. laticeps populations, whose main threats are agriculture activities and hunting pressure. In the Dry Chaco, open areas, such as savannas, grasslands and croplands, are the main habitat of Leptodactylus species. However, the current conservation focus in the Dry Chaco is on forest habitats, for which conservation planning should be redesigned to take into account the natural savannas and grasslands, in order to more comprehensively protect and account for habitat heterogeneity.