Introduction

Lampreys are jawless aquatic vertebrates in the order Petromyzontiformes. They are globally distributed, with two families, Geotriidae and Mordaciidae, in the Southern Hemisphere, and one family, Petromyzontidae, in the Northern Hemisphere (Potter et al. 2015). There is consensus on the placement of taxa within these families (Gill et al. 2003; Renaud 2011), but their systematics at lower levels of classification, such as the relationships among genera and assignment of species within them, are contentious and fluid (Potter et al. 2015; Docker and Hume 2019). This is in part because lampreys lack many of the morphological characteristics used in taxonomic studies of bony fishes, and those few putatively diagnostic characters are often only expressed in the adult life stage (Renaud 2011). At an even finer scale, recent molecular analyses cast doubt on the accuracy of separating some lampreys as distinct species. A recurring pattern among lamprey taxa is the presence of species complexes whose members are often described as paired or stem-satellite species (Docker 2009). The members of a complex are often morphologically indistinguishable during their extended larval stages (Zanandrea 1959; Potter 1980) but develop striking anatomical and behavioral differences as juveniles and adults. Each complex consists of a migratory (anadromous or potamodromous) species with a parasitic, predatory, or necrophagous phase prior to sexual maturity and one or more closely related species that presumably evolved from the former, which as adults either remain in freshwater, do not feed, or both (Docker 2009; Docker and Potter 2019). Although the present taxonomy of lampreys reflects the belief that phenotypic diversity in life history expression is grounds for species designation (Bailey 1980; Potter 1980), this interpretation has not been universally accepted among lamprey researchers (McPhail and Lindsey 1970; Espanhol et al. 2007; Kucheryavyy et al. 2016a, 2016b) and is rarely supported by genetic information. For example, Docker et al. (2012) found that in one species pair—the parasitic silver lamprey Ichthyomyzon unicuspis and non-parasitic northern brook lamprey I. fossor—putative specimens of both species collected from the same rivers could not be differentiated using either mitochondrial (mt) DNA sequences or microsatellite loci (see also April et al. 2011). Similarly, many of the species within the genus Entosphenus lack diagnostic differences in their mtDNA sequences despite different migratory or feeding life histories (Lorion et al. 2000; Lang et al. 2009; Taylor et al. 2012).

There is also evidence that species designations based on phenotype overlook cryptic taxa. For example, Yamazaki et al. (2006) analyzed mtDNA sequences and allozymes of specimens of Lethenteron from across eastern Asia. Although the authors argued that Siberian brook lamprey Lethenteron kessleri should be regarded as a synonym for Far Eastern brook lamprey Lethenteron reissneri, they also observed two groups presumed to be Lethenteron reissneri that were highly genetically divergent; one group appeared to represent a new species within this genus (Lethenteron sp. N), and the other group likely constituted a new species unaffiliated with Lethenteron (albeit labeled Lethenteron sp. S as a placeholder). Likewise, Boguski et al. (2012) found that morphologically similar but geographically disparate populations of Lampetra sp. were genetically divergent and hypothesized that some may represent undescribed species.

Many of the aforementioned issues are evident among species of Lampetra. The type species for this genus—European river lamprey Lampetra fluviatilis—and several congenerics in Europe are distantly related to those in North America (Lang et al. 2009; Boguski et al. 2012), so much so that the lone member from eastern North America has been regarded by some as belonging to a different genus (Okkelbergia; Hubbs and Creaser 1922). Similarly, phylogenetic analyses suggest that Lampetra of western North America may also constitute a separate genus (Docker et al. 1999; Blank et al. 2008; Lang et al. 2009). In western North America, the current taxonomy recognizes four species: (1) western river lamprey Lampetra ayresii, an anadromous, parasitic species with a disjunct distribution restricted to large rivers and their estuaries from Alaska to California; (2) western brook lamprey Lampetra richardsoni, a non-parasitic freshwater species (but see Beamish 1987), considered a species pair with Lampetra ayresii and sharing a similar but more continuous range; (3) Kern brook lamprey Lampetra hubbsi, originally considered to belong to Entosphenus, with the same life history as L. richardsoni but found only in the San Joaquin River basin in central California (Moyle 2002; Potter et al. 2015); and (4) Pacific brook lamprey Lampetra pacifica, also maturing as a non-parasitic species in freshwater that was originally thought to be widely distributed in Oregon and California (Vladykov 1973), then synonymized with L. richardsoni (Robins et al. 1991), and most recently taxonomically resurrected as a distinct species and thought to inhabit an undefined portion of the Columbia River basin (Reid et al. 2011) that apparently does not include Washington (Wydoski and Whitney 2003).

Several approaches could reduce the uncertainty about the taxonomic status of western North American Lampetra within Petromyzontidae. Although phylogenetic analyses have included thorough treatments of particular genera or regions (Yamazaki et al. 2006; April et al. 2011; Boguski et al. 2012; Li 2014; White 2014) and broad overviews of all genera (Lang et al. 2009), no analysis has taken advantage of all publicly available sequences for any gene region to place Lampetra of western North America in their evolutionary context relative to other lampreys. Moreover, formal molecular species delimitation has not been undertaken for any genus of lamprey. Finally, understanding to what degree geography dictates population structure within species is important for guiding management actions.

