Genetic divergence among and within four Italian populations of Aphanius fasciatus (Teleostei, Cyprinodontiformes)

Abstract

Due to the natural fragmentation of its habitats, Aphanius fasciatus, a Cyprininodontid living in coastal brackish‐water habitats, is distributed discontinuously with considerable isolation of its populations, and characterized by a high degree of osteological and genetic differentiation among them. The genetic relationships in four different Italian populations were studied using the control region (D‐loop) mitochondrial DNA sequence and the polymorphisms generated by RAPD‐PCR. The results revealed a well‐defined genetic structuring, that was also confirmed by the high G_st and F_st values, and a close relationship between the populations of Ganzirri and Lesina, as well as the isolation of the other two Sicilian populations. Moreover, with respect to the mitochondrial DNA analyses, the total genetic variance was higher among populations, but lower in the case of the RAPD polymorphisms.

Keywords:

• Killifish
• mtDNA
• RAPD‐PCR
• genetic population structure

Introduction

The genetic structure of a species reflects the geographic distribution pattern of the genotypes and the modality of the gene flow. The levels of gene flow depend not only on the dispersal capacity of individuals and their gametes in relation to geographic, physicochemical and ecological barriers, but also on their fitness, or in other words, on the reproductive capacity that they manifest (Maltagliati 1998). An extensive gene flow can curtail evolution independently of the populations, thereby hindering genetic drift, while a low gene flow generally determines a high level of genetic differentiation between populations (Waples 1987).

The killifish Aphanius fasciatus Nardo 1827 currently lives in coastal brackish‐water habitats; its populations are often isolated and show considerable morphological and genetic differentiation, as pointed out in osteological (Tigano & Ferrito 1985; Parenti & Tigano 1993; Tigano et al. 1999), biochemical (Maltagliati 1998, 1999) and karyological studies (Vitturi et al. 1995; Tigano et al. 2004b). Allozymic and morphological investigations have always shown a significant differentiation in populations, indicating in some cases a genetic divergence in relation to their geographic distribution (Maltagliati 1998, 1999). In other cases (Cimmaruta et al. 2003), this divergence does not seem related to geographic distance, as is shown by the morphological data (Tigano et al. 2001; Ferrito et al. 2003). Recently, analysis of the highly variable D‐loop tract of the mitochondrial DNA (Tigano et al. 2004a; Ferrito et al. 2007) has indicated substantial genetic divergence among populations of this species. Mitochondrial control regions have been regarded as a molecular marker of choice for studying rapidly evolving populations, because they accumulate mutations up to ten times faster than the rest of the mitochondrial genome (Sturmbauer & Meyer 1992; Meyer 1993). Many studies reported the effectiveness of RAPD markers in discriminating between species, subspecies, and populations in a wide range of organisms, including fish (Bardakci & Skibinski 1994; Allegrucci et al. 1995; Gomes et al. 1998; Mamuris et al. 1998). Recently, Maltagliati et al. (2003) detected high levels of genetic divergence between two adjacent populations of A. fasciatus.

The aim of the present study is to use analyses of the mitochondrial DNA sequence and of the polymorphisms generated by RAPD‐PCR to evaluate the level of genetic divergence present in four different Italian populations of A. fasciatus; the work is part of the ongoing studies to analyse the structuring of this species using various investigative approaches: osteology, allozymes, and DNA (Ferrito et al. 2000, 2003, 2007; Tigano et al. 2001, 2006). The geographic distribution of these populations is related to their natural history in order to understand the geographic patterns that may result from divergence.

Materials and methods

Fish samples

A total of 84 specimens of Aphanius fasciatus was examined. They were collected from Ganzirri Lake, the Marina di Modica salt pond, Contrada Manzonara in Sicily and Lesina Lagoon in Apulia. (Figure 1).

Figure 1 Map showing catch localities of individuals used in this study.(LE) Lesina; (GA) Ganzirri; (MA) Manzonara; (MM) Marina di Modica; (LO*) Longarini; (PV*) Pantano Viruca; (FM*) Foce Marcellino. The asterisks indicate the populations from southeastern Sicily (Tigano et al. 2004a).

Tissue samples and DNA extraction

