Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Evolutionary History of the Live-Bearing Endemic Allotoca diazi Species Complex (Actinopterygii, Goodeinae): Evidence of Founder Effect Events in the Mexican Pre-Hispanic Period

  • Diushi Keri Corona-Santiago ,

    diushi.keri.c@gmail.com

    Affiliations Programa Institucional de Maestría en Ciencias Biológicas, Facultad de Biología, Universidad Michoacana de San Nicolás de Hidalgo, Morelia, Michoacán, México, Laboratorio de Biología Acuática, Facultad de Biología, Universidad Michoacana de San Nicolás de Hidalgo, Morelia, Michoacán, México

  • Ignacio Doadrio,

    Affiliation Departamento de Biodiversidad y Biología Evolutiva, Museo Nacional de Ciencias Naturales, CSIC, Madrid, España

  • Omar Domínguez-Domínguez

    Affiliation Laboratorio de Biología Acuática, Facultad de Biología, Universidad Michoacana de San Nicolás de Hidalgo, Morelia, Michoacán, México

Abstract

The evolutionary history of Mexican ichthyofauna has been strongly linked to natural events, and the impact of pre-Hispanic cultures is little known. The live-bearing fish species Allotoca diazi, Allotoca meeki and Allotoca catarinae occur in areas of biological, cultural and economic importance in central Mexico: Pátzcuaro basin, Zirahuén basin, and the Cupatitzio River, respectively. The species are closely related genetically and morphologically, and hypotheses have attempted to explain their systematics and biogeography. Mitochondrial DNA and microsatellite markers were used to investigate the evolutionary history of the complex. The species complex shows minimal genetic differentiation. The separation of A. diazi and A. meeki was dated to 400–7000 years ago, explained by geological and climate events. A bottleneck and reduction of genetic diversity in Allotoca diazi was detected, attributed to recent climate fluctuations and anthropogenic activity. The isolation of A. catarinae occurred ~1900 years ago. No geological events are documented in the area during this period, but the date is contemporary with P’urhépecha culture settlements. This founder effect represents the first evidence of fish species translocation by a pre-Hispanic culture of Mexico. The response of the complex to climate fluctuation, geological changes and human activity in the past and the future according to the ecological niches predictions indicates areas of vulnerability and important information for conservation. The new genetic information showed that the Allotoca diazi complex consist of two genetic groups with an incomplete lineage sorting pattern: Pátzcuaro and Zirahuén lakes, and an introduced population in the Cupatitzio River.

Introduction

The diversity and distribution of freshwater ichthyofauna has been strongly linked to natural events, mainly geological process and historical climatic fluctuations. It is therefore that, taking into account the biogeographic history of each region in the world; the freshwater fish have been used as a model for the study of biogeography, speciation, paleohydrology and evolution. However, translocation and introduction of fish species is well documented in many ancient cultures having an important role in their current distribution, misunderstanding the evolutionary history of the species, In America, particularly in Mexico, the impact of pre-Hispanic cultures in the current distribution of the biota, is little known and mostly documented for birds [13]. In fact, the scientific community has long assumed that translocations and introduction of species in America occurred only after European colonization [4].

The region in West-central Mexico, particularly between the Lerma-Santiago and Balsas rivers, is characterized by varied and rugged physiography, a product of intense tectonic and volcanic activity in the area since the Miocene [5]. Also, this region is an area with a long cultural history, with the P’urhépecha settlement developing one of the most stable empires at least 3000 years ago [6], particularly due to the richness and abundance of lacustrine resources. The main resources currently exploited by settlements around the lakes are fish, as well as turtles, salamanders, clams, and wetland plants [7,8]. Even the area was formerly known by the Náhuatl culture as Mechuacan, “a place of fishermen,” and is currently called Michoacán.

Approximately 100 fish species are in the West-central Mexico, 70% of which are endemic and one of the most representative endemic fish groups of this region is the subfamily Goodeinae [5]. The group presents high species richness in a small area (c. 41 species) and unique characteristics associated with breeding strategies and embryo development, such as internal fertilization, matrotrophy, and sexual selection [913]. The family presents a high diversification apparently influenced by viviparity, with vicariance and adaptive radiation being the most important factors. The speciation rate of goodeids has not been constant and has been impacted by multiple extinctions, estimating that at least 25 000 years is required to establish an evolutionary lineage [14].

The genus Allotoca is the most diverse within the Goodeinae, and is the subject of taxonomic controversy, with the need for taxonomic revision to validate species [9,10,1517]. The A. diazi complex comprises three recognized species [5]: Allotoca diazi [18], endemic to Pátzcuaro basin, A. meeki [16], endemic in the Zirahuén basin, and Allotoca catarinae [15], restricted to the Cupatitzio River, an upper tributary of the Balsas River basin (Fig 1). Allotoca diazi was extirpated from Pátzcuaro Lake and is currently restricted to the Chapultepec Spring, a tributary of Pátzcuaro Lake, while A. meeki was extirpated from Zirahuén Lake and currently restricted to Opopeo Spring, a tributary of the Zirahuén Lake [5]. The lakes are located in the Michoacán-Guanajuato Corridor, an area of over 1000 volcanic cones active from the Pliocene to the present. To explain the evolution of these three sister species and other co-distributed taxa, two biogeographic hypotheses have been proposed. The first argues the existence of a tributary connecting the Lerma River with Cuitzeo, Pátzcuaro, and Zirahuén lakes [19] around 700 000 years ago [20]. The second hypothesis indicates the existence of a tributary connecting the Cupatitzio River with Zirahuén and Pátzcuaro lakes, reaching Zacapu [21].

thumbnail
Fig 1. Geographical distribution of Allotoca diazi complex.

https://doi.org/10.1371/journal.pone.0124138.g001

Additionally, the application of modern biogeographic, genetic, systematic, and bioinformatics methodologies in the study of evolutionary processes within the Goodeinae has questioned recent connections of the Lerma and Cupatitzio Rivers with the Zirahuén, Pátzcuaro, Cuitzeo, and Zacapu basins [9,12,22,23]. The taxonomic status of the three species has also been questioned [13,24], and the possibility of pre-Hispanic translocation of exploited species has been proposed [12].

In order to elucidate and analyze the two biogeographic hypotheses about the connections or disconnections between the Zirahuén, Pátzcuaro and Cupatitzio basins, and the possibility of a species translocation, we examined the evolutionary history of Allotoca diazi complex using two different molecular markers described for the Goodeinae Subfamily [9,25,26]: a conserved molecular marker of mitochondrial DNA (Cytochrome b, Cytb gene) and a less conserved microsatellite nuclear markers. We obtained these molecular markers to explore the genetic differentiation at intraspecific and interspecific level applying phylogeographic and phylogenetic approaches. The results of the genetic study were analyzed based on geological, paleoclimatic and anthropogenic records from central Mexico to access how natural and artificial historical processes are involved in the evolutionary history of A. diazi complex. In addition, we used all the genetic information and ecological niche predictions for identify operational conservation units.

Materials and Methods

Ethics statements

This research was a part of the multidisciplinary work for restoration and conservation of the Cupatitzio River, and the complete work, field and laboratory protocols were approved by the Mexican Government, SEMARNAT (The Ministry of Environment and Natural Resources for Mexico) and CONACyT (National Council of Sciences And Technology) (Permit No. FOMIX-MICH-2009-C02-115897). All individuals were anesthetized using Tricaine mesylate (MS-222) to minimize suffering.

