The identification of watermilfoil, discovery of hybrid watermilfoil, and their implications for aquatic plant management in the Clark Fork River, Western MT, USA

Abstract Eurasian watermilfoil (Myriophyllum spicatum) is an invasive submersed macrophyte that has infested many waters in North America since its introduction. Eurasian watermilfoil has the ability to alter the structure and function of the ecosystems that it invades. Eurasian watermilfoil was first surveyed in the Clark Fork River, MT in 2008, alongside the native northern watermilfoil (Myriophyllum sibiricum). Three reservoirs (Cabinet Gorge, Noxon Rapids and Thompson Falls) on the lower Clark Fork River had entire lake surveys conducted using the point intercept method. Morphological data from these surveys showed that only Eurasian watermilfoil and northern watermilfoil were present during the time of the surveys in 2008. The results of the morphological identification were supported by molecular identification at three different laboratories. In 2015, a genetic survey of watermilfoil species was conducted on the Noxon Rapids Reservoir and found that hybrid watermilfoil (Myriophyllum spicatum × sibiricum) was present. This hybridization event poses a number of issues for aquatic plant management. Hybrid watermilfoil is much more difficult to identify morphologically than its parent types. The hybrid is also more invasive and may be differentially susceptible to some herbicides than the parental type Eurasian watermilfoil.

invasive macrophyte in North America that has been introduced multiple times (Smith and Barko 1990;Zuellig and Thum 2012;Moody et al. 2016). This macrophyte causes damage to the ecosystems it infests by affecting water quality, altering nutrient dynamics, and shading out native plants (Smith and Adams 1986;Madsen et al. 1991;Boylen et al. 1999;Madsen et al. 2015). In addition to these effects on the aquatic ecosystem, Eurasian watermilfoil is detrimental to recreation which in turn has additional economic costs (Newroth 1985). Additionally, infestations have reduced lakefront property values in Vermont by up to 16% (Zhang and Boyle 2010).
The spread of Eurasian watermilfoil is prolific and asexual propagation via fragmentation allows Eurasian watermilfoil to be transported to new lakes rapidly (Smith and Barko 1990). In lakes it infests, Eurasian watermilfoil can form dense canopies at the surface of the water. These canopies allow it to outcompete the native submersed macrophytes, which can lower plant diversity in these lakes (Madsen et al. 1991). Eurasian watermilfoil provides poorer quality habitat for wildlife in waters where it is invasive than the native plant communities (Parsons et al. 2011;Alford and Rozas 2019). An infestation of Eurasian watermilfoil can also change the biogeochemical cycles of lakes. Eurasian watermilfoil can serve as a pump to move sediment-bound phosphorous into the water column, which can alter the community structure of that lake (Smith and Adams 1986). Although Eurasian watermilfoil frequently spreads through vegetative fragments, sexual reproduction can also be an important mechanism for its spread (Smith and Barko 1990;Zuellig and Thum 2012). With the ability for sexual reproduction, Eurasian watermilfoil has hybridized with the native northern watermilfoil (Myriophyllum sibiricum Kom.) on multiple occasions Les 2002, 2007). This hybrid watermilfoil (Myriophyllum spicatum Â M. sibiricum) is highly invasive in North America and sexually viable .
The hybridization between Eurasian and northern watermilfoil poses problems with identification. Typically, Eurasian watermilfoil is differentiated from northern watermilfoil using distinct morphological characteristics. Turions are present on northern watermilfoil but are absent on Eurasian watermilfoil; however, the turions are only present in autumn (Aiken et al. 1979). The typical method to differentiate the two species is to count the number of pinnae (leaflets) on a leaf (Aiken et al. 1979). This method however, becomes less reliable when hybrids are present because they often have intermediate pinnae numbers when compared to the parent populations (Moody and Les 2007). Although morphological identification of Eurasian watermilfoil, northern watermilfoil and their hybrids can still be successful, molecular identification methods are much more reliable (Moody and Les 2002;Pashnick and Thum 2020).
The Clark Fork River is an important water source for multiple reservoirs in Montana that then drains into Lake Pend Oreille in the Idaho panhandle. In 2008, three major reservoirs in the lower Clark Fork River system were surveyed for aquatic invasive macrophytes . These three reservoirs were Cabinet Gorge, Noxon Rapids and Thompson Falls. The surveyors found northern watermilfoil at all three sites and Eurasian watermilfoil at Cabinet Gorge and Noxon Rapids reservoir . In both Cabinet Gorge and Noxon Rapids reservoir, there was no evidence of hybridization between Eurasian watermilfoil and northern watermilfoil in 2008 ). The objectives of this article are to quantify the macrophyte community in the Lower Clark Fork system; quantitatively demonstrate that morphological characteristics can be used to identify milfoils; and demonstrate the rapidity at which watermilfoils can hybridize under field conditions, which in turn can lend insights into population shifts and ultimately management of hybrid watermilfoils.

