Abstract

Waterfalls and drainage divides have been found to act as barriers to gene flow in riverine communities, affecting migration and genetic diversity of fish populations. We studied populations of Rhinichthys atratulus, the Eastern Blacknose Dace, in two tributaries of the Mattabesset River, Connecticut, to test hypotheses about genetic diversity in relation to physical barriers to migration and gene flow. We investigated three specific hypotheses: i) That waterfalls or drainage divides should be barriers to gene flow; ii) That the amount of genetic divergence across these barriers should be greater than among populations not separated by physical barriers; iii) That genetic diversity below waterfalls should be greater than the genetic diversity above waterfalls. We tested these hypotheses using mitochondrial and nuclear genetic markers from 191 and 197 individuals, respectively, from seven locations; the populations were separated by three waterfalls and two drainage divides. We found that: i) The waterfalls and drainage divides were associated with genetic divergence among populations, indicating barriers to migration and gene flow; ii) All populations were significantly different from one another for microsatellites and most for nd2, but that the magnitude of divergence across barriers was not greater than divergence between populations not separated by physical barriers; iii) Genetic diversity was higher above falls for private alleles and private haplotypes but was higher below falls for haplotype and nucleotide diversity. Migration rates among populations ranged from 1 to 3.2 individuals per generation. These migration rates were not sufficient to homogenize the populations over the effects of mutation. Although the mitochondrial data were consistent with an isolation by distance model of evolution, the results were due to the number of private haplotypes in the upper part of the Coginchaug River. We discuss the roles of extinction and geological history, as well as ecological habitat preferences and predation risk in determining genetic divergence among populations.

Keywords

Gene Flow, Barriers, Migration, Mutation, Rhinichthys atratulus, Waterfalls, Evolution

Share and Cite:

1. Introduction

Habitat fragmentation is the natural or human-altered process by which continuous habitats are broken into isolated patches. Patch size not only affects ecological processes but also reshapes evolution by altering selection pressures, gene flow, and population dynamics [1]-[7]. As landscapes are divided by waterfalls, agriculture, urbanization, dams, and other infrastructure, many species must adapt rapidly or face extinction [8]-[10]. Fragmentation consistently decreases biodiversity across taxa and ecosystems over time [11] [12].

Habitat fragmentation has ecological and genetic repercussions for populations. Both anthropogenic and natural barriers disrupt the exchange of individuals between populations, thereby reducing gene flow [1]-[7]. Reduction of gene flow among populations can have a variety of effects such as inbreeding depression and genetic drift, which alter the genetic structure of these populations [10] [13]-[20]. Lacking the available genetic variation to adapt to a changing climate, for example, may make these populations more vulnerable to extinction [14] [18] [20]-[24]. While barriers can reduce biodiversity in small populations, they can also lead to increased genetic diversification or even speciation depending on the life history of a species and the population size [25].

At a broad scale, eco-evolutionary feedbacks dominate fragmented landscapes. Fragmentation alters both the context and the tempo of evolution, especially for dispersal and life history [8]. These dynamics apply equally to terrestrial and aquatic ecosystems where dispersal becomes more necessary and more risky [8] [10]. For example, human-altered fragmentation has effectively become a driver of contemporary evolution shaping species distributions, genetics and ecosystem processes [26] [27]. In many river systems there are measurable evolutionary changes among freshwater fish populations within a few decades [14] [28] [29] or even a few generations [30] [31]. Yet such rapid evolution rarely restores ecological connectivity. Without habitat protection and restoration efforts, adaptation may lead to smaller, more specialized, and less resilient populations.

The evolutionary effects of fragmentation are profound in freshwater fishes where river damming and flow alteration create severe barriers to gene flow [31]-[34]. Although the effects of dams are not part of our study, dams represent a controlled experiment about evolutionary divergence and rates of allelic change that are not necessarily available to studies about fragmentation to waterfalls. Dams divide once-continuous river systems into discrete segments, isolating fish populations and disrupting migration routes that are critical for spawning. The isolation of freshwater fish populations can result in genetic divergence from previously connected demes [33] [34]. Hydrological alterations and damming resulted in high levels of genetic divergence among populations of five species of cypriniform fishes in Iran [35], and in salmonid fishes [28] [34]. The amount of divergence is also a function of the degree of isolation of populations on both sides of the barriers. For example, in Oregon, the largest and least penetrable barriers resulted in increased divergence in Redband Trout (Onchorynchus mykiss gairdneri) [36]. The rate of divergence has also been found to be inversely proportional to the effective population size of the fragmented species [31] [33] [37]. A community wide study found that four species of fishes with high effective population sizes had low divergences relative to dam populations below [37], whereas the Smallmouth Bass, Micropterus dolomieui, which had a low effective population size, was highly divergent. Low effective population sizes allow for drift to work more quickly to shift stochastically allele frequencies [31], as has been shown for the Variegate Darter (Etheostoma variatum) [38].

Natural waterfalls are among the most persistent barriers shaping the population genetics of freshwater fishes. Although they occur less frequently than human-made dams, waterfalls can strongly influence gene flow by restricting upstream dispersal while sometimes permitting downstream dispersal [4] [8] [13] [20] [27] [29]-[31] [39]-[42]. The asymmetric nature of passage beyond waterfalls produces predictable patterns of genetic connectivity and diversity [20] [23] [43]-[45]. Genetic diversity is predicted to be higher below waterfalls because above-falls populations are relatively isolated while below-falls populations may freely exchange individuals within and among river basins and receive an occasional migrant over the falls [3] [18] [31] [41] [45]-[47]. The main expectation is that above-waterfall populations function as small, partly isolated demes prone to genetic drift. Without the replenishment of alleles from downstream sources, heterozygosity and allelic richness should decline over time, and differentiation increases relative to downstream populations. These predictions have been confirmed empirically: in guppies in Trinidad [45], and for Rainbow Trout in the Russian River of CA [28].

Similar to waterfalls and other actions that fragment habitats, the evolution of watersheds and drainages can have profound effects upon the evolutionary divergence of fishes, particularly those that are non-migratory [48]-[50]. Headwater, non-migratory fishes may be effectively isolated from conspecifics that live in different tributaries of the same river drainage or in other drainages [50] [51]. The result is that downstream habitats may serve as effective or leaky barriers to geneflow among tributaries or drainages [8] [51] [52]. The result of riverscape heterogeneity is that the larger the distance among populations and the longer the period of their separation, then the greater we should expect their genetic divergence, such as has been witnessed in darters [53], a leuciscid minnow [49], a characiform tetra [54], the Brown Trout [55], and galaxiids [56]. However, this does not necessarily imply an isolation by distance framework because the genetic divergence will depend on the effective population sizes, habitat diversity within the drainage, normal non-migratory movements of the species, etc.

