Genetic structure of Tupinambis teguixin ( Squamata : Teiidae ) , with emphasis on Venezuelan populations

Tupinambis teguixin, the common tegu, is the only species of the genus found in Venezuela. It is distributed in different bioregions in the Neotropics, some of them separated by geographic barriers that may restrict gene flow among populations. Thus, to assess this possibility, we tested the Paleogeographic hypothesis and the Riverine hypothesis for the divergence among populations. To this end, we evaluated the degree of genetic structuring in six populations of T. teguixin from Venezuela, plus one from Brazil and one from Ecuador. We used two molecular datasets, one with the populations from Venezuela (Venezuela dataset, 1 023 bp) and one including the other two (South America dataset, 665 bp), with 93 and 102 concatenated sequences from cytochrome b and ND4, and 38/37 haplotypes. We used three measures of genetic diversity: nucleotide diversity, haplotype diversity and number of polymorphic sites. Gene flow was estimated with the statistic ΦST and paired FST values. We also constructed a haplotype network. We found genetic structuring with (1) ΦST = 0.83; (2) high paired FST estimates (0.54 0.94); (3) haplotype networks with a well-defined geographic pattern; and (4) a single shared haplotype. The genetic structure does not seem to stem from geographic distance (r = 0.282, p = 0.209), but rather the product of an historic biogeographic event with the Mérida Andes and the Orinoco River (71.2 % of the molecular variance) as barriers. We propose the Zulia population as an Evolutionary Significant Unit and that the other populations be temporarily considered Management Units, pending further data. Populations Delta and Guri should form a single Management Unit since they share a haplotype. Rev. Biol. Trop. 63 (4): 1235-1249. Epub 2015 December 01.

The South American subcontinent harbors the greatest biological diversity on Earth and its origin is usually attributed to historical factors involving vicariant mechanisms or allopatric speciation. The main hypotheses are: the Refuge hypothesis, the Gradient hypothesis, the Disturbance-vicariance hypothesis, the Paleogeography hypothesis, and the River or Riverine hypothesis (Haffer, 1997;Moritz, Patton, Schneider, & Smith, 2000), with the latter two being the most frequently tested. The Paleogeography hypothesis considers changes in the relative distribution of land, sea or landscape as a result of tectonic movements (e.g., mountains rising) or sea level fluctuations as major factors in the distribution of species in the Neotropics (Haffer, 1997). For instance, a number of studies have shown the impact of the rise of the Andes Cordillera on the dispersal of individuals of the same or different species, by assessing the genetic structure of populations on both sides of the Andes and/or by describing their restricted distributions due to the barrier (Mones, 1991;Brumfield & Capparella, 1996;Arrivillaga, Norris, Feliciangeli, & Lanzaro, 2002;Romero, 2003;Bernal, Guarnizoi, & Luddecke, 2005;Albert, Lovejoy, & Crampton, 2006;Brumfield & Edwards, 2006;Weir & Price, 2011). The effect of orogeny has also been tested elsewhere, in regions such as the Northwest Africa Atlas Mountains (Brown, Suárez, & Pestano, 2002); the Central Mountain Range in Taiwan (Lee, Jiang, Su, & Tso, 2004); the mountains of Eastern Australia (Smissen, Melville, Sumner, & Jessop, 2013); the Pyrenees and the Cantabrian Mountains in Spain (Milá, Surget-Groba, Heulin, Gosá, & Fitze, 2013).
Lizards of the genus Tupinambis (Teiidae) are one of the largest lizards in the Americas. They are found only in South America, mostly East of the Andes, with only a few populations to the West of the cordillera, in Colombia and Venezuela. The systematics of this genus has changed several times over the years, with T. teguixin as one of the most controversial species (Peters & Donoso-Barros, 1970;Hoogmoed, 1973;Presch, 1973;Cei, 1993;Avila-Pires, 1995;Manzani & Abe, 1997;Manzani & Abe, 2002;Péres, Jr., 2003;Péres, Jr., & Colli, 2004;Harvey, Ugueto, & Gutberlet, Jr., 2012;Pyron, Burbrink, & Wiens, 2013). The common tegu, T. teguixin, is the only Tupinambis species recorded in Venezuela, distributed on lowlands (0 -650 masl) all over the territory (Beebe, 1945;Staton & Dixon, 1977;Rivas & Oliveros, 1997;Péfaur & Rivero, 2000;Markezich, 2002;Molina, Señaris, & Rivas, 2004). Since 1977, Tupinambis are listed in Appendix II of CITES because of widespread trade in its skin, especially in Argentina and Paraguay. In Venezuela, although tegus are not commercially exploited, difficult economic conditions could spark an interest in the near future, so it is important to have information on patterns of genetic variation of its populations, among other biological data, in order to make recommendations for its sustainable management. In this study, we assessed the genetic diversity and population structure of Tupinambis teguixin in Venezuela, with the addition of one population from Ecuador and one from Brazil, by analyzing mitochondrial DNA sequences, and considering the Orinoco River (OR) and the Mérida Andes (MA) as potential geographic barriers. The study of population structure was aimed at determining whether the populations in this study are in fact a single population or a group of somewhat isolated subpopulations with little or no gene flow between them, due to the geographic barriers that separate them.