Fish sampling and DNA isolation

One hundred and twenty specimens were collected from six populations throughout the distribution range of the species complex (S1 Table). Our sampling did not recover any A. diazi or A. meeki from Pátzcuaro Lake and Zirahuén lacustrine basin, respectively. This is consistent with reports that these species were extirpated from these basins as a result of the introduction of exotic species [5]. Capture was by electrofishing between March and August 2011. Pectoral or caudal fin clips, were obtained and preserved in absolute ethanol and frozen at -20°C. Some specimens of the species complex were identified and deposited in the fish collection of the Universidad Michoacana de San Nicolás de Hidalgo CPUM in México; the remaining individuals were returned to the water. Approximately 2 mm2 of tissue was used for DNA extraction. Digestion was performed with ATL QIAGEN Buffer and Proteinase K and purifying by BioSprint DNA Blood Kit QIAGEN according to the manufacturer’s instructions.

Locus amplification, sequencing, and genotyping

Polymerase Chain Reaction (PCR) was performed to amplify the molecular markers. To amplify the Cytb gene, we used the primers GLuDG [27] and H16460 [28]. The PCR consisted of a 12 μl volume reaction with a final concentration of 0.2 μM of each primer, 0.25 mM of each dNTP, 1.5 M of MgCl2, and 1 unit of Taq DNA Polymerase. The PCR procedure consisted of 2 min at 94°C followed by 35 cycles of 45 s at the 94°C for DNA denaturation, 1 min at 46°C for primer alignment, 1.5 min at 72°C for synthesis, and a final extension of 5 min at 72°C. The PCR products were checked by electrophoresis in agarose gel 1.5% and submitted to Macrogen Inc. (Korea) for sequencing.

For microsatellites, we searched polymorphic loci of the A. diazi complex using 20 primers designed for Goodeinae [25,26]. Seven polymorphic microsatellite loci were selected (S2 Table). The 7 loci were amplified with QIAGEN Multiplex PCR Kit according with the manufacturer specifications, and the products were submitted for genotyping (Secugen S. L. Corp, Spain).

mtDNA analysis

Manual alignment was implemented in Mega v5.2 [29] and examined chromatographically. The genetic diversity was estimated by calculating nucleotide diversity (π), haplotype diversity (h) and proportion of segregating sites (Θ) using the software DNAsp v5.0 [30]. Neutrality tests (Tajima’s D and Fu’s FS) to Cytb were implemented using Arlequin v3.5.1.3 [31]. The genetic P-distances separating the three species were obtained (Dp) with Mega v5.2, and a bootstrapping process was implemented with 1000 repetitions. We conducted an analysis of geographic correspondence by network haplotype reconstruction using Network v4.6.1.0 [32], applying the median-joining method.

Tree-building algorithms were used for phylogenetic reconstruction using Allotoca zacapuesis, a sister group of the species complex [9], as outgroup. Neighbor-joining clustering was implemented in Mega v5.2 based on Dp with bootstrapping with 1000 replicates. Maximum likelihood reconstruction was conducted with RAxMLGUI v1.3 [33,34], using the generalized time reversible (GTR) + gamma + proportion of invariable sites model [33] and 1000 bootstrap repetitions. The evolutionary substitution model selected was the transitional model [35] + invariable site proportion + gamma (TMI+I+G) using jModeltest v1.7 [36] and selected by the Akaike information criterion (AIC) (S3 Table). Phylogenetic reconstruction was by Bayesian inference (BI) with the software MrBayes v3.1.2 [37], using the above selected evolutionary substitution model and implementing 2 reactions for 4 Markov Chains Monte Carlo (MCMC) processes with 3 million generations, sampling every 100 generations. We evaluated the convergence of the log-likelihood (-InL) value of the 2 reactions, with 10% of reconstructions discarded as burn-in to construct the consensus tree (σ = 0.0002). The posterior probabilities obtained based on a confidence limit of 95% (highest posterior density-HPD) were used to evaluate the support values of nodes.

Microsatellite analysis

Quality of the samples genotyped and polymorphism in microsatellite loci was analyzed using GeneMapper v 4. Allele number by loci and species (Na), effective number of alleles (Nae), null allele proportion, and linkage disequilibrium were obtained using Arlequin v3.5.1.3 and Genepop v4.1.2 [38]. Observed (Ho) and expected (He) heterozygosity was determined, and the endogamy coefficient was obtained based on Wright’s fixation index (FIS) [39] to calculate the Hardy-Weinberg (HW) equilibrium across each locus and species using the software Arlequin v3.5.1.3, implementing 100 000 permutations.

The genetic structure level among the three species of the A. diazi complex was estimated with paired fixation indices (FST). Significance was assessed by AMOVA with 20 000 permutations implemented in the software Arlequin v 3.5.1.3 using multiple hypothetical arrangements based on phylogenetic and biogeographic information. A Bonferroni correction [40] was applied to each P-value obtained in the paired test, HW analysis and linkage disequilibrium estimation. Bayesian clustering was conducted with the software Structure v2.3.3 [41]. We performed a series of independent runs from values of K = 1–10 populations assuming correlated allele frequencies and an admixture model [42] with 10 runs for each K value to assess convergence of results without geographic information. For each K value, the MCMC was run with a burn-in of 100 000 steps and chain length of 500 000. Following the method described by Evanno et al. 2005 [43], ΔK was calculated to find the K value that best explained the genetic differentiation level.

To determine whether the difference in allele size contributed to genetic structure and if each locus is consistent with a strict Stepwise Mutation model, we estimated the RST values to assess deviation of equilibrium between RST and FST values. RST values greater than FST imply influence of the size differences of alleles on interpopulation genetic differentiation; thus we can assume that loci mutate faster by multi-step changes [4447]. The RST estimation was performed with SPAGeDi v1.3a [48] without geographic information, applying 20 000 permutations to obtain a test of differences among observed and expected RST and FST values.

Isolation-with-migration and demographic history analyses

This approach allows determination of gene flow before and after isolation within a population and assesses phylogenetic relationship scenarios in recent isolation events within closely related populations or species. Divergence times of populations of members of the A. diazi complex were calculated using three analyses in IMa2 software [49], which has the advantage of implementing models of gene flow.

Pairwise analysis was conducted to find the confidence limits of the time to divergence among the three species, combining the molecular markers to substantiate and support the exchange rates within Markov chains. A molecular clock was calibrated using the mutation rate for Cytb estimated for teleosts of 0.76–2.2%/million years (4.33*10-6-1.25*10-5subst/locus/generation) [5052] and the Hasegawa-Kishino-Yano (HKY) evolutionary substitution model [53]. The mutation rate used for calibration with microsatellites was 0.9*10-4-1.5*10–4 locus/gamete/generation, estimation obtained for the Zebra fish Danio rerio [54] and congruent with the mutation rate for vertebrates [55,56]. The SMM evolutionary model was implemented for microsatellites. Ten independent runs were conducted for each scenario using 150 Markov chains with a burn-in of 500 000 for a total of 5 000 000, sampling topologies every 10 000 steps to recover 100 000 topologies from all runs. It was assumed that effective population size and the migration rate was different in each period.

The demographic history of the members of the species complex was obtained constructing Bayesian skyline plots (BSP) using Beast 1.8.0 [57] for Cytb marker. Chains were run for 70 000 000 steps sampling each 1000 and 10% was discarded as burn-in under HKY substitution model, a relaxed molecular clock calibrated with the information of Cytb mentioned above. All operators were optimized automatically and the results were analyses and visualized using Tracer 1.5 [58]. Effective sample size (ESS) was used to evaluate the strength of the analyses, and we determined the marginal distributions and HPD values in both analyses.

