Phenotypic Variation in Brook Trout Salvelinus fontinalis (Mitchill) at Broad Spatial Scales Makes Morphology an Insufficient Basis for Taxonomic Reclassification of the Species

We compare data presented in Stauffer (2020) to: 1) data presented in Stauffer and King (2014) that were used to describe the neotype of S. fontinalis from Long Island, New York (Fig. 1A), 2) independent meristic counts from 38 Brook Trout populations located in the GSMNP described by Weathers et al. (2019; Fig. 1B), and 3) trait data reported in peer-reviewed literature and governmental agency reports. We provide a more detailed summary of these three data sources in the Supplemental Appendix (see Data Accessibility).

We direct readers to Stauffer (2020) for a description of the author's sample collection and data analysis. Briefly, Stauffer (2020) measured 23 morphometric and 10 meristic traits for ten fish from each of Cosby, Indian Camp, and Greenbrier Creeks within the GSMNP. These three collection locations are all located within the same 12-digit Hydrologic Unit Code (HUC 12) watershed and are separated by a waterway distance of no more than 19 km (Fig. 1B).

We contacted the author, who provided the original, individual-level data used in Stauffer (2020) and the individual-level data from Stauffer and King (2014) from populations in New York. However, analysis of the individual-level data failed to reproduce many of the summary statistics reported in the data tables presented in Stauffer (2020), and the number of vertebrae and basihyal teeth were absent. There were also inconsistencies between the individual-level data from populations in Long Island, New York and data presented in Stauffer and King (2014). Accordingly, we report two sets of values when summarizing data from Stauffer (2020) and Stauffer and King (2014)—one that provides the data as reported in the respective manuscript and another set of values calculated from the individual-level data provided by the author.

We also include data from Weathers et al. (2019), which described the morphometric and meristic variation of Brook Trout from 35 streams in the GSMNP (Table S1; see Data Accessibility). Although not included in the original publication, Weathers et al. (2019) also independently collected data from the three Stauffer (2020) sites (Cosby, Indian Camp, and Greenbrier Creeks; Fig. 1B) that we included in our analyses for a total of 38 populations from the GSMNP. Because Weathers et al. (2019) included an independent analysis from more sites, we can compare Weathers et al. (2019) and Stauffer (2020) to infer the amount of within- and between-stream variation in trait values relative to other streams that occur within the same region (i.e., HUC 4 watershed). Notably, though both studies included morphometric and meristic data, we limit our comparisons to meristic data because morphometric data were not indicated by Stauffer (2020) to be diagnostic, and because there is less subjectivity in collection methods for meristic data.

Stauffer and King (2014) and Weathers et al. (2019) did not include vertebral or basihyal tooth counts in their studies. To understand how vertebral counts compare to other Brook Trout species accounts, we compared the summary data presented in Stauffer (2020) to published species descriptions and experimental studies (sources provided in the Table S2; see Data Accessibility). Vertebral counts are infrequently reported in species accounts, and so the spatial extent of our comparisons is limited with most reported values originating from populations that occur at more northern latitudes. We were unable to locate independent records for basihyal tooth counts and discuss the limitation of this variable in more detail below.

To evaluate claims of unique meristic identity in the proposed new species of Salvelinus, we quantified patterns of phenotypic variation among putative species using the ten meristic traits measured in Stauffer and King (2014), Weathers et al. (2019), and Stauffer (2020; see original publications for variable descriptions). Although specimens from all localities outside of Cosby, Greenbrier, and Indian Camp Creeks would currently be classified as S. fontinalis under the delimitation proposed by Stauffer (2020), we expected broad differences in phenotype between New York and GSMNP populations due to the separation of northern Atlantic and southern Atlantic populations into distinct phylogenetic assemblages (Danzmann et al., 1998; Stauffer and King, 2014). Accordingly, we treated samples from New York and the GSMNP as two separate populations of S. fontinalis.