Clarifying these issues for Lampetra in western North America is more than a taxonomic issue. In 2004, the US Fish and Wildlife Service denied a petition for listing Lampetra ayresii and Lampetra richardsoni (as well as Pacific lamprey Entosphenus tridentatus) in the western continental USA under the Endangered Species Act (US Fish and Wildlife Service 2004). In their decision, the agency concluded that there was no evidence for a widespread decline in either Lampetra species but regarded each as having a broad distribution from Alaska to California. If these taxa constitute a species complex, however, it may be worthwhile to evaluate the conservation status of all of its members. To help fill this information gap, we developed a global phylogeny of all members of Petromyzontidae to assess the phylogenetic position of the western North American Lampetra (Objective 1). We used all available sequences of two mitochondrial gene regions in public sequence repositories (that met minimum length and quality thresholds), as well as sequences of new specimens from much of the Willamette River basin in Oregon. Because sequences for one of these gene regions were taxonomically comprehensive, we performed molecular species delimitation among samples of western North American Lampetra using a variety of approaches to generate hypotheses about the number, identity, and distribution of species (Objective 2). We refined our estimates of the distribution of some of these taxa by using an even larger dataset and performing molecular specimen assignment. Finally, we assessed haplotype diversity and evaluated fine-scale genetic structure of Lampetra among adjacent subbasins of the Willamette River basin to understand population boundaries and identify units for management and conservation (Objective 3).

Methods

Study area

Data collection for this study took place across six, 8-digit hydrologic units (subbasins averaging ~1125 km2; US Geological Survey 2004) located in the Willamette River basin, Oregon, USA: the Upper Willamette, North Santiam, South Santiam, McKenzie, Middle Fork Willamette, and Coast Fork Willamette River subbasins (Fig. 1). Thirteen major hydroelectric dams in this basin create barriers to upstream fish passage, preventing migratory species from reaching hundreds of kilometers of spawning and rearing habitat in the upper subbasins. Consequently, migratory fish species, including Entosphenus tridentatus, were largely extirpated from their historical ranges in these subbasins. Information on species identity and distributions of Lampetra throughout the Willamette River basin is limited. Fish collections dating back to 1951 suggest the presence of Lampetra pacifica or Lampetra richardsoni throughout the Willamette River, including our study area (see below; Table S1; Fishnet2 Portal, www.fishnet2.org, accessed on July 23, 2022 and September 23, 2022). Reid et al. (2011) used paired morphological and genetic analysis to confirm the presence of Lampetra richardsoni in the Tualatin subbasin and Lampetra pacifica in the Clackamas subbasin (both tributaries to the lower main stem Willamette River). No other studies have attempted to identify species of Lampetra or their distribution and population boundaries in the Willamette River basin.

Fig. 1
figure 1

Locations of lamprey samples collected in this study. Labels denote the subbasin (4th-level hydrologic unit). Numbers correspond to Map ID in Table S1. Black triangles represent locations of major hydroelectric dams; green shading represents the US National Forest boundaries

New sample collection

New collections were obtained as part of a larger study to determine the occupancy of lampreys throughout subbasins of the Willamette River basin on the Willamette National Forest. Survey locations were selected using a general random tessellation stratified (GRTS) sampling regime structured at the level of 12-digit hydrologic units (sub-watersheds). Surveys targeted three 100-m segments of floodplain channels as well as perennial streams with an average slope less than 2.5%. Only surveys where lampreys were observed were included in this study. To increase the sample size and geographic scope of this study, additional surveys were conducted outside of the Willamette National Forest. These additional sites were selected based on accessibility (e.g., road crossings) and targeted shallow, low-gradient reaches (< 2.5% stream slope) considered suitable for lamprey.

Between April 2015 and October 2017, larval (and occasionally adult) lampreys were captured at 42 locations in the Willamette River basin using a Smith-Root LR-24 backpack electrofisher (Table S1, Fig. 1). Specimens were sedated using tricaine methanesulfonate diluted in stream water, and tissue samples were collected from each individual. Tissues were preserved individually in vials filled with 95% ethanol and stored at room temperature. Collections were made following the Guidelines for the Use of Fishes in Research, Oregon Department of Fish and Wildlife Scientific Taking Permit 4(d) Permit # 22632 following methods outlined in Jenkins et al. (2014).

During field surveys, Lampetra were distinguished from E. tridentatus based on pigmentation of the caudal ridge and caudal fin (Goodman et al. 2009). Although we did not attempt to identify Lampetra to species in the field, we expected that most specimens would be Lampetra richardsoni or Lampetra pacifica (Reid et al. 2011; Boguski et al. 2012; Potter et al. 2015). We did not expect to find putative Lampetra ayresii because this form is migratory, typically found in coastal systems, and has never been documented in our study area.

DNA extraction and sequencing

