DOI: 10.5061/dryad.kh1893228
title: Dual community assembly processes in dryland biocrusts communities

description: 1.    Biocrusts are critical components of drylands
 where they regulate a wide range of ecosystem functions, however, their
 response to the worldwide phenomenon of shrub encroachment and to
 livestock grazing, the most extensive land use in drylands, are not well
 studied. Grazing by livestock and increases in shrub cover could influence
 biocrust communities directly via trampling or shading, or indirectly, by
 altering biotic interactions amongst biocrust taxa. The extent of these
 changes in biocrust cover, diversity and composition are poorly known.
 2.      We used linear models and structural
 equation modelling to examine the direct effects of grazing and shrubs on
 biocrust community composition and the indirect effects mediated by
 changes in species interactions. 3.     
 Biocrust richness and cover increased with increasing shrub cover at the
 site level. This pattern occurred despite the negative response we found
 (lower cover and richness) under shrub patches vs open areas, which was
 consistent irrespective of grazing level. Functional diversity and
 evenness were similar between shrubs and open at low grazing intensity,
 but at high grazing functional diversity was greater in the open.
 Competition between biocrust species was an important driver of their
 community assembly irrespective of shrub cover, grazing intensity or patch
 type. Structural equation models showed that the effects of grazing and
 shrub cover on functional evenness, functional diversity and richness were
 controlled by biotic interactions within the shrub microsites. In the
 open, however, these effects were either direct or mediated by changes in
 cover. 4.      Biocrust cover, species richness
 and functional diversity increase with shrub cover at the site scale,
 despite the negative effects at the microsite level. We demonstrate here
 that drivers of community assembly differ markedly at small spatial
 scales. Though biocrust communities were directly driven by environmental
 filtering in the open, biotic interactions played a fundamental role in
 their assembly when growing beneath shrubs.
