Dissimilarity among Species and Higher Taxa of Amphibians in a Hotspot of Biodiversity and Endemism in the Neotropics

1. Introduction

In Mexico, the class Amphibia has a high global richness level, ranking fifth in the world in terms of diversity with 376 species, of which 67% (252 species) are endemic to the country (Parra-Olea et al., 2014). Mexico has the second most plethodontid species globally, with 132 species, of which 81% are endemic to the country [1]. At the national level, the state of Oaxaca has the greatest diversity of Amphibia, registering 38 genera (32 anurans, 5 caudates and 1 Gymnophiona) and 149 species (106 anurans, 41 caudates and 2 caecilians). Oaxaca contains 40% of the species, 76% of the genera and 80% of the families of amphibians in Mexico. Of this total, 102 are endemic to Mexico (62 anurans, 39 caudates and 1 caecilian), and 60 are endemic to Oaxaca (40.2%: 27 anurans and 33 caudates) (Mata-Silva et al., 2015). In the case of caudates, endemism is extraordinary since 80.48% are endemic to Oaxaca [2].

Unfortunately, amphibians are experiencing the worst extinction crisis in their history [1,3]; for example, in Mexico, more than 41% of amphibian species are critically endangered, endangered or threatened according to the International Union for the Conservation of Nature (IUCN) [3], and 72% of salamander species (95 species) in Mexico are threatened with extinction [4] due to various factors, such as habitat transformation, which is caused by human activities, including agriculture, livestock, logging, and urban development; introduction of exotic species; and the effects of global climate change and emerging infectious diseases [1,4,5]. In addition to this problem, amphibians have limited dispersal capacities, specialized habitats and habits and high rates of endemism [1], which contribute to the fact that at the national scale, they present the highest beta diversity within terrestrial vertebrate species with a high conservation risk [6,7,8]. Given this background, we consider it important to analyze the diversity of amphibians in the state of Oaxaca with novel methodological approaches. In this work, we analyze amphibian beta diversity and its replacement components and the differences in the richness of the amphibians in the physiographic subprovinces of the state of Oaxaca since this type of analysis offers quantitative information that can contribute to a better understanding of the diversity of this group at different scales and that can be reflected in conservation strategies for these species at the state level [9].

2. Materials and Methods

The units of study in this work were the 12 physiographic subprovinces present in the state of Oaxaca that were used for an analysis of the herpetofauna in Oaxaca by Mata-Silva et al. [2]: Balsas Depression (DB); Mountains and Valleys of the West (MVO); Tehuacán Trench (FT); Sierra Madre de Oaxaca (SMO); Gulf Coastal Plain (PCG); Central Valleys of Oaxaca (VCO); Mountains and Valleys of the Center (MVC); Isthmic Depression of Tehuantepec (DIT); Sierra Madre de Chiapas (SMC); Sierra Madre del Sur (SMS); Pacific Coastal Plain (PCP) and Tehuantepec Coastal Plain (PCT) (Figure 1). From the information on amphibians in this work, a database of presence and absence was constructed.

To determine if there is a relationship between the size of the area of the provinces and the richness of amphibians, anurans and caudates, a regression analysis was applied. Beta diversity or dissimilarity was separated into its two components of replacement and differences in richness [10,11,12]. According to this method, the total dissimilarity is 1 minus the similarity coefficient of Jaccard (βcc) or Sorensen (βsor). This dissimilarity is divided into two components: dissimilarity due to turnover (β.sim) and dissimilarity due to differences in richness (βsne) [12,13,14]. This analysis was carried out in the R program using the BAT package [15,16].