Rhinichthys atratulus, the Eastern Blacknose Dace (Figure 1), is a small, non-migratory, riffle-dwelling minnow found commonly in Atlantic draining streams east of the Appalachian Mountains from Nova Scotia to Virginia and in southeastern tributaries of Lake Ontario [57] [58]. It is an obligate freshwater fish that plays a critical role in trophic interactions of stream invertebrates, feeding on the larval stages of aquatic insects [59]. R. atratulus prefers to live in rocky, riffle habitats with moderate flow, a habitat that occurs above and below waterfalls in our study area. They can also be in pools or below overhanging banks immediately adjacent to riffles [58].

Figure 1. Rhinichthys atratulus, Eastern Blacknose Dace, 65 mm standard length, breeding male, from Laurel Brook adjacent to the Coginchaug River, CT. Photo by Barry Chernoff and A. Machado-Allison.

We examine three specific hypotheses about the effects of large barriers on the effects of migration and genetic differentiation among populations of R. atratulus: i) That waterfalls or drainage divides are barriers to gene flow; ii) That the amount of genetic divergence across these barriers should be greater than among populations not separated by physical barriers; iii) That genetic diversity below waterfalls should be greater than the genetic diversity above waterfalls. We hypothesize that populations above waterfalls will exhibit lower genetic diversity than those below, and that populations separated by barriers will be more genetically dissimilar than those not separated by physical barriers.

2. Materials and Methods

2.1. Sites

We sampled seven locations (Table 1) along the Coginchaug River and in Falls Brook (Figure 2). Both are north-flowing streams that are tributaries of the Mattabesset River which drains into the Connecticut River (Figure 2). All sampling localities were riffle habitats or the adjacent pools. Conductivity, pH, dissolved oxygen, and water temperature were measured using a Yellow Springs Instrument (YSI) probe Model 556 MPS. The conductivities in the Coginchaug River ranged from 0.231 - 0.243 ms (mean = 0.234 ms) [60]; in Falls Brook the conductivity ranged from 0.230 - 0.234 ms (mean = 0.232 ms). The dissolved oxygen was completely saturated at all localities, exceeding 9.5 ml/L. The pH in both the Coginchaug River and Falls Brook localities was slightly acidic, ranging from 6.82 to 6.84 [60].

Table 1. Physical descriptions, coordinates and sample sizes for the seven sampling sites in the Coginchaug and Falls Brook drainages. The data for substrate, flow and canopy cover at the Coginchaug localities are from [60] and this study. *The sample sizes are for both nd2 and microsatellites, except for nd2 in CR-A1 (n = 36) and FB-A (n = 22). Canopy refers to canopy cover.

The sampling localities were distributed above and below three waterfalls, two in the Coginchaug River and one in Falls Brook (Figure 2). Wadsworth “big” Falls (Figure 2, W1) is 10 m high; Wadsworth Little Falls is 5 m high (Figure 2, W2); and the waterfall at Falls Brook is 8 m high (Figure 2, W3). Both W1 and W2 empty into the lower Coginchaug River (Figure 2). The two sampling localities in the lower Coginchaug River (CR-B1, CR-B2) are ca.7 km below the outfall of both Wadsworth falls. We surveyed all of Wadsworth Brook above Little Falls (W2) and found only one location (Figure 2, CR-A3) that had a riffle-pool habit and our study species. A drainage divide, 46 m in elevation, separates Wadsworth Brook from the upper Coginchaug River even during periods of flooding. Falls Brook has only two localities, one above and one below (Figure 2, FB-A, FB-B) with habitats that supported our study species. We surveyed the entire upper part of the brook and found no riffle habitats; the brook ca. 20 m above FB-A had been channelized. Similarly, a weir/low dam has been installed in the lower part of Falls Brook near its confluence with the Mattabesset River; the impoundment comes within 350 m of the riffle with our study species below the falls (Figure 2, FB-B). This weir/dam is partially submerged allowing Falls Brook to flow into the Mattabesset River. There are two dams with impoundments on the lower Coginchaug River, the largest, Starr Millpond (41.650˚N, 72.663˚W) constructed ca. 1812, powered a mill with a “waterfall” of 6 m that would prevent any immigration of R. atratulus into the area of the study sites. Nonetheless, we have sampled the Mattabesset River, including the mouth of the Coginchaug for 16 years and have never encountered an individual of R. atratulus out of the more than 25,000 individuals captured in the region between the mouth of the Coginchaug River and near the mouth of Falls Brook (B. Chernoff, pers. comm.)

Figure 2. Map of study sites along the Coginchaug River, CT, USA. Black circles denote sample sites and bars denote waterfalls. CR-A1 and CR-A2 are above Wadsworth “big” Falls, CR-A3 is above Wadsworth Little Falls, CR-B1 and CR-B2 are below the falls, and the two Falls Brook sites are also separated by a waterfall. Map generated in ArcGIS 10.3.1.

The sites were recolonized approximately 18,500 years ago [61]-[63], after the retreat of the Wisconsinan glaciation. The time of isolation of fishes by waterfalls and drainage divides can be estimated from dates of late-Wisconsinan glacial varves [61]-[63]. The closest varve to the north was dated at 17,900 years ago. Thus, populations of R. atratulus were isolated by these barriers for about 18,000 years or approximately 9000 generations.

2.2. Samples

We collected 197 individuals of R. atratulus (Figure 1) from the seven sampling localities. Almost all individuals were collected between 2012 and 2014. A few individuals from CR-B1 and CR-B2, had been collected in 2006 and 2008, respectively. There were no among-year genetic differences found within localities that had been sampled over time (p > 0.10). We collected R. atratulus using a Smith Root backpack electroshocker and dip nets. A portion of the caudal fins were clipped and stored in 95% ethanol. Scissors were cleaned with 95% ethanol and wiped clean with a Kimwipe.

All of the procedures involving the collection and handling of fishes were humane and ethical, and were approved by state and animal care and use committees. Our study operated under scientific collecting permits from the state of Connecticut: Connecticut Scientific Collection Permits SC-13023 and SC-17031. Permissions for handling fish for scientific study were governed by: IACUC 2006-1212, Chernoff IACUC 2008-1212, Chernoff IACUC 2012-1212-Chernoff-A and IACUC 2014-1212-Chernoff-A.

2.3. DNA Extraction

DNA was extracted from fin clips using a QIAGEN DNeasy Blood and Tissue Kit: QIAGEN Sciences, MD, USA. The provided protocol for “Purification of Total DNA from Animal Tissues” was followed except for the last step, in which 200 μl of buffer was added to the membrane, and the samples were incubated at room temperature for 10 min before being centrifuged in order to increase DNA concentration. Final DNA concentrations for all samples were determined on a Thermo Scientific NanoDrop ND-2000 1-position spectrophotometer.

2.4. Amplification and Analysis of nd2

We analyzed 191 individuals for the nd2 mitochondrial gene. Nd2 was amplified using Polymerase Chain Reaction (PCR) with reaction conditions from [63]. Seven μl of PCR product mixed with 1μl of Loading Dye was run on a 1.5% agarose gel at 100 V for 30 min with 5 μl of SYBRsafe (Invitrogen). Successfully amplified samples that displayed a clear band at ~1045 base pairs (bp) were sent to Yale University’s DNA Analysis Facility for sequencing. Forward and reverse sequences were aligned using ClustalW Multiple alignment in Bioedit v7.2.5, and consensus sequences were generated based on chromatograms viewed in FinchTV v1.4.0 [64].

