Multiple components of plant diversity loss determine herbivore phylogenetic diversity in a subtropical forest experiment

Plant diversity loss can alter higher trophic‐level communities via non‐random species interactions, which in turn may cascade to affect key ecosystem functions. These non‐random linkages might be best captured by patterns of phylogenetic diversity, which take into account co‐evolutionary dependencies. However, lack of adequate phylogenetic data of higher trophic levels hampers our mechanistic understanding of biodiversity relationships in species‐rich ecosystems. We used DNA barcoding to generate data on the phylogenetic diversity of lepidopteran caterpillars in a large‐scale forest biodiversity experiment in subtropical China. We analysed how different metrics of lepidopteran phylogenetic diversity (Faith's PD, MPD, MNTD) and taxonomic diversity were influenced by multiple components of tree diversity (taxonomic, functional, phylogenetic). Our data from six sampling periods represent 7,204 mitochondrial cytochrome c oxidase subunit I (COI) sequences of lepidopteran larvae, clustered into 461 molecular operational taxonomic units. Lepidopteran abundance, the effective number of species (irrespective of the focus on rare or common species) and Faith's PD and MPD (reflecting basal evolutionary splits), but not MNTD (reflecting recent evolutionary splits), significantly increased with experimentally manipulated tree species richness. Lepidopteran MNTD decreased with increasing tree MNTD. Path analyses showed that tree phylogenetic and functional diversity explained part, but not all of the effects of tree species richness on lepidopteran diversity. Importantly, tree diversity effects on lepidopteran diversity were to a large extent indirect, operating via changes in lepidopteran abundance. Synthesis. Our study shows that evolutionary dependencies determine the response of herbivore communities to changes in host plant diversity. Incorporating a wider range of diversity metrics both at the level of producers and consumers can thus help to develop a more comprehensive understanding of the functional consequences of biodiversity change across trophic levels. Moreover, the dependence of trophic linkages on herbivore abundances underlines the need to address the consequences of current declines in insect abundances for ecosystem structure and functioning.


| INTRODUC TI ON
Global environmental change leads to increasing biodiversity loss, with consequences for the structure and functioning of ecosystems (Chapin et al., 2000;Kardol, Fanin, & Wardle, 2018). Previous studies have shown the important role of plant diversity and the consequences of the loss of plant diversity for key ecosystem functions, such as primary productivity and nutrient cycling Tilman, Isbell, & Cowles, 2014). These studies have also demonstrated that changes in plant diversity can cascade up the food web to affect the abundance and species richness of higher trophic levels (Giling et al., 2019;Gossner et al., 2016;Haddad, Crutsinger, Gross, Haarstad, & Tilman, 2011;Scherber et al., 2010). This is an important finding because interactions with higher trophic-level organisms, such as herbivores and predators, play important roles in modifying ecosystem functions and can feed back on plant diversity and performance (Bagchi et al., 2014;Finke & Denno, 2005;Schowalter, 2012). In the light of recent declines in the number of insects and other higher trophic-level organisms (Dirzo et al., 2014;Hallmann et al., 2017), adequate knowledge of the way in which changing abundance and diversity at higher trophic levels affect ecosystem functioning has become a globally recognized, interdisciplinary concern (Eisenhauer, Bonn, & Guerra, 2019).
However, our understanding of the extent to which biodiversity loss at different trophic levels is directly linked, and of the mechanism driving such potential linkages, is still limited (Lewinsohn & Roslin, 2008;Soliveres et al., 2016). Some recent studies on plant diversity have shown the importance of functional and phylogenetic diversity components in explaining the causal links between declines in plant species richness and changes in ecosystem functions Schweiger et al., 2018;Srivastava, Cadotte, MacDonald, Marushia, & Mirotchnick, 2012), whereas others failed to do so (e.g. Huang et al., 2018). For higher trophic levels, we are often still struggling to understand patterns of functional and phylogenetic diversity and their relationship with changes in plant diversity (Ebeling et al., 2018;Jorge, Prado, Almeida-Neto, & Lewinsohn, 2014;Scherber et al., 2010). Identifying and measuring functional traits to quantify animal functional diversity and the causal links with ecosystem functioning can be challenging (Brousseau, Gravel, & Handa, 2018), and using phylogenetic diversity as a proxy for functional diversity has repeatedly been discussed as a potential way forward if key functional traits show some degree of phylogenetic conservatism (Cavender-Bares, Kozak, Fine, & Kembel, 2009;Gravel et al., 2012).
As a biodiversity measure, phylogenetic diversity links the distinct evolutionary history and features of species. It can therefore point to the mechanisms driving patterns of distribution and co-occurrence among species (Castagneyrol, Jactel, Vacher, Brockerhoff, & Koricheva, 2014;Pellissier et al., 2013;Singer et al., 2018;Webb, Ackerly, McPeek, & Donoghue, 2002), which can be difficult to identify with studies based only on taxonomic diversity. Previous studies have shown that both plant species richness and plant phylogenetic diversity are relevant drivers of the species richness of higher trophic levels (Cavender-Bares et al., 2009;Staab et al., 2016;Weiblen, Webb, Novotny, Basset, & Miller, 2006;Whitfeld et al., 2012). In addition, recent evidence indicates that plant functional traits and phylogenetic diversity have effects on the community structure and phylogenetic diversity of higher trophic levels (e.g. Pellissier et al., 2013). However, a general understanding of these relationships across trophic levels is hindered by the very limited number of studies that have attempted to analyse the linkages between the structure and diversity of plant communities and the phylogenetic tree species richness. Lepidopteran MNTD decreased with increasing tree MNTD.
Path analyses showed that tree phylogenetic and functional diversity explained part, but not all of the effects of tree species richness on lepidopteran diversity. Importantly, tree diversity effects on lepidopteran diversity were to a large extent indirect, operating via changes in lepidopteran abundance. 4. Synthesis. Our study shows that evolutionary dependencies determine the response of herbivore communities to changes in host plant diversity. Incorporating a wider range of diversity metrics both at the level of producers and consumers can thus help to develop a more comprehensive understanding of the functional consequences of biodiversity change across trophic levels. Moreover, the dependence of trophic linkages on herbivore abundances underlines the need to address the consequences of current declines in insect abundances for ecosystem structure and functioning.