This partition was performed both for the dissimilarity in the composition of species and for that in the taxonomic structure considering the composition of taxa higher than the species level, using the method of Bacaro et al. [17], where the total taxonomic dissimilarity, here βsorT is equal to the dissimilarity of the Jaccard or Sorensen coefficient but taking into account higher taxa. The values of βsorT will range from zero when the taxonomic structure of both communities is identical to 1 when the taxonomic structure is totally different [17]. Taxonomic dissimilarity measures the proportion represented by nonshared taxa to the total number of taxa in the two communities. Therefore, the partitioning of βsorT with the procedure of Carvalho et al. [12] provides a dissimilarity component due to the change in the taxon (β.simT) and a dissimilarity component due to the difference in taxon richness (βsneT). To calculate the total taxonomic dissimilarity and its components, we used the categories of order, family, subfamily, genus and species. The analysis was carried out in the R program [16].

To observe the patterns of beta diversity among species and the beta taxonomy between physiographic subprovinces, we performed nonmetric dimensional scaling analysis (NMDS). To visualize clusters between subprovinces, cluster analyses were performed using the pairwise grouping method with the unweighted pair group method with arithmetic mean (UPGMA). The clusters formed were superimposed in the NMDS and were described using cutoff levels at 0.60 of dissimilarity for the total beta diversity and 0.30 of dissimilarity for the replacement components and differences in richness. These analyses were performed in the PRIMER V7® program [18]. Simple linear regression analyses were also performed to explore the existence of a relationship between the beta diversity of species and the beta diversity of higher amphibian taxa for the total beta diversity, turnover and differences in richness of the 12 physiographic subprovinces present in Oaxaca. All the analyses were performed for all amphibians; however, the analyses were performed separately for anurans and caudates separately as these orders present important differences in richness and environmental requirements and their dissimilarity patterns can vary.

3. Results

3.1. Species Richness and Endemism

Nationally, the state of Oaxaca has the highest amphibian species richness, with 149 species; in addition, 54.3% of these amphibians have been recorded in only one subprovince, as well as 48 of the 60 species endemic to Oaxaca (80%, Table 1). In comparison to caudates, anurans were distributed more broadly, in as many as 10 subprovinces (X = 2.87 subprovinces), and the caudates were distributed in as many as three subprovinces (X = 1.3 subprovinces, Table 2). Notably, the SMO had 88 species (62 anurans and 26 caudates), with 70.45% endemism (62 species: 36 anurans and 26 caudates), while nationally, there was 21.6% endemism (19 species) at the national level (17 anurans and 2 caudates) and 49% endemism (43) at the state level (19 anurans and 24 caudates), of which 86.7% (36) are endemic to this subprovince (16 anurans and 20 caudates). A positive and significant relationship was recorded between the richness of all amphibians, anurans and caudates separately and the area of the subprovinces (r ² = 0.461 p <0.017; r ² = 0.407 p <0.025; and r ² = 0.51 p <0.0137, respectively).

3.2. Total Beta Diversity of Amphibian Species at the State Level

The total beta diversity for Oaxacan amphibian species was 86% (βSOR = 0.86) and was higher than that of higher taxa (βSORT = 0.79). Turnover was more important than the differences in richness at the species level (βSIM = 0.74 and βSNE = 0.12) and higher taxa level (βSIMT = 0.63 and βSNET = 0.16). In comparison to the anurans, the caudates presented a higher total beta diversity of species and taxa (caudates βSOR = 0.92 and βSORT = 0.81 and anurans βSOR = 0.84 and βSORT = 0.76), and turnover (caudates βSIM = 0.81, βSIMT = 0.63 and anurans βSIM = 0.72 βSIMT = 0.61) was more important than the differences in richness (anurans; βSNE = 0.12 βSNET = 0.15 and caudates βSNE = 0.11 and βSNET = 0.18 Figure 2). However, when all the subprovinces are analyzed, including where the caudates are not distributed, in comparison to the other factors, the differences in richness were more important (βSIM = 0.11 and βSNE = 0.69). There was no relationship between total beta diversity or turnover and the area of the physiographic provinces (r ² = 0.031, p <0.1526; r ² = 0.035, p <0.1318, respectively). The area of the physiographic subprovinces showed a positive and significant relationship (r ² = 0.1033, p <0.0083).