The number of individuals sampled and sequenced per locality are detailed in Table I and Figure 1. The fish were preserved whole in 95% ethanol, and total genomic DNA was then extracted from the muscle tissue (100 mg), which was digested overnight at 55°C in 500 µl of extraction buffer (NaCl 400 mM, Tris 10 mM, EDTA 2 mM, SDS 1%, pH 8, proteinase K 20 µg/µl). The DNA was purified by standard chloroform extraction and isopropanol precipitation (Sambrook et al. 1989).

Table I. Populations and relative number of individuals of Aphanius fasciatus examined. Number and diversity of haplotypes were calculated using DNAsp (Rozas & Rozas 1999).

Populations	Number of specimens	Number of haplotypes	Haplotype diversity
Ganzirri lake (GA)	15	7	0.838
Lesina lagoon (LE)	18	10	0.922
Contrada Manzonara (MA)	10	9	0.978
Marina di Modica salt pond (MM)	41	19	0.952

Mitochondrial control region PCR amplification and phylogenetic analyses

Amplification of a 5′ portion of the mitochondrial control region was accomplished using the CR‐A and CR‐E primers described by Lee et al. (1995). All amplifications (25 µl) contained 10 to 100 ng of DNA, 10 mM Tris HCl (pH 8.3), 50 mM KCl, 1.5 mM MgCl₂, 2.5 units of Taq DNA Polymerase (GE Healthcare, UK), 150 mM of each dNTP and 50 pmol of each primer. A cycling profile of 45 s at 94°C, 45 s at 48°C, and 1 min at 72°C for 35 cycles in an Eppendorf Mastercycler Personal was used. The amplification products were then purified by means of a Nucleospin Extract II kit (Macherey‐Nagel, Germany). An ABI 3100 automated sequencer (Applied Biosystems, U.S.A.) was employed to perform automated sequencing with the CR‐A primer used in the amplification.

The computer program Clustal V implemented by Sequence Navigator (Applied Biosystems) was used to align the mitochondrial sequences. Phylogenetic relationships were assessed using the Neighbour‐Joining (NJ) (Kimura 2 parameter) method and the Maximum Parsimony (MP) method implemented by the PAUP software package (Phylogenetic Analyses Using Parsimony, version 4.0, Swofford 1998). Relationships were also assessed by Bayesian analysis with Markov Chain Monte Carlo (MCMC) sampling in MrBayes v.3.0b4 (Huelsenbeck & Ronquist 2001). Bayesian MCMC analyses with four parallel Markov chains were run for 1,000,000 generations with a sampling frequency of one in every million trees. Due to the large number of samples involved, the MP method on a subset of our data, which included all the individuals with differing haplotypes, was used. Topological confidence was evaluated with 1000 bootstrap replicates (Felsenstein 1985). In both the NJ and MP procedures, bootstrapping analysis was performed with an equal weighting of transitions and transversions. Two Aphanius species, A. dispar and A. mento (Genbank accession numbers U06052 and U06054, respectively) and a single specimen of A. fasciatus coming from Genbank (U06053) were used as outgroups. MtDNA sequences from 29 individuals from southeastern Sicily (Tigano et al. 2004a) were also considered in the phylogenetic reconstruction (accession numbers from AJ605322 to AJ605329). Population structure was calculated by analysing the molecular variance (AMOVA, Excoffier et al. 1992) based on pairwise differences as implemented by Arlequin (vers. 1.1, Schneider et al. 1997). Gene flow (F_st) and haplotype diversity were calculated using the software package DnaSP, version 3 (Rozas & Rozas 1999), following Hudson et al. (1992).

A network was constructed by means of the TCS program (Clement et al. 2000), which uses statistical parsimony (Templeton et al. 1992; Crandall 1994) for network estimation. The maximum number of mutational steps that constitute a parsimonious connection between two sequence types was calculated with 95% confidence. Gaps were treated as a fifth base.

RAPD‐PCR analysis

