Background: Plant-pollinator community diversity has been found to decrease under conditions of drought stress, however research into the temporal dimensions of this phenomenon remains limited. In this study, we investigated the effect of seasonal drought on the temporal niche dynamics of entomophilous flowering plants in a water-limited ecosystem. We hypothesized that closely related native and exotic plants would tend to share similar life history, and that peak flowering events would therefore coincide with phylogenetic clustering in plant communities based on expected phenological responses of plant functional types to limitations in soil moisture availability.

Location: Galiano Island, British Columbia, Canada

Methods: Combining methods from pollinator research and phylogenetic community ecology, we tested the influence of environmental filtering over plant community phenology across gradients of landscape disturbance and soil moisture. Floral resource availability and community structure were quantified by counts of flowering shoots. We constructed a robust phylogeny to analyze spatial and temporal variation in phylogenetic patterns across the landscape, testing the significance of the observed patterns against a randomly generated community phylogeny. Phylogenetic metrics were then regressed against factors of disturbance and soil moisture availability.

Results: Critical seasonal fluctuations in floral resources coincided with significant phylogenetic clustering in plant communities, with decreasing plant diversity observed under conditions of increasing drought stress. Exotic plant species in the Asteraceae became increasingly pervasive across the landscape, occupying a late season temporal niche in drought-stressed environments.

Main conclusion: Results suggest that environmental filtering is the dominant assembly process structuring the temporal niche of plant communities in this water-limited ecosystem. Based on these results, and trends seen elsewhere, the overall diversity of plant-pollinator communities may be expected to decline with the increasing drought stress predicted under future climate scenarios.

Sampling

We analyzed the phenological response of plant communities to seasonal drought using a 2x2 factorial study design contrasting four conditions of disturbance and soil moisture availability: 1) dry semi-natural environments (woodlands and rock outcrops); 2) wet semi-natural environments (wetlands); 3) dry modified environments (disturbed upland areas such as clear-cuts); and, 4) wet modified environments (rural areas including gardens, orchards and fields). Sites were stratified based on available terrestrial ecosystem mapping data (Madrone, 2008), with 6 sites selected per condition (for a total of 24 sites). Limitations in existing site conditions and the logistics of site access resulted in an imbalance in the study design, with 4 sites representing the wet semi-natural condition vs 8 sites representing the dry semi-natural condition (and 6 sites representing the other two conditions). This imbalance may potentially limit insight into effects associated with the wet condition in our study. However, our analyses incorporate soil moisture availability as a continuous variable, thereby overcoming some of the limitations associated with the imbalance in categorical analyses. Sites ranged in size from 0.21 to 6.3 ha and varied in their proximity from 650 m to >23 km (Appendix A), due to limitations in site access and the narrow geography of the island. While grouped site conditions circumscribed similar habitat types based on common soil moisture regimes, disturbance conditions were diverse, including forestry, fire, landscaping, and other anthropogenic effects. Commonalities between site conditions are reflected in the similarity of vegetation communities, as shown in ordination plots (Appendix B).

To sample floral resource availability (FRA), we randomly distributed 6–8 (2x15 m) belt transects throughout each site, using balanced acceptance sampling methods to ensure a random yet spatially balanced distribution of transects (van Dam‐Bates et al., 2018). The number of transects was scaled roughly in proportion to the size of each site, with the aim of capturing variability in floral resource availability across the landscape. Floral resource availability was quantified as counts of flowering shoots, recorded for each plant species at 1 m intervals and used to estimate relative abundance, with each interval surveyed comprising a 1 x 2 m area spanning both sides of the transect line. Soil moisture was recorded at 5 m intervals as volumetric water content (%VWC) using a Field Scout Time-Domain Reflectometry probe. Sampling was conducted on a monthly basis from April through August 2018, resulting in five samples per site, with each sample period elapsing over 13 days. This study design was intended to achieve a high degree of spatial and temporal resolution in our samples while ensuring that field work was logistically feasible.