3.3. Beta Diversity of Amphibians between Pairs of Physiographic Subprovinces

The total beta diversity of amphibians presented the highest values between DB/PCG and DB/DIT (βsor = 1 for both pairs, the DB/SMC βsor = 0.92), and FT/PCT and PCG/VCO (βsor = 0.90 for both pairs, Figure 3). The highest turnover occurred between DB/PCG and DB/DIT (βsim = 1.0), and another 10 pairs showed values greater than 80% (βsim = 0.80) with MVO being involved in five pairs. For the differences in richness, the highest dissimilarity value was obtained between FT and SMO (βsne = 0.72) (Figure 3, Amphibians SP, supplementary material).

For the anurans, 19 pairs showed values above 80%, with the highest being between DB/PCG and DB/DIT (βsor = 1 for both pairs) and DB/SMC and FT/PCT (βsor = 0.90 for both pairs). With respect to turnover, the highest values were between DB/PCG and DB/DIT (βsim = 1.0), followed by those between FT/PCP and FT/PCT (βsim = 0.85 for both). The differences in richness had their highest value between FT and SMO (βsne = 0.63 Figure 3, Anurans SP, supplementary material). The caudates did not have records in five subprovinces, so 21 pairs of 66 possible pairs were analyzed because when comparing subprovinces with and without caudates, the maximum possible beta diversity originated from the effect of the differences in richness; however, when they were eliminated from the analysis, fourteen pairs showed the maximum beta diversity (βsor = 1) caused by turnover (βsim = 1 and βsne = 0). The highest values of the differences in richness were found between SMO/VCO and MVO/VCO (βsne = 0.92 and βsne = 0.8; Figure 3, Caudates SP, supplementary material).

3.4. Determination of Beta Diversity of Species and Higher Taxa

For the total beta diversity of amphibian species, three groups were formed with a similarity of 0.40 (SMO, DB, VCO, MVC, MVO and FT; SMS, PCP, and PCT; and SMC, PCG and DIT) and for the total beta taxonomic diversity, two groups were formed (Figure 4). In the case of anurans at the species level, two groups were formed with a similarity of 40% (MVO, MVC, FT, and VCO and PCT, PCP, SMS, SMC, DIT, and PCG). For the caudates, a group formed with a similarity of 60% (DIT, PCG, and SMC), and the remaining species were not grouped (Figure 4). For the total beta taxonomy, two groups were formed (Figure 4). Regarding the turnover of species and higher taxa of all amphibians, two groups were distinguished with a similarity of 0.60 (MVO, MVC, FT, and VCO and PCT, PCP, SMS, SMC, DIT, and PCG). For the replacement of anurans, three groups formed with a similarity of 0.40 (PCP, PCT, and SMS; DIT, PCG, and SMC; and MVO, MVC, FT, VCO, DB, and SMO). For the caudates, two groups formed (SMO, VCO, and MVO and DIT, SMC, and PCG) as occurred for the replacement of higher taxa (Figure 4). For the differences in richness of amphibians, anurans and caudates at the level of species and higher taxa, a group with a similarity of 0.40 was formed.

3.5. Relationship between Beta Diversity of Species and Higher Taxa

The total taxonomic beta diversity followed the same trend as that of the beta diversity of species since both were positively and significantly related for all amphibians (r ² = 0.85 p < 0.0001) as were turnover (r ² = 0.88, p < 0.0001) and differences in richness (r ² = 0.90, p < 0.0001). The same result was obtained for anurans (total beta: r ² = 0.86, p < 0.0001; turnover: r ² = 0.86, p < 0.0001; differences in richness: r ² = 0.87, p < 0.0001); however, in the caudates, this relationship was less robust (total beta: r ² = 0.70, p < 0.0002; turnover: r ² = 0.65, p < 0.0001; differences in richness: r ² = 0.40, p < 0.0065), probably due to the high species richness and the low richness of higher taxa in this order (Figure 5).