MATERIAL AND METHODS
Samples and laboratory procedures: Six localities were visited for two weeks each between 1999 and 2002, in search for tegus: Zulia (Campo Rosario, Zulia State), West of the Mérida Andes; Delta (Boca de Macareo, Delta Amacuro State), Guatopo (Guatopo National Park, Miranda State) and Llanos (Agropecuaria Puerto Miranda, Guárico State), the latter three East of the Andes and North of the Orinoco River; Amazonas (San Carlos de Río Negro, Amazonas State) and Guri (Campamento Mirador del Frío, Bolívar State), these two East of the Andes and South of the Orinoco (Fig. 1). Lizards were captured in Tomahawk live traps. Five or six individuals were sacrificed at each locality to keep as vouchers at the Museum of Natural Sciences of Universidad Simón Bolívar. Muscular tissue samples were obtained inducing autotomy by bending the tip of the tails of the lizards (García-Muñoz, Ceacero, Pedrajas, Kaliontzopoulou, & Carretero, 2011), preserved immediately in absolute alcohol, and a few days later stored at -20 o C.
DNA was extracted using Promega Wizard Genomic DNA Purification kit (A1120), following the mouse tail protocol from the kit's manual. We estimated the concentration of DNA in a spectrophotometer, and samples were then standardized to an approximate 2 ng/ µL, suitable for PCR amplification. Using Promega PCR Master Mix (M7502), we amplified partial segments of the mitochondrial genes: cytochrome b (cytb1 290bp and cytb2 358bp) and subunit 4 of NADH dehydrogenase (ND4 375pb). PCR cycling protocols and primers used can be found in supplementary material (Appendix 1). PCR products were purified with Promega Wizard SV Gel and PCR Clean-Up System (A9282), and sequenced either in the Center for Sequencing and Analysis of Nucleic Acids (CESAAN, Instituto Venezolano de Investigaciones Científicas, IVIC, Venezuela), in the Genetic Ecology Unit (UEG, IVIC, Venezuela) or at MACROGEN Inc. (Seoul, Korea). Samples with ambiguous or unclear results, as well as a number of other nonproblematic samples, were sequenced in both directions and at least twice to corroborate the resulting sequence.