description: At each site, five shrub and
 five open (interspace) patches (hereafter ‘microsites’) were randomly
 selected, and within each of these microsites, a 5 cm-diameter circular
 subplot was sampled under the canopy or in the open (<i>N</i>
 = 320 in total). Each subplot was carefully removed from the ground,
 placed in a paper bag, and stored in a cool place until biocrust cover by
 species was measured. Once in the laboratory, each taxon was identified to
 species level using keys in Filson and Rogers (1979), McCarthy (1991),
 Catcheside (1980), Scott (1985) and Scott and Stone (1976) and
 more recent generic revisions. Nomenclature followed Buck and
 Vitt (2006) for mosses, McCarthy (2006) for liverworts, McCarthy (2015)
 for lichens, and where appropriate, more recent taxonomic revisions. From
 these data, species richness and cover were directly obtained. Changes in
 species composition were assessed using non-metric multidimensional
 scaling ordination (nMDS) and the Bray-Curtis distance measure on
 log-transformed abundance data with the PRIMER/PERMANOVA statistical
 package (PRIMER-E Ltd., Plymouth Marine Laboratory, UK).
 Functional attributes were assigned to
 each of our biocrust taxa based on eight functional traits: root (rhizine)
 length, taxon height, sediment capture, absorptivity, and the activity of
 four enzymes associated with carbon (b-glucosidase, b-D-cellobiosidase),
 nitrogen (N-acetyl-b-glucosaminidase) and phosphorus (phosphatase)
 cycling. These functional traits are linked to competitive ability and
 response to the environment in biocrust taxa, but also to their effects on
 important ecosystem functions such as nutrient cycling, water infiltration
 and resistance to erosion (see Bowker et al., 2010; Mallen-Cooper
 & Eldridge, 2016). From this information, we calculated two
 complementary metrics of functional diversity: RaoQ (sum of pairwise
 distances in multidimensional trait space), and ii) functional evenness (a
 measure of the extent to which different trait strategies are evenly
 distributed in the community; see Pakeman, 2014 for a
 review). Since RaoQ reflects pairwise differences in multidimensional
 trait space, it should be a good measure of functional similarity between
 two given species. Thus, assuming small fitness differences and
 hierarchical competition, species should differ as much as possible (i.e.,
 high RaoQ) to be able to coexist (limiting similarity).
 Functional metrics were calculated using the FD
 package in R (Laliberté et al.,
 2014). Functional dispersion and richness were also considered, but
 removed in further analyses due to their strong correlations with RaoQ (ρ
 = 0.97 [functional dispersion]) and species richness (ρ = 0.75 [functional
 richness]). Functional trait information was available for 23 of the 53
 species we found, and these represented, which represented 56.1 ± 3.08%
 (mean ± SE) of the total cover within the 320 subplots. All functional
 diversity metrics were abundance-weighted. 
 <i>Changes in biocrusts
 interactions</i> To
 calculate the degree of competition intransitivity, the inverse
 Markov-chain approach of Ulrich et al. (2014) was used as implemented in
 the free software <i>Transitivity</i>.  This approach provides
 a metric for intransitivity (<i>I</i>), which ranges between 0
 (total competition hierarchy) and 1 (total competition intransitivity).
 The metric is based upon the number of competition reversals, i.e., how
 many times the competitor dominant(s) lose (see Laird & Schamp,
 2006; Ulrich et al., 2014 for details). With this approach the number of
 species to analyse cannot be exceed one less than the number of sites,
 therefore in our case, four species. Thus, we analysed the degree of
 intransitivity on the four dominant species (as we had five replicates of
 each microsite × site combination). These were three lichen
 (<i>Collema coccophorum, Psora crystallifera and Psora
 decipiens</i>) and one moss (<i>Didymodon torquatus</i>)
 species. Together with the levels of intransitivity, the software
 calculates a measure (<i>Match</i>) of how well the simulated
 Markov chains fit the observed abundances. This <i>Match
 </i>metric is an indicator of the importance of competition for the
 assembly of the studied species (values close to one indicate that
 competition is strongly driving the observed abundances). The
 <i>Match </i>metric also indicates how reliably the level of
 intransitivity has been estimated. Only in three out of the 64 microsite ×
 site combinations was this metric below 0.6 (60% of the observed
 abundances explained by the simulated Markov chains). Therefore, it can be
 assumed that competition was a dominant interaction among the studied
 species and that the level of intransitivity was reliably estimated (see
 Ulrich et al., 2014; Soliveres et al., 2015 for further details on the
 methods). To account for
 changes in composition, we applied this approach to two different
 datasets: one including the four dominant species across all sites
 (therefore removing changes in composition), and one including the four
 dominant species at each site (therefore accounting for effects of changes
 in composition + environment). The intransitivity metrics obtained using
 the overall dominant species or those most dominant at each site rendered
 similar results (ρ = 0.41; <i>P</i> = 0.017;
 <i>N</i> = 33), although in the former case, this metric could
 only be calculated in half of the 64 microsite × site combinations because
 these species were not present in all subplots. Thus, we only report
 results from the metrics calculated for the four dominant species at each
 site. Overall, we analysed
 changes in biocrust communities by examining total cover, species richness
 and composition, functional composition, the importance of competition for
 community assembly, and the degree of competition intransitivity. For the
 taxonomic and cover measures, we used our full database. For the
 functional diversity metrics, we used a reduced database of the 23 taxa
 for which we had functional trait data. Finally, to assess the
 interactions among different biocrust species, we used the four dominant
 species. <i>Statistical
 analyses</i> Changes in
 biocrust cover, species richness, functional diversity (functional
 evenness and diversity [RaoQ]) and the level of intransitivity among the
 four dominant species in each subplot) were analysed by using linear
 models with site level shrub cover (0-50%), microsite (shrub
 <i>vs</i> open) and grazing level (low <i>vs</i>
 high), and their interactions as fixed factors using the ‘lme4’ package
 (Bates et al., 2015) within R statistical software (Version 3.4.4, R core
 Team 2018). To further
 interpret compositional changes, the degree of association of biocrust
 species in relation to patch type and grazing was measured with Indicator
 Species Analysis in R (De Caceres, 2013) using a data matrix of 56 species
 and 64 combinations of microsite by site. Indicator values combine
 information on relative abundance and frequency of species, and the
 indicator value is maximal (IV=100%) when all individuals of a given
 species are restricted to a particular microsite (e.g. shrub), and all
 samples from that particular microsite contain an occurrence of that
 species. Data were randomized among the factors and a Monte Carlo
 randomization procedure performed with 1000 iterations in order to
 determine the statistical significance of the indicator values.
 Finally, we used Structural Equation
 Modelling (SEM) to examine the direct and indirect effects (i.e. mediated
 by biocrust cover and the way species interact) of increasing grazing
 intensity and shrub cover, on biocrust richness, functional diversity, and
 composition. In order to evaluate possible interactions among the
 importance of each driver and their direct and indirect effects for a
 given microsite, we fitted separate models for open and shrub patches.
 Structural equation modelling tests the plausibility of a causal model,
 based on <i>a priori</i> information (Supporting Information Appendix S3: Fig. S1), in explaining the relationships among different variables. In our model we predicted that grazing would have direct effects on individual crust attributes (e.g. functional composition, richness; Eldridge et al. 2013; Bowker et al. 2013; Concostrina-Zubiri et al., 2017) but also indirect effects <i>via</i> changes in intransitivity (Bowker et al., 2010; Soliveres et al., 2015) and crust cover (as larger patches are more likely to support more species which are likely more competitive, therefore modulating biotic interactions). Similarly, the changes in heterogeneity and productivity induced by shrub patches and increasing density at the plot scale should directly affect the composition of biocrust communities (Maestre et al., 2002), but also their cover and how they interact. Overall goodness of fit probability tests (χ<sup>2</sup> and Bollen-Stine) were performed to determine how well our a priori model structure fit the data. The goodness of fit test estimates the match between the variance/covariance matrix of the observed data and that expected under the <i>a priori</i> model structure. Thus, high probability values indicate that these models are highly plausible causal structures underlying the observed correlations. All SEM analysis were conducted using AMOS Software Version 20.