Study sites
The lower Clark Fork River system is in western Montana and contains multiple reservoirs. The river flows west-northwest and drains into Lake Pend Oreille in Idaho. Eurasian watermilfoil was first documented in this section of the river in 2008, cohabiting with northern watermilfoil. For this study, three reservoirs in the lower Clark Fork River were surveyed: Thompson Falls Reservoir, Noxon Rapids Reservoir and Cabinet Gorge Reservoir.

Vegetation sampling
In 2008, entire reservoir surveys of the aquatic macrophyte communities in all three reservoirs in the lower Clark Fork River were conducted. The surveys were conducted using the point intercept method (Madsen et al. 2015;Madsen and Wersal 2017). The points in these reservoirs were arranged in a systematic grid of 150 m, 250 m and 150 m for Cabinet Gorge Reservoir (n ¼ 334), Noxon Rapids Reservoir (n ¼ 497) and Thompson Falls Reservoir (n ¼ 40), respectively (Figures 1-3). Macrophytes were sampled by tossing a rake and pulling the vegetation to the surface. If no plants were found on the first rake toss, then one more rake toss was made to ensure that plants were not present at that site. Any additional species visible at the point were also recorded. The macrophyte species in the sample were identified on site, and presence (1) and absence (0) data were used to determine the frequency of the macrophytes in the reservoirs. The mean species richness for each lake was calculated by averaging the number of all species present at each point. The same was also done for native species and non-native species.
In addition to the quantitative community surveys, six of the survey points from each reservoir were selected and samples of watermilfoils were harvested using the rake method for laboratory analyses (genetic), and additional specimens were pressed for morphological analysis. If possible, only points with both Eurasian and northern watermilfoil were selected and the species were identified visually, on site. In the Thompson Falls reservoir, only northern watermilfoil was present, so there were no Eurasian watermilfoil samples to analyze. The samples harvested for genetic analysis were stored on ice and shipped to three university laboratories (Bowling Green State University, Grand Valley State University and Mississippi State University). Noxon Rapids Reservoir was surveyed again in 2015 and milfoil populations were sampled at 40 locations throughout the reservoir. Collected samples were sent to Montana State University for genetic analysis.

Morphological analysis
Pressed specimens were analyzed at Mississippi State University. For each pressed specimen, the stem color was coded as green or red; the apical meristem was recorded as rounded or flat, and leaf tips were recorded as rounded or flat. For each specimen sample, six node/internode combinations were selected beginning 210 mm below the apex, the internode length was measured, the number of leaflets from one leaf per node counted, the length of the leaf measured, and the length of the leaflet measured. Stem thickness was measured at the middle of the internode for each internode interval. All measurements were made in mm. These metrics were used to produce a framework for the differentiation of Eurasian watermilfoil and northern water milfoil based on morphological characteristics. A Fisher's exact test was used to statistically differentiate morphological characteristics of Eurasian and northern watermilfoil. All tests were conducted using a p 0.05 significance level.