PCR reactions (25 µl) were performed for the same specimens, using 10 ng of DNA, 10 mM Tris HCl (pH 9), 50 mM KCl, 1.5 mM MgCl₂, 1.5 units of Taq DNA Polymerase (Amersham Pharmacia Biotech), 200 mM of each dNTP, and 50 pmol of each primer. We used a cycling profile of 1 min at 94°C, 1 min at 35°C and 2 min at 72°C for 45 cycles. To ensure positive and reproducible results, we selected three primers that reveal polymorphisms: OPA‐01 (CAGGCCCTTC), OPA‐09 (GGGTAACGCC) and OPA‐10 (GTGATCGCAG) (kit A Operon Technologies, Germany). Amplified products were loaded onto a 12.5% polyacrylamide gel (GeneGel Excel 12.5/24 kit, GE Healthcare, UK). Electrophoresis was carried out at 10°C, with maximum settings of 600 V, 25 mA, 15 W for 2–2.5 h on GenePhor System. The gels were silver stained using Hoefer™ Automated Gel Stainer and the Plus One DNA Silver Staining Kit according to the supplier's instructions (GE Healthcare, UK). Negative reactions were also performed to test the effectiveness of the method.

The phenotype of a distinct locus was considered as the presence (1) or absence (0) of each band. POPGENE version 1.31 (Yeh 1999) was used to calculate Shannon's diversity index and the G_st values, so as to provide a relative estimate of the degree of genetic variation between population. In this case, we also used the Arlequin package to perform the AMOVA test.

Results

Amplifications of the mitochondrial control region from 84 individuals resulted in 414 aligned base pairs (accession numbers from EF440653 to EF440688). A single base insertion (position 39) was observed in the Ganzirri and Marina di Modica populations. No other insertions or deletions were observed in our data sets. Of these 414 bases, 113 were variable, and 112 were phylogenetically informative. After PAUP analysis, transitions (trs) proved more frequent than transversions (trv) (ratio trs/trv = 2.75) and saturation was not present (plots not shown). Thus we did not weight transitions and transversions differently. The number of haplotypes and their diversity are given in Table I.

The phylogenetic trees constructed with the NJ method (Kimura two‐parameter model) and the MP method were superimposable and show a close relationship between the populations of Ganzirri, Lesina and Foce Marcellino, rather than between that of Ganzirri and the other Sicilian populations. On the contrary, these appear to be very different from one another; only Longarini and Pantano Viruca cluster together. All phylogenetic relationships were supported by high bootstrap values (> 72%) in both the approaches used (Figure 2a). The UPGMA dendrogram based on the genetic distance matrix for the RAPD markers confirms the affinity between the two populations of Ganzirri and Lesina, as well as the isolation of the other two Sicilian populations (Figure 2b). The three southeastern Sicilian populations were not analysed using RAPD markers. The same results were obtained with the Bayesian approach (figure not shown).

Figure 2 (a) Neighbour‐Joining tree drawn using D‐loop sequence data. The numbers represent the percentage of 1000 bootstrap replication (NJ/MP) in which a given node appeared. (b) UPGMA dendrogram based on the genetic distance matrix for the RAPD‐PCR markers for the populations examined.

Genetic divergence between populations is shown in Table II. It ranges from 1.6 (between the Ganzirri and Lesina populations) to 3.6% (between the Pantano Viruca and Lesina populations). These data are in agreement with other genetic divergences reported for the mitochondrial control region in fishes (McCune & Lovejoy 1998). The pronounced genetic structuring is also confirmed by the high F_st and G_st (Table III) values obtained from analyses with both molecular markers; the parsimony network shows no connection among the haplotypes of the populations, but PV and LO belong to the same group (Figure 3).

Table II. Average pairwise genetic (Kimura‐2) distances between populations. The values are referred to mtDNA (under the diagonal) and RAPD (above the diagonal) analyses.

	GA	LE	MM	MA	LO	PV	FM
GA		6.4	11.7	16.5
LE	1.6		13.1	17.2
MM	2.0	2.3		17.0
MA	2.9	3.1	1.6
LO	2.7	3.5	1.3	1.9
PV	2.7	3.6	1.4	1.8	0.5
FM	1.5	2.4	1.6	2.3	2.7	2.6

Table III. F_st values for D‐loop mtDNA (under the diagonal) and G_st values for RAPD markers (above the diagonal) calculated using DNAsp (Rozas & Rozas 1999) and POPGENE (Yeh 1999), respectively.

	GA	LE	MM	MA	LO	PV	FM