DnaSP v5.10.01 [65] was used to calculate haplotype counts, nucleotide diversity (π), haplotype diversity (Hd), and to generate input files for Arlequin. We used Arlequin v3.5.2.2 [66] for the following tests of neutrality: Tajima’s D, Fu’s F_s, Fu and Li’s D*, and Fu and Li’s F* [67] [68]. Arlequin was used to determine the number of private haplotypes, which we define as population-specific haplotypes that have not been found in any other population (i.e., endemic). We acknowledge that if sample sizes were to get very large, the probability of finding a population-specific or “private” haplotype in another population increases in proportion to the degree of gene flow among populations [69]. Arlequin was also used to generate Analyses of Molecular Variance (AMOVA). An AMOVA partitions molecular variance and estimates genetic structure based on the allelic contents of haplotypes and their frequencies [70]. AMOVAs were generated for populations by sample site, for groups of populations above or below waterfalls, and for groups of populations separated by drainage divides. We did not correct the table of p-values for the mitochondrial and microsatellite AMOVAs (Table 3) using the sequential Bonferroni procedure [71]. Sequential Bonferroni is too conservative and leads to too many Type II errors [72]. Thus, the Bernoulli equation indicated that the probabilities for finding the number of significant results (p < 0.05) out of the number of p-values generated for the nd2 and microsatellite data were 4.13 × 10⁻²³ and 5.08 × 10⁻³⁴, respectively [72].

A Minimum Spanning Tree estimation of connection lengths between haplotypes was generated in Arlequin v3.5.2.2 to estimate genetic relationships among the haplotypes [63] [73]. The branching diagram also shows the percent of individuals displaying that haplotype from each population. This network was rooted at Haplotype 1 because phylogenetic analysis identified Haplotype 1¹ as the ancestral haplotype of populations that recolonized New England after the Wisconsinan glaciation. Thus, all of the other haplotypes are derived from Haplotype 1.

2.5. Amplification and Analysis of Nuclear Microsatellites

We analyzed 197 individuals for microsatellites. The DNA was screened for 15 microsatellite loci. Seven loci were polymorphic and were used in the analyses. Microsatellites were amplified using PCR with primers and protocol from [58]. Successfully amplified fragments displayed a fluorescent band at the expected length that was visualized with gel electrophoresis and then submitted for fragment analysis by the Yale DNA analysis center. All fragments were scored manually for length using Peak Scanner v2.0 (Applied Biosystems) and binned in Microsoft Excel.

Input files were generated from Excel using CONVERT v1.3.1 [74] and CREATE v1.37 [75]. We used POPGENE v1.32 [76] to determine private, or endemic, alleles at each site. A pooled-variance test along with a test of homogenous subsets were conducted to determine the significance between the mean number of alleles found at each site. FSTAT [77] was used to determine the allelic richness (A_r). A_r is a calculation of the mean number of alleles per site that accounts for the variation in sample size by standardizing the sample size to the smallest n collected. A pooled-variance test along with a test of homogenous subsets were also calculated for A_r. Arlequin v3.5.2.2 [66] generated AMOVAs and estimates of Hardy-Weinberg equilibria for the populations.

The Garza-Williamson index, generated in Arlequin, is a statistic to detect recent population bottlenecks [78]. The index uses the ratio M, the mean ratio of the number of alleles to all possible alleles in the range of allelic sizes in order to determine whether populations are in decline. For datasets with seven or more loci, a value of M < 0.68 denotes that a population has gone through a bottleneck, while M > 0.80 indicates that there has been no reduction in the effective population size. The modified Garza-Williamson index [66] adds the value 1 to the denominator of the original index to prevent division by zero in the case of a monomorphic population.