We used the Qiagen DNeasy® Blood and Tissue Kit to extract DNA from up to eight individuals from a given sampling site (Table S1) and stored the DNA at −20 °C until sequenced. A 988-bp segment of the cytochrome b (cyt b) gene and a 1337-bp segment of the cytochrome c oxidase subunit I (COI) gene was amplified from each sample using existing or newly developed primers (Table S2; Boguski et al. 2012). The 40-μl amplification reactions contained 4 μl (~4–20 ng) DNA, 4 μl of 10× PCR buffer, 4 μl MgCl2 (2.5 mM), 1 μM of each primer, 200 μM each dNTP, 25 μg BSA, and 1 Unit Titanium® Taq DNA Polymerase (Takara Bio USA, Inc.), with the remainder consisting of PCR-grade distilled water. Segments of the cyt b gene were amplified under the following conditions: 2 min at 94 °C, 10 × (2 min at 94 °C, 60 s at 60 °C, 2 min at 72 °C), 10 × (2 min at 94 °C, 60 s at 58 °C, 2 min at 72 °C), 10 × (2 min at 94 °C, 60 s at 55 °C, 2 min at 72 °C), 5 min at 72 °C, then held at 12 °C. Segments of the COI gene were amplified under the following conditions: 12 min at 95 °C, 35 × (60 s at 95 °C, 60 s at 55 °C, 90 s at 72 °C), 5 min at 72 °C, then held at 12 °C. The PCR products were cleaned with Exo-SAP-IT (Affymetrix), and DNA sequence data were obtained using the Big Dye kit and the 3700 DNA Analyzer (ABI; Eurofins Genomics, Louisville, KY). Sequences were viewed and aligned with Sequencher (Gene Codes Corp. MI). We used MEGA version 7 (Kumar et al. 2016) to translate all sequences into amino acids to verify there were no stop codons or nonsense bases which would indicate a sequencing error or the presence of nuclear copies of mitochondrial DNA.

Sequence datasets

To understand the phylogenetic position of western North American Lampetra within the Petromyzontidae, we created two sequence datasets. The “cyt b dataset” comprised the 988-base cyt b sequences obtained in this study (n = 134) and all publicly available sequences (in GenBank; n = 388) representing species in Petromyzontidae (Table S3; Table S4). We excluded public sequences that did not completely overlap our sequences or that contained ambiguous or missing bases; sequences longer than 988 bp were trimmed to ensure consistent coverage across all sequences in the dataset. Similarly, the “COI dataset” comprised 585-base COI sequences from our specimens and an additional 393 publicly available sequences, screened for length and quality (Table S3; Table S4). For both datasets, we used all publicly available sequences (cyt b, n = 6; COI, n = 48) of species of Geotria and Mordacia, representing the remaining two families in Petromyzontiformes, as outgroups. All public sequence labels were updated to reflect the presently accepted taxonomy following Fricke et al. (2022; e.g., sequences labeled as Lampetra tridentata were updated to Entosphenus tridentatus). Apparent misidentified specimens were ignored, except one COI sequence identified as Lampetra ayresii (JQ354155), but which was identical to sequences of Entosphenus tridentatus and thus relabeled.

We created two datasets to address the objective related to species delimitation and geographic distribution. For species delimitation within the western North American Lampetra, we created the “species delimitation dataset.” This dataset used the cyt b dataset restricted to representatives of Lampetra from western North America (n = 294) and a single sequence of Entosphenus tridentatus as an outgroup. Finally, to more precisely describe the geographic distribution of delimited species, we created a “specimen assignment” dataset that also included the shorter sequences (295–438 bases; n = 205) excluded from the species delimitation dataset (Table S3).