GA		0.055	0.087	0.119
LE	0.73		0.098	0.127
MM	0.91	0.87		0.109
MA	0.83	0.79	0.77
LO	0.91	0.87	0.70	0.80
PV	0.97	0.93	0.85	0.92	0.40
FM	0.95	0.89	0.88	0.93	0.90	0.98

Figure 3 Parsimony network constructed using the TCS program(Clement et al. 2000).

The three random primers used in RAPD PCR analyses generated a total of ten unambiguous polymorphic loci, ranging in size from 500 to 2000 bp. The mean genetic diversity within populations estimated by Shannon's indices ranged from 0.4542 to 0.5357 (mean value 0.4846, Table IV). The most variable population was Marina di Modica.

Table IV. Genetic diversity estimates within populations of A. fasciatus using Shannon's diversity index (H₀). The RAPD bands are given as approximate size (bp) of the amplified PCR product.

	OPA‐01	OPA‐01	OPA‐01	OPA‐09	OPA‐09	OPA‐09	OPA‐10	OPA‐10	OPA‐10	OPA‐10	H0
	500	700	750	600	800	900	1030	1200	1500	2000
GA	0.5337	0.6245	0.6245	0.4115	0.5830	0.5712	0.6926	0.6563	0.0000	0.0000	0.4697
LE	0.5830	0.5337	0.3372	0.2512	0.2512	0.6811	0.4767	0.5337	0.5712	0.5712	0.4790
MM	0.6811	0.6623	0.6623	0.4230	0.6365	0.3631	0.6915	0.5461	0.6915	0.0000	0.5357
MA	0.0000	0.6576	0.3372	0.6576	0.6886	0.4452	0.5337	0.6886	0.5337	0.0000	0.4542
Mean values	0.4494	0.6195	0.4903	0.4358	0.5398	0.5152	0.5986	0.6062	0.4491	0.1428	0.4846

The analyses of genetic variance were performed with the populations separated into two distinct groups: those from Sicily, on the one hand (Manzonara, Ganzirri and Marina di Modica), and that from the Lesina peninsula on the other. During the AMOVA test, it was discovered that total genetic variance was higher among populations (87.58% considering the populations as belonging to a single group) in the case of the mitochondrial DNA analyses, but lower in that of the RAPD polymorphisms (39.79%) (Tables V and VI). Moreover, no notable variation relatable to the geographic distribution of these populations was observed in the rates of genetic diversity.

Table V. Results from mtDNA analysis of molecular variance (AMOVA) for the populations studied (P<0.001).

Source of variation	Degrees of freedom	Variance component	% of variance	Φ‐statistics
Among groups	6	14.38925	87.58	0.87
Within groups	106	1.36591	12.42

Table VI. Results from RAPD analysis of molecular variance (AMOVA) for the populations studied (P<0.001).

Source of variation	Degrees of freedom	Variance component	% of variance	Φ‐statistics
Among groups	3	1.24510	39.79	0.39
Within groups	52	5.42528	60.21

Discussion

The paleogeographic changes that took place in the Mediterranean basin starting in the Miocene, (Hsü et al. 1977; McKenzie 1999) had a profound impact on the distribution and genetic structure of the cyprinodontiform fish that populated the Mediterranean waters in this period, and whose European fossils date back to the Oligocene‐Miocene (Gaudant 1978, 1979). Studies by Villwock (1999) on cyprinodontid biogeography showed that, during the Middle Miocene marine transgression, the ancient ancestor of the A. fasciatus/A. dispar group was distributed throughout the eastern Mediterranean, including the Mesopotamian basin. As the Tethys Sea began to recede, this ancestor became divided into two groups, giving rise to the species of the A. fasciatus complex in the central and eastern Mediterranean basin, and to that of the A. dispar one in the southeastern part of the basin. The Messinian Salinity Crisis (Upper Miocene) most probably triggered a process of vicariant speciation in the Aphanius genus (Hrbek & Meyer 2003). The subsequent paleoclimatic changes that took place in Quaternary modified the species range considerably (Hewitt 1999, 2000) and certainly affected, in particular, the present genetic structure of the populations of the different Aphanius and A. fasciatus species.

