The Appalbees menu: a multiyear, multilocus metagenetic assessment of pollen foraging by Appalachian Bombus affinis workers

Sampling

Appalachian pollen samples (n = 86) were collected by Metamorphic Ecological Research and Consulting, LLC, spanning June–August of 2021–2023. Sampling in 2021 and 2022 began in late June whereas in 2023 a small number of samples were obtained in early June. Three additional samples from Illinois and Wisconsin were also included to buttress comparisons between Midwestern and central Appalachian forage diversity; these samples were collected in August 2023 by Environmental Solutions & Innovations, Inc. Collections were made under the following permits where required by state and federal agencies: eS61005D (US Fish and Wildlife Service [USFWS]); 2021.085, 2022.113, 2023.196 (West Virginia Department of Natural Resources); 4200 (Monongahela National Forest); 2600, 2670 (George Washington and Jefferson National Forests); ES02373A-16 (USFWS).

Appalachian survey locations were chosen with the goals of documenting populations not previously known as well as for feasibility of bee identification and capture; these considerations favored roadsides, fields, and forest edges and are further detailed in Hepner et al. (2024). B. affinis workers were netted while foraging on a host plant identified at least to genus, with the exception of a single individual opportunistically netted in flight. Bees were chilled on ice, allowing a single corbicular pollen ball to be removed with forceps and placed in a labeled microcentrifuge tube prior to releasing the bee. Forceps were cleaned with alcohol between sampling events. Sampled bees were photographed for verification and vouchering. Appalachian samples were stored dry whereas Midwestern samples were stored in ethanol (storage method has minimal effect on taxonomic characterization (Quaresma et al., 2021)).

The Midwestern samples were strictly opportunistic samples that were provided to us by the USFWS and had been collected during permitted B. affinis surveys unrelated to the Appalachian survey described here. While the three samples do not represent any coherent effort to obtain pollen samples from a specific region or environment, they nonetheless contain information about B. affinis forage behavior. Any comparison beyond descriptive statistics between the two regions is obviously constrained by the large difference in sample size; here we limit our comparison to three topics: (1) the potential effects of sample-level stochasticity on concordance between internal transcribed spacer 1 (ITS1) and internal transcribed spacer 2 (ITS2) detections (i.e., concordance in a small versus a large sample set) (2) calculating the proportional overlap that the current Midwestern and Appalachian pollen samples have with external data sets, particularly the large Midwestern survey of Wolf et al. (2022), and (3) including the Midwestern samples and the Appalachian samples as two of seven independent data sources, including five previously published data sets, from which a ‘consensus’ forage list was generated. This list was based only on positive detection of a genus in at least one data set from each region (discussed further below).

Library preparation and sequencing

Pollen DNA was extracted using the Plant DNeasy kit (Qiagen) with mechanical disruption in lysis buffer. Disruption was performed with MP Biomedicals lysing matrix A in two mL centrifuge tubes, shaken with a FastPrep-24 5G Bead Beating Grinder (MP Biomedicals) for 40 s at 6 m/s. Two primer sets were used to amplify the ITS1 and ITS2 genetic loci via polymerase chain reaction (PCR). The ITS1 primers consisted of primer ITS5A from Stanford, Harden & Parks (2000) and the reverse complement of primer ITS-3p62plF1 from Kolter & Gemeinholzer (2021). The ITS2 primers consisted of ITS-3p62plF1 and ITS-4unR1 from Kolter & Gemeinholzer (2021). These locus designations are approximate, however, as the ITS-3p62plF1 primer is located near the 3′ terminus of the ITS1 sequence. Thus, the “ITS1” amplicon lacks a short terminal region of ITS1 sequence that is included in the “ITS2” amplicon (along with the highly conserved 5.8S rRNA sequence).

The base amplicon was generated from genomic DNA for each locus using a low cycle-number “preamplification” reaction with the primers cited above. PCRs consisted of 3 µl DNA solution, 0.2 µl Promega GoTaq, five µl Promega 5X reaction buffer, 2.5 µl primer mixture (one µM stock concentration), 1.5 µl Promega magnesium chloride solution, 2.5 µl deoxyribonucleotide triphosphates (0.25 mM stock concentration) and 10.3 µl molecular-grade water. Preamplification PCRs were performed in triplicate and pooled to minimize stochastic variation, then cleaned of residual reagents using the Qiagen Qiaquick 96-well PCR Purification kits for vacuum or centrifuge. Reactions were run in 96-well thermocyclers with heated lids, using the following protocol: an initial denaturing step at 95 °C for 2 min, followed by 12 cycles of 95 °C for 30 s, 50 °C for 30 s, and 72 °C for 30 s. A final extension step of 72 °C was implemented for 2 min, followed by an indefinite hold at 4 °C.