STRUCTURE v2.3.4 was used to generate Bayesian analyses to infer population structure and probabilistically assign a posteriorly individuals to genetic populations [79]. The analysis assumes that there is but a single genetic population; the site from which each individual was captured was not used in the analysis. Under the admixture model, 20 runs per K (# of populations) were generated (K = 1 through 7) with 10,000 burn-in and 10,000 generations. The models each have a maximum likelihood. The likelihoods were tested sequentially with the log likelihood ratio test. We selected the model with the largest number of partitions for which the log likelihood ratio test was significant (p < 0.01) for χ² with 1 degree of freedom [80]. This result was checked by using the likelihoods from Structure to calculate the BIC indices for each value of K [81] [82]. The best fit model has the lowest mean value of BIC. The output from STRUCTURE displays a visualization of estimated population genetic structure by site [79]. The program classifies individuals into color groups, each of which represents a different genetic population (unfortunately, the user has no control of the colors that program assigns to groups).

2.6. Fst Analyses of Divergence

Arlequin v3.5.2.2 [66] was used to calculate the pairwise divergences among populations for both nd2 and microsatellites (Table A2). Pairwise Fst tile plots were generated in R version 4.3.3 [83] using the packages data.table version 1.17.8 [84] and ggplot2 version 3.5.2 [85]. The Kruskal-Wallis test evaluated the significance of pairwise divergences (Fst’s) between populations separated by waterfalls or by drainage divides. The Kruskal-Wallis non-parametric test makes assumptions about normality of the distribution of Fst values. Box-whisker plots illustrate the differences between the treatments for waterfalls and drainages. Because the Fst values were used in multiple tests, we would have applied the Bernoulli method [72] as discussed above to determine the probabilities for finding the number of significant results for the number of comparisons for the Kruskal-Wallis tests of the Fst data. However, the results were not significant and the test was not needed.

2.7. Estimations of Migration and Mutation Rates

Migrate v3.6.8 [86] was used to estimate the mutation rate scaled by population size (θ) and the mutation-scaled effective immigration rate (M). M indicates how much more important immigration is than mutation in bringing genetic variants into a population [87]. The number of individual migrants between populations per generation is calculated from the product of M by θ [87]. The estimates of these parameters were calculated separately for nd2 and microsatellites.

The parameters θ, M and the number of migrants per population are calculated between pairs of populations as part of a network of connections among sampled populations [87]. For example, a hypothesis of panmixia includes connections among all populations. We tested a hydrologically constrained model with 13 parameters against panmixia with 42 parameters using the log-likelihood ratio test. The model of connections is shown in Table 7. This model is most reasonable from a biological perspective as well because a species of the size of R. atratulus should not be able to travel upstream over waterfalls [13] [15] [20] [88]. Differences in the means and coefficients of variation of θ between nd2 and microsatellites were tested in Past v4.03 [89] using 9999 bootstraps with a Monte Carlo procedure.

For analyses of microsatellites, PGDSpider v2.1.1.3 [90] was used to convert Arlequin files to Migrate input files. For analyses of nd2 and microsatellite data, Migrate v3.6.8 was run ten times. Marginal likelihood estimations were determined based on the log-probability of the data given the model to find the most likely migration scenario of the runs.

For microsatellites, SPAGeDi software [91] was used to compare Rst and Fst values. Rst measures divergence based upon a within population model rather than Fst which is based upon an immigration model [92]. If Rst = Fst, then there is no contribution of mutation to genetic differentiation. If Rst > Fst, then mutations are the main contributor to genetic differentiation of populations [91].

2.8. Isolation by Distance (IBD)

Distances among populations were river distances calculated from the connections among the rivers. For example, the distance from the population above Wadsworth Little Falls (CR-A3) to other populations includes the stream distance from CR-A3 to the Coginchaug River and then the distances in the Coginchaug to other populations. A reduced major axis (RMA) regression of log-transformed Fst genetic distance on geographic distance was calculated in PAST v4.03 [89]. Mantel tests were also conducted using PAST v4.03 [89] to compare pairwise F_ST values to geographic distances between populations.

2.9. Data Accessibility

All nd2 sequences were registered in GenBank: MH368787-MH368800. All microsatellite alleles were registered in GenBank: MH368801-MH368807.

3. Results

3.1. Relationships among Haplotypes

There were 14 haplotypes (Table A1) identified from the 191 individuals from the seven populations (Figure 3). All haplotypes are derivatives of Haplotype 1 (Figure 3) which has been identified as the ancestral haplotype of populations that recolonized New England after the Wisconsinan glaciation [63]. Of the eight private haplotypes observed, six were found in two sites (CR-A1 and CR-A2) above the large waterfall (W1) in the Coginchaug River. Only two haplotypes were observed below W1 (Table 2, CR-B1 and FB-B). Haplotype 1 was present in all Coginchaug populations and most frequent in R. atratulus above W1 (Figure 3). Haplotype 1 was absent from the Falls Brook localities (Figure 3). Although Haplotype 1 occurred in low frequencies in the two Coginchaug populations below W1 (CR-B1 & CR-B2; Figure 3), these two populations uniquely possessed Haplotypes 9, 12 and 42, which are derivatives of Haplotype 1. Although CR-B2 only had 9 individuals, it still had the same the number of haplotypes, 5, as the largest sample, CR-A1 (n = 36), but lacked private haplotypes which could be due to the small sample size. Haplotype 10 was the most common but was absent from CR-B2. Haplotype 10 has also been found in the Thames drainage basin of Eastern Connecticut [63]. The populations, CR-A3 and FB-A, above the two smaller waterfalls, W2 and W3, respectively, lacked any of the derivative haplotypes (Figure 3).

Figure 3. 95% maximum parsimony haplotype network displaying relationships among nd2 haplotypes. Each circle represents a haplotype, with each color representing a different population. The size of each circle represents the relative frequency. Each line between haplotypes represents one nucleotide difference between haplotypes.

Table 2. Population statistics for mitochondrial nd2, including haplotype diversity (Hd), nucleotide diversity (π), the average number of nucleotide differences between haplotypes (k), and tests of neutrality for Tajima’s D, Fu’s F_S, Fu and Li’s D*, and Fu’s and Li’s F*. None of the tests were statistically significant (p > 0.05).

3.2. Genetic Diversity of nd2

The average number of nucleotide differences between haplotypes above the falls ranged from 0.000 to 1.432 and below the falls from 0.133 to 2.556 (Table 2). In above-falls sites, mean haplotype diversity (Hd) was 0.341 (range: 0.000 ≤ Hd ≤ 0.656), and mean nucleotide diversity (π) was 0.001 (0.000 ≤ π ≤ 0.001). In below-falls sites, the average Hd was 0.553 (0.113 ≤ Hd ≤ 0.861), and the mean π was 0.0013 (0.000 ≤ π ≤ 0.002). Overall, Hd, π, and the average number of mutations between haplotypes were higher in populations below the falls than above the falls. Because the population at FB-B had one common and one rare haplotype, Hd is approximately 0. The same value of Hd was found at FB-A where there was one haplotype.

The results of Tajima’s D, Fu’s F_s, Fu and Li’s D*, and Fu and Li’s F* for any population failed to reject the null hypothesis of neutrality (p > 0.05), suggesting that nd2 was evolving as a neutral genetic marker (Table 2). The results of these tests also indicated that populations had probably not experienced a recent bottleneck or inbreeding depression [93].

Global AMOVAs (Table 3) were computed to compare sites separated by and not separated by barriers to geneflow (i.e., waterfalls and drainage divides). In all cases, the largest percentage of the variation was explained by haplotypic variation among individuals within populations, ranging from 31.5% - 97.7% (Table 3). There were significant differences between populations that were not separated by barriers whether above the falls (Table 3: CR-A1 vs. CR-A2, p < 0.001) or below the falls (Table 3: CR-B1 vs. CR-B2, p < 0.01). The latter comparison could be the result of small sample size in CR-B2. Of the three over-the-falls comparisons, only that across Wadsworth Little Falls (W2) was significant (p < 0.001, Table 3). The comparisons of populations over Wadsworth Big Falls (W1) and over Falls Brook (W3) were not significant (p > 0.05). The comparisons across drainages were all significantly different (p < 0.01, Table 3).

Table 3. Analyses of molecular variance (AMOVA) for nd2 and microsatellites among populations. Populations inside parentheses are pooled for the particular comparison. The percentage variation among populations within groups is only calculated when there is more than one population within a group. The percentage variation within individuals is the variance explained by heterozygosity and only applies to microsatellite data. N/A means not applicable. Significance levels for p-values are: p < 0.05*, p < 0.01**, p < 0.001***.

The pairwise population Fst values (Figure 4) reveal that the populations above waterfalls in the Coginchaug River (CR-A1-A3) are most similar to each other, as are the populations within Falls Brook (FB1-2). With the exception of the population above Wadsworth Little Falls (CR-A3), the populations from Coginchaug River are most highly diverged from those populations in Falls Brook (Figure 5). Otherwise, there is a patchwork of divergences. The Fst values indicated that there were no significant divergences between populations that were separated by waterfalls (KW H = 0.081, p = 0.78, Figure 6) or by drainage divides (KW H = 0.679, p = 0.09, Figure 5).

Figure 4. Heatmap of pairwise population Fst values for nd2 (A) and microsatellites (B). Colors are scaled so that larger Fst values are darker. Abbreviations: CR is for Coginchaug River; and FB is for Falls Brook; see Figure 2 for location of sampling sites.

Figure 5. Boxplots of Fst values (ordinate) for populations not separated (1), and separated (2) by waterfalls ((A), (B)) and by drainages ((C), (D)). The data are for nd2 ((A), (C)) and for microsatellites ((B), (D)). The boxplots have the following components: i) The box represents the proportion of data that lie between the 25^th and 75^th quartiles; (ii) The horizontal bar within the box is the median; (iii) The vertical line spans the maximum and minimum values.

3.3. Genetic Diversity of Microsatellites