Analysis of geographic variation:
To analyze molecular data, we first edited the cytb1, cytb2 and ND4 sequences using the program Chromas Lite version 2.01 (Technelysium Pty Ltd) and by eye. We then proceeded to align the sequences using the ClustalW option in MEGA 5 (Tamura et al., 2011). We can be certain that the sequences are not pseudogenes since no STOP codons, nor insertions, or deletions were found. We organized the sequences into two molecular datasets: 1) the Venezuela dataset (VDS) which included the cytb1, cytb2 and ND4 sequences concatenated from 93 tegu lizards; and the South America dataset (SADS) with the cytb1 and ND4 sequences from 96 tegu lizards from Venezuela, three from Brazil (Fazenda Nova Esperan ça, Roraima) and three from Ecuador (Reserva Faunística Cuyabeno, Sucumbios). The sequences from Brazil and Ecuador were from Fitzgerald et al. (1999) stored in the GenBank database (cytb1 accession numbers AF151183-88; and ND4 accession numbers AF151200-02 and AF151211-13), and sequences generated by this project were also deposited in the Gen-Bank (combined regions cytb1-cytb2 accession numbers KM259634-729; and ND4 accession numbers KM259730-825).
Genetic polymorphism of T. teguixin populations from Venezuela was estimated through three measures of genetic diversity: nucleotide diversity (P), haplotype diversity (h) and number of polymorphic sites (S). These indices were calculated using the program DnaSP version 5.10.1 (Librado & Rozas, 2009) with data from the VDS.
We estimated the pattern of gene flow between populations with the statistic Φ ST , an index analogous to Wright's F ST , using an Analysis of Molecular Variance (AMOVA, Excoffier, Smouse, & Quattro, 1992) with sequences from VDS and the null hypothesis of no genetic difference between populations. We also compared populations from either side of the geographic barriers using AMOVAs: 1) East v. West of the MA; 2) North v. South of the OR; and 3) West of MA v. North of the OR v. South of the OR. The program Arlequin version 3.5.1.2 (Excoffier & Lischer, 2010) was used to conduct the population structure analysis. Also using Arlequin, we performed the Mantel tests to evaluate the possible correlation between genetic distance (F ST estimates) and geographic distance between T. teguixin populations (using VDS). We estimated geographic distances with the website http://jan.ucc.nau.edu/~cvm/latlongdist.html (Michels, 1997).
Finally, the haplotype networks were constructed with the program TCS version 1.21 (Clement, Posada, & Crandall, 2000), which infers genealogic relations between haplotypes using parsimony (Templeton, Crandall, & Sing, 1992); and with the program Network 4.162 (Bandelt, Forster, & Röhl, 1999; fluxus-engineering.com) which uses a distance algorithm. Since there are no cytb2 sequences of T. teguixin from Brazil and Ecuador, we generated a few networks with the SADS (cytb1 + ND4) and others with the VDS (cytb1 + cytb2 + ND4). Thus, we were able to examine whether the distribution of haplotypes changed when the mitochondrial markers changed.

RESULTS
Geographic distribution of genetic diversity and population structure: We captured 102 T. teguixin individuals in six locations in Venezuela ( Fig. 1 and Appendix 2), Values of haplotypic diversity (h) in all locations were generally high, with a total h of 0.96 (Table 1). Conversely, nucleotidic diversity was low (0.014), with the Amazonas and Delta populations displaying the highest values ( Table 1). As to polymorphic sites (S), the Llanos and Delta populations had the most (21 and 19, respectively) while Guatopo and Guri the least (two and four, respectively) as well as the lowest numbers of haplotypes (two and five) and lowest h (0.33 and 0.73) (Appendix 1).
Global Φ ST (0.834, p < 0.0001) indicates no gene flow between populations, i.e., there is population structure. Similarly, the AMOVA shows that 83.2 % of the genetic variance observed is due to differences among populations. Paired F ST estimates confirm T. teguixin's population structure, with the values with Zulia being the greatest when paired with all other populations (> 0.872, Table 2).
With the AMOVAs, we found an effect of the geographic barriers (MA and OR) only when considering the three groups, i.e., both barriers simultaneously. The three groups were East of MA (Zulia), North of OR (Delta, Guatopo and Llanos), and South of OR (Amazonas and Guri). There was a significant difference between those groups of populations (p < 0.05, Table 3), indicated by 71.2 % of the molecular variance. However, we found a different result when considering each barrier separately. When we compared the groups from the North and South of the OR we found no difference (p = 0.22, Table 3) and the molecular variance was low (26.2 %). When we compared the groups East (Amazonas, Delta, Guatopo, Guri and Llanos) and West (Zulia) of the MA, despite the fact that the percent variation was high (75.0 %), groups were not significantly different (p = 0.18, Table 3).
Because of these contradictory results in the comparison of East v. West of the MA, we proceeded to perform further three AMO-VAs eliminating one individual or population, which appeared somehow anomalous or extreme, to assess their effect on the overall result. This was the case of haplotype 33 (individual RN102) from the Amazonas population, which differed from the other Amazonas individuals (haplotypes 31 and 32) in many mutations (13 and 14, respectively). To this end, we performed the analysis (a) eliminating individual RN102 and (b) eliminating the whole Amazonas population. Again, we found no differences between East and West of the MA (without RN102, p = 0.16; and without Amazonas, p = 0.18, Table 3). We would like to stress that we repeated the PCR and the sequencing of all three markers for individual RN102 and we obtained identical sequences every time. Additionally, since the Delta haplotypes appeared widely dispersed across the haplotype network (see the results below), we performed an AMOVA without this population, and found no significant differences between the East and West groups of populations (p = 0.19, Table 3). It is possible that the differences we found among Venezuelan populations of T. teguixin (according to the global Φ ST and paired F ST estimates) were due to geographic distances rather than the physical barriers, so we proceeded to perform a Mantel test to evaluate this hypothesis, but we found no correlation between genetic and geographic distances among the tegu populations of Venezuela (r = 0.282, p = 0.209; Fig. 2). In the plot of figure 2, we can see an outlier point (the lowest point in the graph corresponding to the Delta-Amazonas distances) which may be affecting the correlation. We carried out the analysis again eliminating individual RN102 from the comparison between Amazonas and Delta but the result was also non-significant (r = 0.356, p = 0.154). When we eliminated the whole Amazonas population, the result was still not significant (r = 0.629, p = 0.066), nor was it significant without the Delta population (r = 0.408, p = 0.173). Table 2 shows that the pairwise F ST estimates from Zulia were higher than the rest and very similar to each other (between 0.872 and 0.940) while the geographic distances varied between 578 km and 1 161 km. This suggested that the Mérida Andes constituted an effective barrier for T. teguixin population migration, and making geographic distances less relevant. A similar situation occurred with the Orinoco River as a barrier, since some of the populations separated by the river showed high values of paired F ST (Guri v. Guatopo and Guri v. Llanos, Table 2).

