Gulf Menhaden Brevoortia patronus is one of the most intensively fished organisms in the Gulf of Mexico. Fishery managers and industry operators have historically worked towards a sustainable fishery and have cooperated on stock assessments to estimate feasible reference points for management. These stock assessments are necessarily rooted in a strong definition of the fishery stock, including the estimation of the number of populations that are exploited. Here, genetic population structure and variability were examined in B. patronus, using genetic markers specifically designed for the species. We observed genetic variability that indicates a relatively high effective population size for a marine finfish (Ne ≥ 1,200), and two analytical approaches implied a single genetic population of B. patronus. We compared the latter finding with the population structure in the closely related Yellowfin Menhaden B. smithi, for which two distinct populations from Florida were identified using the same genetic loci (Fst = 0.015, P = 0.027). These contrasting patterns of population structure between sympatric congeneric species are likely driven by differences in distribution and census size and may relate to factors that originally drove speciation in North American Brevoortia species. The finding of a single Gulf-wide population of B. patronus suggests that there is an extensive migration throughout the species range and supports the notion of a single genetic stock.
The purse seine fishery for menhaden Brevoortia spp. (Clupeidae) is one of the oldest fisheries in United States (Smith 1991) and is supported primarily by two species: the Atlantic Menhaden B. tyrannus and the Gulf Menhaden B. patronus. Presently, B. patronus is the most intensively harvested finfish species in the Gulf of Mexico (“Gulf”) and supports one of the largest commercial fisheries in the United States (Vaughan et al. 2007). This species forms numerically immense, shallow, nearshore schools, and this tendency makes it easy to harvest in high numbers using purse boats deploying seines up to 365 m in length (Smith 1991). Despite this intensive harvest, landings of B. patronus appear to be generally below the maximum sustainable yield, and the fishery is thought to be stable (Vaughan et al. 2000).
The maintenance of a long-term sustainable harvest of B. patronus is important given the prominence of this fishery in the Gulf. A critical element for achieving this goal is to develop stock assessments that are anchored by a strong definition of the unit stock (Begg and Waldman 1999). The definition of a fishery “stock” has evolved through the years from one that was primarily driven by distinct features of local fishing economies to one that is more rooted in the biology of the organism and includes assessments of reproduction, demo-graphics, and genetics (Carvalho and Hauser 1994; Begg and Waldman 1999; Waples et al. 2008). Modern stock concepts take into account practical fishery considerations as well as political concerns, but are also heavily rooted by the idea of a genetic stock (Carvalho and Hauser 1994). For B. patronus, the stock definition used in stock assessments has historically been defined by multiple data considerations, including an estimation of the number of genetically distinctive populations occurring in the fished area (SEDAR 2013). However, this component of the stock definition has been limited by the magnitude of genetic sampling in previous studies (Anderson 2006). An increase in the number of genetic markers used for eliciting discreet genetic stocks would greatly improve the stock definition by improving statistical power (SEDAR 2013; Anderson and Karel 2014).
The possibility of multiple genetically distinct stocks in B. patronus may be confounded by a potential for extensive gene flow across the range of the species based on aspects of reproductive biology at multiple life stages. The species undergoes an offshore migration beginning in October, and the fish are caught offshore at various depths in the northern Gulf (Roithmayr and Waller 1963). Based on tag return data, in addition to this offshore movement there alsomay be a tendency for adults of the species to move closer to the center of the species' range (Ahrenholz 1981; although see Pristas et al. 1976). The seasonal migratory period of adults also coincides with the spawning season (Raynie and Shaw 1994), and spawning in B. patronus occurs from near shore to at least 96 km offshore (Christmas and Waller 1975). The offshore pelagic egg—larval stagemay last up to 10 weeks (Shaw et al. 1985; Deegan and Thompson 1987), whereby migration occurs passively and may rely on nearshore currents that potentially move individuals laterally along shorelines before they settle into estuaries. Thus, two stages of migration (adult seasonal active, larval passive) result in a potential for long-distance gene flow and homogenization of allele and genotypic frequencies across a broad scale. This would imply the potential for a single unit stock, although this hypothesis has not been rigorously examined with a large genetic data set.
An additional consideration with regard to management of B. patronus is whether intensive harvest or other processes have impacted the genetic variability or the effective population size (Ne) in Gulf stocks. There is a high estimated natural mortality in B. patronus (Smith 1991). In marine systems, variances in productivity and survivorship among demes and high estimates of natural mortality have historically been linked to wide disparities between population census size and genetic effective size (Turner et al. 2002). This phenomenon can be exacerbated in heavily harvested species; long-term harvest has been related to low Ne in other exploited species (Diaz et al. 2000; Hauser et al. 2002). A recent meta-analysis of marine populations indicated that genetic diversity is generally lower in overharvested populations than in stable ones (Pinsky and Palumbi 2014).
The commercial harvest of B. patronus is generally heaviest in an area between Cameron, Louisiana, and Moss Point, Mississippi, which is roughly the center of the species' distribution. Estimates from landings data suggest that 4–11.5 × 109 fish are harvested each year (Smith 1991; Vaughan et al. 2007), but the stock is thought to be stable and not overfished (Vaughan et al. 2000). Thus, one question that might be answered with genetic data is whether there are observable negative impacts on genetic diversity in the center of the range of B. patronus due to fishing. Bowen and Avise (1990) previously estimated the effective size of the female population in menhaden (Nfe) to be ~250,000 individuals, suggesting very high genetic diversity. However, the method used in that study was based on nucleotide diversity at a single mitochondrial DNA (mtDNA) locus and relied upon rough estimation under a coalescent model. The estimate of Bowen and Avise (1990) thus may simply reflect long-term (evolutionary) effective size among mtDNA lineages only. In comparison, a short-term estimate of Ne would more accurately reflect the effective size of the contemporary population and would be more useful for estimating existing adaptive potential (Hare et al. 2011).
Here, the question of gene flow, population structure, and genetic variation of B. patronus was addressed by applying a suite of novel genetic markers (Anderson and Karel 2014) to Gulf-wide samples of Brevoortia species. The study design strategy involved sampling primarily juveniles from inshore areas, which also incidentally resulted in relatively large samples of Yellowfin Menhaden B. smithi from Florida. The population structure of B. smithi in Florida was compared with that of B. patronus Gulf-wide. This genetic examination of the unit stock in the B. patronus fishery can be used to improve the stock definition for future Gulf Menhaden stock assessments.
Sample collection and laboratory methods.—Menhaden specimens were collected throughout the Gulf during two targeted sampling periods (Table 1; Figure 1). The first (“early”) period was during the years 2002–2003. These samples were collected during the course of a previous study (Anderson 2006), but DNA extracts from these specimens were reanalyzed in this study using additional genetic markers. Sampling methods for early sampling can be found in Anderson (2007).
The second (“late”) group of samples was collected in 2012–2014. These included five sampling locations in Texas, two locations in Louisiana, and three locations in Florida. When possible, sampling locations for late samples were geographically proximate to those for the early samples. Menhaden were collected using various gears, mainly bag seines and inshore trawls. In most cases, samples consisted of whole juvenile specimens (<100 mm TL). In a few cases, adults were sampled by excising a small portion of fin tissue from the caudal fin. Whole genomic DNA was extracted from caudal fin tissues of all specimens using a Gentra Puregene tissue isolation kit (Qiagen, Valencia, California) following the manufacturer's instructions.
Sample sizes of Brevoortia for each sampling locality and each of two temporal (early versus late) sampling periods sorted from west to east longitude. Samples were broken down by species and hybrids using genetic species identification.
Fourteen microsatellite loci described in Anderson and Karel (2014) were used to develop multilocus genotypes for all specimens. Nine of these were novel loci developed specifically for the genus Brevoortia. The remaining five were developed for use in a closely related clupeid, American Shad Alosa sapidissima, and were included here for comparison with previous genetic examinations that used these loci in Brevoortia species (Anderson 2006, 2007). Each primer set was labeled with a WellRed fluorescent dye (primer manufacture and labeling by Sigma-Aldrich) and used to amplify approximately 100 ng of genomic DNA for each specimen through the use of PCR. Amplified products were multiplexed in four independent panels with a 400-bp size standard and separated with a CEQ 8000 sequencer (Beckman Coulter). Fragment analysis was performed using Gene Marker version 2.4 (SoftGenetics LLC, State College, Pennsylvania), and all allelic bins were established based on estimated allele size in base pairs.
Genetic species identification.—Samples consisted primarily of small juveniles, and it was anticipated that some individuals sampled in the western Gulf may actually have been Finescale Menhaden B. gunteri and that some individuals sampled in the eastern Gulf may actually have been B. smithi. The range of B. patronus overlaps the ranges of each of these species. Morphological identification in juvenile specimens is extremely difficult, so multilocus microsatellite genotypes were used to determine species in all specimens. This was accomplished by using the software Structure version 2.3.4 (Pritchard et al. 2000). The analysis consisted of 10,000 burn-in iterations followed by 90,000 iterations for clustering. Populations were constrained to K = 2 to allow for the identification of two general categories: B. patronus and combined B. gunteri—B. smithi (the latter species likely share common lineage and generally cluster together in phylogenetic analyses: Anderson 2007). Genotypes were sorted into these categories, and if individuals shared a cluster lineage >0.20 with both inferred clusters, they were classified as hybrids (as in Anderson and Karel 2007). A second identical Structure analysis (K = 2, 10,000 burn-in iterations and 90,000 clustering iterations) was used to separate putative B. gunteri and B. smithi. Species identification was aided in this case by the knowledge that B. gunteri generally occur in the western Gulf west of the Mississippi River and B. smithi occur in the eastern Gulf (Dahlberg 1970). When the results of this second Structure analysis did not correspond to expectations based on sampling geography, an additional Structure run was completed with K = 3 to provide further resolution (see Results).
Microsatellite diversity and genetic structure.—The observed (Ho) and expected (He) heterozygosity of each marker locus was calculated for each species observed in the sample, and the within-group coancestry coefficient (Fis) was calculated over all samples (species specific) using the software Fstat version 2.9.3 (Goudet 1995). It was expected that He would generally exceed Ho and that most loci would demonstrate elevated Fis relative to the Hardy—Weinberg equilibrium expectations, as a result of combining two temporal sampling periods, with samples spread across the Gulf. The confidence interval for each locus-specific estimate of Fis (with all samples combined) was estimated in Fstat and used to determine which loci had elevated Fis relative to the overall mean Fis. Marker loci that deviated significantly from the group-wise estimate of Fis were included in genetic structure analysis but were excluded in comparisons of genetic diversity and effective population size (Ne, see below). It was expected that the genetic structure analysis would be fairly robust to violations of Hardy—Weinberg expectations but that estimates of genetic diversity and Ne would not.
Genetic structure was assessed in B. patronus and B. smithi independently using Structure version 2.3.4. All samples from across the sampling geography were combined, and three levels of potential group structure were assessed (K = 1, 2, and 3) for each species. This analytical framework assumes that some level of population subdivision may be detected across the species range, but no assumptions were made as to the geographical arrangement of those populations (that is, geographic data were not used to influence a priori probability distributions). Previous assessments of genetic structure in B. patronus suggested no Gulf-wide genetic structure and high rates of migration across the Gulf (Anderson 2006; Anderson and Karel 2014). Therefore, the Structure analysis was completed assuming admixture among populations and also allele frequencies that were correlated among populations. Structure runs consisted of 10,000 sampling iterations followed by 90,000 cluster iterations; the number of iterations was determined after observing convergence of model parameters in trial data runs. At each level of K, 10 replicate runs with identical run length and model parameters were used to determine whether results achieved consistency among runs. This analytical framework was repeated with samples of B. smithi.
Multiple marine taxa have demonstrated genetic structure across the northern Gulf, and in many cases this structure can be related to taxa occurring roughly east and west of the Mississippi River drainage (Dahlberg 1970; Portnoy and Gold 2012). A previous examination of the distribution of morphological characters among the North American forms of Brevoortia suggested that B. gunteri and B. smithi were sister taxa that had diverged according to this general pattern (Dahlberg 1970). Although no such pattern has yet been identified in populations of B. patronus, tagging data suggests that movement across the Mississippi Delta is infrequent (Pristas et al. 1976). The hypothesis that there may be divergent populations of B. patronus in the eastern and western Gulf was tested explicitly using the analysis of molecular variance (AMOVA; Excoffier et al. 1992). The AMOVA framework assessed genetic structure at three levels: among groups (western versus eastern Gulf), among samples, and within samples (among individuals). Divergence was assessed using F-statistics derived from the AMOVA sum of squares, and permutation tests (n = 1,000) were used to assess statistical significance of each F-statistic. A significantly positive F at any level was assumed to indicate significant genetic divergence. The AMOVA and all associated statistical testing were carried out using the software Arlequin version 3.5 (Excoffier et al. 2005).
Genetic diversity and Ne.—The central—marginal hypothesis was tested in B. patronus with the null hypothesis (H0) that genetic variability was distributed evenly throughout the sampled range. The central—marginal hypothesis posits that, across a species range, genetic variability will be highest in populations at the center of the range due to larger effective population sizes and greater rates of migration than what is expected in marginal populations (Eckert et al. 2008). This mechanism may be even greater in B. patronus given that tagging data suggest that adults migrate towards the center of the species distribution over time (Ahrenholz 1981). Estimates of genetic variability were generated for each B. patronus sample using multilocus heterozygosity and the estimator θH (derived from Ohta and Kimura 1973), as implemented in Arlequin 3.5. Estimates were replicated using early (2002–2003) and late (2012–2014) data sets independently, and a t-test was used to determine whether mean values of θH from each sampling period differed significantly. Genetic variability, θH, was plotted against longitude to qualitatively assess whether θH was higher in middle longitudes or consistent across the species range. To quantitatively test for higher θH in the center of the species distribution, the Mississippi River was used to anchor the distribution of samples, and linear shoreline distance from the Mississippi River to each sample locale was approximated using the plot tool in Google Earth (Google, Mountain View, California). The relationship between this distance and θH was tested for statistical significance using least-squares regression. The mean θH over all samples and all loci was used to indirectly estimate effective population size of B. patronus Gulf-wide, using the formula θH = 4Neμ, where μ is the study-wide mutation rate of all microsatellite loci. Estimates of Ne were generated under fast (μ = 1 × 10-3) and slow (μ = 1 × 10-5) assumed mutation rates.
The effective population size (Ne) of B. patronus was also estimated directly from genotypic data using two methods: (1) the bias-corrected linkage disequilibrium (LD) method (Hill 1981; Waples 2006; Waples and Do 2010), and (2) the moment-based temporal method (Waples 1989). Initial estimates generated using the software NeEstimator 2.1 (Do et al. 2014) suggested that point estimates of Ne using the single-sample LD method usually were very large and had confidence intervals that overlapped infinity. Therefore, confidence intervals around each estimate were used to evaluate a minimum value of Ne. In particular, the lower bounds of estimates of Ne have value from a conservation standpoint and relate directly to the long-term adaptability of populations (Nelson and Soulé 1987). The methods used for Ne estimation are generally sensitive to low frequency alleles; in each case, estimates were generated while throwing out individuals with low-frequency alleles using three threshold frequencies: 0.01, 0.02, and 0.05.
Only the five Texas sampling locales that were identical between sampling time periods were used for Ne estimation. Estimates of Ne based on the single-sample LD method were generated independently for both the early and late samples; in each case, samples were pooled into a single population. Pooling samples was justified by a lack of genetic structure found among sample locales (see Results) and resulted in more robust sample sizes for Ne estimation. The temporal method required the use of both data sets simultaneously and an estimate of the number of generations that occurred between samples. Generation time was set at t = 1.5 years to account for the expected “knife-edge” maturity that occurs between ages 1 and 2 (Vaughan et al. 2007), as well as high estimated natural mortality and short life span of B. patronus (Ahrenholz 1981; Vaughan et al. 2000, 2007). The estimator Fs (Jorde and Ryman 2007) was used for the temporal estimate of Ne; preliminary runs using NeEstimator suggested that this Fs was generally not as sensitive as other estimators to rare alleles in the Gulf Menhaden data set (e.g., Nei and Tajima 1981; Pollak 1983).
Genetic Species Identification
All three species of Brevoortia known to be present in the Gulf were observed in the data set based on genetic species identification (Table 1). The most commonly observed species was B. patronus (overall n = 750), which was present in 17 out of 19 samples. The second most commonly observed species was B. smithi (n = 126), occurring in four Florida samples. Brevoortia gunteri was only observed twice, once each in early (2002) and late (2013) samples in the upper Laguna Madre, Texas. The contrasting distributions of B. patronus and B. smithi in Florida suggest that the latter species numerically replaces the former somewhere south of Apalachicola, Florida, although this finding is tenuous given the limited sampling in this study. Five hybrid genotypes were observed, all of which occurred in Florida samples. No attempt was made to resolve the genetic status of hybrids (F1 versus late-generation back-cross), although the genetics of hybrids has been examined previously (Anderson and Karel 2007).
Interestingly, when the small-scaled menhaden, B. gunteri and B. smithi, were analyzed together using Structure, the model representing the most basic level of population structure (K = 2) did not conclusively resolve the species (Figure 2). Specifically, in the K = 2 model, the two individuals sampled in the upper Laguna Madre (presumed to be B. gunteri) were not conclusively assigned to either population cluster based on admixture scores (individual admixture proportions for individual 1 = 0.40/0.60 and individual 2 = 0.32/0.68). This could be interpreted one of two ways: (1) the two individuals identified tentatively as B. gunteri based on location of catch (Laguna Madre) were actually B. smithi that had travelled farther west than what would be expected based on their species range, or (2) shallow lineage divergence between B. gunteri and B. smithi has resulted in multiple common alleles, resulting in difficulty in resolving these species using the given microsatellite markers. This latter interpretation is strengthened by the fact that themarkers used in this study were heterologous (originally developed for B. patronus and A. sapidissima) and had low genetic variation in B. gunteri and B. smithi (Anderson and Karel 2014). An additional Structure run atK= 3, and using sampling location as a model prior, conformed to geographic expectations and easily distinguished the species.
Microsatellite Diversity and Genetic Structure
The expected heterozygosity He exceeded that of observed heterozygosity Ho at each genetic locus in B. patronus (Table 2). This resulted in positive Fis at all loci and an overall estimate of Fis = 0.087 (±0.029). Only a single locus (AF39660) had an exceptionally high estimate of Fis that exceeded the jackknifed estimate over all loci, indicating the possibility of null alleles. This locus was retained for genetic structure analysis but was excluded when estimating samplewise genetic variability (θH) and effective population size (Ne) in B. patronus. In B. smithi, the overall Fis was higher than in B. patronus (Fis = 0.116 ± 0.055), and five loci had significantly high estimates of Fis (BP003, BP039, BP221, BP500, AF39660). Low genetic diversity and deviation of these loci from Hardy—Weinberg expected genotypic frequencies in B. smithi has been demonstrated previously and likely is the result of using heterologous loci (Anderson and Karel 2014). As a result of high estimates of Fis, as well as limited sampling of B. smithi, sample genetic variability (θH) and effective population size (Ne) were not estimated in this species.
Structure analysis of B. patronus samples suggested a single genetic stock that comprised all samples. The single population (K = 1) model had the highest posterior probability of all tested models (K = 1–3) in each independent run. Models that assumed more complex population structure (K = 2, 3) resulted in admixture proportions for each individual that were approximately equal mixtures of each of the genetic clusters retained (equivocal with respect to geography). The AMOVA that separated samples into eastern and western Gulf components also showed no evidence of genetic divergence between groups (Fct < 0.001, P = 0.403), nor was there evidence for significant divergence among samples (within groups; Fsc < 0.001, P = 0.844). A post hoc examination of pairwise estimates of the genetic differentiation index, Fst, (among all samples) suggested no significant genetic divergence between any two samples after adjustment for multiple tests performed simultaneously.
Observed (Ho) and expected (He) heterozygosity, and estimated coancestry coefficient (Fis) of each genetic locus and overall among Brevoortia patronus and B. smithi samples observed in this study. The confidence intervals around individual locus estimates of Fis were compared with the mean estimate of Fis (jackknifed over loci) to determine which loci exceeded the overall estimate (statistically high Fis values are marked with an asterisk).
In contrast to B. patronus, there was evidence for significant genetic structure among samples of B. smithi. Specifically, samples taken in Apalachicola and Cedar Key demonstrated a majority of genetic contribution (based on individual admixture scores) from cluster 1, which contrasted with Tampa Bay and Charlotte Harbor samples, which had a majority of contribution from cluster 2 (Table 3; Figure 2). This result suggests there is a “northern” and “southern” population of B. smithi in Florida. The AMOVA of B. smithi samples included the three largest samples (Cedar Key, Tampa Bay, and Charlotte Harbor) and also suggested significant divergence among samples (Fst = 0.015, P = 0.027). A post hoc examination of pairwise Fst scores suggested that divergence between Cedar Key (from the northern population) and the two southern samples drove the significant Fst in the AMOVA. The Fst between Cedar Key and Tampa Bay was Fst = 0.019, P < 0.001; the Fst between Cedar Key and Charlotte Harbor was Fst = 0.016, P = 0.006.
Genetic diversity and Ne
The estimates of θH among all samples of B. patronus ranged from 1.77 (upper Laguna Madre, 2002) to 1.90 (Galveston Bay, 2013). There was no significant difference in θH between early and late samples (t = 0.61, P = 0.55). Qualitatively, θH appeared to be elevated in northern Gulf samples relative to the edges of the species range, and this effect was replicated across both sampling periods (Figure 3). However, the model of regression of θH against sample distance from the Mississippi River was not significant (F = 0.07, P = 0.79). The mean θH over all loci and all samples was 1.85, which equated to a range of Ne from 462 to 46,125 depending upon usage of slow versus fast mutation rates.
Proportion of membership of each small-scaled menhaden (Brevoortia smithi and B. gunteri) sample in genetic clusters, generated by Structure analysis, and the samples size (n) used for each sample location. The upper Laguna Madre sample was presumed to be two B. gunteri individuals based on the known distribution of each species; all other samples were assumed to be composed of B. smithi.
Point estimates of Ne using the LD method were unstable, ranging from 1,437 individuals to ∞ in the early samples (depending on low allele frequency threshold; Table 4). In late samples, point estimates at all cutoff levels were ~∞. The lowest bounds of the confidence interval for any allele frequency cutoff threshold were 611 (in early samples) and 1,139 (in late samples). The temporal method yielded more stable estimates and confidence intervals of Ne. Point estimates using the temporal method ranged from 1,201 to 1,435, and the lower and upper bounds of confidence intervals ranged from 819 to 1,822.
Estimates of effective population size (Ne) in Texas samples of Brevoortia patronus. Estimates were generated using two methods: the single-sample linkage disequilibrium (LD) method and the two-sample temporal method. In the former case, estimates of Ne in early and late samples were generated independently. In all cases, Ne was estimated under three lowfrequency allele thresholds (0.05, 0.02, and 0.01).
The modern purse seine fishery for menhaden is centered between Cameron, Louisiana, and Moss Point, Mississippi, based on the location of reduction facilities (Smith 1991; SEDAR 2013). Nominal harvest takes place outside of this range, but close to 95% of the harvest takes place in this general geographic distribution (Smith 1991; SEDAR 2013). Previous assumptions about the taxonomic makeup of the menhaden catch have assumed <1% of the catch to be B. gunteri and B. smithi (Ahrenholz 1981). In this study, neither of the congeneric species observed (B. gunteri, B. smithi) were found in the range from Sabine Lake, Texas, to Moss Point, which roughly encompasses the fished area. These data generally support previous assumptions about the composition of the fishery and suggest that B. patronus likely constitute >99% of the menhaden catch in the targeted area. Harvest in areas outside of this range, particularly in Florida, may encounter other species on a more regular basis. However landings in Florida have historically accounted for <1% of the total annual catch of Gulf Menhaden (SEDAR 2013), and therefore the impact of harvest on B. smithi is likely to be minimal, assuming that sampling in this study is representative of species presence—absence.
With regard to population structure, previous studies of other marine clupeids have demonstrated gene flow over vast spatial scales, resulting in basinwide genetic homogeneity at neutral loci in several species (e.g. Ryman et al. 1984; Kinsey et al. 1994; Mariani et al. 2005; Lynch et al. 2010). The commonly observed pattern of basinwide homogeneity in marine fishes in general has generally been attributed to large population sizes and limited barriers to dispersal that typically characterize marine systems (Waples 1998). The biology of Gulf Menhaden, in particular the presence of an offshore egg—larval stage and the potential for vast migration at both the larval and adult life stages, suggests a predisposition for this species to conform to the expectation of genetic homogeneity over the range of the species. Two independent tests of this assumption here suggest a lack of localized genetic divergence in Gulf Menhaden across the sampled range. A naive Structure cluster analysis that made no assumptions about population boundaries suggested that all samples collected from Texas to Florida resided in a single genetic cluster. Additionally, a model-based AMOVA that assumed the Mississippi River as a boundary between independent populations was not supported by elevated levels of divergence. In fact, this latter analysis suggested that there was not a single pairwise group of samples with elevated genetic divergence between them, confirming basinwide genetic homogeneity in the Gulf. This finding was also implied by results from previous studies (Anderson 2006; Anderson and Karel 2014) and strongly supports a single-stock hypothesis for B. patronus in the Gulf.
In contrast to the basinwide homogeneity observed in B. patronus, there was significant genetic divergence among localized populations of B. smithi. This result was elicited both by a Structure analysis that suggested two genetic clusters in B. smithi and also by elevated levels of a traditional genetic divergence metric (Fst) in comparisons between northern and southern Florida samples. Although the observed genetic divergence among samples was low (Fst = 0.015), it was statistically significant. Moreover, pairwise estimates of divergence among various samples of B. smithi resulted in a repeated pattern; specifically, two independent southern samples (Tampa Bay, Charlotte Harbor) displayed significantly elevated divergence with a single northern sample (Cedar Key). These results imply limited gene flow among localized groups of B. smithi in Florida and are similar to the result described in Anderson (2007), although in that case the divergence was higher (Fst = 0.052) and was observed between Gulf and Atlantic populations of B. smithi rather than two Gulf populations. The biology of B. smithi with respect to distribution may play a role in the establishment of genetically distinct populations in this species. For instance, B. smithi generally reside inshore and do not undertake the vast offshore migrations observed in B. patronus (Dahlberg 1970). Additionally, B. smithi may not form the vast schools observed in B. patronus (Hildebrand 1963) and likely has a smaller effective population size relative to B. patronus and B. tyrannus (Anderson 2007). Taken together, these factors (inshore distribution, limited migration, relatively low Ne) could be expected to result in elevated genetic divergence between localized populations. This interpretation is striking when viewed in the context of genetic homogeneity of the congeneric B. patronus across the entire Gulf, as well as B. tyrannus in the Atlantic Ocean (Lynch et al. 2010), and may hold implications towards the mechanisms that initially drove speciation in the various North American species of Brevoortia.
Estimates of Ne in B. patronus using the LD method were very high and sensitive to the presence of rare alleles. This could potentially be the result of effective population sizes that are higher than can be accurately estimated with the available data; a result that can be interpreted in two interrelated ways: (1) the Ne in Gulf Menhaden is significantly large that the LD method is only useful with extremely large sample sizes, and (2) the number of loci or individuals examined here was not high enough to detect the expected LD generated by a large but finite Gulfwide population. Whatever the reason for unstable estimates of Ne, the point estimates generated using both the historical and recent data sets generally approached infinity, and confidence intervals around point estimates using a variety of estimators suggested an absolute lower bound of Ne = 611. Confidence intervals generated using the temporal method implied a lower bound of Ne = 819 but suggested more conservative point estimates (1,201 < Ne < 1,435). In addition to the instability in LD estimates, high gene flow and insufficient sampling could also create downward bias in estimates of Ne using the temporal method (Wang and Whitlock 2003; Waples and Do 2010; Hare et al. 2011). In general, these estimates that suggest Ne > 1,000 are supported by previous examinations of Gulf Menhaden reporting exceptionally high genetic variability based on mtDNA (Bowen and Avise 1990; Anderson 2007; Lynch et al. 2010) and microsatellite (Anderson 2007; Lynch et al. 2010) diversity.
With regard to the spatial distribution of genetic diversity, estimates of θH from Gulf-wide sampling of B. patronus suggest that genetic diversity is highest in the northern Gulf, towards the center of the species range. This result was repeated in both historical and contemporary sampling, suggesting repeatability in this pattern between two data sets. Based on tagging data, Ahrenholz (1981) suggested that adult menhaden migrate towards the center of the species range as they age. This migratory pattern would result in frequent immigrants into the center of the range, with immigrant frequency declining as distance from the center of the range decreases. Such a pattern could be expected to generate the spatial distribution of genetic diversity observed here. It would also potentially result in the level of gene flow that would be necessary to homogenize allele frequencies across the species range, resulting in a lack of discernible genetic structure. Our finding that no genetic structure was observed among spatial samples in this study may imply just such a mechanism.
If the finding that genetic variability is highest in the center of the species range is a biological reality, it implies that the fishery generally targets individuals in an area with the highest potential for genetic adaptability and that the fishery has likely had a minimal impact on this central—marginal pattern. The caveat to this point is that the data failed to show a statistical pattern with respect to this distribution of genetic variability, as θH was not significantly correlated with distance from the Mississippi River. Furthermore, the genetic variability statistic θH is influenced heavily by sample-level genetic heterogeneity and thus is heavily biased with respect to sampling error or allele-size calling error. Ideally, the central—marginal pattern in B. patronus could be tested using sample-based estimates of effective population size, although the current data suggest that estimates of Ne are unreliable even with very large sample sizes. Genomic sampling methods likely will be necessary to generate proper estimates of Ne and test the central—marginal pattern conclusively.
In conclusion, an inclusive data set of microsatellites designed specifically for use in B. patronus has repeated previous findings of genetic homogeneity in Gulf-wide samples and supports single-stock treatment of B. patronus in stock assessments. Given this repeated finding, it is unlikely that the fishery might exploit rare or exceptional localized genetic diversity in areas where it operates. Because of logistical considerations, the fishery targets populations where the species is most commonly found, in the center of the species range (Smith 1991; SEDAR 2013). Based on estimates of genetic diversity, this area may be a genetic sink, where elevated rates of immigration result in high genetic diversity relative to the edges of the species range. Estimates of population-wide genetic diversity (θH) and effective population size (Ne) suggest relatively high genetic diversity for a marine finfish and reflect the presumed enormous census size projected from fishery-dependent estimates (Vaughn et al. 2000). Estimates of effective population size were similarly large between early and late samples, suggesting that long-term fishing pressure has not had an observable influence on Ne, although uncertainty associated with estimates of Ne makes this conclusion tenuous. In any event, using the temporal method point estimate of Ne = 1,201, the effective size of the B. patronus population in the Gulf appears to be higher than mean estimates generated for other exploited wild populations (Palstra and Ruzzante 2008), similar to or higher than estimates generated for other notable marine fish populations (Hauser et al. 2002; Turner et al. 2002; Hoarau et al. 2005), and above the conservation threshold (Ne = 500) that may be expected to result in long-term population adaptability (Nelson and Soulé 1987).
The first draft of this paper was improved by comments from S. Vanderkooy, Gulf States Marine Fisheries Commission, as well as D. McDonald and W. Karel, Texas Parks and Wildlife Department (TPWD). M. Fisher and three anonymous reviewers made useful suggestions during the course of TPWD internal review. Numerous field staff of the Coastal Fisheries Division of TPWD assisted in field sampling of menhaden specimens. This work was supported in part by funding assistance via the Federal Sportfish Restoration Program, U.S. Fish and Wildlife Service.