There were 125 alleles from seven polymorphic microsatellite loci (Table A3). The variability of the number of alleles within populations was quite high (Table 4), the standard deviations ranging from 4.2 to 9.3; the standard deviation of the mean among all populations was 6.5. The means of the two upper Coginchaug populations and the lower Falls Brook populations exceeded nine alleles per population (Table 4) and the mean among all populations was greater than eight. The number of private alleles varied significantly among populations (Table 4), with the greatest being found in the two upper Coginchaug populations and in two populations below the falls, CR-B1 and FB-B (Table 4). The number of alleles and private alleles exceeded those published for other populations of this species [58] [94]. There were no significant differences between the observed and expected heterozygosities within populations (p > 0.05, Table 4). The lack of departure from Hardy-Weinberg equilibrium indicated failure to reject null hypotheses that the microsatellites were evolving under random (neutral) processes. Thus, we considered these microsatellite loci as neutral markers.

Table 4. Population statistics for microsatellite data of R. atratulus, including the mean number of alleles, allelic richness, and the number of private, or population-specific, alleles. Also listed are the observed and expected heterozygosity of the seven polymorphic loci for each of the seven populations. None of the results were statistically significant, indicating no significant departure from Hardy-Weinberg equilibrium.

For all comparisons, the within-individual variation comprised the largest percentage of the total variance (Table 3). The percentage explained ranged from 84.5% - 99.7% of the total variance and was due to heterozygosity. These values were all very highly significant (p < 0.001).

All the among-group comparisons were significant if not highly significant (p < 0.001): i) Between adjacent populations with no barriers; ii) Between above-falls and below-falls populations; iii) Between populations separated by drainage divides (Table 3). Because there were only two comparisons between adjacent populations (Table 3), we note the results from [94]. They studied three populations of R. atratulus in a tributary of the Coginchaug River that were not separated by barriers and were very close to each other (maximum distance spanned 470 m) with the same genetic markers used in this study. The same result obtained: all three populations were significantly different (p < 0.01) from each other (Table 3) [94].

The log likelihood ratio test and the lowest mean BIC indices indicated that the Bayesian a posteriori classification of individuals best fit a model with k = 5 partitions (i.e., genetic groups; Figure 6). Within the mainstem of the Coginchaug River, there were differences between the upper two populations, CR-A1 and CR-A2. Individuals within the uppermost population, CR-A1, belonged primarily to a single genetic population (yellow) that was infrequent in CR-A2 and the Falls Brook populations. The population above Little Falls, CR-A3, was almost fixed for a single genetic signature that is present in ca. 25% of CR-A2. Interestingly, the two populations below the big falls at Wadsworth, W1, comprise a majority of individuals that were genetically distinct from all other populations (green). The Falls Brook populations were classified primarily into a genetic group (red, Figure 6) that had a low probability of grouping with other populations. The populations at FB-A had a number of individuals classified into the yellow group that was characteristic of CR-A1.

Figure 6. Bayesian a posteriori classification of individuals into genetic groups by the program STRUCTURE. The analysis was based upon seven microsatellite loci for 197 individuals. Each vertical bar represents an individual. The height of each color within a vertical bar (i.e., an individual specimen) is the probability of belonging to a particular genetic group that is represented by a color. The individuals are grouped by the geographic populations in which they were collected, see Figure 2 for localities. Abbreviations: CR is Coginchaug River; FB is Falls Brook.

The heatmap of Fst values (Figure 4) illustrates a patchwork of larger and smaller divergences but relative to nd2, the divergences are quite small. The most similar populations are those from Falls Brook. The most divergent Fst values are between CR-A3 and CR-B1 to FB-A and FB-B (Figure 4). There were no significant differences in pairwise Fst values for either populations separated by waterfalls or populations separated by drainage divides (KW H = 1.64, p = 0.20, and H = 0.714, p = 0.40, respectively; Figure 5).

All the populations had gone through recent bottlenecks or reductions in population size. The means of the Garza-Williamson statistic ranged from 0.340 - 0.552; the means of the modified index ranged from 0.219 - 0.394 (Table 5). Recent bottlenecks or reduction in population sizes are indicated when the value of either index is <0.68 [78] [95]. Both versions of the index revealed that the populations from the most restricted habitats, CR-A3, FB-A and FB-B, had the lowest values.

Table 5. Means and standard deviations of the Garza-Williamson index, M, and the Garza-Williamson Modified Index.

3.4. Mutation and Migration

The population-scaled mutation rates, θ, were very similar among populations (Table 6); with most populations differing by only tenths of a percent per generation. The fastest and slowest population-scaled mutation rates were found for nd2 (Table 6) for populations in the Coginchaug River. Although it appears that nd2 had higher values of θ than microsatellites, neither the means nor coefficients of variation of θ differed significantly (p > 0.1).

Table 6. The estimated mutation rates (θ) for nd2 and microsatellites for each of the seven sampled populations (Figure 1) of Rhinichthys atratulus.

The values of M indicate the importance of immigration over mutation for introducing genetic variation into populations scaled by population size. However, higher values of M may also be due to shared ancestral haplotypes and microsatellite alleles. For both data sources the estimated number of migrating individuals was relatively low, ranging from less than 1 to 3.2 individuals per generation (Table 7). That the species lives two to three years [58] [96]-[98] implied that one or fewer individuals were migrating per year.

The model that best fit the data was the geo-hydrologically-based biogeographic hypothesis that prevented individuals to migrate from sites below waterfalls to sites above waterfalls (e.g., from CR-B1 to CR-A2) or to cross directly over a drainage divide (e.g., from CR-A3 to CR-A2). Thus, in the Coginchaug River populations, upstream migration between adjacent populations was much lower than downstream migration on a per generation basis for both nd2 and microsatellites (Table 7; e.g., CR-A1 → CR-A2 > CR-A2 → CR-A1). The estimates of migrations from above to below the three waterfalls were relatively high (Table 7), which could imply the individuals going over the falls may survive.

Table 7. Estimation of migration rates for nd2 and microsatellite data. The first column shows the direction of migration between sites. The statistic M indicates the role of immigration over mutation in bringing in genetic variants to populations. The number of migrating individuals per generation is determined by multiplying the mutation-scaled effective population size, theta (θ), for the corresponding recipient population, times M.

Furthermore for microsatellites, we tested the importance of θ relative to migration as the source of genetic variation within populations, by comparing Rst to Fst values [91]. The result was that Rst > Fst (M.L, χ² = 7.84, p = 0.0037). This result supports the conclusion that mutations contributed significantly to microsatellite differentiation of the populations.

3.5. Isolation by Distance

IBD predicts a positive relationship between genetic and geographic distances (Wright 1943). RMA regression of Fst values from nd2 on geographic distance was highly significant (p < 0.01; Figure 7(A)). The Mantel correlation test on the matrix of nd2 Fst values and geographic distance was positive (R = 0.5640) and highly significant (p < 0.01). The results mean that the genetic differences among populations can be explained by an isolation by distance model. Both the RMA regression (Figure 7(B)) and the Mantel test calculated for Fst values from microsatellite data were not significant (p > 0.05).