Haplotype network:
We found 38 haplotypes among the 93 concatenated sequences (cytb1-cytb2-ND4) from the VDS, of which 47.4 % (18) were unique (present in a single individual; Fig. 3; Appendix 2 and Appendix 3). In the SADS, we found 37 haplotypes among 102 concatenated sequences (cytb1-ND4), of which 43.2 % (16) were unique (Fig.  3, Appendix 2 and Appendix 3). Only one haplotype was shared among all populations of both datasets (haplotype 7), i.e., all other haplotypes were private. The shared haplotype was the most common (Fig. 3, Appendix 3) and it was found in one Delta individual, 8/9 individuals from Guri (VDS/SADS) and 2 from Brazil. Haplotype 24 from SADS has the same frequency but is found only in the Zulia sample (Fig. 3, Appendix 3). The populations with the greatest numbers of total and unique haplotypes were Llanos and Zulia (Fig. 3, Appendix 2), while Amazonas and Guatopo had the lowest numbers (three and two respectively), possibly because of their small sample sizes (five and six, respectively).
The haplotype networks generated with TCS and Network and for both datasets (VDS and SADS) were similar. In general, we observe a well-defined structure since the haplotypes of each population are grouped (Fig. 3). Still, the haplotypes of the Delta population are related to those of other populations (Amazonas, Guatopo, Guri and Llanos), so Delta is the only population with haplotypes scattered all over the network. As we pointed out above, haplotype 33 of Amazonas is quite different from those of other Amazonian haplotypes and is located 8/12 (VDS/SADS) mutational steps from the ancestral haplotype, which is its closest haplotype in the network (Fig. 3).
The haplotype network generated by TCS with 95 % confidence showed the Zulia haplotypes on an independent network, and to connect with the main network (where most haplotypes lie), 19/24 (SADS/VDS) mutational steps were needed. Similarly, the Ecuador haplotypes appeared on a separate network, but they joined the main one in 44 steps. When we added steps in TCS to force the connection between the main network and Zulia and Ecuador, two loops were formed. The first one  Fig. 1). Each circle represents a haplotype and its size is proportional to the haplotype´s frequency while the fill pattern identifies each population. Each line connecting haplotypes corresponds to one mutational step, unless indicated by numbers. Haplotype 7 (H7) is the ancestral and the only one shared between populations (2 specimens from Brazil, 1 specimen from Delta and 9 specimens from Guri). Haplotype codes refer to Appendix 3. connected one haplotype from the Delta population (individual DE005) with the ancestral (haplotype 7), and the second was formed from two Zulia haplotypes (ZUL72 and ZUL79), which connected to the first loop. The main network connected to the Ecuador haplotypes through the second (Zulia) loop. In contrast, the analysis with the program Network resulted in a single loop between the two Zulia haplotypes, which in turn connected to the Ecuador haplotypes (Fig. 3). A list of the haplotype sequences from VDS can be found in Appendixes.

