Statistical analysis
To examine if tree seedling density (number of seedlings per
m2) was higher overall on nurse logs than on the
forest floor, we used a Wilcoxon rank test because the data did not meet
parametric model assumptions of normality and homogeneity of variance
across groups. We used seedling density per area of the plot because
area of the plot varied with nurse log diameter. We then examined
whether substrate type including forest floor and all three decay
classes of nurse logs influenced tree seedling density, bryophyte depth,
percent canopy cover, percent cover of the three most abundant bryophyte
species (Rhizomnium glabrescens , Hylocomium splendens ,Antitrichia curtipendula ) using 1-way ANOVA or Kruskal-Wallis
tests if parametric model assumptions couldn’t be met. Bryophyte depth
was missing from one plot and percent cover of H. splendens had
to be square-root transformed to meet parametric model assumptions.
To determine the factors associated with seedling density, we used a
generalized linear model (GLM) with a Poisson error distribution.
Explanatory variables included substrate type (nurse log or forest
floor), decay class of nurse log, bryophyte depth, % canopy cover, and
bryophyte species with a total percent cover of >300 across
all plots and that were found in at least 10 plots (Table S1). The
bryophyte species included as predictors were: Hylocomium
splendens , Rhizomnium glabrescens , Antitrichia
curtipendula , Rhytidiadelphus loreus , Sphagnum
girgensohnii , and Kindbergia praelonga . We used the dredge
function in the MuMIn package to determine the model with the lowest
Akaike Information Criterion (AIC) value (Barton, 2009). We calculated
Spearman rank correlation coefficients for all pairs of explanatory
variables using the rcorr function in the Hmisc package in R (Harrell &
with contributions from Charles Dupont and many others, 2020). We
included forest floor as decay class 5 in our model. We also calculated
variance inflation factors (vif) to test for collinearity among our
predictor variables using the vif function in the car package; any vif
< 5 was deemed not correlated.
To examine the influence of substrate (forest floor and nurse logs
characterized by decay class) on bryophyte community composition, we
conducted a non-metric multi-dimensional scaling ordination (NMDS) using
the metaMDS function in the vegan package in R (Oksanen et al., 2010).
We used a Bray-Curtis distance matrix on bryophyte percent cover in each
plot. Only plots that had at least two bryophyte species were included
(n = 98). We used ellipses to denote variances in bryophyte species
composition among substrate types. We fit bryophyte depth (BD), tree
seedling density (seedling), and percent canopy cover (CC) as
environmental vectors to the NMDS ordination using the envfit function
in the vegan package in R with 1000 permutations (Oksanen et al., 2010).
To test for significant differences in species composition among
substrate types, we used the adonis function in the vegan package, which
is essentially a multivariate analysis of variance. We then used the
pairwiseAdonis function in the vegan package to compute all pairwise
comparisons among substrate types with adjusted p-values (Martinez
Arbizu, 2020).
To test the light levels under and beside Hylocomium splendensmats, we used a paired t-test.
RESULTS