Finally, to understand genetic diversity and structure of Lampetra among subbasins in our study area, we created a “concatenated” dataset representing only new samples collected in this study by concatenating a 1337-base segment of the COI gene and a 988-base segment of the cyt b gene to create a single 2325-base sequence for each individual. Sequences were concatenated using Geneious 10.0.6 (http://www.geneious.com, Kearse et al. 2012). Prior to analysis, we reduced sequences in the species delimitation and concatenated datasets to representative haplotypes.

Objective 1: phylogenetic position of western North American Lampetra

To better understand the systematics of Lampetra from western North America within Petromyzontidae, we used the cyt b and COI datasets to estimate maximum likelihood phylogenetic trees in IQ-TREE 2.1.1 (Minh et al. 2020) and implemented via the CIPRES gateway (https://www.phylo.org/). Initially, we removed the outgroups from each dataset, assigned three preliminary partitions based on codon position (see Bofkin and Goldman 2007), then selected edge-linked partitions and the TESTMERGE setting to determine the best-fitting evolutionary model (as measured by BIC) for each partition (Table S5). We then added the outgroup sequences and constructed a consensus maximum-likelihood tree for each dataset, assigning support values based on 1000 ultrafast bootstrap replicates. The final trees were visualized and edited using FigTree v1.4.4 (http://tree.bio.ed.ac.uk/software/figtree) and Inkscape v1.0.1 (http://inkscape.org).

Objective 2: species delimitation and specimen assignment in western North American Lampetra

We used five methods of single-locus species delimitation to generate species hypotheses within Lampetra from western North America using the species delimitation dataset. We recognized species under the phylogenetic species concept, which relies on reciprocal monophyly among lineages (Nixon and Wheeler 1990), further modified by a distance threshold often typical of interspecific differences among fishes (Ward 2009). We began by constructing a maximum likelihood tree for this dataset following the approach described for Objective 1. In the first analysis, we input this tree into the online version of maximum-likelihood, multi-rate Poisson tree processes (mPTP; https://mptp.h-its.org/#/tree; Kapli et al. 2017), which examines branch lengths to estimate the shift from intraspecific to interspecific divergence and can provide a liberal estimate of species diversity relative to other methods (Zhang et al. 2013). Second, we constructed statistical parsimony networks (SPN) with a 95% parsimony limit in TCS v.1.21 (Clement et al. 2000) and visualized these using tcsBU (https://cibio.up.pt/software/tcsBU; Santos et al. 2015). When applied to datasets with multiple taxa, the number of unconnected networks constitutes a conservative estimate of species diversity (Hart and Sunday 2007). Third, we used the online version (https://bioinfo.mnhn.fr/abi/public/asap/#) of ASAP (Assemble Species by Automatic Partitioning; Puillandre et al. 2021), for which pairwise genetic distances among haplotypes are used to distinguish between intraspecific variation and interspecific divergence and a scoring system to identify the best-fitting set of partitions (i.e., candidate species) based on both the probability of panmixia in a given partition and genetic distances within and among partitions. The set of partitions with the lowest score is deemed the most likely number of species. We adopted the default values and distances based on the K80 substitution model because of its similarity to the traditional distance metric used in sequence-based analyses (Ratnasingham and Hebert 2013). The analysis was run 10 times with different initial seeds, for which the highest-scoring set of partitions only changed once. Because the approach is similar to ABGD (Puillandre et al. 2012), we assumed that species counts would be conservative relative to other methods (Puillandre et al. 2021). Fourth, we reexamined the maximum-likelihood tree from the first step and interpreted bootstrap values > 85 as strong support (Minh et al. 2013) for a candidate species identified by one of the preceding methods. Finally, we used MEGA 7 (Kumar et al. 2016) to examine pairwise genetic distances, expressed as p-distance, among these candidate species. We contrasted the maximum intraspecific distances with the minimum interspecific distances. If the latter exceeded ~2% and was larger than the former (i.e., a barcode gap was present; Ratnasingham and Hebert 2013), we interpreted this as supporting the designation of a candidate species.

Due to the limited representation of the COI gene in public databases, we did not perform a species delimitation analysis using this gene. However, we did compare placement of Lampetra sp. from western North America in the COI maximum likelihood tree (above) to the cyt b tree to broadly look for any apparent differences in species placement.

Finally, to assign specimens to putative species, we built a cyt b neighbor-joining tree (in MEGA 7 using the number of base pair of differences) using the specimen assignment dataset.

Objective 3: haplotype diversity and genetic structure among lampreys in the Willamette River basin

To assess haplotype diversity among lamprey in the Willamette River basin (Objective 3), we used the concatenated dataset to create statistical parsimony networks following the method described in Objective 2 above. Unlike phylogenetic trees, statistical parsimony networks are unrooted and do not assume bifurcating divergence (Clement et al. 2000). As a result, they are useful for exploring haplotype relatedness at the population level, particularly when the number of variable sites may be low and parsimony methods will be unable to resolve branching patterns with high levels of support (Clement et al. 2000).

To look for genetic substructure in Lampetra spp. in the Willamette River basin, we grouped individual sequences from the concatenated dataset based on subbasin of capture. We then used MEGA 7 to calculate average p-distance within and among individuals grouped by subbasin. We performed an analysis of molecular variance (AMOVA) in R (R Core Team 2020) using the poppr.amova() function in the package “poppr” to assess population differentiation at the subbasin level (Kamvar et al. 2014). Significance for the AMOVA was assessed at α ≤ 0.05.

Results

Objective 1: phylogenetic position of western North American Lampetra

The cyt b maximum likelihood tree (Fig. 2; Fig. S1) robustly supported (100% bootstrap support) two clades representing the Petromyzontinae subfamily (Ichthyomyzon, Petromyzon, and Caspiomyzon) and the Lampetrinae subfamily (Entosphenus, Tetrapleurodon, Lethenteron, Lampetra, and Eudontomyzon). Within Lampetrinae, five clades also had bootstrap support of 100: (1) three specimens identified as Arctic lamprey Lethenteron camtschaticum or L. reissneri (equivalent to the aforementioned Lethenteron sp. S; Yamazaki et al. 2006; Lang et al. 2009; Li 2014); (2) Entosphenus plus Tetrapleurodon; (3) all remaining Lethenteron (save one; see below) plus Korean lamprey Eudontomyzon morii; (4) western North American Lampetra; and (5) all remaining Eudontomyzon and Lampetra and Western Transcaucasian brook lamprey Lethenteron ninae. This last clade was further subdivided into three additional highly supported lineages: (1) Eudontomyzon; (2) least brook lamprey Lampetra aepyptera from eastern North America; and (3) all remaining European and Asian Lampetra plus Lethenteron ninae.

Fig. 2
figure 2

Maximum likelihood trees showing the phylogenetic relationship of Northern Hemisphere lampreys at the cyt b gene and COI genes. Bootstrap values are shown for branches with > 85% support. See Figures S1 and S2 for an expanded version of each tree displaying all sequences in each dataset

The topology of the COI maximum likelihood tree (Fig. 2; Fig. S2) was remarkably similar to that of the cyt b tree but showed even greater support at intermediate branches. Bootstrap support for subfamilies Petromyzontinae and Lampetrinae was 99 and 100%, respectively, and relationships within Petromyzontinae (Petromyzon + Ichthyomyzon sister to Caspiomyzon) were identical to those observed in cyt b. Within Lampetrinae, membership within the same five clades was identical to that in the cyt b phylogeny and well supported. Intermediate branches relating these clades differed slightly from the cyt b tree but also had strong (> 90%) support. Clade membership, in descending order of divergence, was (1) “Lethenteron sp. S”; (2) Entosphenus + Tetrapleurodon; (3) western North American Lampetra; (4) Lethenteron (save one) and Eudontomyzon morii; and (5) remaining members of Lampetra, Eudontomyzon, and Lethenteron ninae. Division of the last clade into the same three additional clades seen with cyt b was also strongly supported (bootstrap support ≥ 99%): (1) Lampetra from Europe and Asia, along with Lethenteron ninae; (2) all remaining Eudontomyzon; and (3) Lampetra aepyptera in eastern North America.

Objective 2: species delimitation and specimen assignment in western North American Lampetra

In the 988-bp species delimitation dataset, we observed 160 variable nucleotides sites resulting in 67 haplotypes (Table S3). The mPTP and statistical parsimony networks delimited seven candidate species (Figs. 3 and 4), the first six of which had limited distributions: (1) Lampetra pacifica from the Willamette River basin below Willamette Falls (cyt b H63 and H64); (2) Lampetra hubbsi from the San Joaquin River basin (cyt b H54) and one individual from Paynes Creek in the Sacramento River basin, California (cyt b H55); (3) specimens from the Siuslaw River, Oregon (cyt b haplotypes H56, H57, and H58); (4) specimens from Fourmile Creek in the Klamath River basin, Oregon (cyt b haplotypes H59, H60, H61, and H62); (5) specimens from Mark West Creek in the Russian River basin, California (cyt b haplotypes H65, H66, and H67); and (6) specimens from Kelsey Creek, a tributary to Clear Lake, California (cyt b H53). These six candidate species were supported in the maximum likelihood tree with 97–100% bootstrap support, except for those from Kelsey Creek, which did not receive an estimate because it was represented by only one haplotype. The seventh candidate species included 52 haplotypes (cyt b H01–H52) variously identified as Lampetra ayresii, Lampetra richardsoni, and Lampetra sp. and was supported with 89% bootstrap support in the maximum likelihood tree. Within this lineage, haplotypes identified as Lampetra richardsoni and Lampetra ayresii were polyphyletic (Fig. 3). This candidate species contained all 21 Lampetra haplotypes from our study area in the Willamette River basin, 20 of which grouped into a single subtree. The 21st haplotype (cyt b H39) represented three individuals from Muddy Creek in the Upper Willamette River basin, along with an additional 36 individuals from 20 waterbodies across Alaska, British Columbia, Washington, and Oregon.

Fig. 3
figure 3

Maximum likelihood tree displaying relationships among Lampetra spp. of western North America at the cyt b gene. Tip labels correspond to cyt b haplotypes in Table S4. Haplotypes in bold include at least one sequence identified to the species level. Bootstrap values are shown for branches with ≥ 85% support. Black dots denote haplotypes observed in the headwater subbasins of the Willamette River (this project’s study area). Bars to the right of the phylogenetic trees denote species groups identified by ASAP (Assemble Species by Automatic Partitioning; n = 2 species), SPN (statistical parsimony networks; n = 7 species), and mPTP (multi-rate Poisson tree processes; n = 7 species). Labels to the right correspond to the seven species groups identified through species delimitation analysis. All haplotypes identified as L. richardsoni are included in the L. ayresii species group

Fig. 4
figure 4

Statistical parsimony network (SPN) showing relatedness among cyt b haplotypes of Lampetra sp. in western North America. Nodes represent unique haplotypes and are color-coded based on putative species groups. Nodes marked with a cross denote haplotypes observed in headwater subbasins of the Willamette River. Small white nodes show transitionary haplotypes that were not observed in the current dataset. Segments between nodes represent a single base pair mismatch from each connected node. All haplotypes representing specimens labeled as L. richardsoni are included in the L. ayresii species group

The ASAP analysis was more conservative than the others, with support for two candidate species (Fig. 3). Here, only Lampetra from Kelsey Creek, California, was delimited as a single candidate species; all others were pooled into a single candidate taxon. The next best fitting model in 9 of 10 analyses also delimited Lampetra in Mark West Creek, California, as a candidate species.

Divergence among all Lampetra haplotypes in the species delimitation dataset ranged from 0.1 to 6.4% (Table 1). When grouped into the consensus seven candidate species, six of these exhibited intraspecific divergence ≤ 1.21% and a barcode gap (sensu Ratnasingham and Hebert 2013) exceeding 1.70% (minimum interspecific distance observed was 1.72% between Lampetra ayresii and the Siuslaw River lineage). The maximum intraspecific divergence (2.33%) within the highly diverse and widespread seventh candidate species (which included specimens identified as Lampetra ayresii and Lampetra richardsoni) exceeded its minimum interspecific divergence for all but two of the candidate species.

Table 1 Genetic differences (measured as p-distance and expressed as a percentage) within and among candidate species of western North American Lampetra. Values on the diagonal (in bold) are maximum intraspecific differences; values below the diagonal are minimum interspecific differences. Specimens from Kelsey Creek were represented by one haplotype, so the intraspecific difference was not estimated. The L. ayresii species group includes both anadromous/parasitic and resident/non-parasitic ecotypes, commonly known as western river lamprey and western brook lamprey, respectively

The COI maximum likelihood tree (Fig. S2) had limited coverage of Lampetra from western North America. However, it also resolved Lampetra hubbsi and Lampetra pacifica as highly supported (bootstrap support ≥ 99%) clades. These had an uncertain relation to other specimens of Lampetra from Oregon and Washington that were likely to represent the widespread, diverse L. ayresii complex from the species delimitation analysis.

The sequences included for specimen assignment (n = 205) represented 33 river drainages ranging from British Columbia to California, with most in Washington (Table S3). All specimens unambiguously assigned to one of the delimited taxa. Specifically, one sequence (AF177958) represented specimens collected from the Merced and Kings Rivers in California and assigned to Lampetra hubbsi; all other sequences assigned to Lampetra ayresii (Fig. S3; Table S3).

Objective 3: haplotype diversity and genetic structure among lampreys in the Willamette River basin

In the concatenated dataset, divergence among Lampetra captured in headwater subbasins of the Willamette River basin ranged from 0.1 to 1.6% across 71 variable nucleotide sites resulting in 34 haplotypes (Table S3, Table S6). The statistical parsimony network estimated from the concatenated dataset revealed substantial geographic structure (Fig. 5). Only one haplotype (W) was observed in more than one subbasin. Three haplotypes in the Upper Willamette River basin (AE, AF, and AH) were so divergent (1.1–1.6%) that they formed a separate network. The Middle Fork Willamette River basin was the most heavily sampled (n = 59 Lampetra individuals) and the most diverse, representing 11 of the 34 Lampetra haplotypes in the concatenated dataset (Table S3). Average divergence among Lampetra was 0.13 to 0.79%, and within-basin divergence ranged from 0.01 to 0.66% (Table 2). Except for the Upper Willamette River basin, within-basin divergence was less than that observed among basins. The AMOVA comparing genetic variation among versus within subbasins was significant (p < 0.001, d.f. = 133), with 36.0% of variation partitioned among subbasins and 64% within a subbasin.

Fig. 5
figure 5

Statistical parsimony network showing relatedness among concatenated COI–cyt b haplotypes of Lampetra sp. in headwater subbasins of the Willamette River. Nodes are scaled to the number of individuals observed with that haplotype and are color-coded based on subbasin of capture. Nodes that appear as pie charts show the proportion of individuals that originated from a particular subbasin. Letters identify the concatenated haplotype represented by each node as listed in Table S4. Note that haplotypes AE, AF, and AH represent three individuals from Muddy Creek represented by the geographically widespread haplotype cyt b H39

Table 2 Average pairwise genetic distance (p-distance, expressed as a percentage) of Lampetra within and among subbasins in the Willamette River basin. Average within-group genetic distance is listed on the diagonal in bold

Discussion

Our phylogenetic analysis comparing the relationships among Petromyzontidae is among the most robust molecular evaluations of this family to date. Our analysis represented every genus and most recognized species of lampreys in the family, included all publicly available sequences, and was repeated for two mitochondrial genes. The analysis reaffirmed Vladykov’s (1972) concept of the Petromyzontinae and Lampetrinae subfamilies (see Potter et al. 2015), as well as shortcomings in the present taxonomy also identified by other studies with more limited datasets (Docker et al. 1999; Blank et al. 2008; Lang et al. 2009; Li 2014; White 2014). That three of the genera within Lampetrinae (Lethenteron, Eudontomyzon, and Lampetra) are not monophyletic is thoroughly discussed by others (Lang et al. 2009; Li 2014; Pereira et al. 2021) and could be resolved with reclassification of species among genera, and, in the case of Lethenteron sp. S, naming of a new genus (see Yamazaki et al. 2006).

Our study focused on Lampetra and the placement of specimens from western North America relative to other congeneric specimens. Like others, using cyt b sequences, we observed that Lampetra was polyphyletic (Blank et al. 2008; Lang et al. 2009; Li 2014), and these patterns were corroborated by our phylogenetic analysis of COI (Fig. 3). Monophyly could be achieved by revising four strongly supported clades: (1) restricting Lampetra to members of this genus from Europe and Asia, which includes the type species for this genus, Lampetra fluviatilis, and reassigning Lethenteron ninae to this genus; (2) resurrecting Okkelbergia as the genus for Lampetra aepyptera in eastern North America, as was suggested by Lang et al. (2009); see also Hubbs and Potter (1971) and Vladykov and Kott (1976); (3) retaining the genus Eudontomyzon for its current members (save Eudontomyzon morii which should be reassigned to Lethenteron); and (4) assigning species of Lampetra in western North America to a new genus, given that this clade is strongly supported as sister to these three clades in both analyses (as well as Lethenteron sensu stricto in the COI phylogeny). Although proposing a new genus name according to the Rules of Zoological Nomenclature (International Commission on Zoological Nomenclature 1999) is beyond the scope of this study, formal taxonomic revision is needed to better represent the distinctiveness of western North American Lampetra.

Our analyses also favor revisiting species designations within the forms of Lampetra from western North America. Four of five tests provided strong support for seven candidate species. These lineages were represented by Lampetra hubbsi, Lampetra pacifica, four undescribed lineages of Lampetra, and one lineage that includes all specimens identified as Lampetra ayresii and Lampetra richardsoni, as well as a number of unidentified Lampetra specimens. Although Lampetra ayresii and Lampetra richardsoni are considered paired species that can be differentiated based on phenotypic and life history characteristics that develop at metamorphosis (Docker 2009), the species delimitation analysis recognized these as a single species with occasional haplotype sharing and a lack of reciprocal monophyly. Because Lampetra ayresii was first described by Guenther in 1870 (see Vladykov and Follet 1958), 95 years prior to the description of Lampetra richardsoni (Vladykov and Follet 1965), Lampetra ayresii has precedence as the name for members of this clade, in compliance with the principle of priority (Article 23.1) in the International Code of Zoological Nomenclature (International Commission on Zoological Nomenclature 1999). We regard specimens identified as Lampetra richardsoni as a synonymous with Lampetra ayresii and consider them to be a life history variant within the Lampetra ayresii species complex.

These results support a growing body of evidence that the paired or stem-satellite species concepts which designate species based on feeding type do not accurately characterize evolutionary lineages warranting recognition as distinct species. Little research has explored mechanisms underlying speciation and evolution of life history diversity in Lampetra of western North America, although parasitic individuals consistently appear within a western brook lamprey population on Vancouver Island in British Columbia (Beamish 1985; Beamish et al. 2016), and the report of a migratory Lampetra that presumably emerged from a resident population in the Columbia River basin (Jolley et al. 2016) suggests that other west coast Lampetra populations may also be polymorphic for feeding type. Research on the paired Lampetra planeri (resident, non-parasitic) and Lampetra fluviatilis (anadromous, parasitic), both of which are of conservation concern in Europe and at a more local level in some countries (Maitland et al. 2015), may provide insights. Numerous studies across western Europe have demonstrated a polyphyletic relationship between these paired species (Espanhol et al. 2007; Blank et al. 2008; Pereira et al. 2010), and interbreeding is common where they are sympatric (Rougemont et al. 2015; Mateus et al. 2016; Hume et al. 2018). In sympatric populations, Lampetra planeri and Lampetra fluviatilis more accurately represent a single species with mixed life histories in the early phases of speciation (Rougemont et al. 2017; Kucheryavyi et al. 2016; Mateus et al. 2016; Hume et al. 2018; see also review by Docker and Potter 2019). Genetic divergence between these species is greater where the migratory life history is absent, allowing populations to evolve in relative isolation (Pereira et al. 2010; Bracken et al. 2015; De Cahsan et al. 2020; Rougemount et al. 2021). These studies and our own demonstrate that it is inaccurate to regard these groups as distinct species. Instead, they more accurately represent life history variants of a single species. Similar variation in divergence, hybridization, and life history diversity have been observed across fish taxa including forms of Oncorhynchus mykiss (e.g., anadromous steelhead versus freshwater-resident rainbow trout; Docker and Heath 2003), Oncorhynchus nerka (e.g., anadromous sockeye salmon versus lake-spawning kokanee versus stream-spawning kokanee; Tigano and Russello 2022), western forms of Rhinichthys cataractae (longnose and Nooksack dace; Taylor et al. 2015), and many forms of Gasterosteus aculeatus (three-spined stickleback; see Gow et al. 2006; Makinen et al. 2006) and Cottus asper (prickly sculpin; Dennenmoser et al. 2015). In these examples, different common names are sometimes used to distinguish among life history variants of the same species. Similarly, continuing to refer to the migratory/parasitic and resident/non-parasitic forms of Lampetra ayresii as western river lamprey and western brook lamprey, respectively, will emphasize the different life histories and ecological requirements of these ecotypes. Recognition of phenotypic and ecological diversity below the species level is important for conservation even if both forms fall under the umbrella of a widespread and abundant species (Docker and Hume 2019).

The four lineages of undescribed Lampetra in our species delimitations analyses were originally identified by Boguski et al. (2012). Although they did not perform a formal species delimitation analysis, the authors hypothesized that these lineages from Oregon (Fourmile Creek and Siuslaw River) and California (Mark West Creek and Kelsey Creek) were distinct species based on their high levels of divergence. Furthermore, the authors noted that within Northern Hemisphere lampreys, species diversity, incidence of endemism, and intraspecific genetic variation for broader-ranging taxa are greater at the southern end of their distributions (see Yamazaki et al. 2006; Mateus et al. 2011; Spice et al. 2012), as it is for other amphidromous species in this region (e.g., Cottus asper; Young et al. 2022). These patterns are broadly associated with the absence of continental glaciation and its associated climate in the southern portion of their range during the Pleistocene (Reyjol et al. 2009; Griffiths 2010; Su et al. 2022), which has allowed a longer period for divergence of lineages. Based on these patterns, Boguski et al. (2012) hypothesized that the Columbia River may represent the northern boundary of the most genetically divergent lampreys in western North America (see Fig. 6). Furthermore, all four of these candidate species were observed in areas of California and coastal Oregon known for unusually high levels of aquatic endemism (Howard et al. 2015; Markle 2019). Because of limited sampling, it is uncertain whether these four lineages are more widely distributed or whether other cryptic taxa with similarly limited ranges are present (the latter a pattern typical among species of Entosphenus in California; Moyle 2002). Future efforts aimed at recovering additional lineages of Lampetra in western North America will likely be the most successful in these regions.

Fig. 6
figure 6

Location and taxonomic classification of specimens based on species delimitation and specimen assignment analyses. The L. ayresii species group includes both anadromous/parasitic and resident/non-parasitic ecotypes, commonly known as western river lamprey and western brook lamprey, respectively

Surprisingly, the distribution of Lampetra pacifica also appears extremely limited despite a potential connection to many stream networks. Initially, this species was thought to occupy a number of river basins from Oregon to California (Vladykov 1973) and then later restricted to an undescribed portion of the Columbia River basin (Reid et al. 2011). Specimens of putative Lampetra pacifica from the Willamette River basin, including the Upper Willamette and McKenzie subbasins, are represented in museum collections (Table S1). Yet the genetic evaluation of specimens of Lampetra from across the Columbia River basin, including the many new samples from the Willamette River basin from the present study, indicates that Lampetra pacifica has only been detected in the Clackamas River (the type of location, a tributary to the Willamette River downstream from Willamette Falls) and in Crystal Springs Creek (a tributary 12 km downstream from the mouth of the Clackamas River within the city limits of Portland, Oregon). Given that Reid et al. (2011) did not observe Lampetra pacifica in the Tualatin River, whose confluence with the Willamette River is only 5 km upstream from the mouth of the Clackamas River but upstream of Willamette Falls, we hypothesize that the falls constitute the upstream limit of the distribution of Lampetra pacifica. The downstream limit remains uncertain. Five museum specimens (four collected in 1963 and one in 2012; Table S1) identified as Lampetra pacifica on the basis of morphology were collected in tributaries to the Columbia River downstream of the confluence of the Willamette River. Whether these represent members of the Lampetra pacifica lineage is unknown, but further evidence supporting that hypothesis is lacking because individuals genotyped from these and nearby basins have all assigned to Lampetra ayresii. Interestingly, the current distribution of this species would have been inundated during the multiple floods associated with ice dam failure of Glacial Lake Missoula and Lake Bonneville draining (O’Connor et al. 2020), suggesting the presence of refugia for this species somewhere beyond the high- water mark for those events. Where these refugia might have been, whether they remain occupied, and whether Lampetra pacifica is found in any other basin warrant further investigation.

Although the species delimitation analysis indicates that Lampetra ayresii is the only Lampetra species present in our study area in the Willamette River, our extensive sampling effort expanded the known diversity of this species. We observed a total of 21 cyt b haplotypes in our study area. Three of these (cyt b H03, H04, and H39) were previously observed by Boguski et al. (2012); the other 18 haplotypes were previously unobserved (Fig. 3, Table S3). Based on our analysis, cyt b H03 and H04 appear localized to Owens Creek in the Upper Willamette River subbasin. In contrast, cyt b H39 was the only haplotype observed in our study area that was also observed in other river basins of western North America. Furthermore, this haplotype was the most widely distributed in the species delimitation dataset, observed in over three dozen individuals ranging from the Farragut River in southeastern Alaska to the Sacramento Delta in central California. Despite its broad distribution, this haplotype was represented by only three individuals captured in Muddy Creek in the Upper Willamette subbasin. Not only was this haplotype rare in our study area, but it was also genetically distinct from all others in the Willamette River basin (Fig. 3, Fig. 5). A similar pattern was observed among members of the Cottus asper species complex. Forms of Cottus asper in the Marys River (to which Muddy Creek is a tributary) were identical either to forms present in two coastal basins (the Alsea and Yaquina rivers) or to those restricted to the Willamette River basin (Young et al. 2022). These forms of Cottus asper and Lampetra ayresii that are diverged from others in the Willamette may have arrived via headwater capture at some point in the mid-to-late Pleistocene (Markle 2019), but more recent introductions by anglers are also plausible.

Apart from the divergent Muddy Creek fish, we observed a high level of genetic structure at the subbasin level, with little or no haplotype divergence among specimens from the same basin (Table 2; Table S6; Fig. 5). Similar patterns of genetic structure have been observed among populations of European Lampetra, particularly freshwater resident Lampetra planeri (Pereira et al. 2010; Bracken et al. 2015; Rougemont et al. 2021). Our results indicate that Lampetra ayresii populations in the Willamette River basin show little evidence of recent dispersal to other subbasins, suggesting that subbasins constitute reasonable units for management and conservation. We acknowledge, however, that the appropriate scale for management in other areas may vary depending on the evolutionary and geological histories of taxa present. For example, both geological and anthropogenic isolation have played a role in structuring Lampetra populations in Europe (Pereira et al. 2010; Bracken et al. 2015; De Chasan et al. 2020; Rougemont et al. 2021). Management units for Lampetra may encompass several adjacent subbasins in areas with reduced species richness and intraspecific diversity. This may be the case in portions of Washington, British Columbia, and Alaska that were more recently glaciated and likely have less species diversity (Fig. 6; see Boguski et al. 2012). In contrast, smaller conservation units may be required in areas with greater species richness, higher rates of endemism, and in locations where lack of gene flow has led to evolution of divergent lineages (e.g., Pereira et al. 2010; Wang et al. 2020).

In summary, we used mitochondrial gene sequences to address questions related to taxonomic classification, species delimitation, and population structure. Although such data is widely and robustly used in this context, it is not without shortcomings, foremost of which is the separate evolutionary history sometimes exhibited by the nuclear genome (Dellicour and Flot 2018). Further explorations of the patterns we observed, both by consideration of nuclear data (Hess et al. 2013, 2020) and much more extensive and intensive field sampling (Bergsten et al. 2012), are prudent first steps for evaluating our conclusions, in large measure because genetic criteria form part of the basis for designating units of conservation under the US Endangered Species Act (Waples 1991, 1995), Canada’s Species at Risk Act (SC 2002), and the United Nations’ Convention on Biological Diversity (CBD 2010). Once these units of conservation are identified, tools such as environmental DNA sampling could rapidly and comprehensively resolve issues associated with their distributions (e.g., Young et al. 2022). We believe that a similar stepwise approach applied to other members of Petromyzontidae across the Northern Hemisphere, many of which are also of uncertain taxonomic standing and ambiguous distribution, could prove useful.