Genetic and morphometric differences between yellowtail snapper (Ocyurus chrysurus, Lutjanidae) populations of the tropical West Atlantic

Populations of Ocyurus chrysurus were compared genetically and morphometrically along the West Atlantic coast to test the null hypothesis of population homogeneity in the area. Brazilian populations were found to be differentiated in shape (canonical variates analysis; F[48,515] = 10.84, p < 0.0001). Analyses of mitochondrial DNA sequences (663 bp of the control region) did not show any differences between Brazilian populations but could detect differences between Brazilian and Caribbean (Belize) populations. The samples from Pernambuco differed significantly from the other Brazilian populations in allozyme frequencies (11 loci; FST = 0.167; p < 0.05), but this may have resulted from the small number of samples analysed for that population. Sequence variation of Belize samples departed from neutral expectations (Fu’s FS = -8.88; p < 0.001). A mismatch distribution analysis points to an ancient population expansion in that area. We conclude that the genetic data do not allow the rejection of the null hypothesis of panmixia for Brazilian yellowtail snapper populations which should be treated as a single genetic stock, with a latitudinal gradient on their morphology which probably results from phenotypic plasticity. On the other hand, there is a severe restriction to gene flow between O. chrysurus populations from the Caribbean and from the southwestern Atlantic.


Introduction
The yellowtail snapper Ocyurus chrysurus (Bloch, 1791) is a Lutjanid reef fish found on the east coast of the Americas from New England (USA) to southeastern Brazil (Menezes and Figueiredo, 1980;Allen, 1985). It is dioecious and its spawning times vary between geographic regions (Thompson and Munro, 1983;Claro, 1983;Grimes, 1987). O. chrysurus has opportunist and generalist alimentary habits and, hence, it is an important predator in some reef zones (Fallows, 1984;Parrish, 1987). They are longlived (up to 17 years Allman et al., 2005) and their larvae are pelagic (Riley et al., 1995).
Yellowtail snappers are important fisheries resource where they occur. The production of the species in Northeast Brazil is high (Costa et al., 2003) by 2001 averaging 1683 tons.year -1 (Resende et al., 2003). In Southeast Brazil yellowtail snappers accounted for 16% of the total fish biomass captured between 1986-1989(Paiva and Andrade-Tubino, 1998. The formation of spawning aggregations in this species (Claro, 1983) makes it particularly vulnerable to overfishing (Costa et al., 2003), so that it is important that the stock structure of this species be known in order to better preserve it. This is also important because of the reported post-recruitment site fidelity of this species (Watson et al., 2002) which might lead to large stock differentiation and fragility.
In recent years there has been a worldwide decrease in fishing stocks (Garcia and Grainger, 2005;Pauly et al., 2005): the overexploited, fully exploited and exhausted fisheries that were 69% in 1995 increased to 75% in 2002 and only 1% of the stocks are in a state of recuperation (FAO, 2002;Ormerod, 2003). Similar problems are faced by the Brazilian fisheries stocks (Vasconcellos and Gasalla, 2001). For this reason efficient management policies based on unambiguous scientific data are necessary both to protect the fishing stocks and to maximize their exploitation without compromising their integrity. The correct stock as-sessment and identification of the species is fundamental to establish the maximum sustainable effort of a given marine resource (Ryman and Utter, 1987;Ryman, 1991). Many stock concepts can be found in the literature (Booke, 1981;Ovenden, 1990;Smith et al., 1990;Carvalho and Hauser, 1995), but one of the most accepted and used is that "a stock is an intraspecific group of randomly mating individuals with temporal and spatial integrity" (Ihssen et al., 1981) which covers most of the definitions given by other authors.
The identification of fish stocks can be made efficiently by the use of highly polymorphic molecular markers (Blaber et al., 2005;Caddy and Seijo, 2005;Carmen and Ablan, 2006). In this study nuclear (allozymes) and mitochondrial (control region sequences) markers are used for the first time to analyse the stock structure of populations of O. chrysurus along 2380 km of Brazilian coast. A significant latitudinal gradient was observed for geometric morphometric data. The lack of observable heterogeneity in allozymes or in the highly variable mitochondrial DNA sequences indicates that this variation results from phenotypic plasticity. Additionally, samples from the Caribbean differed genetically from the Brazilian ones indicating restriction in gene flow between the two areas for this fish species.

Material and Methods
Sampling Specimens of Ocyurus chrysurus were collected between August and October 2001 by fishing vessels from four locations along the Brazilian coast (Table 1; Figure 1). Samples of muscle and liver tissues from each fish were transported frozen to the laboratory where they were stored in liquid nitrogen until analyses. Muscle samples from O. chrysurus were also collected in Belize (17°29' 34" N, 88°06' 00" W; N = 20) and transported in 90% ethanol to the laboratory where they were used in DNA sequencing.

Morphometric analyses
Morphometric variation among populations was assessed through geometric methods (Zelditch et al., 2004). For this purpose all 192 individuals collected at the four localities were measured. Measurements were taken with a Vernier manual calliper (0.05 mm of precision) between 10 homologous landmarks using a truss network protocol (Strauss and Bookstein, 1982) which provided 21 interlandmark distances (Figure 2).
To retain the relative geometry of landmarks, interlandmark distances were converted to Cartesian landmark coordinates for each specimen using the import truss feature of the Morpheus program (Slice, 1998). Landmark coordinates were then adjusted by the generalized least squares Proscrustes superimposition (GLS) method (Rohlf and Slice, 1990), centering all specimens, and scaling them to an equal centroid size. The GLS-adjusted landmark coor- Vasconcellos et al. 309   dinates were used for thin-plate spline (TPS) analysis (Bookstein, 1989) which provided shape parameters for subsequent multivariate analyses. Both uniform and nonuniform (= partial warps) components of shape were submitted to a canonical variate analysis (CVA) in order to assess shape variation among populations (Cavalcanti et al., 1999;Zelditch et al., 2004). The PAST program (Hammer et al., 2001) was used for thin-plate spline analysis and partial warps computations. Multivariate statistical analysis (CVA) was performed with SYSTAT® package, v. 10.

Allozyme electrophoresis
Muscle and liver samples were tested for thirty allozymes and three buffer systems, using 12.5% starch gel electrophoresis. Based on preliminary results, liver tissue and twelve allozyme systems were chosen for the analyses (Table 2).
Gel slices were stained following standard conditions (Manchenko, 1994;Batista and Solé-Cava, 2005;Lima et al., 2005). Samples that could not be reliably scored for at least half of the allozyme loci were excluded from the analyses. This was the case for most samples from Pernambuco, leading to a severe reduction in allozymes sample size for that locality. Allozyme genotypes were used to estimate allele frequencies, from which heterozygosity levels and inbreeding indices (Global and pairwise F ST ) were calculated, using the Genetix 4.01 program (Belkhir et al., 2002). An exploratory factorial correspondence analysis (Excoffier et al., 1992) was also used, to check for hidden partitions in the allelic hyperspace.

DNA extraction
Total genomic DNA was extracted from muscle tissues of 77 individuals, using a CTAB extraction protocol (CTAB 2%, EDTA 20 mM, 2-mercaptoethanol 0.2% v/v, NaCl 1.4 M, Proteinase K 30 micrograms in 500 μL of Tris 100 mM). The extracted DNA was precipitated with ethanol and 3 M sodium acetate, re-suspended in 20 μL of ultrapure water, and stored at -20°C.

PCR amplification and sequencing
The PCR reactions of the control-region of mitochondrial DNA (D-loop) were done using the external primers THR (5'AGCTCAGCGCCAGAGCGCCGGTCTTGTA AA3') (Lee et al., 1995) and 12S (5'ATAGTGGGGTAT CTAATCCCAGTT3') (Palumbi et al., 1991). PCR reactions were set up using 1 unit of Taq polymerase, 0.2 mM of each dNTP, 0.5 μM of each primer, 1.5 mM of MgCl 2 , 30 μg of BSA (bovine serum albumin) and 1 μL of DNA (approx. 50 ng) as a template, in a final volume of 30 μL of PCR buffer (Promega). Thermocycling conditions were: one initial cycle of 5 min at 94°C, followed by 35 cycles of 1 min each at 94°C, 48°C and 72°C, and one final extension step of 5 min at 72°C. Negative controls (without DNA template) were used in all PCR reactions to check for contamination. PCR products were purified, and both strands were sequenced in an ABI 3730 XL automatic sequencer, using the internal PRO (5'CCCAAAGCTAAAA TTCTAA3') (Kocher et al., 1989) and PHE (5'GCTTTAG TTAAGCTACG3') (Hedgecock and Strong, 1994) primers. The different haplotypes obtained were deposited in GenBank under accession numbers DQ666423-DQ666447 and EF624354-EF624386.

Sequence analyses
Sequences were edited using the SeqmanTM II program (DNAstar Inc.) and aligned using ClustalW (Thompson et al., 1994). For phylogeographic nested clade analysis (NCA), a haplotype network was generated using the TCS 1.21 program (Clement et al., 2000) using the 95% parsi-

310
O. chrysurus fisheries genetics mony criterion (Templeton, 1998). Clades were nested from haplotypes (0-step) to the highest level, each separated by one substitution step (Templeton et al., 1995;Templeton, 1998), and the statistical significance of each group was assessed, jointly with geographical coordinates for each sample, using the program GeoDis 2.2 . Selective neutrality of the control region sequences was tested with Tajima's (Tajima, 1996) and Fu's (Fu and Li, 1993) tests. Exact tests of population differentiation (Raymond and Rousset, 1995) were performed with 1,000 dememorisation steps and 10,000 Markov Chain steps. Selective neutrality tests, mismatch distribution analyses (Rogers and Harpending, 1992), Analyses of Molecular Variance (AMOVA; Excoffier et al., 1992) and exact tests of population differentiation were done using the Arlequin 3.0 program (Excoffier et al., 2005). Haplotype and nucleotide diversity indices were estimated using the DNAsp 4 program (Rozas et al., 2003).

Morphometry
Morphometric variation among populations was visualized via scatter plot of the scores of the first two canonical variables of the CVA applied to the shape variables of the TPS analysis (Figure 3). The first variable accounted for 76.8% of the among-population relative to within-population variation and discriminated the southern-most populations (Espírito Santo, Bahia) from the northern ones (Pernambuco, Ceará). A geographical gradient was also observed when considering the intermediate position of the Bahia specimens between Espírito Santo and the northern populations. The second variable accounted for only 26% of total variation, being related to differences between Ceará and Pernambuco and, to a lesser extent, between Bahia and Espírito Santo. Results of CVA indicated significantly among-populations differences (Wilks lambda = 0.1255; F [48,515] = 10.84; p < 0.0001). Percentage of individuals correctly classified to the four original populations ranged from 80% (Bahia) to 84% (Ceará and Espírito Santo).

Allozyme electrophoresis
Eleven loci were scored for all samples. Another four loci were scored for all populations except that of Pernambuco (Table 3). Observed heterozygosities (average Hardy-Weinberg expected heterozygosity, He = 0.1917; Table 3) were high (Smith and Fujio, 1982;Ward et al., 1994). Overall, the populations were found to be genetically structured (F ST = 0.0277; 95% confidence interval (bootstrapping over loci) = 0.00078-0.05648; standard error (jacknifing across loci) = 0.0195 reject the null hypothesis of population homogeneity). However, this resulted, mostly, from the large gene frequency differences found in two loci (G6PD-2 and GOT, Table 3) of the population from Pernambuco, which was the only one significantly different from the others (Table 4). When the samples from Pernambuco were excluded from the analyses, the F ST became non-significant (F ST = 0.0098; p > 0.30). The factorial correspondence analysis of unclassified samples showed a wide non-significant scatter of the data with only 27.1% of the total inertia explained by the first two factors (Figure 4). When the gravity centroids of each of the four sampling sites were compared, the two first factors accounted for 84.9% of the total inertia, showing significant differentiation between Pernambuco and the other sites.

Sequence variation
From the 663 base pairs of the control region of mitochondrial DNA sequenced, 105 were polymorphic, producing 58 haplotypes in the 77 individuals analysed. Genetic variation over all sampling sites was high (haplotype diversity, h = 0.9614; nucleotide diversity, π = 0.0185). The highest levels of variation were observed in Belize (h = 1.000; π = 0.02378). The null hypothesis of selective neutrality was not rejected for the studied sequences (Tajima's D test; p > 0.05; Fu's Fs test; p > 0.05) from most sampling sites. An exception was the Belize population, where D was not significantly different from zero, but Fu's Fs was significantly negative (Fs = -8.88; p < 0.001). A negative Fs results from an excess of rare alleles (as can be observed in the haplotype network in Figure 5), and it is usually interpreted as an indication of population expansion (Harpending, 1994). Hence, a mismatch distribution analysis (Rogers and Harpending, 1992) was performed on that population (Figure 6), and the observed distribution did not deviate significantly from the null hypothesis of population expansion (SSD = 0.004; p > 0.95). The Harpending's Raggedness index was low (r = 0.010), indicating a smooth distribution, also consistent with the hypothesis of population expansion (Harpending, 1994). Vasconcellos et al. 311 Levels of population differentiation, as indicated by pairwise F ST values, were significant only between the population from Belize and the others (Tables 4 and 5). No differentiation was observed among populations from Brazil with the exact tests of population differentiation (Markov chain test; exact P value = 0.72255 ± 0.02941). This was confirmed by AMOVA, where none of the Φ ST values among Brazilian populations was significant, and less than 7% of the total variance was explained by differences between geographical locations (Table 6). The nested clade analysis did not reveal significant (p > 0.05 after Bonferroni correction) groupings at any hierarchical level. Seven haplotypes, five of which from Belize, were excluded from the nested-clade analysis because they did not meet the 95% confidence level on the TCS haplotype network. This was caused by the very large divergence of those sequences in relation to the others, as can be seen in the neighbor-joining tree of all haplotypes (Figure 7).

Discussion
The analysis of both nuclear and mitochondrial genetic data of populations of yellowtail snappers, Ocyurus chrysurus, did not reveal genetic differences that could justify rejecting the null hypothesis of panmixia along the 2380 km of Brazilian coast, but significant differences were observed between Brazilian and Caribbean populations. In contrast, we observed significant differences in the mor-

312
O. chrysurus fisheries genetics  Figure 4). However, since only 5 individuals were analysed from that population, it is likely that the differences encountered simply result from sampling error. Indeed, when we excluded Pernambuco from the analysis, the mean F ST values became non-significant (F ST = 0.0098; p > 0.30). Therefore, we cautiously conclude that there is no support from the allozyme data to reject the null hypothesis of panmixia in this species along the studied area on the Brazilian coast.
The control region of O. chysurus is extremely variable (haplotype diversity = 0.9614). However, this variability is much higher within than between sampling locations in Brazil (Tables 5 and 6), and both the haplotype network generated by nested-clade analysis ( Figure 5) and the complete haplotype tree (Figure 7) show that the haplotypes from the Brazilian samples are distributed haphazardly along the coast. The Caribbean samples, on the other hand, differ significantly from the southwestern Atlantic ones (F ST = 0.17, p < 0.05; Table 5). This can also be seen in the haplotype network and the haplotype tree ( Figures 5 and 7), where most samples from Belize cluster outside of the southwestern Atlantic ones. The haplotype network has a Vasconcellos et al. 313   signal of recent population expansion (star-shaped topology around haplotype α in Figure 5), but that grouping was not significant in the nested-clade analysis, and none of the Brazilian samples presented deviations from neutrality in both Fu's and Tajima's tests. On the other hand, the sequences of the population from Belize diverged from neutral expectation ; p < 0.001), and the mismatch distribution analysis ( Figure 6) had a mode of 15 nucleotide differences, which is a clear signal of an old population expansion. The analysis of mtDNA thus indicates that the population from the Caribbean is old and currently stable. It shows also that the Caribbean population is different from the Brazilian one which is homogeneous along its distribution.
The results of the morphometric geometric analyses show what seems to be a north-south gradient of body shape, ranging from fusiform individuals, with a long caudal fin in specimens from the north to more rounded individuals, with shorter caudal fin in the south. Gradients in morphometric characters are often associated to phenotypic plasticity (O'Reilly and Horn, 2004). The morphometric differences observed are not incompatible with the existence of discrete fish stocks, since patterns of reef fish community change along the Brazilian coast often show a north-south differentiation (Floeter et al., 2001), but this result contrasts with the homogeneity observed with both nuclear and mitochondrial data. It has been argued that morphometric analyses can sometimes detect subtle differences in stock structure which are undetectable by genetic data, particularly when stock structure is the result of very recent population subdivision (Cadrin, 2000). However, the differences in morphometric data on a background of genetic homogeneity could also indicate phenotypic plasticity. This phenomenon may occur as response to gradients in environmental cues (Via et al., 1995). Such environmental gradient is likely to occur along a geographic range of more than 15°of latitude and ca. 2,400 km of coastline as the one inhabited by the studied Brazilian populations of Ocyurus chrysurus. Differences in oceanographic conditions, mainly sea-water temperature, were reported to this area (Castro and Miranda, 1998) being regarded as an important factor in the genetic differentiation of Macrodon ancylodon populations along the Brazilian coast (Santos et al., 2003). Cases where morphometric and genetic data indicate different scenarios of population structuring are not uncommon (Salini et al., 2004;Levi et al., 2004). In the case of Ocyurus chrysurus, the fact that the morphological variation observed formed a gradient is a clear indication that it has likely resulted from phenotypic plasticity. We thus conclude that the genetic data do not allow the rejection of the null hypothesis of panmixia for Brazilian yellowtail snapper populations which could be treated as a single genetic stock, with a latitudinal gradient on their morphology. On the other hand, there is a severe restriction to gene flow between O. chrysurus populations from the Caribbean and from the southwestern Atlantic.