Figure 7. Reduced major axis regressions of log-transformed pairwise F_ST values on geographic distances between populations. (A): nd2 haplotypes. (B): Microsatellite data.

4. Discussion

That specific habitats are patchy in nature is true for almost all ecosystems. The discontinuous distribution of habitats is especially characteristic of freshwater streams [99]. Our study raises questions about the relative roles of patchy habitat and major barriers, such as waterfalls or drainage divides, in shaping the genetic structures of populations. Many studies have shown that waterfalls serve as barriers that have large influences on population genetic structure [3] [40] [51] [100]-[104]. To the extent that geneflow occurs across barriers, it should occur in a downstream direction as predicted by riverscape theory [20] [40] [43] [44] [47] [104]-[106]. Many but not all of these studies involve salmonids (trouts, salmons, and their relatives) of which some species have the potential to migrate upstream over waterfalls (e.g., Steelhead Trout, Onchorynchus mykiss, [103]). For low-mobility species, the evidence is inconsistent whether waterfalls, rapids, or dams result in genetic differentiation among populations [6] [7] [14] [31] [33] [43].

4.1. Genetic Diversity over the Falls

Genetic diversity is expected to be greater in populations below falls than above falls [3] [18] [31] [41] [47]. This is because populations below waterfalls can receive immigrants from downstream and upstream populations; it is well documented that some fishes can survive the descent over a waterfall [47] [103] [107]-[109]. Our results provided some support for and against the hypothesis of greater genetic diversity in below-falls populations. The conflicting results are both a function of the data and how one measures “diversity”.

Support for the theory was found in both the Coginchaug and Falls Brook populations. The Coginchaug River drainage below-falls populations (CR-B1, CR-B2) had higher Hd and π than above falls populations (CR-A1, CR-A2, CR-A3; Table 2). The population below Falls Brook had more private microsatellite alleles than the above-falls population (Table 4). Hd was greater in the Falls Brook population below than above the falls. This was weak support at best, however, because FB-A is fixed for Haplotype 10 and FB-B had only two haplotypes (Figure 3, Table A1).

Contrary to theory, the above-falls populations in the main-stem of the Coginchaug River (CR-A1, A2) had more private haplotypes and a greater number of private microsatellite alleles than below-falls populations (Table 2, Table 4). The relatively small sample size of CR-B2 may have affected this result, though we note that the number of total haplotypes in (CR-B2), the smallest sample, is the same as the number of haplotypes in the largest sample (CR-A1). Furthermore, the number of alleles and allelic richness did not differ significantly (p > 0.05) between above- and below-falls populations in the Coginchaug and Falls Brook drainages (Table 2).

Not all studies have supported the theory that populations below significant barriers such as waterfalls have higher genetic diversity. For example, there was greater diversity of microsatellites in some above-barrier populations of Semotilus atromaculatus (Creek Chub) in the Mastigouche Wildlife Reserve, Quebec, Canada [20]. They theorized that because greater genetic diversity above impassable waterfalls could not be explained by migration, the diversity might be a result of older influences such as postglacial colonization or introductions. Alternatively, others have argued that waterfalls can create upstream habitat diversity that can, eventually, lead to genetic differentiation in above-falls populations [25] [43]. These suppositions identify three assumptions made by the theory that differentials in genetic diversity are not due to: (i) Extinction or reduction of ancestral genetic variation; (ii) In situ evolution; (iii) Anthropogenic introductions.

Our results lead us to argue that the genetic variation among populations in the Mattabesset River drainage has been largely shaped by in situ mutations rather than by current rates of migration or by maintenance of the genetic diversity of the recolonizing ancestral populations. The effects of both migration rates and the maintenance of ancestral genetic diversity are insufficient to explain the current population genetic structure among our study populations. A single population of R. atratulus recolonized modern-day Connecticut River drainage from a single refugium in which the founding population went through a subsequent bottleneck (e.g., Haplotype 1 was the ancestral nd2 sequence of the founding population) [63]. The legacy of the bottleneck would be low genetic diversity, not high endemic genetic diversity. In the Mattabesset drainage, 14 of the 20 haplotypes are endemic to the drainage (includes Allyn Brook [94]; and of the 14, 12 are private haplotypes found only in a single locality. Furthermore, the number of private microsatellite alleles is also high (e.g., = 40%). Thus, the diversity of haplotypes in the Mattabesset River drainage have evolved since recolonization approximately 18,500 years before present [61]-[63].

The very high percentages of private haplotypes and microsatellite alleles among the individual sites in the Mattabesset drainage are due to in situ mutation rather than migration. Our inference of in situ evolution is also supported by the fact that analyses failed to reject null hypotheses of neutrality for each of 51 nd2 haplotypes and 238 microsatellite alleles among all studies [58] [63] [94] (Table 2). And further, because Rst > Fst for microsatellite alleles implied that mutations are a greater source of genetic differentiation within the populations [91] [92].

The role of in situ evolution due to mutation has not generally been acknowledged [3] [20] [41]. Genetic diversity produced by mutation can produce results that can support or reject the theory based solely upon migration. Our results (discussed below) and those of [94] demonstrated that mutation rates have been sufficient to override homogenization effects of migration or ancestral lineage effects. Furthermore, high rates of mutation were calculated for mitochondrial DNA among the re-colonized populations in New England for R. atratulus [63].

In general, the introduction of organisms by humans to new localities is critical to consider in studies of population genetics [20] [33]. This is especially true for small fishes that may be used as bait. We do not believe that bait bucket introductions formed the basis for the genetic diversity among the populations we sampled. While R. atratulus can be used as bait for sportfishes such as freshwater basses, genus Micropterus, those predators are rare in the Coginchaug and Falls Brook sampling sites; in 22 years of surveys in these rivers we have never encountered people fishing (B. Chernoff, pers. obs.).

4.2. Waterfalls and Drainage Divides as Barriers to Gene Flow

Waterfalls and drainage divides are not only expected to serve as barriers to geneflow but also are expected to increase genetic diversity within river systems [25] [33] [51] [104] [106]. The expectations are that the nature and magnitude of genetic divergence across instream barriers should be greater than among populations not separated by such barriers [4] [7] [20] [40] [44] [103]. Our study on a non-migratory species tests this idea by examining populations above and below three waterfalls, and between two tributaries of the Mattabesset River.

Our results demonstrate that there is significant genetic divergence among all populations for microsatellites and most populations of nd2 (Table 3). The results mostly provide evidence that waterfalls and drainage divides serve as barriers to multi-directional geneflow. For nd2, the only over-the-falls comparison that was significant was that of the populations above (CR-A3) and below (CR-B1, B2) Wadsworth Little Falls (W2) (Table 3). In all cases the variation of nd2 within populations, whether pooled, was high and significant (Table 3). High levels of within-population variation generally lead to lower estimates of among-group variation [22] [110].

For microsatellite data, all comparisons of above versus below falls populations were significantly different (p’s < 0.01, Table 3). Furthermore, the genetic structures of microsatellites were significantly different (p’s < 0.01) for the populations separated by drainage divides (Table 3). These results are confirmed by the Bayesian analysis of the microsatellite data that illustrated the differences in genetic signatures across the barriers (Figure 4).

We also discovered that populations not separated by barriers in the upper and lower Coginchaug River were also significantly different in genetic structure for both data sources (Table 3, Figure 4). These results are identical to the findings of [94] for R. atratulus in Allyn Brook, a tributary of the upper Coginchaug River: populations at distances ranging from 80 - 470 m were genetically distinct from one another for the same neutral markers. Genetic structure among populations in other species has also been found in sections of rivers or streams unimpeded by barriers (e.g., Gadopsis marmoratus [14] [49] [110]-[112]. Similarly, the amphipod, Gammarus fossarum, expressed strong genetic structure among populations at spatial scales of a few kilometers [113]. All of these results point to the lack of geneflow among populations whether or not separated by barriers. The roles of mutation and random processes, such as genetic drift, should not be underestimated in the acquisition of neutral genetic novelties [114] [115].

The theoretical expectations governing genetic differentiation across major barriers in lotic systems [40] [41] [48] [104] [106] [112] [116] make no predictions about the genetic structure of populations that are not separated by barriers. That all populations separated by barriers in this study differed significantly in their genetic structures for at least one data source is consistent with expectations. Furthermore, the degree of genetic divergence among populations separated by major barriers is expected to be greater than the genetic divergence among populations that are unimpeded by barriers. The degree of divergence for neutral markers should increase with the age of separation by a barrier. The sampling localities of our study were opened to recolonization approximately 18,500 years ago based upon the dating of late-Wisconsinan glacial varves [61] [62]. Deglaciation proceeded from the coast to the interior. The closest varve to the north of the study area was dated at 17,900 years ago [61] [62]. Thus, the populations above waterfalls and between drainage divides have been isolated for approximately 18,000 years or approximately 9000 generations. This degree of isolation should result in genetic differences that should far exceed those of unimpeded populations. Our results for the nuclear and mitochondrial data were contrary to this expectation. F_ST values among populations were not significantly greater for populations separated by waterfalls or drainage divides than for those lacking barriers (ANOVA p’s > 0.10; Kruskal-Wallis p’s ≥ 0.2). All populations differed to the same degree indicating that geneflow among all populations is relatively low.

4.3. The Effects of Mutation and Migration

The genetic variation within Connecticut River populations is due to four main factors: i) The genetic variation of the post-glacial recolonizing ancestral population; ii) Mutation of novel haplotypes and alleles; iii) Migration bringing additional alleles and haplotypes; iv) Extinction. The populations that recolonized New England entered the deglaciation of the Connecticut River via an eastward flowing river in what is now Long Island Sound [61]-[63]. The founding population went through a bottleneck that reduced genetic diversity [63]. For example, these researchers demonstrated that Haplotype 1 is ancestral for New England. Though we are unable to estimate the diversity of the microsatellite alleles present in the re-colonizing ancestors, both the allelic richness and heterozygosity would have been decreased [117]. Thus, the processes of mutation and migration have produced a large majority of the genetic variation in the Mattabesset River drainage.