The results obtained from analysing the sequences of mitochondrial DNA control regions and those derived from studying the RAPD polymorphisms are completely superimposable. This, in fact, infers a closer relationship between the populations of Ganzirri and Lesina, rather than between that of Ganzirri and the other Sicilian populations. The haplotype diversity found in the present study (0.838–0.978) is rather high, which, according to Grant and Bowen (1998), is indicative either of a long and stable evolutionary history or of secondary contact among differentiated lineages.

The considerable molecular divergence between the Sicilian populations examined in the present study had already appeared in a previous study, conducted on other Sicilian populations from different areas (Tigano et al. 2004b), that analysed the D‐loop of mitochondrial DNA and various osteological characteristics. Those results revealed a morphological and molecular differentiation between these populations, with genetic distance values similar to those obtained in the present study. Research conducted on a morphological basis (Tigano et al. 2001; Ferrito et al. 2003, 2007) into the same populations as those investigated in the present study corroborates the molecular data obtained here, thus confirming the good overall agreement between morphological and molecular studies (Stepien & Kocher 1997). A low molecular differentiation was found in populations that are geographically remote from one another, such as the Maltese and Tunisian ones (the rate of divergence ranges from 0.2% to 1.3%) or those of western Sicily and Malta (0.0%–0.5%) (Tigano et al. 2006). Consequently, it is still unclear which events might serve as the basis for a model of genetic structuring for A. fasciatus. Taking the paleogeographic events that affected the Mediterranean area into account, the results we obtained can be explained in the light of the climatic and geological events that took place during the periods between the Pliocene and the Pleistocene; indeed, during the Late Pliocene, most of Italy was totally submerged. In the period that followed, with the advent of the last ice age (the coldest ever), the sea level dropped, allowing many of the submerged regions to reemerge. Very likely this caused the gene flow to cease abruptly, in particular among the populations of Aphanius living in the brackish waters of the lagoons near the coast, but, in some cases, it allowed temporary contacts along shorelines (Tigano et al. 2006). The populations of the southeastern part of Sicily, Pantano Viruca, Longarini and Marina di Modica belong to the Hyblean region that was isolated, from the rest of Sicily until Pleistocene, and the high genetic divergence rates of these populations vs the others could be due to the long period of isolation of this region; instead, the closer relationship between Ganzirri and Lesina is difficult to explain because in a recent work (Tigano et al. 2006), it was shown that A. fasciatus living in the Ganzirri Lake could be considered an introduced species.

The seven populations of A. fasciatus from the central and eastern Mediterranean that Hrbek and Meyer (2003) studied (not the same as those investigated by us) showed a genetic differentiation (approximately 1.6% sequence divergence), with one exception: the population from Lake Bafa in Turkey, which had an average divergence from the other populations of 6.86%.

The genetic structure of marine organisms could be the result of habitat connectivity patterns, ecological conditions, adult ecology and dispersal potential. The influence of environmental factors plays a decisive role in the differentiation of populations, and it has, in fact, been hypothesized that individuals subjected to similar environmental conditions, such as the presence of predators and a limited availability of food, can respond in the same way to selection pressures (Mittelbach et al. 1999). Moreover, according to Carvalho (1993) and Bernardi et al. (2001), if localized populations of fish inhabit similar environments or remain interconnected through migration and gene flow, they may display more or less homogenous arrays of phenotypic or genetic traits. However, if they are exposed to contrasting environmental conditions and/or only exchange a few migrants, appreciable population differentiation may arise because of the accumulation of unique mutations in each population due to genetic drift or natural selection.

The high levels of genetic divergence in the populations that we studied might well be related to the fragmented nature of the brackish water habitat, the low dispersal power of A. fasciatus (Maltagliati 1998) and a remarkable degree of adjustment to local environmental conditions. Indeed, particularly stressed habitats such as the brackish waters in which A. fasciatus lives can influence the rate of evolution, reducing the size of the populations and consequently causing the loss of a random number of alleles; this would result in a reduction of genetic variability and heterozygosity (a phenomenon known as genetic erosion) (Cimmaruta et al. 2003). Thus, the selection pressures that the individuals of different populations are exposed to might in part explain the considerable variability that exists even within populations, as demonstrated by the RAPD polymorphism analyses and by the higher Shannon indices, in particular.

Finally, the results obtained from the molecular research provide an effective assessment of the level of genetic differentiation among the four populations of Aphanius fasciatus investigated.