Community composition changes
To analyse the differences in sample composition at MOTU level among the
treatments, we performed Non-metric Multidimensional Scaling (NMDS)
ordinations of the trees (N = 72) based on the arthropod datasets of
each marker. We used Bray–Curtis clustering to measure the
dissimilarity between the trees (function metaMDS , package
‘vegan’). Posteriorly, we ran a Permutation Analysis of Variance
(PERMANOVA) (functions adonis at package ‘vegan’) to analyse the
compositional differences between samples from different treatments. The
effects of Treatment and Vegetation structure on sample composition were
explored by correlating them with the NMDS ordination (functionenvfit (McCune & Grace, 2002), package ‘vegan’), evaluating the
strength of those relationships through the squared correlation
coefficient (r2). Beta-diversity summaries can be
found in Table S11.
We calculated the zeta diversity of each group of samples, that is, the
number of species shared by a given number of trees (Hui & McGeoch,
2014), and also how the values of zeta diversity changed with increasing
number of trees (McGeoch et al., 2019). We used the ratio of the zeta
diversity decline to visualize the species retention rate (Figure S2):
the probability of finding species shared by all trees as the number of
trees increases. We first grouped all trees of the same treatment
together (Figure S2A), and later we separated them by Treatment and Farm
(Figure S2B). In both cases we followed the ALL model (the addition of
trees does not follow any directionality, nor is dependent on proximity
to one another), performing a Monte-Carlo chain sampling with 5,000
samples (function Zeta.decline.mc , package ‘zetadiv’ (Latombe et
al., 2020)). The better fit to a power-law or exponential regression of
the zeta diversity decline was decided on the basis of the Akaike
information criterion.
We also visualized the most abundant MOTUs (those whose total abundance
across all samples was above the third quartile of the abundance
distribution) using a semiquantitative scale in which MOTUs with an
abundance an order of magnitude greater appeared as double the size in
the graph (Turon et al., 2022) (see Supplemental Information: Figures).
To identify the MOTUs most responsible for the differences between
treatments we calculated the Dufrene-Legendre indicator value d(Dufrene & Legendre, 1997) of the MOTUs for each treatment (functionindval , package ‘labdsv’ (Roberts, 2016)). Summary of the
indicator species analysis can be found in Table S12. We visualized a
subset of the datasets containing only these MOTUs using the
semiquantitative scale mentioned above.