DISCUSSION
In this study, we found evidence of genetic structure in populations of T. teguixin in Venezuela and two populations of Brazil and Ecuador. The evidence includes: (1) high global Φ ST and paired F ST estimates; (2) haplotype networks with a well-defined geographic pattern, with each population with its haplotypes grouped (except Delta); Zulia and Ecuador with separate networks; and (3) a single haplotype shared among all populations. This structure does not appear to be a product of geographic distance between populations, but rather a historical biogeographic effect of the Mérida Andes and the Orinoco River. The F ST values suggest that the distribution pattern of the genetic diversity is consistent with the presence of geographic barriers, which therefore have restricted the gene flow between populations.
The geological history of the barriers shows that the main elevation of the Mérida Andes (MA) happened toward the end of the Miocene, between 12 and 8 My ago (figure S5 in Hoorn et al. 2010). Additionally, there are several competing hypotheses on the location of the Proto-Orinoco original course and how it changed to its current location. Hoorn, Guerrero, Sarmiento & Lorente (1995) proposed that in the early Miocene a Proto-Orinoco ran into the Caribbean Sea in what is now Maracaibo lake basin. Díaz de Gamero (1996) hypothesized that the Proto-Orinoco ran into the sea in what is now Falcon state in Northwestern Venezuela, at the same time; Rincón, Solórzano, Benammi, Vignaud & McDonald (2014) suggested that by that time it already drained in central Venezuela (Guarico state). In general, researchers agree that the present Orinoco River course with a delta draining into the Atlantic Ocean exists since the late Miocene (5-6 My; Hoorn et al., 1995;Díaz de Gamero, 1996;Iturralde-Vinent & MacPhee, 1999;Rincón et al., 2014). Tupinambis probably originated in the Oligocene (32 to 33 My) according to Péres Jr. (2003), based on a molecular clock analysis of cytochrome b and by Giugliano, García & Colli (2007) based on 12S and 16S mitochondrial genes. According to the fossil record, the earliest Tupinambislike fossil was found in La Victoria formation (Honda Group of the Magdalena Valley in Colombia) from the late Miocene (13.5 to 12.9 My) (Sullivan & Estes, 1995). This scenario, based on molecular and fossil evidence then suggests that tegu lizards of the genus Tupinambis inhabited Northern South America before the rise of the Andes and before the present location of the Orinoco River (OR), so populations may have in fact be left isolated on either sides of these barriers.
Few studies have looked at the genetic structure of vertebrate populations or species in relation to the possible effect of the MA or -even less -the OR as geographic barriers. Stephen et al. (2013) investigated the genetic structure of the common iguana (Iguana iguana), in relation to several geographic barriers, including the MA and found profound divergence consistent with the barriers. Further, they suggest the presence of cryptic species, although further sampling would be required to confirm this possibility. Also, Miralles et al. (2009) studying the phylogeny of the Caribbean skinks (Mabuya), a genus with a wide altitudinal distribution, found a new species (Mabuya zuliae) on the lowlands of the Maracaibo lake basin, West of the MA, in the general area where we sampled tegu lizards. On the other hand, Vargas-Ramírez et al. (2013) found no molecular differences on populations of the turtle Rhinoclemmys on either side of the Andes but the taxonomy remains unchanged because of the allopatric and parapatric distribution of the genus and the species Rhinoclemmys melanosterna, as well as the caryological and morphological differences. Albert et al. (2006) suggest that the tectonic activity during the Miocene fragmented the ichtyofauna of Northwesten South America. Similarly, the diversification of the Cracidae (Aves) appears to be due to vicariant processes caused by marine transgressions, the rise of the Andes and concomitant changes in the courses and basins of South American rivers (Pereira & Baker, 2004).
AMOVAs in this study have produced contradictory results regarding the effect of geographic barriers on the genetic structure of T. teguixin. On the one hand, they indicate that both barriers (MA and OR) do influence the geographic distribution of the genetic diversity of the species, while on the other hand, when each barrier is considered separately, this is not the case. However, paired F ST values and haplotype networks confirm that the MA has an effect on the populations of these lizards, and that this effect is greater than that of the OR. Specifically the paired F ST estimates of Zulia (West of the MA) are the highest, and the haplotypes of this population form a network that is separate from the rest (main network). Conversely, the OR seems to allow a certain amount of gene flow, since the only shared haplotype is found in the Delta population (North of the OR), in Guri and Brazil (South of the OR). This is expected since the Guri and Brazil populations are on the Guiana Shield, which was never covered by the sea, while other regions of Venezuela where fragmented first by sea incursions and then by geographic barriers (Hoorn et al., 2010). Further, the stretch of the Orinoco where the tegus have a better chance of crossing is precisely the river delta, since in that area the current is milder and there are many islands to facilitate crossing. Moreover, it is known that tegu lizards are good swimmers (Achaval & Langguth, 1973).
The AMOVA comparing populations East and West of the MA (two groups), which showed a high percentage of molecular variance but a non-significant Φ ST , is possibly the reflection of a technical problem. When two groups of populations are compared, there is only one degree of freedom, which reduces the statistical power of the analysis, i.e., there is a lower probability of detecting a difference between the groups. The fact that the differences did turn out to be significant when comparing the three groups (taking both barriers into account) may be due to the greater power of the analysis, plus the addition of an additional factor (the other barrier) which may reduce the error.
Haplotypes in the population from Ecuador are quite different from the rest, with a distance of at least 44 mutational steps. This population is located 1 002 km from the Amazon population (its nearest one), which is a distance comparable to that between the Amazonas and Delta populations (952 km) and between Amazonas and Zulia (1 014 km), but the genetic distance between Ecuador and Amazonas is much greater than the distance between the other populations in this study: for instance, Zulia, with the Mérida Andes barrier separating it from the rest, is only 19 mutational steps from the main network (the network with the rest of the Venezuela and the Brazil haplotypes). Something similar was found by Fitzgerald et al. (1999) who reported a genetic divergence between T. teguixin from Brazil and Ecuador close to that between two species, namely, T. merianae and T. rufescens. It would be important to extend the sampling to populations located towards the center and South of Brazil to assess the genetic variation of T. teguixin in its full geographic range in order to understand such great divergence between Ecuador and the rest. It is possible that the so-called Pebas system (a great wetland composed of shallow lakes and swamps which existed between 23-10 My) may have affected (isolated) the Ecuadorian populations (Hoorn et al., 2010).
The very high pairwise F ST estimates of Zulia irrespective of geographic distance constitute evidence of the effect of the Merida Andes on the genetic structure of the populations of T. teguixin rather than geographic distance. However, the Mantel test could be underestimating the degree of gene flow between populations, since an apparent but not statistically significant trend can be seen in the graph of genetic distance v. geographic distance (Fig. 2), while the only shared haplotype is found in the populations nearest to each other (Appendix 2). Thus, it would be advisable to carry out further sampling in locations located halfway among those used in this study in order to corroborate (or not) the lack of gene flown between the populations. Likewise, it is recommended to include analyses with nuclear markers, to verify that the pattern of genetic variation found with mitochondrial markers is similar to the pattern of nuclear markers.
Finally, the results of this study show that there is a significant population genetic structure in Tupinambis teguixin from Venezuela and that the geographic distribution of the genetic diversity is modulated by the presence of an important geographic barrier, namely the Mérida Andes and, to a lesser degree, the Orinoco River. Thus, both the Paleogeographic and the Riverine hypothesis gain support from the genetic structure we found.
We would like to conclude with brief comments on the conservation implications of our study. Because of the widespread and excessive exploitation of Tupinambis populations, mainly in Argentina and Paraguay, CITES decided to include the genus in the Appendix II. In other countries, exploitation is moderate to non-existent, as in Venezuela, so those populations do not appear endangered. However, economic difficulties, which tend to recur in countries such as Venezuela, may lead to an interest in the species as a source of income. The results of this study can contribute to the decision-making process in the case of the implementation of a management plan for T. teguixin in the wild, which would help preserve the genetic diversity and viability of its populations.
The genetic structure found, where the Zulia population is clearly differentiated and isolated from the rest, in addition to being morphologically different as well (Gols-Ripoll et al. unpublished), leads us to propose that the Zulia population be treated as an Evolutionary Significant Unit (ESU, sensu Moritz, 1994). Further studies are required to confirm that the other populations in this study should also be considered ESUs, so we provisionally suggest that these should be considered separate Management Units (MU, sensu Moritz, 1994), until the genetic-morphological taxonomic status would be assessed. The Delta and Guri populations would be a single unit since they share a haplotype. Thus, for now, only Delta and Guri populations should be allowed to intermix, -in case this was proposed as a management strategy-, pending further data.