K E Y W O R D S
BEF-China, biodiversity and ecosystem functioning, Hill numbers, Lepidoptera, phylogenetic diversity, plant species richness diversity of higher trophic levels (e.g. Lamarre et al., 2016;Peralta, Frost, Didham, Varsani, & Tylianakis, 2015). And yet, better insights into the phylogenetic structure and diversity of higher trophic levels and their relationship with plant diversity may be required to mechanistically explain the consequences of biodiversity loss across trophic levels. These insights can also be used to devise adequate conservation strategies that take into account potentially non-random biodiversity loss (Peralta et al., 2015).
Here, we analyse how taxonomic and phylogenetic diversity of species-rich assemblages of key herbivores -lepidopteran caterpillars -in subtropical forests are affected by changes in multiple components of plant diversity. Herbivores are particularly important to consider in this context, because they are well known to show phylogenetically structured host use that reflects evolutionary adaptations to plant palatability and defense traits (Lamarre et al., 2016;Volf et al., 2018). Specifically, we explore the effects of tree species richness, tree functional trait composition and functional diversity, and tree phylogenetic diversity on abundance, species richness and phylogenetic diversity of lepidopteran larvae in a subtropical forest biodiversity experiment in a highly diverse region of south-east China . We considered Faith's phylogenetic diversity (PD), mean phylogenetic pairwise distance (MPD), and mean nearest taxon distance (MNTD) as metrics of phylogenetic diversity. These metrics are differentially sensitive to changes in taxonomic diversity (Srivastava et al., 2012) and they reflect different degrees of dependence on ancient (MPD) or recent (MNTD) splits in the phylogeny (Webb, 2000), which can help to understand how evolutionary dependencies determine the response of herbivore assemblages to changes in host plant diversity (Lamarre et al., 2016). Previous results from the same study sites at an earlier stage of forest development have shown that tree species richness and tree functional traits can promote caterpillar abundance and species richness (Zhang et al., 2017). However, the phylogenetic diversity of the tree and caterpillar communities were not considered, which might show deviating patterns compared with taxonomic diversity (Pellissier et al., 2013).
We expected that plant species loss will affect both taxonomic and phylogenetic diversity of lepidopteran communities largely via changes in plant functional and phylogenetic diversity. Specifically, we hypothesized that (a) lepidopteran phylogenetic diversity will be best explained by plant phylogenetic diversity, as these metrics account for evolutionary adaptations to specific plant lineages.
However, (b) the strength of these relationships will depend on the metric of phylogenetic diversity and their dependence on the timing and overall number of evolutionary splits in the lepidopteran phylogeny.

| Study area and experimental design
The study was conducted on the 'BEF-China' tree diversity experimental sites, which form the largest tree diversity experiment in the world at present. The experiment is located in Jiangxi province, south-east China (29°08′-29°11′N, 117°90′-117°93′E), a region that is characterized by typical seasonal monsoon climate. The mean annual temperature is 16.7°C and mean annual precipitation 1,821 mm (Yang et al., 2013). In total, 566 plots (25.8 × 25.8 m 2 ) were established on two sites (Site A and B; ~20 ha each and c. 4 km apart). On each plot, 400 saplings were planted in regularly arranged 20 rows and 20 columns (planting distance 1.29 m). Planting took place in 2009 (Site A) and 2010 (Site B) . For our study, we focused on a subset of 64 randomly distributed, intensively studied plots on the two sites (32 plots per site).
The selected plots span a tree diversity gradient from monocultures to 24 species-mixtures (sixteen monocultures, and eight, four, two, one, and one mixtures of 2, 4, 8, 16 and 24 species, respectively, per study site), which allows testing for the effects of plant diversity loss on ecosystem structure and functioning. Altogether, 40 locally common tree species (Table S1) were planted across the entire experiment. Tree species composition at the two sites differed, with two separate species pools of 16 broadleaved species at each site (and an overlap of eight tree species in the 24-species mixtures that were planted on both sites). The species of the less diverse mixtures were selected by randomly subdividing the species of the 16-species plots into two 8-species mixtures, four 4-species mixtures, and eight 2-species mixtures with non-overlapping species composition . Tree species were randomly assigned to the planting positions within each plot, and the total number of individuals per plot was divided equally among the planted species.