The population-scaled mutation rates, θ’s, were very similar among populations for both nd2 and microsatellites (Table 6). These θ’s can easily account for the genetic variations among the sites in the 18,000 years or 9000 generations of residence. While many widespread haplotypes (e.g., haplotypes 10) or microsatellite alleles probably evolved earlier in the history of these drainages, the private haplotypes and alleles are probably more recent in origin.

Our estimates of migration in the best fit model found general levels of concordance as predicted by mitochondrial or nuclear data (Table 7). The number of migrants between sites was predicted to vary from one to three per generation (Table 7). For neutral genetic markers, these values would have been sufficient to swamp genetic variation at the connected sites [118]-[120]. Thus, the large number of private alleles and haplotypes among our sampled populations indicates that the rate of mutation is able to overcome the homogenizing effects of migration. A conclusion that is supported by the fact that the Rst > Fst for microsatellites. This is also the result and conclusion reached for Allyn Brook populations of R. atratulus [94]. While it is possible that extinction could explain the pattern of genetic diversity for these populations, we point out that: i) The null hypotheses of neutrality were not refuted for haplotypes and microsatellite alleles; ii) The genetic effects due to migration or recolonization can be quite insensitive to the effects of extinction under certain theoretical models [120].

Given the geometry of the rivers, one would expect a stepping stone model and IBD to inhere [7] [25] [44] [106] [121]. The mitochondrial data were consistent with IBD (Figure 7(A)), and both the regression and Mantel tests were highly significant (p < 0.01). However, neither test was significant for the microsatellite data. IBD assumes that geographically proximal populations will exchange more genes from migrants and, therefore, will be more genetically similar than populations that are further apart [121] [122]. Potential evidence for geneflow among populations can be found in some of the haplotypes that are shared by a limited number of populations. For example, haplotypes 9, 12 and 44 were found only in the CR-B1 and CR-B2 below Wadsworth Big Falls, W1, (Figure 3, Table A1). Although the presence of Haplotype 14 in CR-A2 and CR-B2 would seem to provide evidence for geneflow over W1, it could also be explained by extinction in CR-A1 and CR-B1 (Figure 3, Table A1). If we failed to sample Haplotype 14 in CR-B1, then geneflow would be more parsimonious.

Our data indicate that migration between and among populations is very low, ranging from 1 - 3 individuals per generation (Table 7). Similar results were found in a small tributary of the Coginchaug River [94]. Mark and recapture studies in the Still River basin of western Connecticut demonstrated no movement of R. atratulus among riffles that were ≥600 m apart (K. Anatone, pers. obs.). Therefore, the positive relationship between genetic and geographic distance of nd2 (Figure 7(A)) was most likely due to the pattern of private haplotypes: relatively many and few at the two ends of the distribution (Coginchaug and Falls Brook, respectively). Unlike the mitochondrial data, the microsatellite data did not support a hypothesis of IBD (Figure 7(B)). The lack of congruence between nd2 and the microsatellites could be due to manner of inheritance and the reduced effective population size for mitochondrial data [123]. Nonetheless, the populations having gone through bottlenecks as indicated by the microsatellite data (Table 5), would also reduce the effective population size and increase the effects of intra-population mutation.

IBD can provide evidence for dispersal and the amount of gene flow between populations [7] [52] [69] [95] [121] [122]. However, IBD assumes that migration and geneflow are the major sources of genetic variation among populations. The effects of mutation, however, may differentially affect populations of species with differing life histories and ecologies. This may be especially true for habitat specialists or species with strong preference for a habitat. R. atratulus strongly prefers riffle habitats and adjacent pools. This results in limited migration among riffle habitats [94].

The predicted migration rate per generation (1 - 3.2 individuals, Table 7) could maintain heterozygosity and minimize the loss of polymorphisms in populations [124], but has not been sufficient to homogenize populations. The estimation of individual migrants between populations is partly based upon shared haplotypes and alleles. Thus, the actual number of migrants may be lower than predicted by our analysis because the analysis does not distinguish between contemporary and longer-ago historic processes.