We used a principal components analysis (PCA) to make multivariate comparisons of meristic traits among populations and calculated pairwise Mahalanobis distances, which are unitless measurements that summarize the differences between group means along all axes of the PCA. To evaluate the ability of meristic traits to correctly classify species of Salvelinus into distinct categories, we conducted a cross-validation linear discriminant analysis (LDA) using the MASS package (Venables and Ripley, 2002) in R v. 4.0.2 (R Core Team, 2020). For these analyses, we used a classification model based on five categories, which included each of three populations surveyed by Stauffer (2020) from Cosby, Greenbrier, and Indian Camp Creeks, S. fontinalis from other watersheds in the GSMNP, and S. fontinalis from New York. For this analysis, we combined data collected by Weathers et al. (2019) with data from Stauffer (2020) for Cosby, Indian Camp, and Greenbrier Creeks. Our assumption was that if fish from the three study streams did represent new species as proposed by Stauffer (2020), then Weathers et al. (2019) would have inadvertently collected data from the three new species of Salvelinus, rather than Brook Trout. If so, adding the additional data from Weathers et al. (2019) would have either no effect, or potentially increase the power to detect significant differences among populations through increased sample size.

To account for the potential confounding effects of regional population structure, which could have occurred when comparing the individual streams from Stauffer (2020) to the collective GSMNP group that included the 35 other sites surveyed by Weathers et al. (2019), we reran this analysis and treated each stream as an independent factor (i.e., 41 total streams including 38 from the GSMNP and 3 from Long Island, New York). Inferences from the two analyses were congruent, and, for concision, we report the results of the individual stream analysis in the supplement (Table S3; see Data Accessibility).

Principal components analysis of meristic traits showed a clear distinction between S. fontinalis in New York and all southern Appalachian populations from the GSMNP on the first axis of the PCA (Fig. 2). This axis, which explained 22.3% of the total variation, was most closely associated with parr marks (variable loading = –0.31), number of teeth on lower left jaw (0.53), and pored scales posterior to the lateral-line scales (0.44). There was considerable overlap on PC1 and PC2 among all GSMNP populations, including the three sites surveyed by Stauffer (2020). In accordance with this finding, pairwise Mahalanobis distances were smaller among GSMNP populations than were Mahalanobis distances between GSMNP populations and populations from New York (Table 1).

Results from the LDA also suggested that the three populations described in Stauffer (2020) were not phenotypically distinct. While the LDA was able to correctly identify Brook Trout from New York, it correctly classified only 50% of individuals from Cosby Creek, 60% from Greenbrier Creek, and 55% from Indian Camp Creek (Table 2).

Stauffer (2020) proposed that the specimens from Cosby, Greenbrier, and Indian Camp Creeks were diagnosable using three traits: pored lateral-line scales, vertebral counts, and number of basihyal teeth. However, when data presented in Stauffer (2020) are compared to additional populations, the three meristic traits do not appear to be diagnostic at regional and range-wide scales and rather support a pattern of gradual variation in phenotype across the distribution of Brook Trout. For example, individual pored lateral-line scale counts documented in Weathers et al. (2019) for 38 streams in the GSMNP differed by up to 39 scales, and this range generally encompassed the variation reported by Stauffer (2020; Fig. 3A, B). Moreover, direct comparison of pored lateral-line scale counts from the three populations surveyed by both studies suggests that there may be significant individual and/or interannual variation within streams, in which case the range of values reported by Stauffer (2020) is not a full representation of the trait's distribution within the populations studied. For instance, data for Cosby Creek are similar between Stauffer (2020) and Weathers et al. (2019); however, counts from Indian Camp and Greenbrier Creeks, which were noted as having diagnostically low and high counts, respectively, in Stauffer (2020; but see Fig. 3B for discrepancies between the provided and published datasets), are close to the median value in the dataset collected by Weathers et al. (2019). Findings from Weathers et al. (2019) are consistent with those of Harris et al. (1978), who indicated that counts of total lateral-line scales—which are presumably correlated to counts of pored lateral-line scales—for Indian Camp Creek were approximately average relative to other populations within the GSMNP and across the range of Brook Trout (Cosby and Greenbrier Creeks were not evaluated by Harris et al. [1978]). Discrepancies between the datasets demonstrate the spatial and temporal variability in morphological and meristic traits and highlight the need for larger sample sizes across space and time in order to characterize intra- and inter-specific variation.

