Animal dispersal of plant propagules fundamentally alters the success of dispersal events, and thus shapes plant community composition through time. While this is well-documented in seed plants, spore-bearing plants have received little attention with regard to this phenomenon. Birds are particularly attractive as a potential bryophyte dispersal vector given their highly motile nature as well as their association with bryophytes when foraging and building nests. Despite this, species-specific dispersal relationships between birds and bryophytes have never been examined. We captured birds in Gifford Pinchot National Forest in the Pacific Northwest of the United States to sample their legs and tails for bryophyte spores. We found 24 bryophyte species across 34 species of bird. We examined the level of specialization 1) within the overall interaction network to assess community-level patterns and 2) at the plant species level to determine the effect of bird behavioral type on the plant-animal interaction. Our results suggest that associations within the network are more constrained (specialized) than expected by chance. Additionally, we found that avian foraging guild impacted the variety of bryophytes found on an individual bird. Foliage gleaners and ground foragers had particularly specialized associations within the overall disperser-bryophyte network. Our findings suggest that diffuse bird-bryophyte dispersal networks are likely to be common in habitats where birds readily encounter bryophytes and that further work aimed at understanding individual bird-bryophyte species relationships may prove valuable in determining nuance within this newly described dispersal mechanism.

(a) Mist netting and sampling birds for bryophytes

To sample bird surfaces for bryophyte spores, we captured birds along a recreational trail adjacent to the Wind River Experimental Forest within Gifford Pinchot National Forest, WA, Pacific Northwest, USA. We deployed ten 12-x-3m, 30 mm mesh mist nets from dawn to midafternoon throughout an Oregon Ash (Fraxinus latifolia Bentham) forest, surrounding a banding station at 45 48’40” N, 121 56’35” W. We checked nets at least every 30 minutes, retaining birds for banding and sampling prior to release. Our site is adjacent to mixed Western Hemlock (Tsuga heterophylla Sargent) and Douglas-Fir (Pseudotsuga menziesii (Mirbel)) forest, and is bounded on one side by a small patch of meadow. The variety of surrounding habitat at our field site generated a diversity of avifauna with different habitat and foraging predilections. We sampled the legs and tails of each captured individual for bryophytes using cotton swabs (see [54] for details). Birds were additionally fitted with individually numbered USGS leg bands to account for recaptures. Avian capture data, along with age, sex, and morphometrics were submitted to the USGS Bird Banding Laboratory. In order to better contextualize the bryophytes sampled from bird surfaces relative to the abundance of bryophytes in our field site, we measured the cover of bryophytes on both tree trunks and the ground by quadrat sampling every 10 meters along multiple transects. We sampling along transects both in the Oregon Ash dominant forest that contained our mist nets, as well as adjacent Douglas Fir dominant forest.

(b) Sample processing and molecular methods

Samples were vortexed twice each for one minute in filtered tap water and vacuum filtered onto gridded 0.45 µm mixed cellulose ester membranes (EMD Millipore). We placed filters onto 60-x-15 mm petri plates containing BCD nutrient agar [55] and grew them under a 12:12 L:D light cycle at approximately 500 lux at room temperature (22-25°C). Plates that germinated over the next year were stored and later sampled for tissue. DNA was extracted using the manufacturer’s protocol (Sigma Aldrich Extract-N-Amp PCR kit) and frozen for storage. Samples were later thawed and the trnF-L region [56–59] of the chloroplast genome was amplified via PCR using the manufacturer’s protocol. PCR products were Sanger sequenced on ABI 3730xl instruments (Functional Biosciences, Inc.). Additionally, we sequenced known species collected at our field site to improve our ability to confirm species identification.   

(c) Sequence processing and tree building

We aligned forward and reverse Sanger reads to generate a consensus sequence for each sample, which we then trimmed to remove primer annealing sites. Samples were then preliminarily identified by comparison with NCBI accessioned trnF-L sequences via BLAST search [60]. Both bird swab and known field sample sequences were aligned in Geneious (v. 8.0.5) and used to build an UPGMA tree using the ape and phangorn packages in the R statistical computing platform version 3.3.3 [61–63]. Known samples from the field were used to confirm avian-derived sample identities before being removed from the tree. One Orthotrichum grouping of sequences failed to align appropriately within our tree and was removed from our analysis of phylogenetic distance. Finally, the single Marchantia sample in our dataset was removed prior to analyzing phylogenetic because it proved to be an overly influential outgroup, resulting in our final tree. Both avian species and family were substituted for bryophyte species names in separate trees for use in analyzing phylogenetic distance, and individual species associations plotted across the bryophyte phylogeny.

(d) Data analysis

Avian-bryophyte dispersal networks were assessed by examining bipartite (two level) ecological network structure. While initially developed for assessing the interaction between two trophic levels, these methods have since been applied to a variety of interaction networks, including dispersal networks [64,65]. The level of specialization (H2’) of the overall network was assessed using an information theoretic approach derived from Shannon’s Diversity Index. This global network index was then compared to a null distribution of specialization developed by iterating species interactions over 1000 runs while maintaining species abundances present in the original dataset [66]. Significance (p value) of a given network-level degree of specialization when compared to a random distribution of associations can be derived as the proportion of values of the random distribution of values that exceed or are equal to the observed H2’ [67]. We conducted a similar analysis on subsets of the interaction network based on avian behavioral group, allowing for a comparison of specialization across behavioral groups. We additionally calculated species-level specialization (d’) of the bryophytes in the overall avian-association network, as well as subset networks defined by avian behavioral guild [68]. We assessed how avian foraging guild impacts bryophyte specialization by building linear models comparing d’ for each bryophyte species in the total dataset to d’ within networks constructed from individual behavioral types.

We analyzed Faith’s Phylogenetic Distance (PD) by treating each avian species as an aggregate “site” type, calculating the minimal path of connectivity across the bryophyte tree of all species found on each individual bird species (Figure S3). For each species, a sample-size controlled null distribution was generated by permuting bryophyte species identities randomly 1000 times across the tree and calculating minimal path of connectivity [69]. Observed values of PD were compared with the mean ± 95% CI of the null distributions to determine whether avian-vectored bryophyte species were more clumped than random. Analyses and visualizations were constructed in R using the bipartite, picante, and ggplot2 packages [70–72].