4. Discussion

Among the states in Mexico, Oaxaca has the highest amphibian richness and endemism, with 149 species, of which 60 are endemic to the state (27 anurans and 33 caudates) [2], surpassing the endemism of amphibians recorded in 13 countries on the same continent [19]. According to our results, this amphibian richness in Oaxaca is related to the size of the physiographic subprovinces, given the species-area relationships initially proposed for island biogeography [20], which is recognized as one of the patterns of global diversity [21,22]. This same species-area relationship has been recorded for the herpetofauna of the Caribbean Islands [20] and the Marías Islands in the Mexican Pacific [23] and for the reptiles of the islands of Mexico [24], as well as for amphibians and reptiles in tropical forest fragments [25,26,27]. In addition to the size of the area, factors such as temperature, precipitation, humidity and environmental heterogeneity are factors that influence amphibian richness [8,9,26,27,28].

The total beta diversity was high for all amphibians and anurans and caudates separately, with dissimilarity values greater than 85% in all cases; however, dissimilarity values were higher for the caudates than for the anurans. This result is consistent with the higher beta diversity values of 80% recorded for amphibians at the continental [9] and national scales [6,7,8]. The same trend was found for taxonomic beta diversity, although with lower values, indicating that the amphibian fauna presents in the 12 subprovinces in Oaxaca drastically change their composition both at the level of species and at the level of higher taxa. The total beta diversity of amphibian species between the pairs of subprovinces also showed generally high values (> 80%), indicating a high rate of change between the composition of the amphibian communities between these units. Similar results were obtained by Mata-Silva et al. [2], who used the biogeographic similarity coefficient for the entire state for the herpetofauna in these same subprovinces, and they found that 48% of the possible combinations presented dissimilarities greater than 80%. Similarly, the highest dissimilarities for amphibians in this study were found between DB and other subprovinces, such as FT, PCT, PCG, SMO and VCO [2]. In addition, Calderón-Patrón et al. [29,30] obtained similar dissimilarity values for the beta diversity of species and higher taxa of amphibians between pairs of ecoregions in Hidalgo and between pairs of quadrants at three scales on the Isthmus de Tehuantepec, and Juárez-Ramírez et al. [31] obtained similar values at the local level among the amphibians of three natural protected areas (NPAs) in the state of Veracruz, recording total beta diversity in the three parks ranging from 60 to100%.

In comparison to the differences in richness, turnover contributed more to the total beta diversity of species and taxonomic beta diversity, both for amphibians and anurans and caudates, with the diversity being higher in the latter order, indicating that species substitution between subprovinces is the mechanism that causes variation in the species composition of the analyzed regions. These results are consistent with those recorded for amphibians in some regions of America [9], those recorded at the state level among ecoregions within the state of Hidalgo [30], and those recorded at the local level among the amphibians from three NPAs in the state of Veracruz, with a total beta diversity ranging from 60 and 100% [31,32]. In two pairs, turnover was more important than the differences in richness, and the opposite occurred in one pair [31], and in the three private NPAs in Veracruz, amphibians presented a relatively low dissimilarity values in species composition (values of βcc between 0.27 and 0.48), mainly due to species replacement, with values between 18 and 28% [32]. However, on the Isthmus of Tehuantepec at fine and medium scales, both components contributed similar amounts [29]. In the case of the differences in richness, their contribution to the total beta diversity did not exceed 13% in any group analyzed, which reflects that there is a very high substitution of species among the subprovinces. These results seem to be influenced by historical and spatial factors [9], such as the Pleistocene glaciations, which caused amphibians to use shelters to survive at high latitudes where the temperature decreased; subsequently, as the planet warmed, the surviving species dispersed and repopulated this region, while in the areas near the equator, this phenomenon was less common, which caused the beta diversity of amphibians to be dominated by turnover below the 37th parallel and by differences in richness about this parallel [9].