A test to identify recent bottlenecks was conducted using Bottleneck v1.2.0.2 [59] with three evolutionary models for microsatellites: Infinite Allele Model (IAM) [60], Two-Phase Model [61], and Stepwise Mutation Model (SMM) [62]. The test was implemented with 10 000 iterations based on a confidence limit of 95%. A mode-shift test was conducted to detect significant deviations in the allelic frequencies in populations that have experienced a recent bottleneck [63].

Identification of operational conservation units and predictions of future distribution

The genetic, biological, and ecological information obtained in this investigation, along with relevant published information [64,65], was used for identification of conservation units. In order to predict a future scenario for range of distribution of the natural members of the species complex, a discriminant analysis of presence-only data through maximum entropy using the program Maxent v3.3.3 [6667] was conducted for the present and for two periods in the future, years 2041 to 2060 and 2061 to 2080. The 6.0 Representative Concentration Pathways (RCP6.0) [68] was selected from downscaling Intergovernmental Panel on Climate Change 5 scenario (IPPC5) (Coupled Model Intercomparison Project 5) [69] implementing the Community Climate System Model 4 (CCSM4) [70]. We implemented independent runs of 5000 iterations with 19 climatic variables (http://www.worldclim.org/bioclim), bootstrapping of 100 replicates to obtain 95% confidence intervals, regularization multiplier (RM) values from 0.01 to 10, evaluating area under the ROC curve (AUC) values and using the omission rates to discard identical AUC of different RM values. Collinear climatic variables were discarded (S4 Table).

Results

mtDNA analysis

One hundred nine sequences of the mitochondrial gene Cytb (1037pb) were obtained for the sampled populations of the three species (GenBank accession numbers: KJ776467-KJ776575) (S1 Table). Overall, A. diazi and A. meeki nucleotide diversity was 0.0005±0.0002, haplotype diversity was 0.65±0.13, and mean proportion of segregating sites per locus was 1.4±0.3. Allotoca catarinae showed null genetic diversity in four sampled populations with a unique haplotype. Allotoca diazi was the most genetically diverse with five haplotypes (n = 35), followed by A. meeki with three haplotypes (n = 22), of their respective and unique populations (Table 1). No significant deviation from neutrality was found for the A. diazi complex. The mean of the genetic distances of A. catarinae with respect to A. diazi and A. meeki was 0.6% and between A. diazi and A. meeki was 0.3%. The three phylogenetic approximations were congruent. The topologies were similar with all DNA sequences and haplotypes. The tree based on haplotypes is presented in Fig 2. The consensus tree was formed by two haplogroups: a non-monophyletic group of A. diazi/A. meeki and a group with the complete of A. catarinae. The haplotype network recovered the same two genetic groups, with four mutational steps separating them (Fig 3). Spatial congruence was not found in haplogroup A. diazi/A. meeki. Two haplotypes of A. diazi were shared with A. meeki, and the unique haplotype of A. catarinae is not shared with the other members of the species complex.

thumbnail
Fig 2. Phylogenetic inference based on haplotypes using Neighbor-joining,

Maximum likelihood and Bayesian inference of Cytb gene. Red asterisk: Allotoca diazi (n = 35 Hn = 5), Blue asterisk: A. meeki (n = 22 Hn = 3) and Green asterisk: A. catarinae (n = 52). Support values are represented by numbers above the nodes for analyses in the order mentioned above.

https://doi.org/10.1371/journal.pone.0124138.g002

thumbnail
Fig 3. Haplotype network based on Cytb gene: Red asterisk: A. diazi (n = 35 Hn = 5), Blue asterisk: A. meeki (n = 22 Hn = 3) and Green asterisk: A. catarinae (n = 52).

Values represent the number of individuals per haplotype. The yellow circles represent the median-vectors.

https://doi.org/10.1371/journal.pone.0124138.g003

thumbnail
Table 1. Genetic diversity parameters and neutrality test based on Cytb gene sequences.

https://doi.org/10.1371/journal.pone.0124138.t001

Microsatellite analysis

For the 7 microsatellites loci amplified in populations of the A. diazi complex, we obtained a total of 96 genotypes. The average number of alleles per locus was 10.4±3.5; the number of genetic variants contributing to heterozygosity in each locus of the species complex was 5.8±2.2; and the proportion of null alleles was 0.06±0.03 (P > 0.05). Significant deviation from linkage equilibrium was not found at most loci, LD = 0.06±0.03 (P > 0.05). Significant linkage disequilibrium was observed in two loci of A. diazi (S5 Table). Significant deviation from HW equilibrium associated with lower than expected heterozygosity was found in A. diazi and A. catarinae (P < 0.05) (Table 2).

thumbnail
Table 2. Genetic diversity, HW equilibrium, and bottleneck test for microsatellite loci.

https://doi.org/10.1371/journal.pone.0124138.t002

The observed RST values for all loci were not significantly different from expected (Ho:RST = pRST) or from the FST values (S6 Table), suggesting that allele size difference did not contribute to the genetic structure of the A. diazi complex. Low but significant genetic differentiation was detected (FST = 0.113 P < 0.05) via AMOVA analysis (S7 Table). Genetic structure by BI recovered two genetic groups (A. diazi+A.meeki—A. meeki) with a high proportion of individuals assigned to other groups (0.3–0.7) (Fig 4). However, significant genetic differentiation was detected between the three members of the species complex. The lowest FST value of genetic differentiation was found between Allotoca diazi and A. meeki (0.084), while genetic differentiation of A. catarinae from A. diazi and A. meeki was 0.133 and 0.109 (P<0.05), respectively (Table 3). The number of migrants value between members of the species complex were not significant.

thumbnail
Fig 4. Genetic structure by Bayesian inference based on microsatellite loci.

https://doi.org/10.1371/journal.pone.0124138.g004

thumbnail
Table 3. Estimated genetic differentiation based on microsatellites.

https://doi.org/10.1371/journal.pone.0124138.t003

Genetic structure by BI recovered two genetic groups (K = 2 A. catarinae + A. diazi/A.meeki) with high proportion of individuals assigned to another group (0.3±0.7) (Figs 3 and S1).

Estimate of divergence times and demographic history

The A. diazi complex divergence from its common ancestor A. zacapuensis was estimated to have occurred in the past million years. The isolation event of A. diazi from A. meeki occurred in the past 7000–400 years (HPD≥95%) with the peak of probability estimated as 1700 years ago (Table 4 and S2 Fig). Divergence time and isolation of A. catarinae from A. diazi and A. meeki was calculated at 8000–700 years ago, with the highest probability being 1900 years ago for A. diazi and 2500 years ago for A. meeki.

thumbnail
Table 4. Divergence times based on mitochondrial Cytb gene and microsatellite loci in an Isolation-with-Migration model.

https://doi.org/10.1371/journal.pone.0124138.t004

The migration rate estimated among the members of the A. diazi complex was fewer than 5.5 individuals/generation and was not significant (Table 3). BSP analysis using Cytb revealed demographic decline in the last 10 000 years for A. catarinae (HPD≥95%) (S4 Fig), and recent bottlenecks for the three members of the species complex were detected out of the HPD confidence range. However, recent bottlenecks were identified in the three species based on different evolutionary models using microsatellites. The bottlenecks detected in A. diazi and A. catarinae via SMM were supported by the FIS coefficients: 0.141 and 0.151, respectively. The bottleneck detected in A. meeki was associated with higher than expected heterozygosity and allele size differences in the IAM analysis.

Operational conservation units and predictions of future distribution

We included in the analysis the species considered as natural members of the species complex: Allotoca diazi and A. meeki. The species Allotoca catarinae was excluded of the analysis because could be result of translocation. The predictions were performed taking into account that the gene flow between species will be null.

The ecological niches area predicted with present environment conditions for Allotoca diazi and A. meeki is bigger than their real current distribution. In the case of the ecological niche models for two periods in the future, high AUC values for were obtained for both species (AUC = 0.940–0.995 RM = 0.01–0.1, see S4 Table). In the first period between the years 2041 to 2060, the model suggest that the ecological niches available for the two species will be located at higher altitudes than they are currently distributed (Figs 5 and 6). In the second period of the future (2061–2080), the area of ecological niches predicted will increase for A. diazi but the probability of occurrence in the Pátzcuaro Lake is low. In the same period, the ecological niches predicted to A. meeki could decrease drastically. According with the genetic and future predictions, we propose the existence of two OCUs corresponding to natural populations in Pátzcuaro and Zirahuén basins.

thumbnail
Fig 5. Ecological niche projections for Allotoca diazi in present and two future periods.

https://doi.org/10.1371/journal.pone.0124138.g005

thumbnail
Fig 6. Ecological niche projections for Allotoca meeki in present and two future periods.

https://doi.org/10.1371/journal.pone.0124138.g006

Discussion

Allotoca diazi complex

The isolation of the Allotoca diazi complex from its common ancestor A. zacapuensis, does not contradict the hypothesis of Álvarez (1972) [21] with respect to the connection between the Zacapu and Pátzcuaro basins during the Pleistocene, resulted from the formation of the El Zirate mountain and the northern Pátzcuaro Lake shoreline during the Late Pleistocene [71].

The Allotoca diazi complex consists of a genetic group with low genetic structure and its isolation within the members was more recent, during the Holocene. The isolation of A. diazi, A. meeki, and A. catarinae is reflected in the low genetic differences, non-monophyletic patterns, shared haplotypes, and genetic groups assignment, and we considered that the species complex consist in an incomplete lineage sorting pattern [7275]. The genetic distances among the three species in Cytb are smaller than reported to sister Goodeinae species (1.7–11%) [9] and sister species, congeneric species, and cofamilial genera within and across the major vertebrate taxonomic classes (~2%) [7679]. The genetic variation estimated in the species complex based on the Cytb gene is attributed to genetic drift and mutation (DT and FS values did not show significant deviation from neutrality). Low genetic variation was detected as a genetic pool arrangement through AMOVA analysis. However, BI detected this minimum mark of genetic structure as two genetic groups. The results obtained are supported by the minimal morphological differences previously reported [15,17,8084].

Allotoca diazi and A. meeki isolation

The hypothesis of connection and disconnection of Pátzcuaro and Zirahuén lakes during the last 700 000 years [19] is supported by stratigraphic [85]; limnologic data [86,87]; shared native icthyofauna such as Chirostoma estor, Chirostoma pátzcuaro (Atherinopsidae), Algansea lacustris (Cypronidae), Alloophorus robustus, Skiffia lermae, Allotoca dugesii, and Goodea atripinnis (Goodeidae) [88]; and genetic data [12,89].

Shared haplotypes of A. diazi and A. meeki and the high proportion of individuals assigned to other species (0.5 to 0.7) by microsatellite data, suggest recent gene flow. Microsatellite information showed low level of genetic differentiation (FST = 0.084). Low genetic diversity was estimated for A. diazi in the Cytb gene (h = 0.78, π = 0.0030) and even lower for A. meeki (h = 0.52, π = 0.0007). The migration rate between the two species is not significant any time.

We estimated that the isolation of A. diazi from A. meeki occurred in the past 400–7000 years, which is consistent with geological, climatic, and anthropogenic events involving Pátzcuaro and Zirahuén lakes. Tecto-volcanic events giving rise to the separation of Pátzcuaro and Zirahuén lakes began 8000 years ago with the formation and activity of the La Tasa volcano southwest of Lake Pátzcuaro and climate fluctuations causing decline in water level and drying of streams [21,9093]. Fluctuations in the level of Lake Pátzcuaro during the Holocene associated with climate change, human activity [90,91,9499], and occurrence of tsunamis [100] could be involved to the isolation and demographics changes of A. diazi and A. meeki. These demographic changes were reflected in genetic diversity including the positive DT value, which was found to be low in A. diazi, and associated with recent bottlenecks that promoted inbreeding (FIS = 0.141), reduced heterozygosity, and the loss of allelic richness. In contrast, a bottleneck identified in A. meeki is attributed to the rapid loss of allelic diversity; however, the RST value and the BSP analysis (in the HPD range) indicates that the bottleneck detected by IAM were overestimated. The population dynamics of A. meeki can be explained by the variations in nucleotide diversity and proportion of segregating sites (π < Θs), which is interpreted as a recent population reduction followed by expansion, as demonstrated by the estimated value of DT = -0.86 and non-significant positive values of Fu’s FS estimates [101103].

Based in our results, we conclude that natural processes, such as geologic events and climate fluctuations, are the main factors responsible for the isolation of A. diazi from A. meeki during the Holocene, not enough time to diverge and showing an incomplete lineages sorting.

Allotoca catarinae

Divergence times and isolation of A. catarinae from A. diazi and A. meeki was calculated to be within the past 8000–700 years and 8000–800 years, respectively. The peak of probability for isolation of A. caratinae was 1900 years ago. The low number of native species fish shared [12,88], and the geological and climatic records not support the connection and disconnection of the region in ancient periods: as the formation of Cupatitzio River ~30 mya [104,105], or in recent periods: as the activity and collapse of Tancítaro volcano 792 000±22 000 years ago [106,107], and the origin of Zirahuén Lake 17 000 years ago [92,93,108], contrasting the Álvarez hypothesis (1972) [21].

We propose that the existence of this taxon in the Cupatitzio River could be the result of a pre-Hispanic introduction, causing a founder effect, which resulted in rapid fixation of a single Cytb haplotype in populations, along with a loss of genetic diversity. The lack of shared haplotypes between Cupatitzio populations and Pátzcuaro-Zirahuén populations cloud this interpretation, but the specimens sampled in the Pátzcuaro basin were from a small spring 10 km east of the lake and represent the only known local population of the species, which has been extirpated form the main basin [5]. We speculate that the ancestral haplotype introduced into Cupatitzio came from the main basin, and ancestral haplotype is lost. The haplotypes collected from the small spring (50 m2 approximately) were isolated from the main basin, and the genetic drift fixed rare haplotypes. In addition, null genetic diversity was found in A. catarinae, associated with recent bottlenecks. Allotoca catarinae is the only one in which bottlenecks were detected in several models of evolution, and also supported with the BSP analysis. These bottlenecks drive in all cases to low heterozygosity (FIS = 0.151).

Relevant anthropogenic events in the region of Lakes Pátzcuaro and Zirahuén and the Cupatitzio River are congruent with our estimated divergence times. There is evidence that fishing was a major activity among the P'urhépecha. Fish, specifically goodeines and atherinopsids, had importance as food, even today goodeids are important part of the food supply for natives in the region, as well as culturally, in all socioeconomic strata, especially for Cazonci, the supreme ruler of the P’urhépecha state. Offerings of fish found in tombs at the archaeological site Uricho southwest of Lake Pátzcuaro, dating to the Late and Epi-Classic periods (~500–900 A.D.) [7,8].

The “islander”, “Chichimeca” and the “Coringuaro”, considered as first human groups to settle in the region of Uruapan and south of Michoacán, probably had influence on the management of aquatic resources and translocation of species may have occurred [109,110]. However, the information about their biological impact is less well documented. The rise of the P'urhépecha culture occurred around 1350, during the formation of the empire, with establishment of new villages, and species translocations possibly took place [109112]. The expansion of the P'urhépecha Empire occurred during the second period with establishment around 1400 (second period) and 1450–1530 (third period) (S3 Fig). During these periods, settlements ranged from the region of the La Palma River, to the Marqués River region, in the Balsas Depression. Evidence from the P'urhépechas settlements in the Cupatitzio River subbasin and the Tepalcatepec-Balsas River suggests that, in this period, A.diazi translocation may have occurred. We conclude that the biogeographic history of A. catarinae was possibly determined by a founder effect mediated by pre-Hispanics in the past 3000 years, and taking into account this inference we cannot rule out species translocations in Pátzcuaro and Zirahuén lakes.

Conservation implications

The factors involved in the genetic diversity and demographic history of A. diazi and A. meeki, were not only the founder effect (A. meeki) and climate fluctuations (A. diazi), but also anthropogenic activity [9399].

The two OCU’s proposals (Pátzcuaro and Zirahuén lakes) are areas with a high degree of alteration, including eutrophication, introduction of exotic species, drought, and overfishing [11,12,22,24,113115], which led to the extirpation of A. diazi from the main Pátzcuaro basin and A. meeki from Zirahuén Lake. Neither Pátzcuaro nor Zirahuén are designated protected areas. In addition, the results of the future ecological niche modeling for A. diazi and A. meeki show a reduction of potential habitats questioning the survival of the species, taking into account its genetic vulnerability to climate fluctuations and human activity, which indicates the necessity of the implementation of effective conservation and management strategies of this micro-endemic species with a long history of economic and cultural importance to the lacustrine settlements.

Supporting Information

S1 Fig. Number of groups (K) to explain the genetic structure by ΔK criteria.

https://doi.org/10.1371/journal.pone.0124138.s001

(EPS)

S2 Fig. Marginal distribution with Isolation-with-Migration model of three scenarios of ancestry based on mitochondrial Cytb gene and microsatellite loci (ESS>100 000).

https://doi.org/10.1371/journal.pone.0124138.s002

(EPS)

S3 Fig. Extension of the P’urhépecha Empire during the period ca. 1450–1520 A.D.

Modified from [1].

https://doi.org/10.1371/journal.pone.0124138.s003

(TIF)

S4 Fig. Demographic history of the members Allotoca diazi complex using Bayesian skyline plots from Cytb sequences.

Dotted lines represent the location of the upper bound (Max), the mean () and lower bound of the HPD≥95%.

https://doi.org/10.1371/journal.pone.0124138.s004

(EPS)

S2 Table. Selected polymorphic microsatellite loci.

https://doi.org/10.1371/journal.pone.0124138.s006

(DOC)

S3 Table. Evolutionary substitution model from mitochondrial Cytb gene by Akaike Information Criterion (AIC).

-lnL = log likelihood, Ti = Transitions, Tv = Transversions.

https://doi.org/10.1371/journal.pone.0124138.s007

(DOC)

S4 Table. Variables included in ecological niche models.

RM = Regularization multiplier value selected according to the of AUC value.

https://doi.org/10.1371/journal.pone.0124138.s008

(DOC)

S5 Table. Estimated allelic diversity by locus and species.

n = sample size, Na = mean of alleles, Nae = mean number of effective alleles, Pan = null allele proportion, LD = linkage disequilibrium test with significant value after Bonferroni correction (* = P > 0.05).

https://doi.org/10.1371/journal.pone.0124138.s009

(DOC)

S6 Table. Genetic structure based on RST values compared with FST of permutation analysis.

Ho RST = pRST, ns = non-significant, *significant (P < 0.05) after Bonferroni correction, RST = observed values, pRST = expected values, IC = 95% confidence interval.

https://doi.org/10.1371/journal.pone.0124138.s010

(DOC)

S7 Table. Genetic structure inferred via Analysis of Molecular Variance AMOVA based on microsatellite data.

https://doi.org/10.1371/journal.pone.0124138.s011

(DOC)

Acknowledgments

Thank to Rodolfo Pérez, Fernando Mar, Eloisa Torres and Silvia Perea for their help in sampling trips and valuable comments to this work. Also, thank to the two anonymous reviewers and editor for their observations and contributions to this work.

Author Contributions

Conceived and designed the experiments: DKCS ID ODD. Performed the experiments: DKCS. Analyzed the data: DKCS. Contributed reagents/materials/analysis tools: DKCS ID ODD. Wrote the paper: DKCS. Acquisition sample: DKCS ID ODD. Final approval of the version: ID ODD.

References

  1. 1. Christensen AF. The fifteenth- and twentieth-century colonization of the basin of Mexico by the great-tailed grackle (Quiscalus mexicanus). Global Ecology and Biogeogrpahy. 2000;9, 415–420.
  2. 2. Peer BD. Invasion of the emperor's grackle. Adeola. 2011;58, 405–409.
  3. 3. Haemig PD. Aztec introduction of the great-tailed grackle in ancient Mesoamerica: Formal defense of the Sahaguntine historical account. NeoBiota. 2014;22, 59–75.
  4. 4. Telles JL. The unknown extent of ancient bird introductions. Ardeola. 2011;58, 399–404.
  5. 5. Miller RR, Minckley WL, Norris SM. Freshwater fishes of México. University Of Chicago Press. 2005.
  6. 6. Bravo-Espinosa M, Mendoza Cantú ME, Medina Orozco LE. Escenarios de erosión bajo diferentes manejos agrícolas en la cuenca del lago de Zirahuén, Michoacán, México. Investigaciones geográficas. 2009; 73–84.
  7. 7. Williams E, Weigand PC. Las cuencas del occidente de México. Colegio de Michoacán. 1996.
  8. 8. Guzmán AF, Polaco ÓJ, Pollard HP. Ofrendas de peces asociadas a entierros del Clásico-Epiclásico en Urichu, Michoacán, México. Archaeofauna. 2001; 149–162.
  9. 9. Doadrio I, Domínguez O. Phylogenetic relationships within the fish family Goodeidae based on Cytochrome b sequence data. Molecular Phylogenetics and Evolution. 2004;31, 416–430. pmid:15062784
  10. 10. Webb SA, Graves JA, Macias-Garcia C, Magurran AE, Foighil DÓ, Ritchie MG. Molecular phylogeny of the livebearing Goodeidae (Cyprinodontiformes). Molecular Phylogenetics and Evolution. 2004;30, 527–544. pmid:15012937
  11. 11. Domínguez-Domínguez O, Doadrio I, Perez-Ponce de Leon G. Historical biogeography of some river basins in central Mexico evidenced by their goodeine freshwater fishes: A preliminary hypothesis using secondary Brooks parsimony analysis. Journal of Biogeography. 2006;33, 1437–1447.
  12. 12. Domínguez-Domínguez O. Filogeografía de Zoogoneticus quitzeoensis, Xenotoca variata y Alloophorus robustus (Cyprinodontiformes: Goodeidae) en el centro de México: implicaciones taxonómicas y de conservación. Doctoral Dissertation. Instituto de Biología, UNAM. Distrito Federal, México. 2008a.
  13. 13. Corona-Santiago DK, Domínguez-Domínguez O. Goodeids phylogeny. In Kempkes M, Köck M, Stawikouski R (eds): Contributions to the biology and conservation of splitfins (Goodeidae). Westarp Wissenschaften, Germany. 2013; 681.
  14. 14. Ritchie MG, Hamil RM, Graves JA, Magurran AE, Webb SA, Macías García C. Sexual differentiation: population genetic divergence and sexual dimorphism in mexican goodeid fish. Journal of Evolutionary Biology. 2007;20, 2048–2055. pmid:17714321
  15. 15. De Buen F. Una nueva especie del Neoophorus diazi (Meek). Anales del Instituto de Biología. 1942;13, 342–249. pmid:19707616
  16. 16. Álvarez J. Contribución al conocimiento del género Neoophorus (Pisc., Goodeidae). Ciencia: Revista Hispano-Americana de Ciencias Puras y Aplicadas. 1959;19, 13–22.
  17. 17. Meyer MK, Radda AC, Domínguez-Domínguez O. Notes on the genera Neoophorus Hubbs & Turner, 1937, with a description of a new species of Allotoca from Laguna de Zacapu, Michoacán, Mexico (Teleostei, Cyprinodontiformes: Goodeidae). Annalen des Naturhistorischen Museums in Wien. 2001;103 B, 453–460.
  18. 18. Meek SE. A contribution to the ichthyology of Mexico. Zoological Series, Field Columbian Museum. 1902;3, 63–128.
  19. 19. De Buen F. Los lagos michoacanos. II. Pátzcuaro. Revista de la Sociedad Mexicana de Historia Natural. 1943;4, 211–232.
  20. 20. Israde-Alcántara I, Garduño-Monroy VH. Lacustrine record in a volcanic intra-arc setting: The evolution of the Late Neogene Cuitzeo basin system (central-western Mexico, Michoacan). Palaeogeography, Palaeoclimatology, Palaeoecology. 1999;151, 209–227.
  21. 21. Álvarez J. Ictiología Michoacana V. Origen y distribución de la ictiofauna dulceacuícola michoacana. Anales de la Escuela Nacional de CIencias Biológicas. 1972;19, 155–161.
  22. 22. Domínguez-Domínguez O, Alda F, Pérez-Ponce de León G, García-Garitagoitia JL, Doadrio I. Evolutionary history of the endangered fish Zoogoneticus quitzeoensis (Bean, 1898) (Cyprinodontiformes: Goodeidae) using a sequential approach to phylogeography based on mitochondrial and nuclear DNA data. BMC Evolutionary Biology. 2008b;8, 161. pmid:18503717
  23. 23. Pérez-Rodríguez R, de León GP-P, Domínguez-Domínguez O, Doadrio I. A new species of Algansea (Actinopterygii: Cyprinidae) from the Ameca River basin, in Central Mexico. Revista mexicana de biodiversidad. 2009;80, 483–490.
  24. 24. Domínguez-Domínguez O, Pedraza-Lara C, Gurrola-Sánchez N, Perea S, Pérez-Rodríguez R, Israde-Alcántara I, et al. Historical Biogeography of the Goodeinae (Cyprinodontiforms). In Uribe-Aranzabal MC, Grier H (eds): Viviparous fishes II. New Life Publications, Florida. 2010; 13–30.
  25. 25. Boto L, Doadrio I. Polymorphic microsatellites in two different species of the genus Zoogoneticus Meek, 1902 (Goodeidae, Actynopterygii). Molecular Ecology Notes. 2003;3, 70–72.
  26. 26. Hamill RM, Webb SA, Macías Garcia C, Graces JA, Magurran AE, Ritchie MG. Comparison of genetic diversity at microsatellite loci in near-extinct and non-endangered species of Mexican goodeine fishes and prediction of cross-amplification within the family. Journal of Fish Biology. 2007;70, 16–32.
  27. 27. Palumbi S, Martin A, Romano S, McMillan WO, Stice L, Grabowski G. Simple fool's guide to PCR. University of Hawaii, Honolulu. 1991.
  28. 28. Perdices A, Bermingham E, Montilla A, Doadrio I. Evolutionary history of the genus Rhamdia (Teleostei: Pimelodidae) in Central America. Molecular Phylogenetics and Evolution. 2002;25, 172–189. pmid:12383759
  29. 29. Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Molecular Biology and Evolution. 2011;28, 2731–2739. pmid:21546353
  30. 30. Librado P, Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics, 2009; 25, 1451–1452. pmid:19346325
  31. 31. Excoffier L, Lischer HEL. Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Molecular Ecology Resources. 2009;10, 564–567.
  32. 32. Bandelt H-J, Forster P, Röhl A. Median-Joining networks for inferring intraspecific phylogenies. Molecular Biology and Evolution. 2000;16, 37–48.
  33. 33. Stamatakis A. RAxML-VI-HPC: maximum likelihood-based phylogenetic analyses with thousands of taxa and mixed models. Bioinformatics. 2006;22, 2688–2690. pmid:16928733
  34. 34. Silvestro D, Michalak I. raxmlGUI: a graphical front-end for RAxML. Organisms Diversity & Evolution. 2012;12, 335–337.
  35. 35. Posada D. Selecting models of evolution. The phylogenetic handbook. A practical approach to DNA and protein phylogeny. Cambridge University Press, Cambridge. 2003; 256–282.
  36. 36. Posada D. jModelTest: Phylogenetic Model Averaging. Molecular Biology and Evolution. 2008;25, 1253–1256. pmid:18397919
  37. 37. Huelsenbeck JP, Ronquist F. MRBAYES: Bayesian inference of phylogenetic trees. Bioinformatics. 2001;17, 754–755. pmid:11524383
  38. 38. Rousset F. genepop'007: a complete re-implementation of the genepop software for Windows and Linux. Molecular Ecology Resources. 2008;8, 103–106. pmid:21585727
  39. 39. Wright S. The genetical structure of populations. Annals of Eugenics. 1951;15, 323–354. pmid:24540312
  40. 40. Rice WR. Analyzing tables of statistical tests. Evolution. 1989; 223–225.
  41. 41. Pritchard JK, Sthepens M, Donnely P. Inference of populations structure using multilocus genotype data. Genetics. 2000;155, 945–959. pmid:10835412
  42. 42. Hubisz MJ, Falush D, Stephens M, Pritchard JK. Inferring weak population structure with the assistance of sample group information. Molecular Ecology Resources. 2009;9, 1322–1332. pmid:21564903
  43. 43. Evanno G, Regnaut S, Goudet J. Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Molecular Ecology, 2005;8, 2611–2520. pmid:15969739
  44. 44. Hardy OJ, Charbonnel N, Fréville H, Heuertz M. Microsatellite allele sizes: A simple test to assess their significance on genetic differentiation. Genetics. 2003;163, 1467–1482. pmid:12702690
  45. 45. Cutter AD. Molecular evolution inferences from the C. elegans genome. In The C. elegans Research Community (eds): WormBook-The Online Review of C. elegans Biology. Pasadena (CA). 2005.
  46. 46. Seyfert AL, Cristescu MEA, Frisse L, Schaack S, Thomas WK, Lynch M. The rate and spectrum of microsatellites mutation in Caenorhabditis elegans and Daphnia pulex. Genetics. 2008;178, 2113–2121. pmid:18430937
  47. 47. Haasl RJ, Payseur BA. The number of alleles at a microsatellites defines the allele frequency spectrum and facilitates fast accurate estimation of θ. Molecular Biology and Evolution. 2010;12, 2702–2715. pmid:20605970
  48. 48. Hardy OJ, Vekemans X. SPAGeDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. Molecular Ecology Notes. 2002;2, 618–620.
  49. 49. Hey J. Isolation with migration models for more than two populations. Molecular Biology and Evolution. 2010;27, 905–920. pmid:19955477
  50. 50. Zardoya R, Doadrio I. Molecular evidence on the evolutionary and biogeographical patterns of european Cyprinids. Journal of Molecular Evolution. 1999;49, 227–237. pmid:10441674
  51. 51. Machordom A, Doadrio I. Evidence of a cenozoic Betic-Kabilian connection based on freshwater fish phylogeography (Luciobarbus, Cyprinidae). Molecular Phylogenetics and Evolution. 2001;18, 252–263. pmid:11161760
  52. 52. Near TJ, Benard MF. Rapid allopatric speciation in logperch darters (Percidae: Percina). Evolution. 2004;58, 2798–2808. pmid:15696757
  53. 53. Hasegawa M, Kishino H, Yano T. Dating of the human-ape splitting by a molecular clock of mitochondrial DNA. Journal of Molecular Evolution. 1985;22, 160–174. pmid:3934395
  54. 54. Shimoda N, Knapik EW, Ziniti J, Sim C, Yamada E, Kaplan S, et al. Zebrafish genetic map with 2000 microsatellite markers. Genome. 1999;58, 219–232.
  55. 55. Weber JL, Wong C. Mutation of human short tandem repeats. Human molecular genetics. 1993;2, 1123–1128. pmid:8401493
  56. 56. Yue GH, David L, Orban L. Mutation rate and pattern of microsatellites in common carp (Cyprinus carpio l.). Genetica. 2007;129, 329–331. pmid:16897452
  57. 57. Drummond AJ, Suchard MA, Xie D, Rambaut A. Bayesian phylogenetics with BEAUti and the BEAST 1.7. Molecuar Biology and Evolution. 2012;29, 1969–1973. pmid:22367748
  58. 58. Rambaut A, Drummond A. Tracer v1.5. Available: http://tree.bio.ed.ac.uk/software/tracer. 2007.
  59. 59. Cornuet JM, Luikart G. Description and power analysis of two tests for detecting recent population bottlenecks from allele frequency data. Genetics. 1996;144, 2001–2014. pmid:8978083
  60. 60. Kimura M, Crow JF. The number of alleles that can be maintained in a finite population. Genetics. 1964;49, 725–738. pmid:14156929
  61. 61. Di Rienzo A, Peterson AC, Garza JC, Valdes AM, Slatkin M, Freimer NB. Mutational processes of simple-sequence repeat loci in human populations. PNAS. 1994;91, 3166–3170. pmid:8159720
  62. 62. Kimura M, Ohta T. Stepwise mutation model and distribution of allelic frequencies in a finite population. PNAS. 1978;75, 2868–2872. pmid:275857
  63. 63. Luikart G, Allendorf FW, Cornuet JM, Sherwin WB. Distortion of allele frequency distributions provides a test for recent population bottlenecks. The Journal of heredity. 1998;89, 238–247. pmid:9656466
  64. 64. Doadrio I, Perdices A, Machordom A. Allozymic variation of the endangered killifish Aphanius iberus and its application to conservation. Environmental Biology of Fishes. 1996;45, 259–271.
  65. 65. Domingues VS, Santos RS, Brito A, Almada VC. Historical population dynamics and demography of the eastern Atlantic pomacentrid Chromis limbata (Valenciennes, 1833). Molecular Phylogenetics and Evolution. 2006;40, 139–147. pmid:16603388
  66. 66. Phillips SJ, Anderson RP, Schapire RE. Maximum entropy modeling of species geographic distributions. Ecological modeling. 2006;190, 231–259.
  67. 67. Phillips SJ, Dudík M. Modeling of species distributions with Maxent: new extensions and a comprehensive evaluation. Ecography. 2008;31, 161–175.
  68. 68. Fujino J, Nair R, Kainuma M, Masui T, Matsuoka Y. Multi-gas mitigation analysis on stabilization scenarios using AIM global model. Multigas Mitigation and Climate Policy. The Energy Journal, 2006;3 (Special Issue).
  69. 69. Taylor KE, Stouffer RJ, Meehl GA. An overview of CMIP5 and the experiment design. B. Am. Meteorol. Soc. 2012;93, 485–498.
  70. 70. Moss RH, Edmonds JA, Hibbard KA, Manning MR, Rose SK, van Vuuren DP, et al. The next generation of scenarios for climate change research and assessment. Nature. 2010;463, 747–756. pmid:20148028
  71. 71. Robles-Camacho J, Corona-Chávez P, Morales-Gámez M, Guzmán AF, Polaco ÓJ, Domínguez-Vázquez G, et al. Estratigrafía y paleoambiente asociados a un Gomphoteriidae (Cuvieronius hyodon) en Tzintzuntzan, Michoacán, México. Revista Mexicana de Ciencias Geológicas. 2010;27, 530–544.
  72. 72. Maddison WP. Gene Trees in Species Trees. Systematic Biology. 1997;46, 523–536.
  73. 73. Arbogast BS, Edwards SV, Wakeley J, Beerli P, Slowinski JB. Estimating divergence times from molecular data on phylogenetic and population genetic timescales. Annual Review of Ecology and Systematics. 2002;33, 707–740.
  74. 74. Yu Y, Than C, Degnan JH, Nakhleh L. Coalescent histories on phylogenetic networks and detection of hybridization despite incomplete lineage sorting. Systematic Biology. 2011;60, 138–149. pmid:21248369
  75. 75. Rosenfeld JA, Payne A, DeSalle R. Random roots and lineage sorting. Molecular Phylogenetics and Evolution. 2012;64, 12–20. pmid:22445448
  76. 76. Brown WM, George M, Wilson AC. Rapid evolution of animal mitochondrial DNA. PNAS. 1979;76, 1967–1971. pmid:109836
  77. 77. Avise JC, Arnold J, Ball RM, Bermingham E, Lamb T, Neigel JE, et al. Intraspecific phylogeography: The mitochondrial DNA bridge between population genetics and systematics. Annual Review of Ecology and Systematics. 1987;18, 489–522.
  78. 78. Avise JC. The history and purview of phylogeography: A personal reflection. Molecular Ecology. 1998;7, 371–379.
  79. 79. Bradley RD, Baker RJ. A test of the genetic species concept: cytochrome-b sequences and mammals. Journal of Mammalogy, 2001;82, 960–973.
  80. 80. Regan CT. Pisces. Biología Centrali-Americana. 1906–1908;8, 1–203.
  81. 81. Hubbs CL. Studies of the fishes of the Order Cyprinodontes. 1926;VI. 3–87.
  82. 82. Hubbs CL. Studies of the fishes of the order Cyprinodontes, XII: A new genus related to Empetrichthys. 1932; 1–5.
  83. 83. Hubbs CL, Turner CL. Studies of the fishes of the Order Cyprinodontes, XVI: A revision of the Goodeidae.1939.
  84. 84. Hieronimus H. "Leven baren" wat is dat wigenlijk?. Poecilia Nieuws. 1995;2, 22–27.
  85. 85. Israde-Alcántara I, Wade M, Garduño-Monroy VH, Barron J. Estratigrafía y encuadramiento geodinámico de las cuencas lacustres del centro de México. Unión Mexicana de Estudios del Cuaternario. Fondo de Cultura Económica, México, DF. 2008.
  86. 86. Bradbury JP. Limnologic history of the Lago de Pátzcuaro, Michoacán, Mexico from the past 48,000 years: impacts of the climate and man. Paleogeography, Palaeoclimatology, Palaeoecology. 2000;163, 69–95.
  87. 87. Israde-Alcántara I, Garduño-Monroy VH, Ortega MR. Paleoambiente lacustre del cuaternario tardio en el centro del lago de Cuitzeo. Hidrobiología. 2002;12, 61–78.
  88. 88. Villaseñor-Gómez LE. La biodiversidad en Michoacán: Estudio de estado. 2005.
  89. 89. Bloom DD, Piller KR, Lyons J, Mercado-Silva N, Medina-Nava M. Systematics and biogeography of the silverside tribe Menidiini (Teleostomi: Atherinopsidae) based on the mitochondrial ND2 gene. Copeia. 2009;2, 408–417.
  90. 90. Israde-Alcántara I, Garduño-Monroy VH, Fisher CT, Pollard HP, Rodríguez-Pascua MA. Lake level change, climate, and the impact of natural events: The role of seismic and volcanic events in the formation of the Lake Patzcuaro Basin, Michoacan, Mexico. Quaternary International. 2005;135, 35–46.
  91. 91. Garduño-Monroy VH, Chávez-Hernández J, Aguirre-González J, Vázquez-Rosas R, Mijares-Arellano H, Israde-Alcántara I, et al. Zonificación de los periodos naturales de oscilación superficial en la ciudad de Pátzcuaro, Mich., México, con base en microtremores y estudios de paleosismología. Revista Mexicana de Ciencias Geológicas. 2009;26, 623–637.
  92. 92. Ortega B, Vázquez G, Caballero M, Israde I. Late Pleistocene: Holocene record of environmental changes in Lake Zirahuen, Central Mexico—Springer. Journal of Paleolimnology. 2010;44, 745–760. pmid:20964053
  93. 93. Vázquez G, Ortega B, Davies SJ, Aston BJ. Registro sedimentario de los últimos ca. 17000 años del lago de Zirahuén, Michoacán, México. Boletín de la Sociedad Geológica Mexicana. 2010;62, 325–343.
  94. 94. Street-Perrott FA, Perrott RA, Harkness DD. Anthropogenic soil erosion around Lake Patzcuaro, Michoacan, Mexico, during the Preclassic and late Postclassic-Hispanic periods. American Antiquity. 1989;54, 759–765.
  95. 95. O'Hara SL. Historical evidence of fluctuations in the level of Lake Pátzcuaro, Michoacán, Mexico over the last 600 years. Geographical Journal. 1993;159, 51–62.
  96. 96. O'Hara SL, Metcalfe SE, Street-Perrott FA. On the arid margin: The relationship between climate, humans and the environment. A review of evidence from the highlands of central Mexico. Chemosphere, 1994;29, 965–981. pmid:7953466
  97. 97. O'Hara SL, Metcalfe SE. Reconstructing the climate of Mexico from historical records. The Holocene. 1995;5, 485–490.
  98. 98. Fisher CT, Pollard HP, Israde-Alcántara I, Garduño-Monroy VH, Banerjee SK. A reexamination of human-induced environmental change within the Lake Pátzcuaro Basin, Michoacán, México. PNAS. 2003;100, 4957–4962. pmid:12671066
  99. 99. Fisher CT. Demographic and Landscape Change in the Lake Patzcuaro Basin, Mexico: Abandoning the Garden. American Anthropologist. 2005;107, 87–95.
  100. 100. Garduño-Monroy VH, Soria-Caballero DC, Israde-Alcántara I, Hernández VMM, Rodríguez-Ramírez A, Ostroumov M, et al. Evidence of tsunami events in the Paleolimnological record of Lake Pátzcuaro, Michoacán, Mexico. Geofísica Internacional. 2011;50–2, 147–161. pmid:22915666
  101. 101. Hudson RR, Kreitman M, Aguadé M. A test of neutral molecular evolution based on nucleotide data. Genetics. 1987;116, 153–159. pmid:3110004
  102. 102. Tajima F. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics. 1989;123, 585–595. pmid:2513255
  103. 103. Fu YX. Statistical Tests of Neutrality of Mutations Against Population Growth, Hitchhiking and Background Selection. Genetics. 1997;147, 915–925. pmid:9335623
  104. 104. López Ramos E. Paleogeografia y tectonica del Mesozoico de Mexico. Revista Mexicana de Ciencias Geológicas. 1981;5, 158–177.
  105. 105. Gómez-Tuena A, Orozco-Esquivel MT, Ferrari L. Petrogénesis ígnea de la Faja Volcánica Transmexicana. Boletín de la Sociedad Geológica Mexicana. 2005;57, 227–283.
  106. 106. Ownby S, Delgado Granados H, Lange RA, Hall CM. Volcán Tancítaro, Michoacán, Mexico, 40Ar/39Ar constraints on its history of sector collapse. Journal of Volcanology and Geothermal Research. 2007;161, 1–14.
  107. 107. Morelli S, Monroy V, Gigli G, Falorni G. The Tancítaro debris avalanche: characterization, propagation and modeling. Journal of Volcanology. 2010;193, 93–105.
  108. 108. Tamayo LT, West RC. The hidrology of middle America. In Wanchope R, West RC (eds): Handbook of Middle America Indians. University of Texas Press, US. 1964.
  109. 109. De Alcalá FJ. Relación de Michoacán (F de Cultura Económica, Ed,). Fondo 2000, 1540;88.
  110. 110. Roskamp H. La historiografía indígena de Michoacán. Research School Cnws. 1998.
  111. 111. Goggin JM. An archaeological survey of the Río Tepalcatepec basin, Michoacán, México. American Antiquity. 1943;9, 44–58.
  112. 112. Nichols DL, Pool CA. The Oxford Handbook of Mesoamerican Archaeology—Google Books. Oxford Univesity Press. 2012.
  113. 113. Lyons J, González-Hernández G, Soto-Galera E, Guzmán-Arroyo M. Decline of freswater fishes and fisheries in selected drainages of west-central Mexico. Fisherie. 1998;23, 10–18.
  114. 114. Zambrano L, Cordova-Tapia F, Ramírez-Herrejón JP, Mar-Silva V, BustamanteL , Camargo T, et al. Las especies exóticas en el lago de Pátzcuaro, Michoacán, México. In Huerto-Delgadillo R, Vargas-Velázquez S, Ortiz-Paniagua CF (eds): Estudio Ecosistémico del lago de Pátzcuaro. Aportes en gestión ambiental para el fomento del desarrollo sustentable. Instituto Mexicano de Tecnología del Agua, México. 2011; 133–156.
  115. 115. Ramírez-Herrejón JP, Zambrano L, Mercado-Silva N, Torres-Téllez A, Pineda-García F, Caraveo-Patiño J, et al. Long term changes in the fish fauna of lago de Pátzcuaro in central México. Latin American Journal of Aquatic Research. 2014;42, 137–149.