Weathers et al. (2019) did not collect vertebral or basihyal tooth counts, and we do not have individual-level data for these traits from Stauffer (2020) or Stauffer and King (2014). However, the three populations of Salvelinus surveyed by Stauffer (2020) did not differ significantly from one another in either trait. (Fig. 4). Additionally, vertebral counts reported in Stauffer (2020) were comparable to the range of values reported for Brook Trout in other studies (Fig. 4A; Supplemental Table S2; see Data Accessibility). While we were unable to locate independent records of basihyal tooth counts in Brook Trout, the range of values reported by Stauffer (2020) was small with significant overlap among each of the populations surveyed (Fig. 4B). Overall, the limited data available for vertebrae and basihyal tooth counts in Brook Trout reflect the inherent challenges in quantifying these traits in a species of conservation concern, where researchers seek to minimize individual mortality when conducting field studies.

Given data insufficiencies, the use of vertebral and basihyal tooth counts for diagnosing species of Salvelinus is problematic as it is not possible to compare the distribution of these traits across many populations. Moreover, vertebral counts have also been shown to be a plastic trait in Brook Trout, with up to a three vertebrae difference in populations incubated at difference temperatures and dissolved oxygen concentrations (Garside, 1966; Fowler, 1970). Less is known about variation in the extent and number of basihyal teeth, but, given that it is a trait that is related to foraging success, and others have documented rapid changes in traits associated with foraging kinematics (Adams et al., 2003), we would expect the variability in habitat across streams and the fine-scale partitioning of within-stream habitat that occurs from the realization of social dominance hierarchies to create similar levels of variability in Brook Trout. Accordingly, basihyal tooth counts are unlikely to be diagnostic across thousands of populations of Brook Trout across eastern North America, and at a very minimum would require data collected from populations at broader spatial scales before the trait could be considered diagnostic of a putative species.

Although Stauffer (2020) relied primarily on meristic traits for species delimitation, the author noted qualitative differences in mouth shape when compared to the Brook Trout neotype—specifically the presence of a more retrognathous versus isognathous mouth position. However, no morphological differences are evident in the sheared PCA in Stauffer (2020), and it is difficult to evaluate the potential for sex-based differences in the trait as the sex of the neotype is not noted by Stauffer and King (2014). Mouth morphology is difficult to compare across populations as jaw shapes are rarely reported in Brook Trout species accounts, are highly subjective, and often vary considerably within a single population (Zastavniouk et al., 2017; Fig. 5). In addition to the significant effects of sexual dimorphism and ontogeny, jaw morphology can change within a few generations in response to competition and the relative abundance of autochthonous and allochthonous prey sources in other species of Salvelinus (Nakano et al., 2020). Overall, given the potential for substantial phenotypic plasticity from multiple, interactive factors, the use of morphological traits associated with foraging has been cautioned against when conducting phylogenetic studies in Salvelinus (Reist et al., 2013; Muir et al., 2016).

As previously discussed, small sample sizes are inadequate for controlling for the effects of phenotypic plasticity and sexual dimorphism, and the minimal comparisons to regional and/or range-wide data risk a spurious conclusion of diagnosability when traits display widespread continua across the Brook Trout range. Moreover, the available data suggest broad overlap in the morphometrics and meristics of the three populations surveyed by Stauffer (2020) from the GSMNP. This shortcoming is apparent in the statistical analyses of Stauffer (2020) where, despite the author's conclusions of significant differentiation in meristic characters, the three populations broadly overlapped one another in multivariate ordination space. This result is unsurprising given that trait variation can be large in small, isolated populations (Zastavniouk et al., 2017), which highlights the need for more robust sample sizes.