Regarding the turnover between subprovinces, in comparison to the other subprovinces, MVO, DB and FT presented the highest values, which means that in addition to changes in the composition of amphibian communities, the communities also present important environmental differences and are geographically distant. These results may be due to differences in temperature, precipitation and environmental heterogeneity and the geographical distance between them and/or the size of the area of the subprovinces, as well as the number of records and collections [8,9,28,29,30,31,32,33,34,35] (Figure 1, supplementary materials). For the differences in richness, the highest values were obtained between DB and SMO and the other subprovinces, which was the product of the low species richness (seven) recorded in the DB and the high number of species in the SMO (88).

The highest values of total beta diversity and turnover were higher in caudates than in anurans, which reflects the more restrictive environmental requirements and smaller distribution areas of anurans [28]; for example, Dalmolin et al. [36] evaluated the beta diversity of amphibians in different pools in Brazil and determined that total beta diversity and turnover were related to the depth and vegetation of the pools. In Mexico, the beta diversity of amphibians is associated with the heterogeneity of precipitation [8], which is related to their reproduction water requirements since adults need environmental humidity and colder temperatures to survive [28]. Thus, a greater spatial variation in precipitation within a region indicates a greater variation in local conditions and, therefore, the possibility that numerous species of amphibians inhabit a restricted subset of the region [8]. In addition, some caudates are particularly difficult to study because their detectability is biased by researchers' knowledge of their microhabitats [29]. In addition, salamanders have a lower dispersal capacity [28], as observed in Oaxaca with the amplitude of their distribution.

Amphibians are the most threatened class of terrestrial vertebrates in terms of conservation [1], and according to Alroy [37], at the global level, the extinction of amphibians has accelerated since the 1980s and 1990s, with between 28 and 201 species going extinct and a current extinction probability for frogs of 66.1% [38]. Unfortunately, globally, Mesoamerica has the second highest probability of frogs going extinct at 17.74%, with the main causes being deforestation, the introduction and expansion of the distribution of exotic species and the presence of the fungus Batrachochytrium dendrobatidis [1,37,39] estimated the environmental vulnerability score for the amphibians of Mexico and determined that salamanders are the most threatened order with 87% of them having a high environmental vulnerability; at the same time, 42.6% of anurans have a high environmental vulnerability. Within the Plethodontidae family, some genera in Mexico face critical conservation problems, as indicated by the high proportion of threatened species: Thorius (96%), Chiropterotriton (83%) and Pseudoeurycea (78%) [4].

Regarding the evaluation of amphibian conservation according to the Secretariat of Environment and Natural Resources SEMARNAT [40], of the 149 amphibian species recorded in Oaxaca, only 58 are considered in NOM-059, and 91 species have not been evaluated. In addition, 56.5% of the amphibians endemic to the area are not under any conservation category [2]. Despite the complex scenario for amphibians, it is noted that in Oaxaca, the SMO subprovince implements successful conservation actions in most of its territory based on the efforts of local communities [34]. This information is important because it is the most important subprovince for the conservation of amphibians at the state level as it has the greatest richness (88 species) and endemism (43 species) levels in Mexico; this subprovince contains 95% of the caudates endemic to Mexico and 80.5% of those endemics to Oaxaca. The same observation was found by Mata-Silva et al. [2] for all herpetofauna since the SMO contains the greatest total number of species at 216, of which 88 are amphibians. High richness levels of vascular plants, gymnosperms, oaks, birds and mammals have also been found in this subprovince [6,35,41,42]. Although the SMO is the best preserved region of Oaxaca, 10 species of anurans are experiencing declines in their populations, and 11 anurans in MVC and SMS are also experiencing declining populations [43]; in addition, salamanders of the genus Thorius have been evaluated as critically endangered, endangered and threatened according to the IUCN [3,44] and the high mountain Pseudoeurycea species of Cerro San Felipe have not been observed recently despite search efforts [45].