Genetic analysis
Each university lab (Bowling Green State University, Grand Valley State University and Mississippi State University) received a total of 30 plants samples (18 samples of northern watermilfoil, and 12 of Eurasian watermilfoil) that were collected from Thompson Falls (northern watermilfoil only), Noxon Rapids and Cabinet Gorge Reservoirs in 2008. At the Bowling Green University lab, total genomic DNA was extracted from these plant tissues using DNeasy kit (Qiagen, CA). Measurement of DNA concentration was performed by a spectrophotometer (NanoDrop ND-1000, Thermoscientific). PCR amplification of the nuclear ribosomal DNA of the 3 0 end of the 18S-like gene to the 5 0 end of the 28S-like gene was performed using primers ITS4: 5 0 TCCTCCGCTTATTGATATGC3 0 and ITS5: 5 0 GGAAGTAAAAGTCGTAACAAGG3 0 (White et al. 1990), 50-100 ng genomic DNA and Phusion TM High-Fidelity DNA Polymerase (NEB, MA). The amplified products were subjected to electrophoresis on 1.5% agarose gels containing 5.0 ug/ml of ethidium bromide. The amplicons were purified using a QIAquick PCR purification kit (QIAGEN, CA). The sample DNA was sequenced using Basic Local Alignment Search Tool (BLAST).
At Grand Valley State University, phylogenetic analysis of DNA sequence data was used to identify both Eurasian and northern watermilfoil, and to identify potential hybrids. Specifically, the internal transcribed spacers 1 and 2 (ITS) were sequenced and then compared to samples of known identity (Moody and Les 2002). In addition, PCR was performed on a chloroplast gene, trnL-F, which was compared to known sequences from GenBank as well as sequences from collections maintained at Grand Valley State University (Taberlet et al. 1991).
The lab at Mississippi State University used the Polymerase chain reaction restriction fragment length polymorphisms (PCR-RFLPs) method (Grafe et al. 2015). A PCR-RFLP analysis involved the PCR amplification of DNA that was variable across individuals of interest. The PCR product was then digested with a restriction enzyme that targeted polymorphisms. That is, genetic divergence between target taxa would have resulted in a gain or loss of a restriction enzyme cut site. The digested PCR product was run on agarose gels, and different restriction enzyme profiles became apparent (Grafe et al. 2015). Ribosomal ITS sequences for both Eurasian and northern watermilfoil were downloaded from NCBI's Genbank (Accession #'s DQ786012-DQ786027). Sequences were then analyzed and compared. Restriction enzyme cut sites that were present in representative sequences of one species, and absent in representative sequences of the other were targeted for marker development. One restriction enzyme presented itself as particularly useful for the purpose of PCR-RFLPs; HhaI was found to cut specific segments that differentiate between Eurasian and northern watermilfoil. HhaI cut ribosomal ITS from Eurasian watermilfoil at three restriction cut sites that were absent in northern watermilfoil. A 239 bp fragment was formed when digesting northern watermilfoil ribosomal ITS with HhaI. This 239 bp fragment was absent after digesting Eurasian watermilfoil ribosomal ITS with HhaI. Ribosomal ITS from samples collected from Montana lakes were PCR amplified, and the amplicons were digested with HhaI. Digests were then run on 4% metaphor agarose gels, and restriction fragment profiles were scored visually.

Vegetation sampling
During the 2008 survey, Eurasian watermilfoil was found at 6.3 and 3.3% of the survey points in Cabinet Gorge and Noxon Rapids Reservoirs, respectively (Tables 1 and 2; Figures 4 and 5). Eurasian watermilfoil was not found in Thompson Falls Reservoir during the 2008 survey (Table 3). Additionally, northern water milfoil was observed at 5.  The mean frequency was calculated as the average presence of an individual species across the sample points. The mean richness values were calculated as the average number of individuals across the sample points. The mean frequency was calculated as the average presence of an individual species across the sample points. The mean richness values were calculated as the average number of individuals across the sample points.