4.4. Habitat and Predatory Effects on Migration

There are listed six potential reasons why R. atratulus was not migrating or migrating minimally in Allyn Brook [94]. Of these, we believe that two are relevant to this study: i) ecological habitat preference and the patchiness of that habitat; and ii) the effects of predation. Some fish populations have been found to differ genetically in the absence of physical barriers [13] [52] [94]. The lack of migration between populations may be due to subtle riverine landscape features. For instance, the depths and lengths of corridors between pool habitats, or riffles, have been found to affect rates of fish migration [1] [125] [126]. As R. atratulus is a riffle dwelling specialist, its movement may be restricted by long stream runs and deep pools that separate riffle habitats because the runs and pools may contain predators and may lack preferred food resources.

Inter-riffle movements may be low for riffle-dwelling specialists because of predation. Based upon 16 years of sampling, there are few large predatory fish at our study locations (B. Chernoff, pers. obs.); but the sites have a relatively open canopy, allowing for predation by birds. Since R. atratulus is often the most common species in our streams, it is preyed upon by many fish-eating birds [98]. Great Blue Herons and Belted Kingfishers have been found to alter the abundance of common prey and affect aquatic food webs [127]. Great Blue Herons, Green Herons, Kingfishers and Ospreys have been seen frequenting the study areas. In addition to bird predation of mature dace, larvae may also be preyed upon. R. atratulus bury their eggs during spawning [96]. If the eggs are swept downstream by the current or some perturbation, then fish that feed by filtering zooplankton may consume the eggs or larvae by chance. Thus, the probability of risk from predation would also diminish migration and could further affect their population genetic structure.

5. Conclusions

Ecosystems are characterized by patchworks of habitats. Yet patchiness and discontinuities have important implications for evolution and maintenance of populations within species [33] [128]-[131]. The more impenetrable the barriers that exist between habitat patches, the larger the degree of separation among populations on both sides of the barriers. Over time, populations separated by barriers such as waterfalls should manifest in genetic divergence [1] [4]-[6] [22] [25] [40] [44] [109].

Our work examined three specific hypotheses about the effects of large barriers on the effects of migration and genetic differentiation among populations of a riffle-dwelling minnow, R. atratulus in the northeastern United States and found that: i) The waterfalls and drainage divides were associated with genetic divergence among populations, indicating barriers to migration and gene flow; ii) All populations were significantly different from one another for microsatellites and almost all for nd2, but that the magnitude of divergence across barriers was not greater than between adjacent populations that were not separated by physical barriers; iii) Genetic diversity was higher above falls for private alleles and haplotypes, but genetic diversity was higher below falls for haplotype (Hd) and nucleotide diversity (π).

The latter two findings indicate that migration of individuals among populations with subsequent gene flow may not be the prevalent factor determining genetic variation at the level of the population. Few studies consider that the rate of in situ mutation may be driving the differentiation among populations [14] [33] [94] [113]. The roles of mutation and random processes, such as genetic drift, should not be underestimated in the acquisition of neutral genetic novelties [43] [95] [109] [114] [115] [120] [132].

Other confounding effects are time of isolation across barriers and extinction rates. The isolation across waterfalls and drainage divides in our study area is at least 18,000 years old [61]-[63]. Is the presence of shared alleles and haplotypes across the major boundaries due to contemporary migration and gene flow or due to more historic influences such as the post-glacial colonization of the rivers by genetically polytypic ancestors? Based upon the morphology of the Coginchaug River, one might expect a steppingstone model of evolution with isolation by distance [7] [41] [44] [106] [109] [133]. Some of our results were in conflict, because IBD was demonstrated for nd2 but not for microsatellites. We concluded that our populations were probably not evolving by IBD because the pattern was driven by the larger number of private haplotypes in the upper portion of the Coginchaug River.

The divergence of all populations with low migration, 1 - 3 individuals per generation, is most likely due to mutations of neutral genetic markers [43] [95] [109] [114] [115] [120] [132]. The effects of mutation will certainly be greater in populations with low effective population sizes or those that have recently gone through bottlenecks, such as the populations in the Mattabesset River drainage. In conclusion, mutation, migration, extinction, and time of colonization all play important roles in the patterns of genetic divergence among populations within species. The strength of influence of these factors will depend upon the geological history of the region, the frequency of habitat patches, the presence of major barriers to gene flow, and in the case of aquatic systems, the pattern of connectivity in the river basin. Thus, habitat fragmentation and degradation remain a threat to populations, as all these factors are critical for the evolution and maintenance of genetic diversity in wild populations.

Acknowledgements

We are grateful to Joel LaBella, Laurie Kenney, and Suzanne Bussolari for all the help with field gear, arranging and accounting for supplies, shipping samples, etc., that we so depend upon. CN, SL, LB, AM, DM, and NN received funding from undergraduate internship programs provided by the Essel and Menakka Bailey Experience Funds of the Bailey College of the Environment, and from the College of Integrated Sciences. BC received Schumann Funds from Biology Department and College of the Environment, Smith Funds from Earth and Environmental Sciences Department, and project grants from Academic Affairs, Wesleyan University, to undertake this study. We are grateful to F. Cohan for discussions, and to A. Machado-Allison and, particularly, to an anonymous reviewer for comments on the manuscript.

Authors’ Contributions

Kayla Anatone: Conceptualization (supporting); project administration (equal); fieldwork (equal); data curation (equal); formal analysis (equal); methodology (equal); writing—original draft (equal); writing—review & editing (supporting). Laura Bither: fieldwork (equal); data curation (equal); formal analysis (equal); methodology (equal); writing—original draft (equal); writing—review & editing (supporting). Nola Neri: fieldwork (equal); data curation (equal); formal analysis (supporting); methodology (supporting). Sage Loomis: fieldwork (equal); data curation (equal); formal analysis (supporting); methodology (supporting). Michelle L. Kraczkowski: Conceptualization (supporting); project administration (equal); fieldwork (equal); data curation (equal); formal analysis (equal); methodology (equal); writing—review and editing (supporting). Timothy S. Earley: methodology (supporting); formal analysis (supporting); graphics (equal); writing—review and editing (supporting). Barry Chernoff: Conceptualization (lead); funding acquisition (lead); project administration (equal); fieldwork (equal); data curation (equal); methodology (equal); formal analysis (equal); methodology (equal); writing—original draft (equal); writing—review and editing (lead).

Appendices

Table A1. Haplotype frequencies found in the sampled populations of Rhinichthys atratulus. Haplotypes 1, 9-12 and 14 correspond to haplotypes A, I, J, K, L, and N of Tipton et al. (2011) [63].

Table A2. Genetic distances, F_ST values, between localities for nd2, above diagonal, and microsatellites, below diagonal. See Table 1 and Figure 2 for locality information.

Table A3. The allele frequencies for the seven polymorphic microsatellite loci found in the seven study locations.

NOTES

¹Labelled as Haplotype A by Tipton et al. (2011).

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References