Similar mechanisms underlie beta diversity of bryophytes in two archipelagos with different isolation time

Community assembly of island plants is a central topic in island biogeography. Most previous studies have focused on seed plants, while our understanding of bryophytes is limited. Specifically, how dispersal limitation and habitat heterogeneity shape beta diversity of bryophyte communities on islands remains unknown. We used two datasets of bryophyte communities from an artificial (~60 yr of isolation) and natural archipelago (~8500 yr of isolation) to understand their beta diversity patterns. The null model, in which species could disperse randomly and colonize successfully between islands, but the number of species on each island and the frequency of each species occurrence across all islands remain unchanged, was used to calculate expected beta diversity. Although we found significant differences in species composition between the two archipelagos, the difference in observed beta diversity was negligible. Further, there was no significant difference between observed and expected beta diversity in each archipelago. The difference in island area, rather than isolation, was significantly correlated with the beta diversity of bryophytes in both archipelagos. Our results suggest that the dispersal limitation is not the major ecological process driving the assembly of bryophyte communities in both archipelagos, but the habitat heterogeneity related to island area determines the beta diversity of bryophytes.


INTRODUCTION
Understanding mechanisms of community assembly on islands has long been of interest to researchers, who frequently find that dispersal limitation and/or habitat heterogeneity shape community assembly in fragmented systems (e.g., Ernest et al. 2008, Rominger et al. 2016, Liu et al. 2018, 2020. However, most such studies in fragmented landscapes focus on the assembly mechanisms of woody or herbaceous plants (e.g., Burns et al. 2010, Keppel et al. 2010, Carvajal-Endara et al. 2017, and we know little about other organisms, such as bryophytes (but see Tiselius et al. 2019), which are broadly distributed (Heden€ as 2007, Frahm 2008 and play important functional roles in terrestrial ecosystems (Lindo andGonzalez 2010, Deane-Coe andStanton 2017). As bryophytes require specific microhabitats, establish on spatially discrete habitats rapidly (Zartman 2003), and have strong dispersal ability due to their small unicellular spores (Hutsemekers et al. 2008, Vanderpoorten et al. 2019, they can give us the opportunity to examine the role of habitat heterogeneity and dispersal limitation in community assembly processes.
Beta diversity (i.e., community dissimilarity) describes the site-to-site variability in species composition (Condit et al. 2002). It can be driven by community assembly mechanisms involving habitat heterogeneity (e.g., when different species with different niches favors varied habitat conditions across sites), dispersal limitation, and the interaction between the two processes (Condit et al. 2002, Myers et al. 2015, Liu et al. 2016, Mori et al. 2018. By analyzing beta diversity across fragments, we can better understand community assembly mechanisms in fragmented habitats (Si et al. 2016, Liu et al. 2018, Mori et al. 2018. To disentangle the mechanisms of community assembly, null models can be used to simulate expected beta diversity where the interested ecological process does not exist (Chase et al. 2011, Segre et al. 2014. For example, when studying dispersal limitation, the widely used null model proposed by Chase et al. (2011) generates specieslevel stochasticity by randomly drawing species based on their occurrence frequency in the relevant set of sites, and species could disperse randomly and colonize successfully between islands in the null model. As such, if the observed beta diversity and the expected beta diversity show no significant difference, it can be interpreted as that dispersal limitation does not limit the colonization of species. Besides, when sampling effort across sites is unequal, null models can be used to test that the observed pattern is not a spurious artifact due to differences in sampling effort (e.g., Myers et al. 2013Myers et al. , 2015. The main driver shaping beta diversity of bryophytes is still controversial. For example, Mota de Oliveira and ter Steege (2015) indicated that species composition of bryophytes was mainly regulated by environmental conditions (e.g., vertical microenvironmental gradient related to tree height) rather than by geographical structure (see also Mota de Oliveira et al. 2009). Cacciatori et al. (2019) suggested the relative importance of environmental and geographical distance on beta diversity in epilithic bryophyte communities varies with spatial scale. Patiño et al. (2018) found that the age of host trees affected bryophyte beta diversity. Silva et al. (2014) concluded that the distribution of bryophytes on rocky outcrops is driven by stochastic rather than deterministic processes. Few studies, however, have assessed community assembly mechanisms of bryophyte communities across fragmented landscapes, which is a crucial part to understand community assembly under global land use change.
Because of varied fragmentation histories (e.g., isolation time) and species traits (e.g., dispersalrelated traits), debate remains as to which ecological process is the main driver of community assembly in fragmented landscapes (e.g., Negoita et al. 2016, Rominger et al. 2016. Even for the same taxon, community assembly mechanisms might be different in different fragmented systems (Horsak et al. 2012). For fragments with different isolation time, Badano et al. (2005) found that on newly created islands, ant species showed a random pattern of co-occurrence, but ant species on older islands co-occurred less often than expected because of stronger competition than on newly formed patches. Furthermore, taxa with weak dispersal ability would need more time to colonize a habitat and be poorly represented in more newly isolated islands (MacArthur and Wilson 1967, Kadmon and Pulliam 1995, Goosem et al. 2016, resulting in nonrandom distribution patterns of species across fragments. As such, isolation time of fragments and dispersal ability of species will affect community beta diversity across fragments. Such interaction between species dispersal ability and fragments isolation time is clear for species with weak dispersal ability, but we do not know for species with strong dispersal ability, such as bryophytes, and the co-occurrence patterns will also be affected by the isolation time of fragments. A large number of artificial archipelagos have recently formed worldwide, with a short isolation history compared to some oceanic islands (Jones et al. 2016). In this study, we compared bryophyte beta diversity and inferred community assembly mechanisms between islands of Thousand Island Lake (TIL, a human-made reservoir, isolated~60 yr ago) and Zhoushan (ZS, a natural oceanic archipelago, isolated 8,500 yr ago) in eastern China (Fig. 1). These two archipelagos are geographically close to each other but vary greatly in their isolation time, which makes them ideal systems for a comparative study on species composition and community assembly (e.g., Yuan et al. 2015). Here, we aimed to answer which factor drivers beta diversity of bryophyte assemblages on islands? Bryophytes have strong dispersal ability (Frahm 2008), so we tested the hypothesis that isolation among islands in these two archipelagos might be not far enough to affect bryophytes colonization, while the habitat heterogeneity on islands could affect habitats requirements of bryophytes, and maybe the main driver of the beta diversity of bryophyte assemblages. Specifically, we test the following predictions in this study: (1) There is no significant difference between the expected beta diversity in which species could disperse randomly and colonize successfully between islands and the observed beta diversity within each archipelago; (2) the difference in island area, which often correlates with habitat heterogeneity, rather than the difference in isolation, affects the observed beta diversity of bryophytes in both archipelagos.

Study area
The Zhoushan (ZS) archipelago (29°31 0 -31°04 0 N, 121°30 0 -123°25 0 E) is located in the East China Sea, northeastern Zhejiang Province, China, and the Thousand Island Lake (TIL, 29°22 0 -29°50 0 N, 118°34 0 -119°15 0 E) is located~310 km near the ZS archipelago ( Fig. 1). ZS archipelago was formed naturally 7000-9000 yr ago by the marine transgression during the Holocene (Xu et al. 2017), while TIL was created in 1959 by the construction of the Xin'anjiang Dam for hydroelectric power production (Liu et al. 2020). The area of TIL is about 580 km 2 , and there are 1078 islands larger than 0.25 km 2 . In the ZS archipelago, the land area of the ZS archipelago is about 1240 km 2 and there are 1339 islands larger than 0.005 km 2 (Xu et al. 2017). Natural vegetation in this region is dominated by subtropical evergreen broad-leaf forest (Wu 1980). The climate is a typical subtropical monsoon climate with average annual temperature of 15.6-17.0°C and annual precipitation of 930-1430 mm, with a rainy season between April and September in this region (Xu et al. 2017, Liu et al. 2018).

Data collection
We compiled the presence-absence data of bryophytes in Thousand Island Lake and Zhoushan archipelagos from the peer-reviewed published papers of Zang (2018) and Wang et al. (2011). The bryophytes were surveyed using the same sampling method on 10 islands in TIL and 12 islands in ZS. The principle of choosing sample islands is to consider island size distribution and isolation gradients across islands in an archipelago. An overview of characteristics of selected islands are shown in Table 1. To ensure the sufficiency of sampling, specimens of stony, arboreal, and native bryophytes were collected on islands. During the surveys, the incidence of bryophyte species was determined through multiple visits to entire islands until no additional bryophyte species were found to ensure a complete species list of bryophytes for each island . The species accumulation curves in both archipelagos showed that the surveys on the largest island were complete (Wang et al. 2011, Zang 2018. All specimens (1102 and 5570 specimens in TIL and ZS separately) were dried naturally and dissected into slides and identified under a microscope (Olympus BH-2 and Nikon Ds-Fi) in the laboratory. Although the two datasets are from different sources, the comparison of beta diversity between the two archipelagos will not be affected due to the complete collection of species composition on each island. Detailed species composition information and island attributes of each island could be found in Data S1.

Statistical analysis
All statistical analyses were performed in the software R 3.5.1 (R Core Team 2018). We used indicator species analysis (function multipatt, package indicspecies) to detect indicator species in TIL and ZS. To examine differences in community composition between the TIL and ZS archipelagos, we used permutational analysis of variance (PERMANOVA test, with Bray-Curtis distances) using 999 permutations in the function adonis of the R package vegan (Anderson 2001(Anderson , 2017. PERMANOVA tests whether there is a significant difference in species composition between two or more groups (Anderson and Walsh 2013). We then used non-metric multidimensional scaling (NMDS) with Bray-Curtis v www.esajournals.org distances for visualization (function metaMDS, package vegan).
We measured observed beta diversity as the dissimilarity between each pair of islands within TIL and ZS using an incidence-based (Bray-Curtis distance with binary = TRUE) metric. As a result, we got 45 and 66 beta diversity values in TIL and ZS separately. To test which mechanism drives beta diversity within each archipelago, we used a null model method to evaluate differences in species composition between observed and expected communities (Chase et al. 2011, Myers et al. 2015. We used the null model that retained differences among species in the number of sites occupied (row sums) and differences among sites in the number of species they contained (column sums) to calculate the expected beta diversity (Ulrich and Gotelli 2012). Thus, in the null model, species could disperse randomly and colonize successfully between islands, but the  number of species on each island and the frequency of each species occurrence across all islands remain unchanged. This fixed-fixed algorithm performs well when confronted with heterogeneous but random matrices and it is a good choice for analyzing assembly rules (Gotelli 2000, Ulrich and Gotelli 2007. We then calculated standardized effect size (beta deviation) as the difference between the observed and mean expected beta diversity from 2000 iterations of the null model, divided by the standard deviation of expected values (Myers et al. 2013). Zero beta deviation indicates that observed beta diversity does not differ from the expected beta diversity, and a positive or negative beta deviation indicates higher or lower beta diversity than expected. We also tested for a difference between observed and expected beta diversity between ZS and TIL using nonparametric analysis of homogeneity of multivariate dispersions (a straightforward multivariate extension of Levene's test) based on the distance to centroid values (PERMDISP test, Anderson 2006, Anderson et al. 2011. We used the nonparametric Mann-Whitney test to examine the difference in beta deviation between TIL and ZS. We also plotted the relationship between the observed beta diversity and expected beta diversity to examine the difference between the expected and observed beta diversity within each archipelago, and the diagonal dashed line (1:1 relationship) in the plot represents identical observed and expected beta diversity. Because there was a significant correlation between island size and elevation, thus, to test the main factors that affected the observed beta diversity of bryophytes, we only selected island area (ha) and isolation (i.e., distance to the nearest mainland, m) as the explanatory variables to avoid the collinearity, and used multiple regression on dissimilarity matrices (MRM) to examine the relationships between beta diversity matrices and the euclidean distance matrices of island area (i.e., the difference in island area) and isolation (the difference in isolation; Lichstein 2007). We evaluated the significance of the regression coefficients using 2000 permutations of the response variable (function MRM, package ecodist; Lichstein 2007). A combination of mantel correlation and multiple regression, MRM allows a regression-type analysis of two or more dissimilarity matrices, using permutations to determine the significance of the coefficients of determination and has been used as a method to disentangle the influence of space and environmental factors in ecological data (Lichstein 2007, Legendre andLegendre 2012). We used linear regression models to test the relationships between observed beta diversity and difference in island area as well as isolation in each archipelago.

Species composition of the two archipelagos
According to the indicator species analysis, the species Jungermannia truncate, Pylaisiadelpha yokohamae, Bazzania tridens, Weissia crispa, and Brotherella erythrocaulis were indicator species in TIL, and the species Hypnum cupressiforme, Tortula muralis, Weissia longifolia, etc., were indicator species in ZS. The similarities analysis and NMDS analysis further revealed clear distinctions in bryophyte species composition between TIL and ZS archipelagos (Table 2, Fig. 2; PERMANOVA F = 9.648, P < 0.001).

Comparison of beta diversity between archipelagos
There was no significant difference in observed beta diversity between TIL and ZS (PERMDISP test F = 0.534, P = 0.478; Table 2, Fig. 3a). A similar result was found for the expected beta diversity (F = 0.911, P = 0.3513; Fig. 3b). Beta deviation also showed no significant difference between TIL and ZS (Mann-Whitney test W = 1439; P = 0.785; Fig. 3c). However, the relationship between observed and expected beta diversity showed a 1:1 regression line in TIL and ZS (Fig. 4), meaning that the mean values of the beta deviations were equal to zero, reflecting random co-occurrence of species across islands.

Relationship between island attributes and beta diversity
The difference in island area was significantly correlated with observed beta diversity in TIL and ZS (Table 3). Further, the pairwise difference in observed beta diversity between islands increased with an increase in the difference in island area in TIL (mantel test: r = 0.398, P = 0.032; Fig. 5a) and ZS (r = 0.363, P = 0.057; Fig. 5a). However, the difference in isolation was v www.esajournals.org not significantly correlated with observed beta diversity (TIL r = À0.057, P = 0.601; ZS r = 0.033, P = 0.391; Fig. 5b).

DISCUSSION
Generally, dispersal limitation and habitat heterogeneity are regarded as the main mechanisms driving community assembly on islands and fragmented habitats (Myers et al. 2015, Mori et al. 2018. However, bryophytes, with high dispersal ability, have been rarely studied in the context of habitat fragmentation. Here, we compared the beta diversity of bryophyte communities between two archipelagos with different isolation time. Consistent with our predictions, our results suggested that species cooccurrence patterns of bryophytes were driven by similar ecological processes in both archipelagos.

Difference in species composition between archipelagos
We found significantly different species composition between TIL and ZS. For example, in TIL, Jungermannia truncate, Pylaisiadelpha yokohamae, Bazzania tridens, and Weissia crispa were indicator species, while Hypnum cupressiforme, Tortula muralis, and Weissia longifolia etc. were indicator species in ZS. The two archipelagos have different climates: subtropical monsoon climate in TIL and subtropical ocean monsoon climate in ZS. Such climate differences may be the main factor contributing to the different species composition of the two archipelagos. In addition, because of the longer history of succession in the forests of the ZS archipelago, the main vegetation is dominated by Pinus thunbergiana, Pinus massoniana, and some remnants of secondary evergreen broad-leaved forests on large islands (Yuan et al. 2015, whereas the vegetation is still at the early succession stage, dominated by Castanopsis sclerophylla in southeast TIL (Wang et al. 2019). Different forest types with different dominant angiosperm taxa are likely to be another factor driving bryophyte species composition in the TIL and ZS archipelagos. Similarly, Jagodzi nski et al. (2018) found that different tree species provided different habitats for bryophytes.

The effect of isolation on beta diversity
To disentangle the community assembly mechanisms of bryophyte communities in TIL and ZS, we used a null model that species could disperse Notes: R 2 is the proportion of variation in community composition explained by location. Statistical significance levels were calculated with 999 permutations. Fig. 2. Non-metric multidimensional scaling (NMDS) ordination of bryophyte community composition of all islands from two archipelagos together (stress = 0.96). Ellipses show 75% confidence limits within each archipelago. Two dimensions and Bray-Curtis distance were applied in the analysis. The circle size refers to the relative size of islands.
v www.esajournals.org randomly among islands within an archipelago (Chase et al. 2011). We found the observed beta diversity showed no significant difference between TIL and ZS (Fig. 3a) and that the observed beta diversity was not significantly different from the expected beta diversity within each archipelago (Figs. 3c, 4), indicating that species co-occurrence patterns across islands were random (i.e., no dispersal limitation) in their distributions as expected by the null model in TIL and ZS. Considering the high dispersal ability of bryophytes (Frahm 2008), our results suggest that the bryophytes in our study systems are free to colonize across islands. Fig. 3. Beta diversity (Bray-Curtis dissimilarity) in the Thousand Island Lake (TIL) and the Zhoushan (ZS) archipelago. (a) Observed beta diversity. (b) Expected beta diversity from a null model that retains differences among species in the number of sites they occupy and retains differences among sites in the number of species they harbor within each landscape. (c) Standardized effect size of beta diversity (beta deviation). Boxes represent the median and 25th/75th percentile for all pairwise comparisons between plots, and whiskers extend to 1.5 times the interquartile range. Fig. 4. Beta deviation in TIL (a) and ZS (b), measured as the relationship between the observed beta diversity and expected beta diversity from the null model. The diagonal dashed line (1:1 relationship) represents identical observed and expected beta diversity. Table 3. Results (regression slopes and R 2 ) of multiple regression on distance matrices (MRM) using a permutation method with the matrix of observed beta diversity in bryophyte species composition among islands as the dependent matrix and matrices of differences in island area (log10-transformed) and differences in isolation (log10-transformed) among islands as independent matrices. In addition, we found that variation in isolation among islands had no significant effect on the beta diversity of bryophytes (Table 3). Several studies have suggested that many bryophytes are highly mobile and easily overcome dispersal problems posed by habitat fragmentation (Moen andJonsson 2003, Pharo et al. 2004). The dispersal ability of bryophytes is much greater than vascular plants in general because the diaspores of bryophytes are usually small unicellular spores between 7 and 100 lm that can be easily dispersed by water, wind, or animals (Frahm 2008). The rate of distance decay in similarity is 1.5-1.9 times higher for vascular plants than for bryophytes; that is, vascular plants have low dispersal ability than bryophytes (Nekola and White 1999). Highlighting this ability, Hutsemekers et al. (2008) found that for 56% of bryophyte species, recruits dispersed at least 6-86 km from source populations in <50 yr. In previous studies, bryophyte communities have been found to restore successfully within~30 yr after catastrophic disturbance (e.g., wildfire; Smith et al. 2014, Jagodzi nski et al. 2018). In contrast, Pharo et al. (2004) found evidence for lower richness and shifts in species composition in newly formed fragments (<25 yr) compared to near contiguous habitats, suggesting that some species with weaker dispersal ability may exhibit greater time lag between colonization events. Moreover, Qian (2009) showed that the beta diversity of vascular plants was negatively related to their dispersal ability. However, in our sites, the longest distance from island to mainland was 30.38 km in ZS archipelago and all islands have been isolated for more than 60 yr in both archipelagos. Given that bryophytes appear to be able to successfully colonize appropriate microhabitats of every island in a very short period of time (e.g., Hutsemekers et al. 2008), fragmentation time (~60 yr) in our study sites may be long enough for bryophytes to disperse to all studied islands. Therefore, isolation is unlikely to be the main factor structuring the beta diversity of bryophyte communities within either archipelago.

The effect of island area on beta diversity
The difference in island area, rather than isolation, affected the beta diversity of bryophytes (Table 3), with beta diversity increasing with the difference in island area within each archipelago (Fig. 5). Our findings match previous findings for spiders (Wu et al. 2017), birds, and lizards (Si et al. 2015) in TIL, and some flowering plants on oceanic islands (Carvajal-Endara et al. 2017), while other studies found isolation mainly affected beta diversity (e.g., Wang et al. 2017, Liu et al. 2019. The emergence of these inconsistent results is likely related to the species dispersal ability. Isolation is not the main factor limiting bryophytes' distribution, but it has been noticed that the presence of bryophyte species is likely to depend on substrate availability, habitat quality, For example, varied habitat quality among host types would increase dissimilarities among bryophyte communities (Kraichak 2014), and numerous studies have found that community composition (of epiphytic bryophytes) was mainly regulated by habitat heterogeneity Zartman 2007, Kraichak 2014, Mota de Oliveira and ter Steege 2015). Liu et al. (2018) and Hu et al. (2019) also found that the habitat heterogeneity associated with island area affected beta diversity of woody plants. Similarly, we observed a significant island area effects on beta diversity of bryophyte communities in TIL, as well as in ZS. Therefore, the habitat heterogeneity related to island area mainly determines the beta diversity of bryophytes in both archipelagos.

Uncertainties and limitations
To compare the two archipelagos, we used species presence-absence data on 10 and 12 islands in TIL and ZS separately to calculate beta diversity, without considering species abundance. The number of islands sampled may be relatively small, and the presence-absence data may increase the influence of rare species when analyzing community assembly mechanisms. Thus, it is still necessary to expand the number of sampled islands and to assess whether species abundance will affect the beta diversity of bryophytes in archipelagos in future researches. Besides, we assumed that island area is a reliable indicator of habitat heterogeneity, which has been widely used to test ecological process on islands (Liu et al. 2018). However, we are unable to test direct environmental drivers of bryophytes beta diversity due to a lack of environmental factors across studied islands. To better understand community assembly mechanisms of bryophyte communities, it is crucial to include both environmental variables and species abundance data in future studies.
We found similar beta diversity patterns in these two archipelagos, indicating that the time since isolation is longer enough than the time needed for the successful colonization of bryophytes, that is,~30 yr (Smith et al. 2014, Jagodzi nski et al. 2018). Thus, we cannot conclude that isolation time has no effect on species composition of bryophyte communities, unless assessing more natural archipelagos with a larger range of isolation time (especially those newly formed archipelagos), or conducting experiments to manually control the effect of isolation time.

CONCLUSIONS
Overall, although the species composition of bryophyte assemblages showed significant differences between the naturally formed archipelago and the newly formed human-made archipelago, patterns of beta diversity within each archipelago were consistent across both, and island area rather than isolation was the main driver of beta diversity. Further, dispersal limitation did not appear to be the main driver of community assembly in both archipelagos, supported by a lack of difference between observed and expected beta diversity. Our results suggest that for taxa with strong dispersal abilities, such as bryophytes, habitat heterogeneity plays an important role in community assembly within fragmented landscapes.