The absence of sex data in Stauffer (2020) is problematic given the commonality of sexual dimorphism, including the potential for stronger phenotype-by-environment interactions in females and stronger selection for primary sex characters in males (Zastavniouk et al., 2017). Given that specimens used by Stauffer (2020) were sacrificed, it seems reasonable that fish could have been sexed directly (via dissection or using molecular sex markers [Schill et al. 2016]) or indirectly using secondary sexual characteristics (Kazyak et al., 2013). In the absence of reliable sex data, the robust literature on sexual dimorphism in Salvelinus could have been used to evaluate the potential influence of ontogeny and sexual dimorphism on inferred patterns of variation (Proulx and Magnan, 2004; Nitychoruk et al., 2013; Zastavniouk et al., 2017).

While Stauffer (2020) asserts that all individuals were sexually mature, the range of individual sizes from each stream makes this assertion speculative, and sexually dimorphic characters tend to become exaggerated with size (Kazyak et al., 2013; Nitychoruk et al., 2013). Acknowledging that food and habitat limitations often limit the size structure of Brook Trout populations within the species' southern range (Ensign et al., 1990), and that Stauffer (2020) measured standard length (SL) rather than a total length (TL), a threshold of 100 to 120 mm TL is commonly applied for distinguishing adults and subadults/juveniles (Etnier and Starnes, 1993; Whiteley et al., 2012; Davis et al., 2015). The size range of individuals included in Stauffer (2020), particularly individuals from one stream (which we believe to be Greenbrier Creek, but inconsistent data reporting makes it unclear which stream used individuals that ranged in size from 77.0 to 107.7 mm), raises concerns that a significant proportion of individuals in the analysis were not sexually mature. Moreover, the sex for the holotype from Greenbrier Creek was not indicated, but the specimen was only 97.2 mm suggesting the holotype may have been described from a sexually immature individual. The inclusion of subadults would likely have limited effects on meristics, but it could impart significant biases on morphometric relationships, particularly those related to mouth and jaw shape/size and body depth.

The potential for familial relationships within the samples considered by Stauffer (2020) was also not considered. In small headwater streams, full-sibling families are often restricted to short reaches of habitat (Hudy et al., 2010), which suggests the small sample sizes may represent a limited number of families. If morphology and meristics are heritable, which is an implicit assumption of using these traits to delimit species, then undocumented family structure could bias assessments of taxonomy based on phenotypic characters.

The potential for intraspecific variation and plasticity highlights the importance of molecular data for a thorough taxonomic investigation in this clade (Muir et al., 2016; Kuparinen and Hutchings, 2019). Stauffer (2020) did reference previous neutral microsatellite genetic analyses by Stauffer and King (2014) and suggested significant population structure and isolation. While microsatellites provide valuable insights on recent population demographic parameters (Wan et al., 2004; Selkoe and Toonen, 2006), they are not informative for inferring relationships at deeper evolutionary time scales. When microsatellites are appropriately applied, it is not uncommon to find significant spatial structuring at scales commensurate to, or sometimes finer than, those detected by Stauffer (2020). Accordingly, the finding of population genetic structure is unsurprising given the ecology of Brook Trout and the limited number of populations analyzed (Kelson et al., 2015; Kazyak et al., 2016; Weathers et al., 2019).

While there are currently not sufficient data to support the distinctiveness of additional species currently classified as Salvelinus fontinalis, we acknowledge that populations of Brook Trout, particularly those in the southern extent of the species' range, are not interchangeable. Genetic and phenotypic differences do exist (Kazyak et al., 2015; Weathers et al., 2019), sometimes at fine scales, and identifying the causes and consequences of this variation is critical to the ecology and conservation of the species. Taxonomic resolution may prove challenging as Brook Trout are frequently characterized by small, isolated populations. In these scenarios, phenotypic variation can be a weak indicator of natural selection due to the stronger influences of drift, founder effects, or the decoupling of phenotypic correlations due to interspecific competition (Lowe et al., 2012; Zastavniouk et al., 2017). Conversely, small, isolated populations may be under stronger selective pressures than larger populations and can harbor unique phenotypic and genotypic variation that is important for species conservation at larger, range-wide scales (Stelkens et al., 2012; Fraser et al., 2014). Ultimately, it is important to consider how neutral and adaptive variation can confound inferences about speciation. Stauffer (2020) suggested that reduced offspring fitness observed by Kulp et al. (2017) after artificially crossing individuals from the three GSMNP populations was suggestive of speciation; however, the potential for outbreeding depression due to the breakup of coadapted gene complexes that evolved in differentially locally adapted populations was not considered (Edmands and Timmerman, 2003).

