Received Date: March 21, 2016;Accepted Date: May 17, 2016; Published Date: May 25, 2016
Citation: Ruiz-García M, Luengas-Villamil K, Pinedo-Castro M, Leal L, Bernal-Parra LM, et al. (2016) Continuous Miocene, Pliocene and Pleistocene Influences on Mitochondrial Diversification of the Capybara (Hydrochoerus Hydrochoeris; Hydrochoeridae, Rodentia): Incapacity to Determine Exclusive Hypotheses on the Origins of the Amazon and Orinoco Diversity for This Species. J Phylogen Evolution Biol 4:166. doi:10.4172/2329-9002.1000166
Copyright: © 2016 Ruiz-García M, et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Visit for more related articles at Journal of Phylogenetics & Evolutionary Biology
The capybara (Hydrochoerus hydrochoeris) is the largest rodent on the world and it is strongly linked to the river systems of a large fraction of South America and part of Central America (Panama). Thus, it is an interesting species to test hypotheses about the origin of the high biodiversity within the Amazon Basin and in a sizeable fraction of the Neotropics. To test these hypotheses, we sequenced two mitochondrial genes (control region and Cytochrome b) of 78 wild capybaras sampled in Colombia, Peru, Ecuador and Brazil. At least, five different “populations” or ESUs were detected in well delimited geographical areas. However, our results do not support the more recent view that two different species of capybara are present (H. hydrochoeris and H. isthmus), unless chromosomal speciation (stasipatric or parapatric) can be demonstrated between these two groups. A Bayesian tree with the aforementioned two mitochondrial genes, and another Bayesian tree with a subset of 25 capybaras for 10 mitochondrial genes, showed that the initial diversification of the mitochondrial haplotype in capybaras was initiated in the Late Miocene. The trees also showed that the other haplotype diversification processes extended into the Pliocene and Pleistocene. We also detected population expansion events during different moments of the Pleistocene. Although some authors strongly suggest that the Miocene diversification explains the extreme biodiversity in the Amazon Basin and in surroundings areas (for instance, the Paleogeography hypothesis), others consider it the result of available forest refugia (Refuge hypothesis) during the Pleistocene. However, our results suggest that both hypotheses (and others, such as the Riverrefuge, the Recent Lagoon and the Hydrogeological Recent Change hypotheses) could have affected the evolution of the capybara to generate the current mitochondrial diversity. Thus, it is difficult to generalize a unique Amazon biota diversification hypothesis because each species or taxon could be affected by different processes and because the temporal antiquity of each taxon in South America is also different. Many mammalian taxa, and others, migrated into South America during the Great American Biotic Interchange (GABI) and their diversification processes in South America were mainly driven by Pleistocene events as those proposed by the Refuge hypothesis. Older taxa within this continent could have begun their current genetics diversification processes earlier, such as in the case of the capybara.
Hydrochoerus hydrochoeris; mitochondrial genes; Amazon biodiversity; Paleogeography; Riverine-refuge; Hydrogeological Recent Change; Recent lagoon and Refugia hypotheses; Bayesian analyses
The great and amazing biodiversity, at all organismic levels, of the Amazon Basin has been a topic, which has drawn the attention of numerous naturalists since the 19th century. To help explain this noteworthy diversity, a series of hypotheses have been proposed. The most important ones are the Stability-time Hypothesis (STH), Paleogeography Hypothesis (PH), Riverine Barrier Hypothesis (RBH), Recent Lagoon Hypothesis (RLH), Refugia Hypothesis (RH), River Refuge Hypothesis (RRH), Canopy-density Hypothesis (CDH), Disturbance-Vicariance Hypothesis (DVH), Museum Hypothesis (MH), Gradient Hypothesis (GH), and Hydrogeological Recent Change Hypothesis (HRCH) (all hypotheses are explained in Table 1). Excellent reviews on the question can be found in the literature [1-12].
|Name||Main arguments of each hypothesis||Authors|
|Stability-time hypothesis, STH||Amazonian forests have existed for an extremely long time, since the Cenozoic.In such stable environment, the rate of extinction would be less than the speciation rate, and this led to large levels of biodiversity||[231-233]|
|Paleogeographyhypothesis, PH||Three mechanisms for vicariance events: •Uplift of the Northern Andes in the Late Miocene and the formation of the current Amazon River’s drainage
•Island model by recurrent Miocene marine introgressions caused by fluctuations of world sea-level and the appearance of internal Amazonian lakes
•Appearance of geological arches, which changed the geomorphology of the Amazonian basin
|RiverineBarrierhypothesis, RBH||The largest rivers of the Amazonian Basin formed in the late Miocene as a barrier to the dispersion of organisms||[235,236]|
|Recentlagoonhypothesis, RLH||Most of Amazonia was covered by a huge lake or lagoon during the Pliocene, and successively smaller portions of Amazonia were covered during a series of assumed high sea-level stands during the Pleistocene||[144,145]|
|Refugia hypothesis, RH||There are dry/humid cycles caused by the Milankovitch cycles. During the dry periods rainforest communitiessplit into isolated refugia separated by savannah or arid pampas.This hypothesis was originally established for the Pleistocene, but later expanded to the Miocene-Pliocene periods as well||[4,5]|
|RiverRefugehypothesis, RRH||A combination of the RBH and RH||[237-239]|
|Canopy-densityhypothesis, CDH||Changes in vegetation type and structure during the Last Glacial Maximum, around 18,000 YA, occurred near the southern margins of the Amazonian basin. The results also suggest repeated reductions and increases in forest canopy density over large areas during glacial-interglacial cycles. Under this hypothesis, biological separation of taxa does not necessarily require forest fragmentation, as in the RH, but only a change from wet to dry forests. Although this hypothesis was initially originated for the Pleistocene, it could be expanded to prior periods|||
|Disturbance-Vicariancehypothesis, DVH||Cold/warm cycles during the Pleistocene period||[241,242]|
|Museumhypothesis, MH||This hypothesis claimed, in a first step, speciation in very localized stable habitat pockets in mountainous regions around the periphery of Amazonia due to climatic fluctuations without major vegetational changes, followed, in a second step, by range expansion of new species into the Amazonian lowlands where they accumulate and are preserved as in a museum||[140,141]|
|Gradienthypothesis, GH||This model predicts parapatric speciation across steep environmental gradients, boundaries, while all the other hypotheses are based on allopatric speciation|||
|Hydrogeological Recent Change hypothesis, HRCH||This hypothesis considers the role of Pliocene-Pleistocene geological changes in the formation of the Neotropical rivers to be more important than the role of the Miocene North Andean uplift in the origin of the Amazon drainage.This ‘‘Young Amazon’’ model suggests that the current Amazonian drainage system was established between 2-3 Millions of years ago, MYA, and the major Amazonian river tributaries originated within the last 2 MYA. Watershed breakdowns and headwater captures influence the diversity of both aquatic and non-aquatic organisms||[11,111,154,205,210]|
Table 1: Eleven different hypotheses on the origin of the Amazonian biodiversity tested throughout molecular population genetics in the capybara (Hydrochoerus hydrochoeris).
Capybara (Hydrochoerus hydrochoeris), also named poncho (Panama), chigüiro (Colombia), chigüire (Venezuela), ronsoco (Peru), carpincho (Argentina, Uruguay, Paraguay), is the world’s largest living rodent. On average, adult capybaras reach a body length of 1-1.5 m, a height of 0.5-0.7 m, and a weight of 48.9 kg . However, their weight has a clinal distribution, with smaller values (35-50 kg) to the north (Colombian and Venezuelan Llanos) and larger values (up to 91 kg) to the south (Southern Brazil and Argentina) [13,14]. This is a semi-aquatic species tightly linked to Neotropical rivers and flooded savannas, with water and temperature the main factors in its distribution. The maximum altitude registered for this species is around 1,300 meters above sea level.
Capybara is a very interesting species to test some of these Amazon (and Neotropical) biodiversity hypotheses because it inhabits the basins of the main South American rivers: Orinoco, Amazon, Sao Francisco and Parana-La Plata rivers. The Hydrochoeridae family first appeared in the fossil record of South America about 10 million of year ago (MYA). Capybara probably evolved from an unknown species of Cardiatherium about 2 MYA in Southern South America . It is closely related with the genera Neocherus and Hydrochoeropsis (Upper Pliocene to Recent; ). All the known fossil records of capybara are from the Pleistocene in Curacao, Uruguay, Argentina, Colombia and Bolivia which suggests that only Pleistocene biodiversity hypotheses are germane to mitochondrial diversification in capybara.
Capybara’s geographical distribution is divided into two populations. A trans-Andean population extends from middle Panama (from the Panama channel to the Darien forests) to Northern and Pacific Colombia (Atlantic area, lower valleys of the Sinu, Atrato and Cauca rivers, lower and middle valleys of the Magdalena and Cesar rivers and some populations on the Pacific coast, Choco and Valle Departments) and Northeastern Venezuela (around Maracaibo Lake). The cis- Andean population is distributed from Eastern Colombia and the rest of Venezuela, Guyana, Surinam, French Guiana, the Amazonian areas of Ecuador, Peru, Bolivia, the majority of Brazil, Paraguay, Uruguay and Northern Argentina . The first population was traditionally considered as a differentiated sub-species (H. isthmus), whilst the second was considered another subspecies, H. hydrochaeris . However, Goldman  claimed that the Panama population could be a different species (H. isthmus). The individuals of this population have a lower weight (26-28 kg)  and a shorter gestation period (104- 111 days compared to 115-125 days for the other capybara population; ). Wilson and Reeder  recognized both species on the basis that H. hydrochoeris was larger in all external and cranial characters. In contrast, H. isthmus had a wider frontal bone to total skull length proportion, with the lower diastema proportionally longer as well as the pterygoid bones shorter and thicker than found in the other capybara population. However, within H. hydrochoeris, no subspecies has been defined .
Mitochondrial markers have been used in a large variety of studies with mammals to resolve evolutionary significant units, taxonomic conflicts, and to determine phylogeographic patterns and biogeographic history [21-26]. The control region of mt DNA (D-loop) is a portion that evolves particularly rapidly and thus allows fine-scale resolution of population structure and micro-evolutionary divisions [27,28]. The mt Cytochrome b gene (Cyt-b) contains both slowly and rapidly evolving codon positions, as well as conservative and variable regions; it constitutes a suitable marker for phylogenetic purposes at various divergence levels . A subset (25 specimens) of the overall sample analyzed (78 specimens) and sequenced at the mt control region and mt Cyt-b, was also sequenced at other 10 mitochondrial genes.
The main aim of this work are as follows: 1- to provide molecular genetics evidence to test if only Pleistocene biodiversity hypotheses could be applied to this species or, in contrast, if the mitochondrial diversification process preceded the Pleistocene; 2- to provide molecular genetics evidence about of the number of species or subspecies within the genus Hydrochoerus; and 3- to compare if the molecular results obtained with different mitochondrial markers provide identical information.
Materials and Methods
The geographical distribution of the 78 capybara samples analyzed was as follows: Northern and Western Colombia (Córdoba and Valle del Cauca Department, n=24), Casanare Department (n=19), Meta Department (n=7), Guainia Department (n=4), Colombian-Peruvian- Brazilian Amazon frontier (n=2), Negro River at the Central Brazilian Amazon (n=1), Northern Peruvian Amazon, including the Napo River (Loreto region, n=8), Peruvian San Martin Department (n=3), Medium and Southern Ucayali River at the Peruvian Amazon (n=7) and Ecuadorian Amazon (n=3) (Figure 1). The samples were obtained from animals hunted by indigene communities for food (the animals were hunted near their respective communities; tissues were muscle or integument, which were stored in a solution of 20% dimethyl sulfoxide (DMSO), supersaturated with NaCl) or live animals maintained in these communities as pets (hairs with follicles stored in ethanol 96°).
Molecular markers employed
The DNA of the samples of muscle and pieces of skin was extracted using the phenol-chloroform procedure , while DNA samples from hair with follicles were extracted with 10% Chelex resin . For Cyt-b gene amplification via polymerase chain reaction (PCR), we used the forward primer L15579 (5’-CCT AGT TTA TTT GGA ATG GAT CGT AG -3’) and the reverse primer H15049 (5’- GCC TGT ACA TCC ACA TCG GAC GAC G-3’)  under the following touchdown PCR profile: 95°C for 5 min, followed by three cycles of 95°C for 30 s, 50-57°C for 30 s, 72°C for 75 s, followed by 27 cycles of 95°C for 30 s, 38-45°C for 30 s, 72°C for 75 s and a final extension of 72°C for 7 min. For the amplification of the control region (D-loop), the forward and reverse primers were, respectively, CAPIF3 (5’-CAG GAA ACA GCT ATG ACC CAA TTA TTC TAY YTG ACA TAA GAC-3’) and CAPIR3 (5’-TGT AAA ACG ACG GCC AGT GAG CGG GTA TAA YRT TAT GG-3´) being the PCR profile of 95°C for 3 min, followed by 30 cycles of 94°C for 30 s, 52°C for 30 s, 72°C for 1 min, and a final extension of 72°C for 5 min . For both genes, the PCRs were performed in a 25 μl volume with reaction mixtures including 2.5 μl of 10 x buffer, 2 μl of 2 mM MgCl2, 0.625 μl of 1 mM dNTPs, 1.5 μl of 0.6 mM of each primer, one unit of Taq DNA polymerase, 13.8 μl of H2O and 2 μl (20-80 ng/μl) of DNA. PCR reactions were carried out in a BioRad thermocycler. All amplifications, including positive and negative controls, were checked in 2% agarose gels, using the molecular weight marker μ1 74 DNA digested with
A subset of the samples analyzed were sequenced for other additional mitochondrial genes. These 10 genes were 12S rRNA (945 bp), 16S rRNA (1,560 bp), ND1 (960 bp), ND2 (1,060 bp), COI (1,540 bp), COII (683 bp), COIII (780 bp), ATP8 (190 bp), ATP6 (670 bp) and ND4 (1,370 bp), with a total of 9,758 bp (their amplification conditions is shown elsewhere). These genes were selected because we are doing a mitogenomics project with different mammals. Many mitogenomic analyses have shown to be extremely useful in determining phylogenetic relationships and in estimating divergence times in many different groups of organisms, such as birds [34,35] and mammals [36-45]. Mitochondrial gene trees are more precise in reconstructing the divergence history among species than other molecular markers . Cummings et al.,  showed that mitochondrial genomes have higher information content per base than nuclear DNA. Additionally, many of the samples we analyzed (hairs) could provide a good amount of mitochondrial DNA but not of nuclear DNA. The sequences were concatenated by means of the SequenceMatrix v. 1.7.6 Software . Overlapping regions were examined for irregularities, such as frameshift mutations and premature stop codons. The lack of such irregularities agrees quite well with the absence of numt sequences in our mitochondrial sequences.
Population genetics and Phylogenetic analyses
The sequences were edited and aligned with the BioEdit Sequence Alignment Editor  and the DNA Alignment Software (Fluxus Technology Ltd). The Modeltest  and the Mega 6.05 Software  were applied to determine the best nucleotide substitution model. We used the Ln L criteria.
For all these population genetic analyses we employed four geographical regions [1-Western Amazon (Colombian, Peruvian and Brazilian Amazon); 2-Meta, Casanare, in the Colombian Eastern Llanos, and Guainia; 3-Northern Colombia; 4-Napo River in Ecuador and Peru]. These four geographical haplogroups were determined by the phylogenetic analyses. In some analyses, only three geographic regions were employed due to the small sample size of the Napo River.
To determine the genetic diversity of the capybara, we used the number of polymorphic sites (S), the number of haplotypes (H), the haplotypic diversity (Hd), the nucleotide diversity (π), the average number of nucleotide differences (k) and the θ statistic by sequence. These gene diversity statistics were carried out with the DNAsp 5.1  and Arlequin 220.127.116.11  Software.
We relied on the HST , KST , KST *, Z, Z* tests [54,55], Snn test  and the chi-square test on the haplotypic frequencies with permutation tests using 10,000 replicates to measure genetic heterogeneity. We also estimated the GST statistic from the haplotypic frequencies and the γST, NST and FST statistics  from the nucleotide sequences.
We estimated the possible historical female effective population numbers among the capybara populations studied-using the Migrate 3.6 Software . This program calculated θ (=Neμ) for each population considered where Ne is the effective female population size and μ is the mutation rate per generation and M (=m/μ) among the population pairs analyzed, where m is the migration rate per generation. The average μ values were 2.575×10-7 for the D-loop and 6.18×10-8 for the Cyt-b. The first average value was obtained from different studies with different mammalian species [58-66]. The second average value was estimated from other different studies but also with mammals [67-71]. Two strategies were used to calculate the effective female population size and gene flow estimates. For the first strategy, we used a maximum likelihood procedure [72-74] with 10 short Markov chains, with 500 recorded steps, 100 increments and a total of 50,000 sampled genealogies with a burn-in (number of discard trees per chain) of 10,000, and one long Markov chain, with 5,000 recorded steps, 100 increments and a total of 500,000 sampled genealogies with a burn-in of 10,000. The second strategy was a Bayesian procedure [74,75] with one long Markov chain, with 5,000 recorded steps, 100 increments, one concurrent chain and a total of 500,000 sampled genealogies with a burn-in of 10,000. This analysis is interesting to carry out because the populations with the highest historically effective female numbers are probably the original ones.
We used several methods to determine possible demographic changes across the natural history of the capybara: 1-A mismatch distribution (pairwise sequence differences) was obtained following the method of Rogers and Harpending  and Rogers et al., . The raggedness rg statistic [78,79] was used to determine the similarity between the observed and the theoretical curves. This procedure let us estimate th76e time of the beginning of a demographic change and the initial and the final population sizes . 2-We also used the Fu & Li D and F tests , the Fu Fs statistic , the Tajima D test  and the R2 statistic , originally created to detect natural selection affecting DNA sequences to determine possible changes in populati80on size [83,84]. All of these statistics and tests were obtained by means of the DNAsp 5.10 and Arlequin 3.0 Software. 3-A Bayesian skyline plot (BSP) was obtained by means of the BEAST v. 1.6.2 and Tracer v1.5 Software. The Coalescent-Bayesian skyline option in the tree priors was selected with five steps and a piecewise-constant skyline model with 40,000,000 generations (the first 4,000,000 discarded as burn-in). In the Tracer v1.5 Software, the marginal densities of temporal splits were analyzed and the Bayesian Skyline reconstruction option was selected for the trees log file. A stepwise (constant) Bayesian skyline variant was selected with the maximum time as the upper 95 % High Posterior Density (HPD) and the trace of the root height as the treeModel.rootHeight.
Two phylogenetic trees were obtained. The first one was a Bayesian tree BI; [85-87] by using the concatenated Cyt-b and D-loop genes. We concatenated both genes because we carried out two analyses which showed that both trees with each individual gene showed the same topology. First, we undertake a likelihood ratio test LRT; . We performed this test separately for mt Cyt-b and mt D-loop and for combined both mitochondrial genes. The last case showed similar likelihood to the individual likelihood of each separated tree. Second, we performed a posteriori significance tests Swofford-Olsen-Waddell- Hillis test, SOWH test; [89-91] to compare a maximum parsimony (MP) tree with the mt Cyt-b gene with MP, maximum likelihood and Bayesian trees with mt D-loop. The MP with mt Cyt-b gene was employed for generating 100 replicate data sets in the software Seq-Gen 1.2.5  which presented a uniform base composition. Goldman et al.,  demonstrated that this procedure can increase power in rejecting the null hypothesis and is better than typical non parametric tests for comparisons of a posteriori hypotheses. The differences between the log likelihood of the mt Cyt-b gene tree and the other trees herein obtained were compared with the distribution of the differences between each parametric replicate and the mt Cyt-b gene tree. As there were not significant differences among all these trees, we decided to concatenate the sequences of both mitochondrial genes.
We used the evolutionary substitution model GTR+G+I (General Time Reversible model with the gamma distributed rate varying among sites, and variable rate categories) because this was the best model estimated with Modeltest Software (Posada and Crandall, 1998). This Bayesian analysis was completed with the BEAST v1.6.2 Software . Two separate sets of analyses were run, assuming a Yule speciation model and a relaxed molecular clock with an uncorrelated log-normal rate of distribution . Results from the two independent runs (60,000,000 generations with the first 6,000,000 discarded as burn in and parameter values sampled every 1,000 generations) were combined with LogCombiner v1.6.2 Software . The effective sample size (ESS) for the parameter estimates and convergence were checked using Tracer version 1.5. The lower and upper 95% HPD of these parameters as well as the means, geometric means, medians, marginal densities and traces were also estimated with the Tracer v 1.5 Software. Posterior probability values provide an assessment of the degree of support of each node on the tree. The final tree was estimated using the TreeAnnotator v1.6.2 Software  and visualized in the FigTree v1.3.1 Software . Additionally, this program was run to estimate the time to most recent common ancestor (TMRCA) for the different haplotype lineages found. We used the value of 12.28 ± 2.3 MYA as a prior for the treeModel. rootHeight in the Bayesian analysis with a normal distribution and a standard deviation. This value was based on an estimated split between the ancestor of the capybaras and Kerodon to have occurred 16-12 MYA [95-97]. Another prior included the split of the ancestors of paca and capybara (25.25 + 2.4 MYA, average from : 26.5 MYA, and from : 24 MYA).
The second tree was also a Bayesian one with a subset of samples sequenced for another 10 concatenated mitochondrial genes. All the conditions and priors were the same as in the previous tree. To estimate possible divergence times among the haplotypes found in the capybara studied, a Median Joining Network (MJ)  was constructed by means of the Network 4.6 Software (Fluxus Technology Ltd). Additionally, the ρ statistic  was estimated and it was transformed into years. The standard deviation of ρ was also calculated . The ρ statistics is unbiased and highly independent of past demographic events. We used two different mutation rates for this task. The first was that employed in the previous analysis for obtaining the effective female size (2.575×10-7; medium mutation rate for mitochondrial markers) and the second was 7.23×10-6 (fast mutation rate for mitochondrial markers). This translates into one mutation every 200,000 and around 30,000 years, respectively.
For the Ln L criteria, the best substitution evolutionary model was GTR+G+I (-1296.276). We located 43 haplotypes and the genetic diversity statistics were Hd=0.977 + 0.019, π=0.026 + 0.009, k=13.580 + 6.314 and θ=22.246 + 7.439 at the two concatenated mitochondrial genes.
The genetic heterogeneity was considered for the capybara at the two concatenated mitochondrial genes (Table 2). Out of seven genetic heterogeneity tests, seven were significant and the gene flow estimates (Nm) were very low (around 0.4-0.6). Table 3 shows the FST values for all the population pairs. All them were highly significant, being the highest value that for Western Amazon vs. Napo River (FST=0.761). Therefore, significant amounts of genetic heterogeneity were found among capybara populations from diverse geographical areas throughout Colombia, Ecuador, Brazil and Peru. A fraction of this genetic heterogeneity is related with the geographic distance among the populations considered. Four regression models (genetic distances vs. geographic distances) were highly significant for the two mitochondrial genes we used. The geographic distance explained 78.1 to 88.3 % of the genetic distances found. The best regression equation was: log (genetic distance)=2.07 (+ 0.71) log (geographic distance)–7.94 (+ 2.23) with r=0.94 and p=0.03.
Table 2: Genetic heterogeneity for four populations of capybara analyzed at the concatenated mt D-loop and mt Cyt-b genes [1-Western Amazon (Colombian- Peruvian-Brazilian Amazon); 2-Meta-Casanare-Guainia Departments (Eastern Llanos; Colombia); 3-Cordoba and Valle del Cauca Department (Northern Colombia); 4-Napo River]. P=Probability; **P <0.001; Nm=Gene flow; df=Degree of freedom.
Table 3: FST population pairwise tests among four geographical capybara regions. 1=Western Amazon; 2=Colombian Eastern Llanos; 3=Northern Colombia; 4=Napo River in Ecuador and Peru. *P <0.01; **P <0.001.
However, the Kimura 2P genetic distances among the capybara populations studied are low or moderate (Table 4). At the mt D-loop, the Western Amazon population and the Northern Colombian population showed the highest value (3.0 %), whereas the Northern Colombian- Napo River population pair yielded the lowest (1.5%). At the mt Cyt-b, again, the Western Amazon population presented the highest genetic distances with reference to the Northern Colombian and the Napo River populations (3.6%), whilst the lowest genetic distance was among the Napo River and both the Northern Colombian and Eastern Llanos-Guainia capybara populations (1.0%).
Table 4: Kimura 2P genetic distances at the mitochondrial D-loop gene (A) and at the mitochondrial Cyt-b gene (B) among four different capybara populations. Genetic distance values (below) are in percentages (%). Standard errors (above) are in percentages. 1=Western Amazon; 2=Eastern Llanos (Guainia, Meta, Casanare in Colombia); 3=Northern Colombia (Córdoba and Valle del Cauca in Colombia); 4=Napo River in Peru and Ecuador.
The effective female population sizes were estimated by means of the Migrate 3.6 Software. Both procedures (maximum likelihood and Bayesian) offered the same results. For brevity, only the Bayesian results are commented upon here. At both mitochondrial genes (Table 5), we considered three populations. The Western Amazon population showed the highest value (Ne=852,104; 97.5 % confidence interval: 235,113-1,590,129). The mean value for these three populations was 582,794 + 422.950.
|Western Amazon||0.0527||0.0145; 0.0983||852,104||235,113-1,590,129|
|Colombian Eastern Llanos-Guainia||0.0059||0; 0.0176||95,307||0-284,790|
|Northern Colombia||0.0495||0.0081; 0.0930||800,971||130,583-1,504,854|
Table 5: Bayesian estimations of the female effective numbers by means of the Migrate Software for different capybara populations with the mt D-loop region and the concatenated mt Cyt-b gene. θ=Neμ, Ne=Female effective number and u=Mutation rate. CI=97.5% confidence interval.
Both mitochondrial gene showed, with the mismatch distribution, a clear population expansion for capybaras. The values of τ and θ0 were 7.824 and 32.615, respectively. Assuming that a generation in capybara could be 1 or 2 years, this could mean than this population expansion began for this sample around 43,111 to 86,222 YA with an initial female population size of around 1 million individuals and with a final female population of around 30.9 million of individuals (Figure 2).
The statistics used to determine population changes are in Table 6. At the global level, in the Colombian Eastern Llanos-Guainia and in the Northern Colombian, only one out of five statistics showed evidence of population expansion (Fu Fs). For the Western Amazon sample, no statistics showed evidence of population changes.
|Tajima D||Fu & Li D*||Fu & Li F*||Fu’sFs||raggednessrg||R2|
|NorthernColombiancapybarasample||P[D <0.396]=0.3661||P[D*<0.707]=0.2142||P[F*<0.716]=0.2147||P[Fs<5.111]=0.0081*||P[rg<0.0281]=0.0412+||P[R2 <0.109]=0.3173|
|Eastern Colombian Llanos-Guainiacapybarasample||P[D <0.100]=0.5143||P[D*<0.199]=0.3965||P[F*<0.198]=0.3902||P[Fs<5.581]=0.0140+||P[rg<0.0175]=0.0161+||P[R2 <0.142]=0.7422|
|Western Amazon capybarasample||P[D <0.091]=0.5156||P[D*<0.459]=0.3022||P[F*<0.408]=0.3401||P[Fs<1.736]=0.2067||P[rg<0.0177]=0.0451+||P[R2 <0.133]=0.5710|
Table 6: Demographic statistics applied to the overall, Northern Colombian, Eastern Colombian Llanos-Guainia, and Western Amazon Basin capybara samples studied at the mt D-loop region and at the concatenated mt Cyt-b gene. +P<0.05; *P<0.01, significant population expansions.
We conducted one BSP analysis for the capybara (Figure 3). Based upon both mitochondrial genes, for the overall sample, the female population slightly underwent growth in the last 1 MY, although this population expansion hardly increased in the last 0.5-0.3 MYA. Thus, some evidence of population expansions were obtained for the capybaras.
Figure 3: Bayesian skyline plot analysis (BSP) to determine possible demographic changes across the natural history of the overall capybara sample sequenced at two concatenated mitochondrial genes (D-loop and Cyt-b). On the x-axis, time in millions of years; on the y-axis, size of the effective female population.
The first BI tree (Figure 4), with both concatenated mitochondrial genes, showed that the temporal split between pacas and the clade Kerodon + capybara was around 24.57 MYA and the split between Kerodon and capybara was around 11.55 MYA. Within the capybaras, Western Amazon was the first clade to diverge (p=1). The other large cluster was comprised of the Colombian Eastern Llanos group (in this case, including the specimens of the Guainia Department; p=0.78) and the Northern Colombian group (p=0.67) plus the Amazonian animals from the Napo River in Ecuador and Peru (p=1). The split between the main Western Amazon clade and the others occurred around 6.6 MYA (95% HPD: 6.2-10.18 MYA). The temporal separation of the central Brazilian Amazon haplotype from the other western Amazon haplotypes of this clade was around 4.6 MYA. The temporal split process within the Peruvian individuals (excepting those from the Napo River) began around 1.75 MYA (95% HPD: 0.25-3.67 MYA). The split between the Colombian Eastern Llanos population and the trans-Andean North Colombian population was around 2.8 MYA (95% HPD: 1.05-6.61 MYA), whereas the haplotype diversification within the first quoted group was around 2.6 MYA (95% HPD: 0.42-4.44 MYA). The temporal split between the trans-Andean and Napo River capybara populations occurred approximately 1.9 MYA. The haplotype diversification process within the trans-Andean population began around 1 MYA.
Figure 5 shows the BI tree for a subset of capybaras sequenced at 10 concatenated mitochondrial genes. Some results, such as the separation between capybara and the sister group (Kerodon; p=1), agree quite well with that observed at the previous tree. Also, as in the previous case, the Western Amazon clade and the Colombian Eastern Llanos-Northern Colombian clade were clearly separated (p=1). In this case, different individuals from the Meta and Guainia Departments were intermixed with the trans-Andean individuals sampled in Northern Colombia. The temporal split between the capybara and Kerodon occurred about 11.88 MYA (95 % HPD: 9.99-13.88 MYA). The separation of the capybara ancestors from the Western Amazon and the Colombian Eastern Llanos- Northern Colombia clade happened about 9.09 MYA (95% HPD: 5.68- 11.8 MYA). The diversification within the Western Amazonian clade began around 7.75 MYA (95% HPD: 3.28-10.38 MYA). Within the Colombian Eastern Llanos-Northern Colombia clade, it began around 7.41 MYA (95% HPD: 5.01-8.63 MYA). These temporal estimates are higher than those obtained in the first tree, although 95% HPD of both trees are overlapped.
Figure 5: Bayesian tree (BI) for 25 capybara samples sequenced at 10 concatenated mitochondrial genes (12S rRNA, 16S rRNA, ND1, ND2, COI, COII, COIII, ATP8, ATP6 and ND4) studied. In nodes, the numbers are the average divergence times and, after slash, the 95% HPD divergence times (range) among clades in millions of years. In blue, the Western Amazon individuals; in green, the Guainia and Eastern Colombian Llanos; in red, the trans-Andean individuals and in lilac, one Napo River individual.
Figure 6 shows the results of the MNJ analysis. With the first mutation rate, the main Amazon haplotype and the Northern Colombian haplotypes split around 3.0 + 0.659 MYA, while this main Amazon haplotype diverged from the Casanare and Inirida haplotypes, 2.4 + 0.606 MYA and 1.0 + 0.304 MYA, respectively. The temporal splits between the Eastern Llanos and the Northern Colombian haplotypes and the last ones and the Napo haplotypes were around 1.574 ± 0.419 MYA and 0.675 ± 0.248 YA, respectively. Although the Bayesian split estimates are somewhat higher than those obtained with this first MNJ analysis, the last values are within the range of the Bayesian ones. However, if we employed the second mutation rate in this analysis, all the temporal splits in capybara were during the Pleistocene. For example, the main Amazon haplotype and the Northern Colombian haplotypes diverged around 0.450 ± 0.098 YA, while the temporal splits between the first haplotype with reference to the Casanare and Inirida haplotypes were 0.392 ± 0.094 MYA and 0.150 ± 0.045 MYA, respectively. The temporal splits between the Eastern Llanos and the Northern Colombian haplotypes and the last ones and the Napo haplotypes were around 0.236 ± 0.062 MYA and 0.101 ± 0.037 YA, respectively. This last procedure was unique and refers to the molecular evolution of capybara exclusively within the Pleistocene.
Figure 6: Median Joining Networks (MJN) for the two concatenated mitochondrial genes sequenced (D-loop and Cyt-b) in capybara; blue circles=Western Amazon (D), Guainia & Eastern Colombian Llanos (A) haplotypes (cis-Andean); yellow circles=Northern Colombia (Cordoba and Valle del Cauca) haplotypes (B) (trans-Andean); related to these haplotypes, there are three Napo River haplotypes (C) (in blue). Red circles are intermediate haplotypes not found.
Strong genetic heterogeneity among capybara populations of different river basins and possible Western Amazon origins
Capybara showed a striking genetic heterogeneity between populations linked to different river basins. Related with this, there was a highly significant positive relationship between geographical and genetic distances between these populations placed in different river basins. This agrees extremely well with recent findings by other authors  for populations of capybara in Argentina, Paraguay and Venezuela. In Paraguay, the rivers were totally disconnected from each other and there were noteworthy differences in genetic structure among them. In contrast there is no striking contrast in genetic structure among the Argentinian rivers or among the Venezuelan rivers because periodic flooding connected these rivers within these respective countries. However, similar to the present work, there was a significant difference in heterogeneity among the major river basins of the three countries. Therefore, the evolution of the capybara seems to be related to the climatic and geological changes which affects Neotropical rivers. This allows us to exclude some of the hypotheses on the origins of the Amazon biodiversity quoted in Table 1: STH, RBH, CDH, DVH and GH. In contrast, PH, RLH, RH, RRH and HRCH could be important to our understanding the evolution of capybara.
Many of our results support the Western Amazon population (less the Napo River population) as the original population which derived Guainia and the Eastern Llanos populations. This in turn originated the trans-Andean Northern Colombia population (Córdoba and Valle Departments). The ancestor of this population originated another population which returned to the Northern Ecuadorian and the Peruvian Amazon (Napo River) (Figure 7). Many data support the Western Amazon population as the original population. The Western Amazon population showed the greatest effective female numbers (and the Napo River population, the lowest ones) and there was no clear expansion for this population. Also, the ancestors of this population were the first to diverge in all of the trees. The Western Amazon area was also the focus for the dispersion of many other mammalian taxa such as the lowland tapir [102-105], the jaguar , and diverse Neotropical primate genera [107-110].
This results for some mammals disagrees with the fact that, for other species, the focus of dispersion is from the Eastern Amazon. Some authors [111,112] claimed that the chronology of the vicariant events within the Amazon supports an Eastern-Western establishment of the Amazon tributaries. This could agree with some geological data (but not others), which previously showed that the South Western Amazon began to take its current shape only in the last 3 MYA  while the Brazilian and Guyana Shields have been stable during at least the last 10 MYA [114,115]. The analysis of the Curimatidae, Steindachnerina, supported that an old splitting event gave rise to the species from the Tocantins River . Similar results were reported for a ro115dent, Proechimys, [117,118] and the Neotropical bird, Xiphorhynhus  which is correlated with the hypothesis that the Western Amazon arches were the atest to uplift.
Miocene and Pliocene Amazon biodiversity hypotheses for the mitochondrial diversification in capybaras
The two BI trees showed the initial differentiation of the ancestor of the Western Amazon population relative to the ancestor of all the other capybara populations around 6-9 MYA, during the late Miocene. Primordial events did occur during the Miocene (20-5 MYA). 1- The large Pebas wetland of shallow lakes and swamps (and/or the hypothetical Pebasian Sea) developed in the Western Amazon, 17- 11 MYA. This system of lakes was formed by the Early Miocene and grew to 1,000,000 km2 in the Middle Miocene (around 15 MYA) as the Andean uplift delivered more water to the developing Amazon Basin and it was bounded on the east by the Purus Arch. This aquatic environment would have precluded extensive and long-term terra firme habitats until the Late Miocene, when the Purus Arch was breached and a transcontinental Amazon drainage was created around 11.8- 10.0 MYA [113,115the present, which could be increased sedimentation during Plio-Pleistocene glacial cycles. 2- There were two Miocene marine ingressions, one 15-13 MYA (21-16 MYA for some authors ) tectonically and eustatically controlled, and the other, younger than 11-10 MYA, and tectonically controlled. 3-The Magdalena Drainage Basin was isolated (12-11 MYA) from the Amazon basin [113,125]. This was provoked by the intensified uplift in the Andes around 12 MYA, especially the tectonic uplift of the Eastern Cordillera in Colombia, the Santander Massif and the Merida Andes which caused the closing of the Amazon–Caribbean connection. 1254- Lastly, the Amazon River was completed formed around 11-7 MYA [113,115,120-123,126-128]. The final breakthrough of the Amazon River towards its modern course occurred with the uplift of the Central Andean Cordillera  and related to the rise of the Purus arch . Therefore, a large fraction of the current Western Amazon was under the water. This ecological scenario could favor the adaptation of some rodents to the aquatic systems. Indeed, family Hydrochoeridae appeared in the fossil record of South America about 10 MYA , which agrees quite well with these marine introgressions during the Miocene. However, Mones and Ojasti  claimed that the capybara probably evolved from an unknown species of Cardiatherium about 2 MYA in Southern South America, with all the fossil record obtained during the Pleistocene. If we consider the Bayesian temporal splits we obtained, the fossil record is clearly incomplete and the notion of the capybara appearing in Southern South America is incorrect. Nevertheless, we have not analyzed capybaras from the southern range of the current distribution (Bolivia, Paraguay and Argentina). Therefore, we still need to determine if the original ones are from the Western Amazon or the Southern South America. Additionally, both the Amazon and Orinoco Basins were isolated around 8-10 MYA (Miocene) after the uplift of the Vaupes Arch in the foreland of the Andes [113,120-122]. It’s likely that the formation of several arches and the orogeny of the Merida Andes could have helped to isolate Amazon and Orinoco capybara populations as well as capybara populations within of the Orinoco drainage Basin. El Baul Arch (although the uplift structures created by this arch are transient and migrated across time; ) and the Arauca Arch (in the Apuré area, which may be the northern extension of a large arch located in the Colombian Llanos Orinoco Basin, which has not yet reported) could have uplifted around 8-5 MYA (Late Miocene-Pliocene). Also the Merida Arch (related with the isolation of the Maracaibo drainage Basin around 10-8 MYA [115,131]), could play an important part in creating barriers within the Orinoco Basin responsible for the molecular heterogeneity of the capybara in this area. Our Bayesian data agree quite well with these temporal geological changes. If our estimates are right, then, the PH hypothesis can explain the origins of the molecular diversification within the capybara.
In an identical sense, the divergence between Pseudoplatystoma magdaleniatum (Magdalena River) and P. corruscans (Parana River) from the other cat fishes was found to be synchronous at around 11.8 MYA. The ancestor of the Orinoco’s species (P. orinocoense and P. metaense) diverged from the remaining species at around 8.2 MYA . These estimates are consistent with Hoorn et al.,  and Lundberg et al., , who dated the establishment of the Magdalena, Parana and Amazon Basins during the Late Miocene between 11.8 and 10 MYA and the primary isolation of the Orinoco between 8.0 and 5.0 MYA.
The temporal split of the Negro River (Central Brazilian Amazon) individual with regard to the Peruvian Amazon individuals was around 4.6 MYA, which coincides with the initiation of the Pliocene. During the Pliocene, another PH event was very important for the generation of biodiversity in the Amazon Basin. Espurt et al., [3,133] demonstrated that the Nazca Ridge subduction imprint had a significant influence on the eastern side of the Andes by means of the Fitzcarrald Arch. This uplift is responsible for the atypical three-dimensional shape of the Amazon foreland Basin. During the Miocene, the Fitzcarrald Arch did not exist. Related with the Nazca Ridge subduction, the arc volcanism in the Peruvian Andes ceased around 4 MYA , thus the older time estimation of the Fitzcarrald Arch uplift is not older than 4 MYA (Pliocene), which is relatively similar to the split between the Western Amazon capybara population and the central Brazilian Amazon capybara population. To ratify the influence of this PH event on the molecular diversification of capybara, it is necessary to analyze animals from Southern Peru and Bolivia, mainly. However, some authors [135,136] suggested that many arches could have been buried under subsequent sedimentary deposits and thus not be responsible for the differentiation among taxa.
Another possibility is also in the Pliocene and very related to the RLH. Many phylogenetic studies of Amazon organisms emphasize that populations and species from the lowlands of the Upper Amazon Basin originated more recently than other organisms in the Brazilian and Guyana shields [116-119,137-139]. This finding was attributed to the Miocene-Pliocene marine incursion that occurred around 5 MYA and provided considerable species range contractions for terrestrial and freshwater biota [8,140,141]. This isostatic change prompted a rise of sea level (100 m) that considerably reduced the area of suitable habitat for freshwater organisms, such as the capybara. Consequently, the colonization of the lowlands from the Upper Amazon Basin occurred after the marine retreat around 4 MYA.
The split between the Eastern Llanos-Guainia capybara population and the trans-Andean Northern Colombia capybara population and the diversification within the Eastern Llanos-Guainia population began in the Middle and the Late Pliocene, following our Bayesian estimates. Also, the initial mitochondrial capybara diversification, by using the MJN procedure with the first mutation rate estimate, between the Western Amazon and the Eastern Llanos and the trans-Andean Northern Colombia capybara populations coincided with the middle and the late Pliocene. Although, these results could be explained as a consequence of the PH, some Pliocene events following other Amazon biodiversity hypotheses could help to understand the molecular diversification of the capybara if these temporal estimations are right.
The existence of a relatively recent Amazon Lake in the Late Pliocene until the Middle Pleistocene (by eustatic marine changes, RLH) could be of great importance for the diversification of the capybara. This event connected the Orinoco and Amazon Basins after its separation during the Miocene [122,123]. This agrees quite well with the views of many authors [8,9,135,142-145]. A point supporting RLH is the existence of a confirmed Pliocene (5 MYA) sea rise of 100 m for a minimal duration of 0.8 MY . Thus, marine transgressions could have had a great influence on biotic diversification within the Amazon long after the Miocene. Our data agree quite well with some recent studies which indicate that the wetland system may have persisted until the Late Pliocene or Pleistocene [135,147,148]. Aleixo and Rossetti  claimed that the Midwestern Amazon may be the most dynamic area for new colonization, due to a relatively recent reduction of the Amazon Lake and subsequent floodplains in that area. Although, some authors have resisted this view [5,150], the mitochondrial capybara diversification could perfectly agree with this hypothesis. When this recent Amazon Lake began to break apart or dry out (Late Pliocene), the capybara populations were isolated. Thus, the interaction of the RLH during the Pliocene and the RRH when the lake dried (late Pliocene and Upper Pleistocene) could generate an important fraction of the diversification in the capybara.
All of these hypotheses regarding capybara diversification could be supported by molecular data. This confirms that an important diversification stage for Amazon organisms occurred between 10 MYA and 3 MYA [151-154], and therefore challenged the Pleistocene changes as a potential explanation for the origin of many South American aquatic biotas.
Mitochondrial diversification occurred within the Western Amazon population during the Upper and the Middle Pleistocene (both for the Bayesian trees and the MNJ with the first mutation rate). It also occurred within the trans-Andean Northern Colombia population and the split between this last population and that of the Napo River in Ecuador and Peru. At least 20 major glacial period happened (with dry periods) during the Pleistocene (beginning around 1.6-2.4 MYA). The mean global temperature was 4°C lower than today . The combination of the RH and RRH could have some importance for the most recent clade splits during the Pleistocene. During the dry periods, some lagoons and rivers may have decreased their channels and surface waters. The social and reproductive system in the capybara could have helped to isolate populations in correlation with the RH and RRH. The capybara herds have an average group size of 5.6 individuals in July (peak rains) to 15.9 in March (driest period) in the Venezuelan Llanos . Capybara spend the majority of their time on a small core area of less than 1 ha in the Venezuelan Llanos. During drought, capybara herds congregate around water holes conforming temporal aggregations in the hundreds. This behavior can prevent gene flow between populations of different rivers when they decreased their size during dry periods. Thus, by using Bayesian and MNJ with the first mutation rate estimates, we found possible evidence of genetic diversification in the capybara during the Miocene, Pliocene and Pleistocene epochs, with possible noteworthy contribution of the PH, RLH, RH and RRH.
Another possibility: very recent Pleistocene diversification in capybara?
In contrast, if we used some MNJ results with the first mutation rate (slow) estimate and all the results with the second mutation rate (fast) estimate, all the capybara time splits were within the Middle and Late Pleistocene. It is interesting to note the discrepancy between the BI (fossil-calibrated DNA phylogeny) and borrowed molecular clock (MNJ). We found this to be true in many of our other research projects involving a variety of Neotropical mammalian species [107-110,156- 161]. If these last temporal splits are accurate, the PH and RLH could be excluded as main hypotheses in the generation of genetic capybara diversity. The RH and RRH are very important along with the HRCH. It’s interesting to note that for the capybara, the Western Amazon showed the highest levels of gene diversity. This coincides with the Napo and the Inambari Refugia stated by the RH [162,163]. The Napo and Inambari Refugia have been determined to be of importance for butterflies [137,164], characiform fishes , amphibians and lizards  and birds [165-167]. Eight aquatic freshwater refuges have also been detected for characiform fishes .
The HRCH claims that important changes occurred in Neotropical rivers during the Pleistocene. The PH sustains that the Amazonas River became transcontinental around 11.8-11.3 MYA  and other hypotheses sustain that during the Pliocene several changes controlled by tectonic processes, such as flexural and dynamic uplift, subsidence and ridge subduction [3,168], were decisive in the drainage of the current Amazon River. In particular, the uplift of the Iquitos Forebulge [169,170], the Fitzcarrald Arch  and the Purus Arch  shaped the Pliocene Amazon drainage network. However, mounting evidence suggests that many changes in the Amazon River occurred during the Pleistocene. Neo-tectonic reactivations occurred during this period along the Negro and Amazon Rivers [172-174]. Also, studies of Neogene sediments from the Ecuadorian  and Amazon  Basins, together with the fourfold increase of sedimentation rates in the Amazon Fan 2.4 MYA , also indicate that important and complex changes in the Amazon Basin drainage occurred during the Pleistocene. Rossetti et al.,  estimated a maximum age of 43,000 YA for the Quaternary terrace deposits along the Amazon River. Coimbra-Horbe et al.,  also determined important changes in the composition of the sediment source in the Amazon River, which correlates with a significant reorganization of the architecture of the Amazon River since the Miocene uplift of the Andes. Furthermore, some important climatic changes during the Pleistocene have been described in the Amazon Basin [178-180]. Additionally, the Pleistocene was a time of broad oscillations in the sea level relative to the present, most resulting in lower levels but some were higher [181-183]. Particularly large oscillations occurred in the last 1 MYA  and major sea level rises occurred around 50,000 and 400,000 YA when the sea rose by 100 m [184,185]. Such a rise is enough to flood the entire Eastern and Central Amazon Basins, and the eastern portions of the Western Amazon basin . A more moderate temperature increase and associated sea level rise occurred in the last interglacial (around 130,000-115,000 YA) with sea levels rising between 6 and 9 m above present levels . Some examples of Pleistocene river changes are as follows. The Bolivian Beni River has experienced anticlockwise movement over the course of the last 0.01 MY. This has led to a change in flow orientation of 45 degrees and is a response to the continuous uplifting of the Andes piedmont . Additionally, the Peruvian rivers are extremely dynamic , with striking shifts in flow recorded over very recent geological time. In favor of the HRCH, such dynamic behavior by river courses and periodic river-switching events could generate fragmentation and a degree of isolation of riparian populations .
If the MNJ estimate with the second mutation rate were accurate, the HRCH should be extremely important in understanding the mitochondrial diversification of the capybaras following the more recent analyses of Amazon stratigraphy and the paleoenvironment [11,147, 148,190]. This suggests a much younger, Plio-Pleistocene history for the transition from a lacustrine system in Western Amazon to the current drainage system. However, both HRCH and RH may simultaneously act. Iriondo and Latrubesse  have offered evidence for a dry Late Glacial climate in the central portion of Lower Amazon and commented on the reduced discharges of the rivers in this area and highly seasonal climatic conditions with enhanced Trade Winds. The oxygen isotopic composition of planktonic foraminifera recovered from a marine sediment core in a region of the Amazon River discharge shows that the Amazon Basin was extremely dry during the very cold Younger Dryas period (13,000-11,600 YA). The discharge was reduced by at least 40 % as compared with that of today . However, recent isotopic evidence suggests that precipitation was lower in the Eastern Amazon, but not in the Western Amazon (where we studied capybaras), during the last glaciation . Unfortunately, knowledge of the recent history of forest cover in the Amazon is limited and contentious .
Some of our genetics results could agree with this second view. The mismatch distribution at the concatenated mt D-loop and Cyt-b detected a clear population expansion around 43,000-86,000 YA. This is during the last phase of the Pleistocene. The BSP also detected a population expansion for the overall capybara sample around 0.5-0.3 MYA, also during the Pleistocene. These results also could agree quite well with the aforementioned remark by Mones and Ojasti . They claimed that the entire capybara fossil record to date is exclusively from the Pleistocene.
Another point which could favor both MNJ temporal split estimates (more recent diversification) is the fact that our data don’t seem to support the existence of different capybara species. The Northern Colombian capybara population could not be considered a different species (H. isthmus) because is clearly descended from the original Amazon and the Orinoco capybara populations. Also, part of their descendants are now in the Amazon (Napo population) and both the original Amazon and the descended Napo populations likely intermixed and there is no evidence of morphological differences between these mitochondrial lineages. This means that reproductive cohesiveness probably exists between the capybara populations, following the Biological Species definition [194,195]. They are the same species. Peceño  showed that the capybara from the Venezuelan Llanos has 66 chromosomes (FN=102), whilst one specimen from the Maracaibo Lake Basin (Northern population) had 64 chromosomes (FN=104). This karyotype could be derived from the first by one pericentric inversion and one Robertsonian change. But polymorphic and variable chromosome numbers can be found within the same species. Additionally, Peceño  analyzed 44 enzymatic blood loci and estimated very small genetic distances (D=0.0056) between two samples representing the two “a priori” capybara species (37 individuals from the Apure State and 16 from the Maracaibo Lake basin). Thus, these enzymatic loci, as well as our mitochondrial sequence results, do not agree with two well differentiated capybara species. Thus, the Northern Colombian population (and probably the Panama one) should be considered as a subspecies (H. isthmus). It may also be preferable for these different capybara populations to be considered ESUs (Evolutionary Significant Units; ). If we follow this line of reasoning, we should consider, at least, four ESUs: (1) Western Amazon, (2) Colombian Eastern Llanos and Guainia, (3) Northern and Pacific Colombia and (4) the Napo River. However, the level of hybridization must be investigated as well as the geographical frontiers of the Napo River population and also the possible existence of chromosomal differences among these ESUs. Other ESUs probably exist especially in the Southern distribution of the species. However, the morphological differences are very limited among the capybara populations and may align better with relatively recent differentiation (in the last 2-3 MYA). This agrees with a drainage of the Western Amazon wetland system approximately 2.7- 2.0 MYA and the formation of the Lower Amazon River . This is a better fit than an older differentiation (6-9 MYA), unless some stasis morphological evolution exists in the capybara. However, the genetic distances we obtained agree better with a unique species of capybara and thus more recent diversification in this species. Bradley and Baker  claimed, for mt Cyt-b, that values <2% would equal intraspecific variation, values between 2% and 11% would merit additional study, and values >11% would be indicative of specific recognition. Avise  determined 5-7% of differences at the mt control region for different species and around 2% for subspecies in mammals. The highest genetic distances at mt Cyt-b and mt D-loop (2.6-3.0%) are between the Western Amazon and the other populations (3.3-3.6% and 2.6-3.0 % respectively). These values are typically of the status of different populations or even subspecies but not of different species. Therefore, our data could support the existence of a unique species with relatively recent differentiation (Late Pliocene and Pleistocene). The HRCH implies that different regions of the Amazon may have undergone distinct rates of landscape change, with the Western area of the Madeira River being the most dynamic (the area sampled here for capybara), Today it is covered by the 200Amazon sedimentary basin [135,149,200]. This area has only recently transitioned from a lacustrine system to floodplain and in turn to lowland forest and retains many patches of seasonally inundated forest. In the HRCH, the formation of the floodplain forest (várzea) habitat as the Mamirauá Reserve is very recent (5,000-15,000 YA; Late Pleistocene and Holocene; ). The process of destruction and reconstruction of floodplain habitat can be rapid, constantly transforming the river ways [202,203], which could highly affect capybara populations. Klammer  suggests that floodplain forest formation has probably played an important role in the distribution of many species. This formation could have similarly influenced capybaras. The HRCH seems to play a fundamental role in the diversification process of numerous avian taxa [11,204-209] and some Primates taxa [110,210]. The temporal splits obtained with the MNJ procedure could also support the RH, because glacially-driven diversification mainly occurred within the last 0.9-0.02 MYA. There was a peaking at the Last Glacial Maximum at 26,000-18,000 YA when the forest was likely at its smallest size. However, as capybara are highly associated with water bodies, the HRCH may be more relevant than RH. Nevertheless, mountain formation, small-scale neotectonic changes at different landscape levels as well as Milankovich cycles can influence river dynamics. Extremely dry climatic events in the past could have significantly reduced the Amazon River’s size and discharge. These changes have the potential to eliminate or dampen the effects of riverine vicariance barriers, especially in the headwaters, which would promote expansion in many species.
Nevertheless, there is one point which doesn’t support this exclusive Pleistocene hypothesis based on the HRCH. The continuous uplift of the Andean cordilleras and the split of the Magdalena and Amazon River Basins is not a very recent process during the Pleistocene and clearly our results showed that the capybaras sampled at the Ecuadorian and Peruvian Amazon Napo River were derived from the Northern Colombian Magdalena capybara population. The divergence must have occurred when the Andean cordilleras were considerably lower than during the Pleistocene and some connection existed between Magdalena and the Northern-Western Amazon Basins. We know that the geological isolation of the Magdalena Basin from the Amazon drainage Basin began with the rise of the Eastern Cordillera of Colombia (12 MYA; ) suggesting a minimum age for the divergence of lineages in cis- and trans-Andean Basins. Indeed, the fishes of La Venta fauna in the Villavieja Formation at the Magdalena Valley include many living forms of fishes that are currently known only211 in the Amazon and Orinoco Basins. Many of these species are indistinguishable from living species (Arapaima gigas and Colossoma macropomum, for instance), or are closely related to current living Amazon species (Brachyplatystoma and Hoplosternum, for instance) . This data doesn’t support a rec216ent connection for the Magdalena and Napo capybara populations. An alternative hypothesis is that human beings transported capybaras from Magdalena valley to Northern-Western Amazon when our species colonized South America (around 10,000-30,000 YA; [213- 216]). However our more recent temporal split between these capybara populations preceded (100,000 YA) the human colonization. Another alternative explanation could be similar to that of Torrico et al.,  for the catfish P. reticulatum. They found a relatively young ages inferred for the divergence of P. reticulatum populations from the Parana and the Amazon Basins (between 0.8 and 1.5 MYA), which was not consistent with the primary establishment of the Parana in the Late Miocene (around 11-10 MYA; [113,115]). They stated that headwater capture events and temporary connections between the headwaters of the Amazon and the Parana Basins promoted speciation by longdistance dispersal and further allopatric divergence [112,152,154]. The presence of river connections during the Pleistocene among basins that had been previously disconnected since the Miocene could also explain the connection between the Magdalena and Napo River Basin capybara populations.
The differentiation of the Amazon and Orinoco capybara populations
There is an interesting point independent of the three different temporal split estimates (6-9 MYA, around 2 MYA, and around 0.3 MYA) between the Amazon and the Orinoco capybara populations. All the estimates precede the current connection, by means of an anastomosing phenomenon named stream capture , between the Upper Orinoco River and the Amazon Basin (Negro River) throughout the Cassiquiare Channel (250 km). Recent elevation changes (since the Holocene 10,000 YA) favor the capture by the Negro River [218, 219], although this elevation increase is only about 20 m between its source at the Upper Orinoco and its mouth at the Negro River . This indirectly shows the previous existence of other river systems in the Western and/or Central Amazon that connected the Orinoco and Amazon Basins around 0.5-2 MYA. Evidence in favor of this comes from geological studies. Almeida-Filho and Miranda  detected the relicts of a large ancient drainage system hidden by the tropical rain forest of this area. It was not visible by optic or synthetic aperture radar images. For the oldest BI estimates, the PH and RLH could explain why the split of the capybara populations was during or slightly after the separation of the Orinoco and the Amazon Basins (8-10 MYA).
Hubert and Renno  determined the characiform fishes of the Upper Negro River to be nested within the Orinoco River cluster because there were a high number of shared species between these two rivers. This result is consistent with extensive fish exchanges between the Orinoco and Negro Rivers throughout the Cassiquiare Channel, as was previously reported . However, the strong ichthyological difference between the Upper Negro River and the Amazon River suggests a recent headwater capture event between the two rivers. The fact that the Negro River capybara individual (sampled at the Middle Negro River) was clustered with the Western Amazon cluster and the animals sampled at the Inirida and Guaviare Rivers (affluent from the Upper Orinoco River) were clustered with the Colombia Eastern Llanos capybara population, agrees quite well with that reported for fishes . But this could not happen with the current Cassiquiare Channel because its recent formation (neither of our estimates coincided within the last 10,000 YA). This provides indirect evidence in favor of other, older, connections between the Amazon and Orinoco drainage Basins.
Which Amazon biodiversity hypothesis is the best for the capybaras?
There is a long history of controversy about which of the Amazon diversity hypotheses is the most important . Here we show, with the example of the capybara, the difficulty in selecting one hypothesis on the generation of molecular differentiation in the Amazon. We have a wide range of temporal splits that could overlap with the majority of the generated hypotheses that explain the high biodiversity of the Amazon. A problem is the selection by the researcher to focus on geological, environmental or climatic changes during a specific epoch (Miocene, Pliocene and Pleistocene) to justify the divergence within the taxa they studied. This study clearly reflects this problem. We can justify the generation of divergence among the capybara populations anytime within the last 10 MYA because we do not have precise fossil data nor exact mutation rates for the markers used. We suggest that the best way to fix this problem is to use exact mutation rates for borrowed molecular clocks since it should be easier to estimate these exact mutation rates for several mitochondrial and nuclear markers than to discover all the fossils (and their correct ages) of a determined lineage. Indeed, some of our findings suggest that, rather than there being a single predominant Amazon biodiversity hypothesis, the evolution of the capybara should be the result of an interaction between Miocene marine incursions, uplift of the paleoarches and mountains (PH), more recent marine incursions and Amazon islands (RLH), relative recent connections allowing cross-drainage dispersal (HRCH), and riverine and forest refugia (RRH and RH).
Another problem in trying to determine a single hypothesis as the key one to generate Amazon biodiversity is that each species represent a different “experiment.” First, the geographical origin of each taxa which generated the current species could be different. Each colonization process should be viewed as a Markovian process, where the immediate anterior step conditions the following step. Most of the current mammalian species in South America are thought to have colonized this part of the world subsequent to the completion of the Isthmus of Panama, about 3 MYA, as part of the Great American Biotic Interchange (GABI, ). These mammals were basically affected “in situ” by Pleistocene events, being this not the case of the capybara. Second, each species (or genome) has a different capacity and adaptability to colonize new areas because the additive effects of the natural selection is different for each taxon. Many researchers have searched for shared phylogeographic (or phylogenetic) relationships among taxa that would suggest a common mechanism of biological diversification [137, 223- 225]. Nevertheless this is not an easy task to generate global identical hypotheses to explain the high levels of genetic variability found in the many Amazon taxa. Therefore we agree with previous studies reporting strikingly different patterns and timings of the origins of Neotropical species. No single hypothesis is sufficient [226-230] because multiple factors on multiple occasions have differentially impacted the genetic structure of the Amazon’s diverse taxa.
Thanks to Dr. Diana Alvarez, Pablo Escobar-Armel, Luisa Fernanda Castellanos-Mora and Nicolás Lichilín for their respective help in obtaining capybara samples during the last 18 years. Many thanks go to the Peruvian Ministry of Environment, to the PRODUCE (Dirección Nacional de Extracción y Procesamiento Pesquero from Peru), Consejo Nacional del Ambiente and the Instituto Nacional de Recursos Naturales (INRENA) and to the Environment Ministry at Coca (Ecuador) for their role in facilitating the obtainment of the collection permits in Peru and Ecuador. Special thanks to the Ticuna, Yucuna, Yaguas, Witoto and Cocama Indian communities in the Colombian Amazon, as well as the Bora, Ocaina, Shipibo-Comibo, Capanahua, Angoteros, Orejón, Yaguas, Cocama, Kishuarana and Alama communities in the Peruvian Amazon, and the Kichwa, Huaorani, Shuar and Achuar communities in the Ecuadorian Amazon for helping to obtain capybara samples.