Cleaned preamplification products were then input into a second PCR to append exogenous adapter sequences for sequencing on the Illumina platform, using primers modified as described in the manufacturer’s protocols (Illumina document 15044223 Rev. B, available at https://support.illumina.com/documents/documentation/chemistry_documentation/16s/16s-metagenomic-library-prep-guide-15044223-b.pdf). Four equimolar variants of the modified primers were used in this reaction, the minimal primer described in the Illumina protocol and three additional primers containing one to three completely degenerate bases inserted between the Illumina adapter and the locus-specific sequence. These additional linker bases generate phase diversity (Wu et al., 2015) that improves the calibration of the sequencing run by moderating signal intensity at highly conserved positions within a barcode. Reactions consisted of three µl cleaned preamplification product, 2.5 µl of 10 µM primer mixtures, 12.5 µl KAPA2G DNA polymerase reaction mix (Roche), 0.25 µl bovine serum albumin solution, and 6.75 µl molecular grade water. Thermocycler conditions were as described above except the cycle number was increased to 24 and the extension time was increased to 40 s. These reactions were performed in duplicate, then pooled and cleaned as described above.

The final amplification step of the amplicon libraries appended sample-specific indexing oligonucleotides using Illumina Nextera CD combinatorial dual index kits. These reactions followed the manufacturer’s protocol without modification and were cleaned of residual reagents as described above. Prior to the indexing PCR, amplicons of the two primer sets were quantified with a Qubit Model 2 fluorometer and Broad-Range DNA Assay (Invitrogen) and then pooled to approximate equal concentrations within samples. Libraries were then diluted as per manufacturer’s protocol and sequenced with Illumina MiSeq v.3 cartridges to generate 250-bp paired-end reads. Sequencing output for biological samples was deposited in the Sequence Read Archive (SRA) of the National Center for Biotechnology Information (NCBI) under BioProject PRJNA1235776.

Reference database generation

A regional reference database was generated to support the assignment of metabarcode sequences to plant taxa of origin. The database was generated by subsetting the PLANiTS ITS sequence database (Banchi et al., 2020) to include only genera listed to occur in the sampling region in the USDA Plants database (available at https://plants.usda.gov/). The sampling region was broadly construed so that the database would be generalizable to similar applications in the Mid-Atlantic and Appalachian regions of the US. Furthermore, because plant distributions are imperfectly known at the state level, overly narrow geographic limits can exclude taxa that are actually present. Therefore, we selected all genera listed as present within any of the following ten US states: Virginia, West Virginia, Maryland, Pennsylvania, North Carolina, Tennessee, Ohio, Kentucky, Delaware, and New Jersey. PLANiTS sequences belonging to these genera were selected regardless of specific epithet, provided the species was recognized in the NCBI taxonomic database (as official taxon or recognized synonym). Geographically unexpected species within geographically expected genera can therefore be included in the database; this is desirable because invasive and horticultural congeners are often absent from geographic checklists and because including many representatives of a genus should improve discrimination at the genus level, the primary level of analysis in this study (discussed further below).

Database sequences were screened for fungal contamination, which can occur with these primers (e.g., Cornman et al., 2015), by classifying them with the UNITE fungal database (Kõljalg et al., 2019), which is pre-formatted for SINTAX (Edgar, 2016) and was accessed on 02/19/25. This contamination screen resulted in the removal of two reference sequences. An additional 31 reference sequences were obtained from NCBI representing taxa that were indicated by the data but absent from the PLANiTS database (i.e., by searching abundant unassigned operational taxonomic units (OTUs) against the comprehensive nucleotide database of NCBI). Two taxonomic names were replaced with accepted synonyms to avoid artificial family-level assignments: Menziesia pilosa was equated to Rhododendron pilosum (Craven, 2011) and Chrysoma pauciflosculosa equated with Solidago pauciflosculosa (Semple et al., 2023). Reference sequences shorter than 250 bp were then removed as insufficiently complete, and identical sequences with the same taxonomy were de-replicated with the “derep_fulllength” command of vsearch v. 2.21.0 (Rognes et al., 2016). The final reference database included 6,107 sequences representing 1,125 genera and 2,751 species (File S1); sequences were converted to uppercase if necessary because lowercase characters can be interpreted as masked by some software. ITSx v. 1.1.3 (Bengtsson-Palme et al., 2013) and the bbmap package (available at https://sourceforge.net/projects/bbmap) were used to extract ITS sequences and expected amplicon sequences, respectively, from the reference database for inspection, such as to determine whether taxa detected at only one primer set were fully represented in the database.

Bioinformatic processing

Read pairs were trimmed of low-quality bases using the bbduk program of the bbmap package, using a quality trim value of 10 (Phred-scaled). Trimmed read pairs were merged with vsearch specifying a minimum overlap of 65 bp and a maximum percent difference of 8% within the overlap. Unmerged pairs and merged pairs less than 300 bp in length were discarded (the mode amplicon lengths were 329 bp and 345 bp for ITS1 and ITS2 primer sets, respectively). OTUs were generated by clustering merged reads at 99% with vsearch using similarity definition “1” of that program. Primer sequences were removed from the OTUs prior to taxonomic assignment with bbduk using a search kmer value of 11.

A taxonomy was assigned to each OTU with the vsearch SINTAX command using the reference database described above. The SINTAX method assigns a confidence score to each rank of the assignment; we assigned the OTU the lowest rank of the SINTAX assignment hierarchy that had a confidence score of at least 0.9 (possible range is 0–1 with higher scores indicating higher confidence). The number of reads attributed to each taxon within each sample was the sum of cluster sizes for all OTUs with the same assignment (multiple OTUs may have the same taxonomic assignment within a single sample). Unfiltered taxon counts and counts filtered, scaled to counts per million (cpm), and aggregated at the genus level are provided in File S2. By policy, the raw data are also contained within a U.S. Geological Survey (USGS) data release (Cornman & Otto, 2025). As the analyses presented here are based on prevalence and Shannon’s diversity, which depend only on proportions, no log-ratio transformation was performed to open the inherently closed, compositional nature of the data (reviewed by Gloor et al., 2017).

A common filtering practice (Alberdi et al., 2018; Drake et al., 2022) is to censor low count values within each sample (i.e., convert those cells of the taxon matrix to missing or zero) rather than count them as detections, because low taxon counts within individual samples may arise through erroneous index sequences (referred to as ‘demultiplex error’ or ‘cross-talk’, e.g., Edgar, 2018; MacConaill et al., 2018). For both loci, sequence counts in the negative control sample were consistent with known rates of demultiplex error for Nextera eight-base dual indexes (MacConaill et al., 2018) for both loci, indicating no systemic contamination occurred (Fig. S1). An additional consideration is that taxa that comprise a small proportion of a sample total could potentially arise from ‘natural’ contamination, i.e., wind- or insect-mediated transfer of heterospecific pollen. The frequency and mechanisms of heterospecific pollen incorporation in B. affinis corbicular pollen loads are unknown and likely contingent on factors such as weather, plant community, and insect community (discussed further below). Given our primary goal of identifying common pollen forage sources relevant to the management of B. affinis habitat and the ambiguity as to whether plants at low abundance in a sample were actively foraged, we limited positive detections in a sample to those plant taxa with at least ten sequence reads and constituting at least 0.1% of the total plant-assigned sequence reads. Counts in samples for taxa with nonzero counts in the negative control were not further adjusted or removed, as the rationale for these thresholds includes mitigating the cross-talk detected in the negative control.

After this censoring, counts were aggregated at the genus level for further analyses. Because of the high number of rare genera detected, some analyses were limited to genera detected in at least three samples (see ‘Results’). Shannon’s diversity index was estimated with PAST v. 4.17 (Hammer, Harper & Ryan, 2001), which was also used to perform nonparametric analyses of variance.

Overlap with genera detected in previous studies of B. affinis foraging was assessed by referencing Table S2 of Simanonok et al. (2021), Table 1 of Wolf et al. (2022), Table 2 and Table 3 of Simanonok et al. (2024), and Figure 2 of Hepner et al. (2024), as summarized in File S3. Taxa above genus rank in those references were excluded. Wood et al. (2019) examined corbicular pollen loads of multiple Bombus species including B. affinis, but did not list these by species and thus could not be included in this comparison.

Evaluation of bioinformatic parameters with simulated data

To evaluate merge parameters and taxonomic assignment accuracy, paired 250-bp sequence reads were simulated for the ITS2 locus with ART (Huang et al., 2012) to generate 1000X coverage of a modified version of the reference database. References were oriented uniformly and trimmed so that the ITS2 forward primer sequence was 5′; the 3′ end was not specifically trimmed to the reverse primer so that partial reference sequences could be included. References with greater than 1% ambiguous characters or lacking an explicit species-level taxonomic assignment were excluded. Illumina MiSeq amplicon error profiles were modeled by selecting settings “MSv3” and “amp” in ART. Simulated reads were then trimmed, merged, and assigned as described above but with key parameters systematically adjusted.

Merge rates for simulated reads were determined for 54 combinations of overlap and maximum percent mismatch values, in which the overlap incremented from 60 to 100 in steps of 5 and the mismatch value incremented from 5 to 10 in steps of 1. Little difference was seen among these combinations, with merge rates greater than 98% for all parameter combinations.

Taxonomic error rates were calculated for three SINTAX scores: 0.8, 0.85, and 0.9 (Fig. S2). At a threshold of 0.9, 92.5% of species assignments were correct and 0.26% incorrect, with the balance unassigned. At a threshold of 0.8, 95.2% of species assignments were correct and 0.37% incorrect, with the balance unassigned. At the genus level, greater than 99% of assignments were correct at all thresholds and error rates were an order of magnitude lower than at the species level (ranging from 0.0259%–0.0538%). Furthermore, genus-level error rates were concentrated in only a few species. For example, at a threshold of 0.9, the uncommon crop species Momordia charantia and two Centaurea species accounted for over 90% of all genus-level errors. On the other hand, the simulation results likely underestimate the true error rate because biological samples often contain valid but unrepresented sequence variation for a taxon, and often include taxa that are missing from or misassigned in the reference database. Nonetheless, the simulated data support the bioinformatic parameters chosen and indicate a high-level of confidence in taxonomic assignments when high-quality reference sequences for a taxon are available, at least for the ITS2 locus. The observation that simulated sequences were much more likely to be unassigned than incorrectly assigned suggests that ITS1 sequences assigned at a conservative (high) threshold are also unlikely to generate significant genus-level false positives.

Adaptive targeting of focal plants recovers high pollen diversity

Appalachian samples were collected from early June to late August over three consecutive years (2021–2023), necessarily clustered around periods of greater survey effort and greater worker abundance. A total of 89 corbicular pollen samples were sequenced, 86 from Appalachian sites and three from comparison sites in the upper Midwest (Fig. 1). The number of sequences assigned to plant taxa ranged from 87–144,221 (mean 44,385) for ITS1 and 34–84,144 (mean 40,321) for ITS2. Although most sequences were assigned at the species level (84.5% for ITS1 and 76.8% for ITS2), taxonomic assignments are analyzed herein at the genus level to maximize both the amount and robustness of data included (see simulated error rates in Materials and Methods). The minimum number of assigned sequences for inclusion of a sample in analyses was 5,000; this threshold excluded two samples at the ITS1 locus and one sample at the ITS2 locus. After these exclusions, the mean genus-level richness detected per sample was 3.40 (2.72 SD) at ITS1 and 6.08 (3.50 SD) at ITS2. While many of the additional ITS2 detections are taxa with low relative abundance and low prevalence, rarefaction to a fixed sampling effort of 10,000 sequence reads further demonstrates the consistently greater genus diversity detected at the ITS2 locus (Fig. 2).

The proportion of sequences consistent with pollen of the host plant from which the worker was captured was highly variable among samples but overall comprised roughly three-quarters of assigned sequence reads (Fig. 3, File S2). The mean proportion of non-host reads was 0.190 for ITS1 and 0.273 for ITS2. Notably, all four samples taken from bees visiting Monarda fistulosa (wild bergamot) failed to yield sequences of that host plant above the threshold detection level for both loci. This is likely a biological rather than methodological phenomenon, as bumblebees are frequently observed robbing Monarda of nectar in a manner that precludes contact with anthers (Wolf et al., 2022) and Monarda metabarcode sequences have been detected previously from honeybee-collected pollen with comparable primers (Smart et al., 2017; Otto et al., 2020).

Comparison of plant diversity detected at each locus

In total, 69 plant genera were detected in Appalachian samples, 45 (65.2%) of which were detected at both loci (Fig. 4A). Seven genera were detected only at ITS1 and 17 genera were detected only at ITS2 (File S2). For the three Midwestern samples, eight of 21 total genera were detected at both loci (38.1%), two were detected only at ITS1 and 11 were detected only at ITS2. The total number of detections summed across samples was 268 at ITS1 and 469 at ITS2.

Many of the differences between loci appear to be due to stochasticity in the detection of rarer taxa per se rather than intrinsic biases between loci. For example, five of the genera detected at only one locus in Midwestern samples (Hypericum, Rudbeckia, Senecio, Impatiens, and Ageratina; Fig. S2) were also detected in Appalachian samples, but in that larger data set three of the five were detected at both loci (Hypericum, Rudbeckia, and Impatiens). Moreover, the single Ageratina reference sequence did not include the region amplified by the ITS1 primer set and therefore that genus would likely have been assigned at a higher rank if amplified from samples. More generally, genera detected at only one locus in Appalachian samples mostly had low prevalence at that locus (Fig. 4A), with the exception of Tilia (Ageratina is uninformative in this regard as it lacked ITS1 reference sequence). While detected at both loci, Daucus was notably more prevalent at ITS1 than ITS2.

The two loci were also concordant with respect to relative abundances of sequences attributed to each genus, again declining somewhat for rarer taxa as expected with sampling stochasticity (Fig. 4B). The taxa that strongly differed in prevalence also strongly differed in relative abundance: ITS2 primers amplified much less Daucus than ITS1 overall, whereas Tilia was abundant in ITS2 sequences but absent in ITS1 sequences. The cpm values for genera were highly correlated when both were detected (Pearson’s r = 0.999, P = 4.77E−61). While the primer sets were broadly consistent in this regard (with the stated exceptions), consistency does not imply that taxon relative abundances are unbiased, but does suggest that metrics based on relative abundance, such as frequency-based diversity indices, would also be reasonably consistent between the two loci. For example, Shannon’s H was well correlated between the two loci (Pearson’s r = 0.753, P = 2.8E−7) after removing very low diversity samples (H < 0.1 at either locus).

Forage patterns in Appalachian samples

Sorting Appalachian samples chronologically (Figs. 5–6) highlighted a forage pattern in which diversity was lower and more consistent in early summer, albeit somewhat confounded by the tendency for more captures during this period, particularly on Hydrangea host plants. Note the size of the ellipses in Fig. 5 are natural-log scaled to increase the visibility of proportionally rare taxa, and only taxa detected in at least three samples are shown. Genus diversity was higher and more variable among samples in later summer (Figs. 5–6), a pattern observed each year. When samples were binned by calendar month across years (Fig. S3), median Shannon’s H differed significantly by month (Kruskal–Wallis test, P = 0.00031 for ITS1 and P = 0.0011 for ITS2). For both loci, August had significantly higher median Shannon’s H than June and July whereas the latter two months did not significantly differ in pairwise comparisons (see Fig. S3 for Bonferroni-corrected P-values).

Figure 5 also illustrates a relative paucity of Fabaceae in this data set compared with some other studies (Simanonok et al., 2021; Simanonok et al., 2023; Simanonok et al., 2024) but which is consistent with the lower visitation of Fabaceae by B. affinis relative to congeners reported by Wood et al. (2019). Wolf et al. (2022) also noted the relative rarity of observed B. affinis visits to several common Fabaceae widespread in Wisconsin, such as Lotus corniculatus, Vicia purpurea, Melilotus sp., and Trifolium pratense. The pea-flower genera Trifolium and Securigera were nonetheless repeatedly detected in Appalachian samples early in the season and in August after the more dominant forage plants were no longer flowering (Fig. 5). Other late-summer forbs identified as forage sources were predominantly within the Asteraceae family, e.g., Ageratina, Carduus, Helianthus, Verbesina, and Solidago (Fig. 5).

Shrubs and trees commonly detected in Appalachian pollen forage included Rubus, Rhododendron, and Sambucus as well as an apparent abundance of Tilia, with the caveat that the latter was detected only with the ITS2 primers as noted above. Tilia is known to be attractive to Bombus (Anderson, 1976; Jacquemart et al., 2018) and flowers abundantly in June and July but not in August (Anderson, 1976, https://www.inaturalist.org/taxa/54854-Tilia-americana), consistent with our data (Fig. 5). Interestingly, a potential toxicity of Tilia nectar to bees has been reported, particularly in Europe, although effects on bee mortality have been questioned (Koch & Stevenson, 2017; Jacquemart et al., 2018; Lande et al., 2019). While the majority of Tilia sequence reads were assigned at the genus level (File S2), those assigned at the species level were almost exclusively T. americana (basswood), which has not been reported to have a negative effect on Bombus pollinators to our knowledge. Sorbus, Oxydendrum, and Kalmia were only infrequently detected (File S2), however the overlap between survey efforts and flowering was likely low for the latter two taxa based on phenology records in the Maryland Biodiversity Project online database (https://www.marylandbiodiversity.com/). Another common genus of shrubs in the region, Vaccinium, potentially overlapped in flowering with early June survey dates but went undetected in our data. Vaccinium is known to be important forage for many bumblebee species (Cane et al., 1985; Jacquemart, 1993; Lichtenberg & Graves, 2023), although Jacquemart (1993) found that most bumblebee visits in similar European habitat (the Ardennes) were by queens.

Comparisons with other B. affinis forage surveys

While the present study does not systematically review historical B. affinis forage plants documented in the literature, the genera identified here (combined across loci) can be compared to several recent surveys of contemporary foraging patterns that employed various observational and pollen-binning methods. These additional studies include two metabarcoding data sets (Simanonok et al., 2021; Simanonok et al., 2024), the palynological component of Simanonok et al. (2024), the observation- and imagery-based compilation of Wolf et al. (2022), and the observational survey of Hepner et al. (2024). Hepner et al. (2024) and the present study represent Appalachian surveys, whereas Midwestern surveys include Wolf et al. (2022) and Simanonok et al. (2024), and the present Midwestern samples. Simanonok et al. (2021) included museum samples from across the northern tier of the historical range and overlaps spatially with the Midwestern data sets but not with the Appalachian data sets.

Overall, 180 plant genera were detected at least once in these studies, but only 77 of these were found in more than one data set (File S3). As with the between-locus comparison, many of the genera unique to individual studies are likely a consequence of rarity per se. For example, of the 28 genera unique to the current Appalachian data set, 23 (82.1%) were detected in three or fewer samples. Furthermore, the number of unique genera detected appeared to scale with the total number of genera detected overall rather than relate to methodology (Fig. S4). Other sampling and biological factors likely to contribute to differences among studies include geography, season, sampled habitats, scale, and ascertainment biases (e.g., visual methods may be biased against trees, as noted by Wolf et al., 2022). Regional differentiation is indicated by the fact that fewer than half of the genera detected in Appalachian samples (29 of 69, or 42.0%) were also listed as Midwestern forage plants by Wolf et al. (2022) from surveys and digital imagery (File S3). In contrast, 15 of 21 genera (71.4%) detected from the three Midwestern samples (with identical methodology) were also identified by Wolf et al. (2022).

While it is difficult to interpret differences in genera detected among studies for the above reasons, it is possible to identify a ‘consensus’ list of forbs identified at least once by both pollen taxonomy and visual observation, and at least once in both Midwestern and Appalachian data sets. These taxa (Table 1) can be viewed as high-confidence forage plants visited by worker B. affinis that could be targeted with seed mixes or other interventions. The frequency of use of those floral resources should then be readily quantifiable by multiple independent means in subsequent monitoring. Of the 36 listed genera, 32 (88.9%) were found in the current Appalachian data set and 33 (91.7%) were found visually by Wolf et al. (2022), indicating the two approaches were comparably effective at monitoring this core set of commonly visited forbs (File S3). The remaining data sets individually detected 30.5–63.8% (mean 50.5%) of the genera in Table 1. Note Rhododendron and Hydrangea are only minimally present in the upper Midwest based on USDA Plants data (https://plants.usda.gov/plant-profile/RHODO, https://plants.usda.gov/plant-profile/HYDRA), but both genera include common horticultural plants likely present in that region and are also widespread throughout the remaining historical range. Importantly, the systematics of some listed genera have been subject to revision, such that genetic and distributional records may not fully align with current treatments, for example genus Rhododendron (Craven, 2011) and tribe Eupatorieae (Robinson & King, 1985), which includes Eutrochium, Eupatorium, and Ageratina.