A future challenge for Brook Trout conservation will be to decouple the relative influence of genotype and phenotypic plasticity in shaping local morphology and behavior. Plasticity to local conditions can have long-term and/or transgenerational effects, as demonstrated by consistent differences in spawning timing in the three populations examined by Stauffer (2020) even after being held in captivity under identical conditions for 1.5 years. While plasticity can be particularly beneficial for surviving in highly stochastic habitats (Miner et al., 2005), it can be misconstrued as evidence of divergent evolution (Price et al., 2003). For example, the phenology of spawning in stream salmonids is often adapted to location conditions, to the extent that spawning occurs across multiple months, and sometimes all year (Rinne, 1980; USFWS, 2009). However, phenological plasticity may not be as strongly correlated with ecological speciation as previously believed (Hendry, 2009), and populations often return to less variable behavior when environmental conditions stabilize.

Phylogenetic uncertainty stemming from vast intraspecific variation is not unique to the discussion of Brook Trout taxonomy and has long been part of conversations related to the conservation of Salvelinus and other stream salmonids (Guinand et al., 2021). Repeated radiation through increasingly diverse environments has led to significant uncertainty in how to appropriately classify the extensive taxonomic and functional diversity in many species of Salvelinus (Reist et al., 2013). Given these challenges, insular efforts to describe microendemic species that ignore the larger body of research and literature are likely to yield greater uncertainty. Ultimately, this highlights the need for holistic phylogenies that incorporate equally the influence of biology, ecology, geography, and genetics for inferring past and present relationships among Brook Trout and its congeners.

Addressing questions about genotypic-by-phenotypic interactions and the appropriate scales of conservation, particularly as it relates to conservation under future climatic and anthropogenic disturbance regimes, will be aided by large-scale studies of Brook Trout phylogenetics and genomics. Range-wide assessments would be particularly beneficial for investigating the adaptive significance of broad patterns of genetic and morphometric differences that have been documented between regions and, if warranted, modifying conservation objectives to meet the unique phylogenetic histories and future conservation concerns at these larger spatial scales. Widespread interest in Brook Trout conservation has already motivated the establishment of a unique, range-wide collaboration among non-profits, universities, and state and federal agencies. Using these partnerships to investigate Brook Trout phylogenetic relationships offers great potential for developing management plans that effectively balance the need to conserve ecological processes across multiple temporal and spatial scales with the challenges of executing conservation efforts given finite resources. Ultimately, while conservation of important sources of genetic and phenotypic diversity is of utmost priority, poorly supported taxonomic revisions can negatively impact conservation initiatives. As such, caution is warranted when applying species concepts that allow for microendemic distributions when characterizing a broadly distributed taxon that is noted for high levels of phenotypic plasticity.

Supplemental material is available at https://www.ichthyologyandherpetology.org/i2020154. Unless an alternative copyright or statement noting that a figure is reprinted from a previous source is noted in a figure caption, the published images and illustrations in this article are licensed by the American Society of Ichthyologists and Herpetologists for use if the use includes a citation to the original source (American Society of Ichthyologists and Herpetologists, the DOI of the Ichthyology & Herpetology article, and any individual image credits listed in the figure caption) in accordance with the Creative Commons Attribution CC BY License.

We thank Paul Super for facilitating the receipt of raw data, Duane Diefenbach and Ken White-Sands for assisting with preliminary data analyses, and Eric Hallerman for providing comments on a previous draft. Mark C. Colaw assisted with the collection of the Brook Trout images used in Figure 5.