| Sampling
Lepidopteran larvae were collected six times in 2017 and 2018 (April, June and September in each year). We collected all caterpillars by beating individual trees and knocking down resident insects.
Lepidopteran caterpillars are much more restricted in their mobility than adult and flying insects and have a high probability of being collected from the trees they actually feed on (see also Wardhaugh, Stork, & Edwards, 2012). We beat the trees with a padded stick over a white sheet (1.5 × 1.5 m 2 ) and collected all caterpillars knocked down from the trees Schuldt, Baruffol, et al., 2014). We sampled all trees in the first rows of each plot for a total 80 living trees in each plot. Our collection completely covered the tree species composition and species richness at the plot level due to the completely random planting design. To avoid contamination of samples, we collected all lepidopteran larvae individually and stored them in separate tubes filled with 99.5% ethanol.
All samples were kept in a −20 ℃ freezer until further processing.
In addition to the plot-level caterpillar data, we collected and taxonomically identified adult moths to construct a reference DNA barcode library. Moths were collected by light trapping near the experimental sites during the season of sampling for caterpillars, and identified by one of the authors (Chun-Sheng Wu).

| DNA extraction, amplification and sequencing
All specimens were sequenced for a region of the mitochondrial cytochrome c oxidase subunit I (COI) gene (Hebert, Ratnasingham, & de Waard Jeremy, 2003), which has been widely used for species delimitation in molecular biology. We used sterile equipment to cut open small-to medium-sized larvae from the head to the abdomen and put them individually into Eppendorf tubes prior to DNA extraction. For larger sized caterpillars, we used 3 or 4 body segments for DNA extraction. Whole genomic DNA was extracted using DNeasy Blood & Tissue Kits (QIAGEN GmbH, Hilden, Germany), following the manufacturer's protocols. COI sequences of samples were amplified using universal primer pairs, LCO1490 and HCO2198 (Folmer, Hoeh, Lutz, & Vrijenhoek, 1994). In cases when a COI sequence was not generated successfully, one of the alternative primer pairs LepF1 and LepR1 (Hebert, Penton, Burns, Janzen, & Hallwachs, 2004) were employed to amplify those samples. Polymerase chain reactions (PCR) were carried out in 96well plates in 30 μl reaction containing 10 μl ddH 2 O, 15 μl Premix PrimeSTAR HS (TaKaRa), 1 ul of each primer at 10 μM, and 3 μl template genomic DNA using a thermo cycling profile of 2 min at 94°C; 29 cycles of 50 s at 94°C, 50 s at 50°C, 1 min at 72°C; followed by a final extension at 72°C for 6 min. All PCRs were performed on an Eppendorf Mastercycler gradient. PCR products were visualized on a 1% agarose gel, and samples with clean single bands were sequenced after PCR purification using BigDye v3.1 on an ABI 3730xl DNA Analyser (Applied Biosystems).
Overall, 7,204 COI sequences were successfully generated, whereas sequencing failed for 1,267 COI sequences of caterpillars.
However, the percentage of sequencing failures was independent of tree species richness, meaning that our analyses use comparable data across the tree diversity gradient of our study (Pearson's r = −0.10, p = .47 for a correlation of percent of failed sequences and tree species richness per plot).

| Sequence alignment and phylogenetic analysis
Molecular Operational Taxonomic Units (MOTU) were inferred from the plot level data, and then a phylogeny was constructed based on the MOTU and the reference data. First, we made a curated reference alignment of the 177 COI barcodes of adult moths collected by light trapping. The sequences were aligned using MAFFT (Misawa, Katoh, Kuma, & Miyata, 2002), and translated into amino acid sequences using MEGA v7.0 (Kumar, Stecher, & Tamura, 2016) to test the presence of stop codons, and manual adjustments were made. Caterpillar sequences were then aligned against the references using a Perl-based DNA barcode aligner (Chesters, 2019). Three methods were used for inferring MOTU: threshold-based hierarchical clustering with BLASTclust, Automatic Barcode Gap Discovery (ABGD), and Poisson Tree Processes model (PTP). The BLASTclust module of the ncbi-blast package (https ://blast.ncbi.nlm.nih.gov/Blast.cgi) was run under a clustering threshold of 97.8% identity (Ratnasingham & Hebert, 2013). MOTU were inferred using ABGD (Puillandre, Lambert, Brouillet, & Achaz, 2012) with the online tool at https ://bioin fo.mnhn.fr/abi/publi c/abgd/, under the Kimura (K80) model with parameters: Pmin = 0.001, Pmax = 0.1, Steps = 50, X = 0.5 and Nb bins = 50. The PTP model (Zhang, Kapli, Pavlidis, & Stamatakis, 2013) was run with default parameters on the maximum likelihood phylogeny of the unique plot-level barcodes. Alternative clusterings produced from the three methods were compared using the Hubert and Arabie-adjusted Rand index using r package clues, and the most consistent clustering selected for further analyses.
A phylogeny was constructed of the finalized MOTU for calculation of phylogenetic indices. For improving phylogenetic structure of a DNA marker with limited information content, we integrated fully identified references to anchor the phylogenetic analysis (e.g. Zhou et al., 2016). For setting the backbone topol-  Figure   S1). In addition to construction of a single maximum likelihood phylogeny, phylogenetic uncertainty was quantified by bootstrapping (Felsenstein, 1985). As branch lengths were required in downstream analyses, we conducted both a maximum likelihood tree-search and branch-length optimization on each of 100 bootstrapped alignments.
We used the Statistical Assignment Package (Munch, Boomsma, Huelsenbeck, Willerslev, & Nielsen, 2008) to taxonomically identify MOTU via reference data. We selected the 'ConstrainedNJ' algorithm with 24 homologues retrieved from GENBANK and our moth adult database for each query sequence (>85% sequence similarity, Ransome et al., 2017). MOTU that could not be confidently assigned by either database were labelled unidentified Lepidoptera (Table S2).
In constructing the species level phylogeny, topological constraints were applied only partially due to a number of factors. In the current instance a high degree of overlap was expected between plot-level MOTU and reference data (as a site specific reference library was constructed). However, limitations possibly still occurred due to the general structure of the taxonomic framework and the patterns in monophyly of the backbone phylogeny. The latter was an apparent issue caused by non-monophyly of families within the Noctuoidea (Heikkilä et al., 2015; Figure S1). Thus, despite the use of topological constraints, some phylogenetic uncertainty was expected to remain, particularly when using single genes. To test the impact of this on the results, we calculated the three phylogenetic diversity indices on Lepidoptera phylogenies from bootstrapped alignments.

| Plant traits, plant diversity, and environmental covariables
We used a range of morphological and chemical leaf traits of the tree species to characterize plot conditions in accordance to nutritional quality and potential defense traits of the trees. For our main analyses, we focused on five general traits that have frequently been shown to affect the composition and performance of insect herbivores (Pérez-Harguindeguy et al., 2003;Schuldt, Baruffol, et al., 2014;Zhang et al., 2017). As morphological traits, we included specific leaf area (SLA), leaf dry matter content (LDMC) and leaf toughness (LT). As for chemical leaf traits, we used leaf carbon (C) concentration and the ratio of leaf carbon to nitrogen (C:N) concentrations. SLA is an important indicator of growth rate potential (Pérez-Harguindeguy et al., 2003) and often positively related to herbivory, as leaves with a high SLA are comparatively soft and easy to attack (Pérez-Harguindeguy et al., 2013). In contrast, LDMC and leaf toughness are expected to decrease herbivory, as tough and structurally robust leaves are more difficult to attack (Pérez-Harguindeguy et al., 2003;Poorter, Plassche, Willems, & Boot, 2004). However, previous research in our study region frequently found strong positive relationships between LDMC and leaf herbivory Schuldt, Baruffol, et al., 2014;Schuldt et al., 2012), possibly because many of the dominant leaf chewers might be adapted to coping with robust leaves in areas where trees with such leaves are common (Pérez-Harguindeguy et al., 2003). Such a positive relationship might reflect compensatory feeding to make up for low nutrient contents (see Schuldt et al., 2012), or correlation with unmeasured traits that affect herbivory, such as secondary plant metabolites (Blonder et al., 2018). Carbon concentrations and the C:N ratio of leaves are related to palatability, and low C:N ratios in particular have commonly been found to promote herbivory (Pérez-Harguindeguy et al., 2003;Poorter et al., 2004).
To test whether consideration of a larger number of leaf traits changed the inferences that could be drawn regarding the effects of functional traits and functional diversity on herbivore diversity, we additionally used five further leaf traits for alternative analyses. These traits were leaf area (LA), leaf potassium (K) content, leaf magnesium (Mg) content, leaf sodium (Na) content, leaf phosphorus (P) content, i.e. traits that might add information particularly on the nutritional quality of the tree species (e.g. Borer, Seabloom, Mitchell, & Cronin, 2014;Poorter et al., 2004). All of the traits were measured on pooled samples of sun-exposed leaves of a minimum of five tree individuals per species following standard protocols (Pérez-Harguindeguy et al., 2013). Because our leaf-sampling design did not allow us to quantify intraspecific variability, we used single mean trait values for each species. Previous studies indicated that variability in trait-environment relationships in our study region is much more pronounced at the interspecific than the intraspecific level (Kröber, Böhnke, Welk, Wirth, & Bruelheide, 2012;Schuldt et al., 2012).