Modelling

We assembled a phylogeny of seed plants known to Galiano Island based on a previously published robust species-level phylogeny of 353,185 seed plants derived from GenBank sequence data (Smith & Brown, 2018; Appendix C). A total of 207 species of entomophilous flowering plants were sampled across sites, 173 of which were represented in the study area phylogeny. This subset of Galiano Island’s study area phylogeny comprises the community phylogeny analyzed in this study (Fig.1). Species omitted from this phylogeny due to lack of sequence data do not demonstrate any taxonomic biases that could affect the results of this study.

We investigated phylogenetic patterns at two spatial scales. First, we pooled transects to assess patterns within each site condition for each sample period. Second, we assessed finer scale temporal patterns across the landscape within each transect for each sample period. To compare phylogenetic patterns between native vs exotic plant assemblages, we compared phylogenetic patterns for each group separately across each community and sample period. Interspecific phylogenetic distance matrices were used to calculate abundance-weighted mean pair-wise distance (MPD) and mean nearest taxon distance (MNTD) metrics for community and transect scale matrices. Whereas MPD is a measure of mean phylogenetic distance between all pairs of species in a sample, MNTD measures the average phylogenetic distance from each species to its closest relative. While MPD is sensitive to deep phylogenetic patterns, MNTD is more sensitive to terminal branching (Tucker et al., 2017). Both metrics are routinely used in community ecology to test hypotheses about the mean relatedness of taxa within communities.

To determine significance in the phylogenetic structure of community data, we compared observed phylogenetic metrics (MPD and MNTD) calculated for each subset of community data against a null model. The null model was generated by randomly shuffling labels (randomizations = 999) representing taxa across the tips of phylogenies generated based on each subset of community data, thereby randomizing assemblages while holding species richness and occurrence frequency constant (Kembel et al., 2010). Standardized effect sizes (z) resulting from null model tests were analysed for each sample period and site condition, to determine the significance of phylogenetic dispersion patterns in community data. Positive z values are associated with phylogenetic over-dispersion; negative values with phylogenetic clustering. The R packages ‘pez’ (Pearse et al., 2015) and ‘picante’ (Kembel et al., 2010) were used to implement community phylogenetic diversity analyses.

Transect-scale phylogenetic metrics were regressed against environmental variables with linear mixed effects models (LMM) using the R package ‘nlme’ (Pinheiro et al., 2019). Random slope-intercept models were fit by incorporating phylogenetic metrics as a response to log-transformed soil VWC and categorical factors of disturbance (fixed effects), with transects nested within sites (random effects). Models failed to converge using mean pairwise distance (MPD z), so mean nearest taxon distance (MNTD z) was adopted for these analyses. To test differences in soil moisture availability between sites, random intercept LMMs were fit using the same nested random effects structure, incorporating log transformed soil VWC as a response to site conditions as fixed effects (Pinheiro et al., 2019).

Finally, we fit negative binomial generalized linear mixed effects models (GLMMs) incorporating counts of flowering shoots as a response to site conditions (fixed effects) and transects nested within sites (random effects) to test differences in floral resource availability (FRA) between site conditions. Similar models were fit for each site condition and sample period to test differences in the contributions of native vs exotic species to FRA, and to test differences in FRA contributed by different plant functional types within clades. These models were implemented using the R package ‘lme4’, and ‘glmmTMB’ to address issues of zero-inflation (Bates, Maechler, Bolker & Walker, 2015, Brooks et al., 2017). The best models were selected based on AIC test scores (R Core Team 2019), assuming ΔAIC of 2.0 as a threshold for model improvement (Burnham & Anderson, 2002). Model effects are reported as Incidence Rate Ratios (IRR), and as marginal R² values (variance explained by fixed effects) where appropriate. Marginal effects were calculated using ‘sjPlot’ (Lüdecke, 2018). Statistical analyses were implemented in R Version 3.6.0 (R Core Team, 2019).

Please consult READMe.txt