Morphological analysis
In comparing Boolean characteristics of Eurasian and northern watermilfoil, there was a difference (p < 0.01) in flattened versus round leaf ends and meristems (Table 4). All but one of the Eurasian watermilfoil specimens had a flattened meristem, while all but four of the northern watermilfoil specimens had rounded apical meristems. Stem color was not different (p ¼ 0.26) between Eurasian watermilfoil and northern watermilfoil and thus would not be a good characteristic for identifying these species. There was no difference (p ¼ 0.44) in stem thickness between Eurasian and northern watermilfoil (Table 5). Leaflet pair number, leaf length, leaflet length and internode length were all different when comparing specimens of Eurasian and northern watermilfoil. On average Eurasian watermilfoil, specimens had 16.3 ± 0.2 leaflet pairs per leaf whereas norther watermilfoil specimens had 8.1 ± 0.1 leaflet pairs per leaf. Based the morphological analysis the three most useful characteristics for differentiating northern watermilfoil and Eurasian watermilfoil were the rounded apical meristem, rounded leaf tip and leaflet pairs. Both leaflet length and leaf length were also different between northern watermilfoil and Eurasian watermilfoil. However, the mean lengths for both leaflet length and leaf length are relatively close in value; and it would be difficult to use these characteristics for differentiating between northern watermilfoil and Eurasian watermilfoil while in the field.

Genetic analysis
All three university labs were successful in identifying and distinguishing between milfoil species using their respective genetic analysis methods. At Bowling Green University, the PCR products of rDNA sequences of both species ranged from 700 to 1000 nucleotidelong. Subsequent sequencing and BLAST analyses revealed that the Clark Fork River sequences exhibited 94-99% similarity to sequences publicly available through NCBI for both Eurasian and northern watermilfoil. Based on this analysis, a DNA methodology was successfully developed to differentiate Eurasian watermilfoil from northern watermilfoil, and no evidence of hybridization was found. Analyses conducted at Grand Valley State University based on both ITS and trnL-F for all field milfoil samples were able to differentiate both species of milfoil and indicated that genetic identifications were identical to morphological identification made in the field. Furthermore, Grand Valley State University found no evidence of hybridization from samples collected during the 2008 survey. Likewise, analyses by Mississippi State University found that in all cases, field identifications of milfoil species were confirmed by genetic analysis, and for all samples, using all three genetic methods, there was no evidence of any hybrid milfoil during the 2008 sampling event.
Although there was no evidence of hybridization in 2008, within seven years extensive hybridization had occurred in Noxon Rapids Reservoir. Of the milfoil samples that were collected in Noxon Rapids Reservoir in 2015, 39% of these samples were hybrids (Myriophyllum spicatum Â sibiricum) (Figure 6). These hybrids exhibited physical characteristics of both northern and Eurasian watermilfoil. Hybrid milfoil exhibiting intermediate characteristics would make field identification more difficult and reliance on genetic testing would be necessary as hybrids begin to spread and cross. The mean frequency was calculated as the average presence of an individual species across the sample points. The mean richness values were calculated as the average number of individuals across the sample points.