| Statistical analysis
All analyses were conducted in r 3.4.3 (http://www.R-proje ct.org) with the packages ade4, inext, lavaan, MuMIn, picante and vegan. Prior to analyses, samples from the six sampling times were pooled at the plot level. Nine plots were excluded in our study because of high tree mortality (with very low numbers of living trees), resulting in 55 plots that were included in the final analysis. We tested the correlations among all predictors ( Figure S4) and checked variance inflation factors (VIF) of our statistical models to ensure that the analyses were not strongly affected by multicollinearity.

| Hill numbers of lepidopteran larvae
Observed samples are often incomplete due to many rare species in an assemblage that may lead to a larger number of undetected species (Chao et al., 2014). As a consequence, species richness values are often underestimated when based on observed species richness. To account for this, we followed the approach by Chao et al. (2014) and used the first three Hill numbers (Chao et al., 2014) to estimate species richness (q = 0), the exponential of Shannon's entropy (q = 1; referring to Shannon diversity) and the inverse of Simpson's concentration (q = 2; referring to Simpson diversity) of the lepidopteran communities. Hill numbers provide statistically rigorous assessments of the effective number of species, weighted by the abundance of rare or common species (Hsieh, Ma, Chao, & McInerny, 2016). The calculation was based on sample coverage (level = 0.6), which is less affected by differences in total sampling effort than other methods (Chao & Jost, 2012).

| Community-weighted mean trait values, functional and phylogenetic diversity
We used the community-weighted mean values (CWMs) of each trait as well as the functional diversity of these traits as predictors of lepidopteran abundance and diversity. CWMs were calculated as the mean value of each trait weighted by tree wood volume.
Tree wood volume was estimated from data on basal area and tree height (Fichtner et al., 2017) measured in October 2016 on trees in the centre of each plot. Values were upscaled to represent plotlevel wood volume of each tree species for our analyses, and the CWM of wood volume was used as an additional predictor in our models.
Tree functional diversity was calculated as the mean pairwise distance in trait values among tree species, again weighted by tree wood volume, and expressed as Rao's Q (Ricotta & Moretti, 2011). We also calculated an alternative tree functional diversity that included additional leaf traits (i.e. LA, K, Mg, Na and P) to test whether results qualitatively differed when more leaf traits were considered. Similar to functional diversity, we quantified the phylogenetic diversity of the tree communities by calculating wood volume-weighted phylogenetic Mean Pairwise Distance (MPD), which in the abundance-weighted case is equivalent to Rao's Q (Tucker et al., 2017). Moreover, we calculated Mean Nearest Taxon Distance (MNTD), which is more sensitive to variations towards the tips of the phylogeny than MPD, because MNTD is a measure that quantifies the distance between each species and the nearest neighbour on the phylogenetic tree (Webb, 2000).
Phylogenetic metrics were calculated based on a maximum likelihood phylogenetic tree available for the tree species of our study region (Purschke, Michalski, Bruelheide, & Durka, 2017). We tested for phylogenetic signal in functional traits across the tree phylogeny using Blomberg's K (Blomberg, Garland, & Ives, 2003), implemented in the r package phylosignal (Table S3).
To characterize the heterogeneous topography of the study plots, plot means of elevation, slope, 'eastness' (sine-transformed radian values of aspect) and 'northness' (cosine-transformed radian values of aspect) were included in our analysis as environmental covariables. Data were obtained from a 5 m digital elevation model based on differential GPS measurements (Scholten et al., 2017).
We quantified the phylogenetic diversity of the lepidopteran communities in three ways. First of all, we calculated Faith's Phylogenetic Diversity (PD) as a simple metric, which, however, is usually strongly correlated with species richness (Tucker et al., 2017). As above for the tree communities, we therefore also calculated abundance-weighted phylogenetic MPD and MNTD. All phylogenetic metrics were calculated using the r package picante.

| Lepidopteran larvae abundance and biodiversity
We used linear models to analyse the effects of tree diversity, plant traits and environmental covariables on caterpillars. As response variables, we modelled abundance, observed species richness, Hill numbers and phylogenetic diversity of the caterpillars. As biotic predictors, we used tree species richness, the CWMs of the five major plant traits and of tree wood volume, tree functional diversity based on Rao's Q of the five major traits, and tree phylogenetic MNTD. We did not include tree phylogenetic MPD in the same models, because it was correlated with tree species richness (Pearson's r = 0.74, p < .001). However, we additionally analysed alternative models where we replaced tree species richness by tree MPD and compared the AICc values of both model variants to assess whether tree species richness or tree MPD more strongly affected lepidopteran diversity. In the same way, we tested in separate models whether replacing functional diversity based on five leaf traits by functional diversity based on ten traits yielded qualitatively different results.
As abiotic predictors, we included in all models sites, elevation, northness, eastness and slope. We also included the interactions between site and tree species richness and site and tree functional diversity as predictors. We simplified the linear models in a stepwise procedure based on values of the Akaike Information Criterion corrected for small sample sizes (AICc), and chose subset models with the lowest AICc.
To improve normality and variance homogeneity of the model residuals, lepidopteran larvae abundance and observed species richness as response variables were log-transformed, and estimated species richness, Shannon diversity and Simpson diversity as response variables were square-root transformed. Likewise, tree species richness as a predictor was log-transformed in all analyses. All continuous predictors were standardized (mean = 0, standard deviation = 1) prior to the analyses to enable direct comparisons of model estimates.
Finally, we assessed whether effects of tree species richness on lepidopteran species richness were driven by the combined effects of different tree species (as hosts for different herbivore species), or were caused by sampling effects due to the inclusion of tree species with particularly high numbers or lepidopteran species. For this, we followed the approach by Hector et al. (1999) and calculated regressions for lepidopteran species richness with plot-specific biomass of each tree species as a predictor, testing whether many or only a few tree species showed significant effects.

| Path analyses
To shed light on potential causal direct and indirect pathways that illustrate the effects of tree diversity on caterpillar abundance and diversity, path analyses (Grace, 2007) were fitted with the r package lavaan. We fitted separate models for lepidopteran species richness, PD, MPD and MNTD and used tree species richness, tree functional diversity, tree phylogenetic MPD (and MNTD, but only in the analysis of lepidopteran MNTD, because it did not affect other metrics of lepidopteran diversity), and lepidopteran abundances as predictors.
We assumed that tree species richness as the experimental treatment influences the phylogenetic and functional diversity of the tree communities, and that functional diversity is influenced by phylogenetic diversity because functional traits might show a phylogenetic signal. For these three components of plant diversity, we hypothesized that they influence caterpillar diversity both directly (direct paths from tree diversity to caterpillar diversity) and indirectly (via effects on caterpillar abundance). We sequentially removed nonsignificant pathways if their removal resulted in increased model fit (lower AIC) (Scherber et al., 2010). Indirect effects via abundance were recently shown to play an important role in modifying arthropod species richness (Schuldt et al., 2019) and might operate by influencing local population persistence or host choice behaviour of a larger number of species (Scheirs & De Bruyn, 2002;Storch, Bohdalková, & Okie, 2018). In addition, we tested two alternative model variants: first, we recalculated the models in the same way as described above, but assuming that lepidopteran diversity influenced lepidopteran abundance (i.e. turning around the pathway between lepidopteran abundance and diversity). Second, we tested for significant residual covariances between the lepidopteran larvae abundance and all of diversity metrics (instead of direct pathways connecting lepidopteran abundance and diversity), because the two might have reciprocal effects on each other (e.g. Storch et al., 2018). In both cases, we simplified models as described above and compared the resulting AICc values to those of the final models of our initial approach.