Discussion
The current study has defined and verified a reliable framework for the morphological differentiation of northern watermilfoil from Eurasian watermilfoil with identifications verified molecularly. The most commonly used characteristic for differentiating Eurasian and northern watermilfoil is the number of leaflets pairs per leaf. Additionally, identification reliability increases when the leaflet pair number is used in tandem with the shape of the leaf ends and apical meristems. However, hybrids often express intermediate phenotypes for leaflet numbers (Moody and Les 2007). This can make the identification of watermilfoil more difficult when hybrid watermilfoil is present with northern and Eurasian watermilfoil.
The absence of hybrid watermilfoil in the Noxon Rapids reservoir in 2008 and its subsequent discovery in 2015 is likely a product of one of two scenarios: (1) two genotypes of Eurasian watermilfoil and northern watermilfoil in the Noxon Rapids Reservoir hybridized and became established in the reservoir, or (2) the hybridization of the watermilfoil species occurred in another system and was transported and established in Noxon Rapids Reservoir. If the first scenario is the most accurate, then that means that Eurasian watermilfoil and northern watermilfoil had hybridized in fewer than seven years, which is shorter than the life of many lake management plans. In either scenario, the presence of hybrid watermilfoil in the Clark Fork River system may affect the structure and function of the aquatic community.
The hybridization of milfoils can cause potential ecological issues. According to the 'Genes to Ecosystem Cascade' model, changes in an invader's genotype from admixture can compound in scale and ultimately result in the alteration of the introduced ecosystem's structure and function (Molofsky et al. 2014). This relationship was initially studied in reed canary grass (Phalaris arundinacea L.) (Lavergne and Molofsky 2007;Molofsky et al. 2014). Like reed canary grass, Eurasian watermilfoil can also grow in extremely high densities where it is introduced. The admixture in reed canary grass occurs between different genotypes of the same species. Conversely, admixture in populations of Myriophyllum occurs predominantly between northern and Eurasian watermilfoil. The interspecies admixture is well-documented, and there is laboratory evidence of intraspecific hybridization between two distinct genotypes of Eurasian watermilfoil (Thum and Mcnair 2018). Intraspecific hybridization has occurred in the lab with Eurasian watermilfoil, but to date has not been documented in field populations. Given current molecular markers it would be difficult to detect this in field populations. Considering lab and field data on hybridization of milfoils, the Genes to Ecosystem Cascade model could also apply to Eurasian watermilfoil.
Genetic variation and admixture may also be driving forces for invasiveness in exotic macrophytes. Previous work on Phalaris arundinacea has shown that admixture between different genotypes has produced more phenotypically plastic individuals with greater rates of vegetative colonization (Lavergne and Molofsky 2007). The introduction of novel genotypes could alleviate genetic bottlenecking, raise adaptive potential, and ultimately, cause greater detriment on the invaded system (Lavergne and Molofsky 2007). Many instances of hybridization between Eurasian and northern watermilfoil have been recorded (Moody and Les 2002;Sturtevant et al. 2009;Zuellig and Thum 2012). Few studies have assessed the ecological fitness of hybrid milfoils. However, preliminary data has found that hybrid watermilfoil grows faster than its parental Eurasian watermilfoil, in a laboratory setting Taylor et al. 2017).
Watermilfoil hybridization also poses issues with aquatic plant management through reduced accuracy of identification using morphological characteristics between the parental strains of milfoil when the hybrid is present (Parks et al. 2016). If hybrids cannot be readily identified in the field, it will necessitate the development of a rapid genetic test to verify field populations. Control of hybrid watermilfoil can also be more difficult. Some genotypes of hybrid watermilfoil show differential susceptibility to some herbicides when compared to Eurasian watermilfoil. The hybrids have expressed differential susceptibility to the common herbicide 2,4-D in laboratory and field experiments even when used in tandem with triclopyr (Glomski and Netherland 2010;Parks et al. 2016). Additionally, hybrid watermilfoil is also differentially susceptible to fluridone, another common herbicide, as shown in the lab and the field Berger et al. 2015). A recent study also found that different hybrid genotypes responded differently to 2,4-D treatments (Taylor et al. 2017). These findings suggest that future management projects may need to identify the genotypes of hybrid in the system in order to maximize chemical control.
Although the morphological framework for watermilfoil identification was verified in the Clark Fork River system, additional work is needed to determine if morphological characteristics can be used to determine if a hybrid has been found. This may include the development of common characteristic classes. Future studies should also examine the ecological factors that mediate hybridization as well as fitness characteristics of different hybrid genotypes compared to parental strains. Finally, more studies are needed to document genotypic responses of hybrid milfoils to chemical management techniques. These later studies will be critical for resource managers in developing effective management plans based on the hybrid genotype present in their particular waterbody.

Disclosure statement
No potential conflict of interest was reported by the author(s).

Data access statement
The data that support the findings of this study are available from the corresponding author, [JDM], upon reasonable request. Mark E. Welch, PhD, is a Professor of Biological Sciences at Mississippi State University in Starkville, Mississippi. His research foci include the evolutionary genomics and conservation genetics of plants and animals. In his recent work, he has shown that cis-regulatory elements may facilitate rapid adaptive change in plants, and that inbreeding depression may be a powerful force influencing the dynamics of small populations.

Notes on contributors
Vipaporn Phuntumart, PhD, is a Professor at Bowling Green State University, Bowling Green, Ohio, USA. Her work explores the molecular interactions between plants and microbes, both pathogenic microbes and beneficial microbes.