| RE SULTS
In total, we collected 8,471 lepidopteran larvae across the six sampling periods, for which 7,204 COI sequences of larvae were grouped into 461, 408 and 524 MOTU by hierarchical clustering, ABGD and PTP, respectively. We selected the hierarchical clustering for further analyses as it was the most consistent when compared to the remaining two methods (pairwise clustering congruence measured by the HA-adjusted Rand index; hierarchical clustering vs. ABGD: 0.995; hierarchical clustering vs. PTP: 0.990; ABGD vs. PTP: 0.985).
Moreover, the lepidopteran species richness was very similar at the plot level for the three methods (Pearson's r > 0.97, p < .001 for a correlation of lepidopteran species richness by all pairs of delimitation methods at plot level). In total, 5,181 sequences were identified to family, 2,743 to genus and 1,178 to species using SAP analysis (Table S2).
The abundance and all taxonomic and phylogenetic metrics of the diversity of lepidopteran larvae (with the exception of MNTD) were significantly positively related to tree species richness (Figures 1 and 2, Table 1). This positive relationship was consistently found for a wide range of the studied tree species when analysed separately in single regressions: more than one third F I G U R E 1 Relationships between tree species richness and (a) abundance of lepidopteran larvae, (b) species richness of lepidopteran larvae, and (c) esitimation of species diversity (q = 0: estimated species richness; q = 1: Shannon diversity; q = 2: Simpson diversity). Regression lines (with 95% confidence bands) show significant (p ≤ .05) relationships. Note that axis values are on a log-scale for tree species richness, lepidopteran abundance, and lepidopteran species richness and a square root-scale for estimation of lepidopteran species diversity F I G U R E 2 Relationships between tree species richness and (a) Faith's phylogenetic diversity (PD) and (b) phylogenetic mean pairwise distance (MPD) of lepidopteran larvae. (c) Relationship between phylogenetic mean nearest taxon distance (MNTD) of trees and lepidopteran larvae (adjusted for covariates of the final regression model). Regression lines (with 95% confidence bands) show significant (p ≤ .05) relationships TA B L E 1 Summary results of linear models for abundance, species richness, Hill numbers and phylogentic indices of lepidopteran larvae across a tree species richness gradient. Standardized parameter estimates (with standard errors, t and p values) are shown for the variables retained in the minimal models (13 out of 38) of all tree species showed a significant or close to significant relationship between their plot-level biomass and lepidopteran species richness, and many more also showed positive, but non-significant relationships (Table S4). Lepidopteran species richness was predicted to more than double (from 30 to more than 60 species) from monocultures to the 24-species mixtures ( Figure 1b). Importantly, the species richness of lepidopteran caterpillars strongly increased with tree species richness irrespective of whether the many rare species were allowed to influence species richness patterns (q = 0) or the focus was on common species (q = 2; Figure 1c). Tree species richness was a stronger predictor of lepidopteran species richness than tree functional diversity. The latter was only kept in the minimal model of the species richness of common lepidopterans (q = 2), where tree functional diversity showed a significant positive effect, but only at one of the study sites (Table 1) (Table S5). Only in the case of the species richness of common lepidopterans (q = 2) and MPD did functional diversity based on ten traits show a slightly stronger effect on both (common species) or one (MPD) study site (Table S5).
The phylogenetic diversity of lepidopteran caterpillars showed very similar results when measured as Faith's PD, because the latter was strongly correlated with taxonomic species richness (Pearson's r = 0.99, p < .001; Figure 2a, Table 1). However, the phylogenetic mean pairwise distance (MPD) of lepidopteran communities (which was less strongly correlated with caterpillar species richness; Pearson's r = 0.43, p = .001) increased with tree species richness as well, although less strongly (Figure 2b, Table 1). Only the phylogenetic mean nearest taxon distance (MNTD) of caterpillars was unrelated to tree species richness (Table 1). However, MNTD of caterpillars significantly decreased with increasing tree phylogenetic MNTD.
Lepidopteran PD and MNTD further showed a significantly positive and negative relationship, respectively, with increasing mean values of leaf toughness of the tree communities (Table 1). In contrast to the taxonomic species diversity metrics, the CWM of SLA had no significant effects on lepidopteran phylogenetic diversity (Table 1).
The direction and magnitude of correlation between tree species richness and Lepidoptera phylogenetic diversity indices were indistinguishable whether using a single Maximum Likelihood phy- Tree phylogenetic MPD had similar, positive effects on lepidopteran diversity as tree species richness when tree species richness as a predictor was replaced by the correlated tree MPD (Table   S6). However, AICc values for the models with tree MPD tended to be higher than AICc values of the models with tree species richness (except for the model of common (q = 2) lepitopteran species richness), indicating that tree phylogenetic MPD did not predict lepidopteran diversity better than tree species richness (Table 1,   Table S6). Variance inflation factors for all linear models were <2.2 in all cases, indicating that multicollinearity could only have little influence on the analyses. Of the ten leaf traits considered in our study, four showed a significant phylogenetic signal (SLA, LDMC, LT, LA; Table S3).
The path analyses -which simultaneously incorporated the potential effects of tree species richness, tree functional diversity, and tree phylogenetic diversity -showed that to some extent the effects of tree species richness on caterpillar taxonomic and phylogenetic diversity operated indirectly via tree functional and phylogenetic diversity. Tree species richness directly influenced tree phylogenetic MPD, which in turn affected tree functional diversity ( Figure 3). The majority of effects of tree diversity on caterpillars operated by influencing caterpillar abundances, which in turn promoted lepidopteran species richness and PD (Figure 3a,b), and decreased lepidopteran MNTD (Figure 3d). Moreover, MNTD of caterpillars was negatively affected by tree phylogenetic MNTD, which in turn was influenced by tree phylogenetic MPD.
Lepidopteran MPD was not directly related to caterpillar abundances, leaving only weaker effects via tree functional diversity on MPD (Figure 3c). Directional effects of lepidopteran abundance on lepidopteran diversity received stronger support (lower AICc values, except for lepidopteran MPD) than alternative models assuming reciprocal effects between the two variables (using a covariance term instead of a direct path; AICc values for directional effects vs. Note: CWM C, Community-weighted mean value of leaf carbon concentration; CWM LDMC, Community-weighted mean value of leaf dry matter content; CWM LT, Community-weighted mean value of leaf toughness; CWM SLA, Community-weighted mean value of specific leaf area; CWM TV, Community-weighted mean value of tree volume; FD, Functional diversity. and effects of tree diversity were weaker (lower standardized path coefficients) in models with pathways going from lepidopteran diversity to abundance ( Figure S3, Tables S6-S8).

| D ISCUSS I ON
Our study provides insights into the linkages between tree diversity and the diversity of insect herbivores in a species-rich subtropical forest ecosystem. It indicates that changes in both the taxonomic and phylogenetic diversity of the tree communities can be tracked in the associated herbivore communities particularly when more basal splits in host and herbivore lineages (as captured by MPD) are taken into account. More recent splits among herbivore species (as captured by MNTD) seemed to be not directly affected by tree diversity.
Importantly, effects of tree diversity on herbivore diversity were to a large part indirect, because tree diversity had strong effects on herbivore abundances, which in turn affected herbivore diversity.
This finding has important implications considering the currently observed decline in insect abundances in many ecosystems .
The increase in the species richness of lepidopteran caterpillars with increasing tree species richness at the plot level is consistent with theoretical expectations, because forest plots with more tree species may support a higher abundance of consumers (which, in turn, can promote diversity, see below) and provide more niche opportunities for herbivore species adapted to different environments or host plants (Peterson et al., 2011). Our finding that many tree species contributed significantly to increasing lepidopteran species richness (Table S4) indicates that complementarity among tree species, and not a sampling effect of having few tree species with a particularly high lepidopteran species richness, was the driving mechanism underlying these relationships (see e.g. Hector et al., 1999;Loreau et al., 2001 for a general discussion of this issue in biodiversity studies). Because many lepidopteran species show preferences for specific host plants (Forister et al., 2015) which often is the result of long co-evolutionary adaptations (Jorge et al., 2014;Peralta et al., 2015;Volf et al., 2018), positive effects of tree species richness on herbivores can also be expected to be reflected by herbivore phylogenetic diversity (Tucker et al., 2017). We found strong positive relationships between tree species richness and herbivore phylogenetic diversity expressed as Faith's PD. However, the latter is known to strongly correlate with taxonomic diversity (Tucker et al., 2017), which explains the high similarity of the results for analyses of herbivore species richness and Faith's PD in our study. Additional metrics of phylogenetic diversity might therefore provide further insight into evolutionary constraints involved in the relationship between tree and herbivore diversity. In our study, the for example, in the monoculture plots of our study. However, this requires closer inspection, as others have emphasized that competition between folivores might be limited (Lawton & Strong, 1981), and ecological processes might be difficult to infer from phylogenetic patterns alone (Godoy, Kraft, & Levine, 2014;Mayfield & Levine, 2010). In comparison to MNTD, MPD reflects more basal splits in the phylogeny which are probably representative of evolutionary adaptations of different lepidopteran lineages to differences in palatability and defense traits of specific plant lineages (Pellissier et al., 2013;Volf et al., 2018). Our results show that such evolutionary adaptations that occurred deep in the phylogeny play an important role in driving the response of herbivore species richness and diversity to changes in tree species richness. This has implications for our understanding of the consequences of ongoing environmental changes for the structuring of multitrophic communities as well as for biodiversity conservation, because changes in producer diversity or composition will have non-random effects on herbivore communities and the diversity of evolutionary information represented by these communities.
The effects of tree species richness on both taxonomic and phylogenetic diversity of lepidopteran caterpillars were partially explained by changes in tree phylogenetic and functional diversity.
This indicates that functional traits -and in particular phyloge- However, significant effects of tree species richness that were not mediated by tree phylogenetic or functional diversity remained, suggesting that tree species richness as a diversity metric can include information that is not easily captured by functional diversity or phylogenetic diversity metrics based on a limited set of traits (Devictor et al., 2010). Indeed, a recent study found that in the same experiment as the one studied here, functional and phylogenetic diversity did not offer better explanations of stand-level productivity and carbon storage than did species richness per se (Huang et al., 2018). We note that our analyses are limited by the This is also reflected in analyses of functional identity and diversity effects on early growth in the BEF-China experiment (Kröber et al., 2015).
It is important to note that effects of plant diversity on the diversity of lepidopteran caterpillars were strongly mediated by lepidopteran abundances, as indicated by the path analyses and their support of direct pathways from lepidopteran abundance to diversity (rather than the other way around or reciprocal effects). This suggests that mechanisms related to the more-individuals hypothesis (i.e. species richness being promoted by factors that increase abundance, as more individuals can support viable populations of more species; see e.g. Storch et al., 2018) play an important role in structuring host plant-consumer relationships at our study site (see also Schuldt et al., 2019). Direct pathways between tree diversity and lepidopteran abundance were generally stronger than the direct pathways between tree diversity and lepidopteran diversity tested in our alternative path model approach. Lepidopteran abundances might influence lepidopteran diversity in multiple ways, e.g. reflecting local extinction dynamics or optimal foraging and ovipoisition behaviour (Scheirs & De Bruyn, 2002;Storch et al., 2018)  Faith's PD were positively related to lepidopteran abundance, MPD was not directly linked to abundance and MNTD declined with increasing lepidopteran abundance in the path models.
The latter is likely due to an increase in the number of common lepidopteran species in species-rich families ( Figure S2), which decreases the mean dissimilarity to the nearest relative. The covariation between MPD and lepidopteran abundances might indicate more complex, reciprocal interactions between the two, rather than a unidirectional effect of abundance on MPD (Storch et al., 2018).
An important general finding of our study is that changes in tree species richness had consistent effects on the taxonomic diversity of lepidopteran caterpillars, irrespective of whether the effective number of lepidopteran species based on Hill numbers was determined largely by rare or only the more common species. This indicates that with increasing tree species richness there is not only an accumulation of less abundant species, which might be functionally less effective (Smith & Knapp, 2003, but see Leitão et al., 2016), but also an increase in very abundant herbivore species which can be assumed to particularly strongly influence the ecological impact of the herbivore communities. An investigation at an earlier stage of forest development at the same study site showed that tree species richness strongly promoted herbivore abundances (Zhang et al., 2017). Our study adds insight into the potential ways in which herbivore communities are structured by tree species richness and how they feed back on producer communities and their functioning by showing that the general increases in herbivore abundances translate into an increased number of common and therefore functionally important herbivore species. These patterns might help explain the finding that herbivore damage at our study sites increased with increasing tree species richness (Schuldt et al., 2015). We note that our study, except for the experimental manipulation of tree species richness, is based on observational data, and further experiments manipulating the species diversity and composition of higher trophic levels, as well as their overall abundance, will be helpful to establish causality beyond that of tree species richness on lepidopteran abundance and diversity.
Overall, our study shows that incorporating a wider range of diversity metrics (from taxonomic to functional and phylogenetic) both at the level of plant communities and herbivore communities can help us to develop a more comprehensive understanding of how changes in biodiversity at different trophic levels will affect ecosystem structure and functions. Moreover, the important role of herbivore abundances in mediating the effects of tree diversity on herbivore diversity underlines the importance of recent calls for a more thorough investigation of the ecosystem-level consequences of current declines in the abundances of insects and other consumers (Dirzo et al., 2014;Eisenhauer et al., 2019). Therefore, we need to take multi-diversity metrics and effects of consumer abundance into account to assess the stability and susceptibility to environmental changes of forest communities and their functioning.