UvA-DARE (Digital Academic Repository) Root traits explain rhizosphere fungal community composition among temperate grassland plant species

(cid:1) While it is known that interactions between plants and soil fungi drive many essential ecosystem functions, considerable uncertainty exists over the drivers of fungal community composition in the rhizosphere. Here, we examined the roles of plant species identity, phylogeny and functional traits in shaping rhizosphere fungal communities and tested the robustness of these relationships to environmental change. (cid:1) We conducted a glasshouse experiment consisting of 21 temperate grassland species grown under three different environmental treatments and characterised the fungal communities within the rhizosphere of these plants. (cid:1) We found that plant species identity, plant phylogenetic relatedness and plant traits all affected rhizosphere fungal community composition. Trait relationships with fungal communities were primarily driven by interactions with arbuscular mycorrhizal fungi, and root traits were stronger predictors of fungal communities than leaf traits. These patterns were independent of the environmental treatments the plants were grown under. (cid:1) Our results showcase the key role of plant root traits, especially root diameter, root nitrogen and speciﬁc root length, in driving rhizosphere fungal community composition, demonstrating the potential for root traits to be used within predictive frameworks of plant – fungal relationships. Furthermore, we highlight how key limitations in our understanding of fungal function may obscure previously unmeasured plant – fungal interactions


Introduction
Soil fungi are highly diverse and drive many critical ecosystem functions (Blackwell, 2011;Frac et al., 2018), including nutrient cycling (Gui et al., 2017) and the decomposition (Zak et al., 2019) and stabilisation of soil organic matter (Clemmensen et al., 2013).Moreover, through their many complex interactions with plant roots, fungi aid plant nutrient acquisition (Averill et al., 2019), pathogen (Marx, 1972) and drought resistance (Jayne & Quigley, 2014), and play a key role in shaping plant productivity and community dynamics (Mommer et al., 2018;Liang et al., 2020).Given this, understanding the key determinants of the diversity and composition of rhizosphere fungal communities is recognised as an important goal in terrestrial ecology (Singh et al., 2004;Francioli et al., 2020).
Many factors influence the composition of rhizosphere fungal communities, both biotic and abiotic.Among them, plant species identity (Burns et al., 2015), plant phylogeny (Barber an et al., 2015) and plant functional traits (Leff et al., 2018) have been proposed as important determinants of rhizosphere fungal communities, although debate remains as to their relative importance.Patterns of plant host preferences (Dickie, 2007) and host specificity (Zhou & Hyde, 2001) have been reported for rhizosphere fungi, suggesting that plant species identity may be an important factor explaining variation in rhizosphere fungal community composition.Indeed, plant species identity at a local scale (Leff et al., 2018), and plant community composition at larger scales (de Vries et al., 2012;Prober et al., 2015), have been shown to be major drivers of soil fungal community composition.However, disagreement remains as to the importance of plant species identity as a driver of fungal community composition (Lekberg & Waller, 2016), with the suggestion that the effects of a given species on fungi are context dependent and driven by environmental or edaphic conditions (Tedersoo et al., 2016).Alternatively, should fungal preferences for plant hosts be conserved between related plant species, plant phylogenetic relatedness may determine rhizosphere fungal community composition.This might be due to the phenotypic and functional similarity that results from the phylogenetic signals found within many plant traits (Flores et al., 2014;Valverde-Barrantes et al., 2015) driving the observed fungal preferences for phylogenetically similar plants (Gilbert & Webb, 2007;Schroeder et al., 2019).However, significant uncertainty remains around the importance of plant phylogeny in explaining rhizosphere fungal community composition, with studies both providing support (Barber an et al., 2015;Schroeder et al., 2019) or not (Leff et al., 2018) for its role.
While plant species identity and phylogeny appear to contribute to patterns of plant-fungal relations, plant functional traits are also likely to mediate these interactions.While the abundance of different groups of fungi, such as saprotrophs, pathotrophs and arbuscular mycorrhizal fungi (AMF), have been shown to be related to plant traits (Eissenstat et al., 2015;Semchenko et al., 2018;Francioli et al., 2020), others have failed to find relationships between plant traits and fungal community composition (Barber an et al., 2015;Leff et al., 2018).This is despite many known interactions between plants and fungi, such as the 'collaboration gradient' (Bergmann et al., 2020), which defines how plant tissue construction strategies influence nutrient foraging via associations with fungal symbionts (Eissenstat et al., 2015;Chen et al., 2018;Zhang et al., 2019), and the growing understanding of the importance of root exudates in shaping fungal communities (Hu et al., 2018).Evidence from studies of soil microbial functioning also suggests that root traits are more important predictors of rhizosphere fungal community composition than aboveground leaf traits (Orwin et al., 2010;Cantarel & Pommier, 2015).However, despite both aboveground (de Vries et al., 2012;Semchenko et al., 2018;Boeddinghaus et al., 2019) and belowground plant traits (Legay et al., 2014;Steinauer et al., 2017;Semchenko et al., 2018) having been linked to measures of rhizosphere fungal community composition, it remains unknown whether easier to measure aboveground traits can be used as proxies for belowground trait-fungal relationships.
While plant traits may help us to understand patterns of rhizosphere fungal community composition, trait values can be highly variable due to intraspecific variability induced by both biotic and abiotic factors (Callaway et al., 2003;Berg & Ellers, 2010;Matesanz et al., 2010).Trait plasticity is an important mechanism by which plant species can both coexist in diverse communities (Abakumova et al., 2016;P erez-Ramos et al., 2019) and persist under environmental change (Henn et al., 2018).It is well known that leaf and root traits display considerable plasticity in response to environmental change (Bardgett et al., 2014;Brunner et al., 2019) such as variation in light (Keenan & Niinemets, 2016) or water availability (Fry et al., 2018).There is also growing awareness that intraspecific trait variation, including that caused by plasticity, can be similar and, in some cases stronger, than trait variation across species, with potential ecological impacts (Albert et al., 2010;Des Roches et al., 2018).Interactions between plant traits and soil fungal communities may therefore be mediated by trait plasticity in response to environmental change (Li et al., 2017;Defrenne et al., 2019;Zhang et al., 2019), however our understanding of the importance of trait plasticity for fungal communities is limited.
Here, we assessed the role of plant traits, plant species identity and plant phylogenetic relatedness in driving the rhizosphere fungal community composition across 21 species of common British grassland plants.We also aimed to quantify how these relationships are modified by changes in environmental conditions and resulting intraspecific trait plasticity.We hypothesised that root traits, rather than leaf traits, are primary determinants of rhizosphere fungal communities.We also hypothesised that trait plasticity in response to changes in environmental conditions will lead to significant shifts in rhizosphere fungal community composition.Finally, we hypothesised that the 'collaboration gradient' of root resource acquisition explains a significant portion of variation in rhizosphere fungal community structure, particularly driven by relationships between plant roots and AMF.

Field soil collection
Field soil was collected from Selside Meadows in the UK Ingleborough National Nature Reserve, Yorkshire Dales (54°10 0 47.9″ N, 2°20 0 11.1″W) in June 2017.The soil at this site has been well characterised in previous studies (Leff et al., 2018;De Long et al., 2019b), and is a brown earth (60% clay, 39% sand, < 1% silt, pH c. 5.7, 4.9% carbon (C), 0.46% nitrogen (N); De Long et al., 2019b) of the Malham series (Eutric Endoleptic Cambisols) (Cranfield University, 2020).The plant community is typical of northern UK meadow grasslands (UK National Vegetation Classification, MG3b -Anthoxanthum odoratum-Geranium sylvaticum grassland, Rodwell, 1992).Sheep graze the site throughout the year; aside from April to July when a hay cut is taken.The site is located at 303 m asl with an approximate average annual precipitation of 1550 mm and average annual minimum and maximum temperatures of 4.3°C and 10.5°C, respectively (De Long et al., 2019b).Soil was passed through a 4mm sieve and homogenised after roots and large debris were removed, and soil water holding capacity (WHC) of sieved soil was determined gravimetrically.

Plant species selection and experimental establishment
In total, 21 plant species were selected for this experiment, consisting of seven species each of forbs, grasses and legumes (Supporting Information Table S1) and representative of the communities occurring in the area from which the soil was harvested.These species exhibit a range of nutrient acquisition strategies, from fast growing, resource acquisitive species to slow growing, resource conservative species (Grime et al., 2014).Therefore, they span a wide range of root and leaf trait values and syndromes that can be present within common UK grassland species (Table S1).Seeds were sourced from Emorsgate Seeds (Norfolk, UK), and germination took place in seed trays under glasshouse conditions.Based on a preliminary trial (data not shown) seed sowing was staggered over a 10-d period so germination occurred approximately at the same time.Seedlings were transplanted into 2 l pots (10 cm diameter 9 c. 25 cm height) once the first true leaves emerged.One seedling was present per pot and those that died within the first week were replaced.Seedlings were grown for 12-wk in total under three environmental treatments.

Plant and soil harvest
At the end of the 12-wk growth period, plants and soil were harvested.In total, 305 plants were harvested, but at least four samples were present per species per treatment.We collected rhizosphere soil, defined as the soil directly adhering to the root system, once the loose soil was gently shaken away from the plant.Rhizosphere soil was immediately frozen at À20°C before it was transferred to À80°C within 8 h for DNA analysis.The remaining soil was passed through a 2-mm sieve and stored at 4°C before analysis.Leaves were sampled, wrapped in tissue paper saturated with Milli-Q water, and stored in 50-ml falcon tubes at 4°C to encourage maximum leaf hydration before weighing.Leaf traits requiring fresh biomass were measured within 48 h.Aboveground and belowground plant biomass was separated and roots were stored at 4°C before washing and trait analyses.Leaf fresh weights were recorded and trait measurements followed methodology described in P erez-Harguindeguy et al. (2013).Leaf thickness was measured using digital callipers by taking an average of three measurements across the leaf.Leaves were pressed gently under a clear Perspex sheet, photographed alongside a ruler and IMAGEJ software (Schindelin et al., 2015) was used to process the photographs and extract leaf area measurements.Leaves were subsequently dried for 72 h at 60°C and the dry weights were recorded.Specific leaf area was calculated as leaf area per unit dry mass, and leaf dry matter content was calculated as the leaf dry mass per unit leaf fresh mass.Dry leaves were ground in a ball mill (MM400 Mixer Mill; Retsch, Haan, Germany) for 40 s and analysed for total leaf C and N concentration using a Vario EL Cube (Elementar UK Ltd, Stockport, UK).Roots were washed of soil and other debris within 2 wk of harvest and stored for up to 4 wk in a 20% (v/v) ethanol solution.Fresh root mass was recorded.
Roots were scanned on a flatbed scanner (Epson 11000 XL) at a resolution of 800 dpi and WINRHIZO software (Regent Instruments, Quebec City, Canada) was used to analyse the images for root diameter, root length, root volume and root surface area, before roots were dried for 72 h at 60°C and dry weights recorded.Specific root length was calculated as root length per unit dry mass, and root dry matter content was calculated as root dry mass per unit root fresh mass.Specific root area was calculated as root surface area per unit dry mass and root tissue density was expressed as root dry mass per unit volume.Dried root samples were analysed for total root C and N concentration as above.

Soil chemical and biological properties
Plant available N was measured by extracting 2.5 g fresh soil in 12.5 ml of 1 M KCl, shaking for 1 h at 200 rpm and filtering through Whatman 42 filter paper (De Long et al., 2019a).Water extractable C and N were measured by extracting 4 g of fresh soil in 28 ml of Milli-Q water, shaking for 1 h at 200 rpm and then filtering through Whatman 42 filter paper (Jones & Willett, 2006).Water extractable N was measured as total N, nitrate and ammonia, but levels of ammonia in the soil were below the limit of detection.An c. 20 g soil sample was dried at 105°C and soil moisture content was determined as percentage of fresh soil weight.Soil total C and N was measured via elemental analysis using a Vario EL Cube (as above), following grinding of the dry soil for 40 s in a ball mill.All N analyses were performed on a Continuous Segmented Flow Analyser (AA3, Seal Analytical, Southampton, UK), while C measures were conducted on a total organic carbon analyser (TOC-L; Shimadzu, Buckinghamshire, UK).Microbial biomass C was measured using the chloroform extraction method (Brookes et al., 1985).Briefly, 5 g fresh soil was fumigated in chloroform for 24 h before extraction in 25 ml of 0.5 M K 2 SO 4 , shaken for 1 h at 200 rpm and then filtered through Whatman 42 filter paper.Microbial biomass C was calculated by subtracting the C in unfumigated K 2 SO 4 extracts from those that were fumigated.Values were corrected with a kEC value of 0.45 to account for extraction efficiency (Vance et al., 1987).
Demultiplexed raw sequence data were processed using the DADA2 pipeline (Callahan et al., 2016).Cutadapt (Martin, 2011) was used to remove primers in all orientations from the dataset while maintaining amplicon length variation.We imposed a minimum length of 50 bp to remove any small fragments at the filtering stage, at which, the error in the maxEE argument was 2.5 as this optimised the retention of reads throughout the pipeline.Error rates were calculated by the DADA2 algorithm before dereplication and merging of paired end sequences.Chimeras were subsequently removed and taxonomy was assigned using the UNITE database (Nilsson et al., 2019).The 'even' mock community from Bakker (2018) used to optimise the pipeline.Data were held within a Phyloseq (McMurdie & Holmes, 2013) object and filtered to exclude any nonfungal associated sequences variants (ASVs).Identified ASVs were filtered such that the mean relative abundance of an ASV had to be greater than 1 9 10 À5 .Sequence data were not rarefied, New Phytologist (2021) 229: 1492-1507 Ó2020 The Authors New Phytologist Ó2020 New Phytologist Foundation www.newphytologist.comas rarefaction curves indicated that all samples had reached the asymptotes of the species accumulation curves.

Plant phylogeny construction
We constructed a phylogeny of our 21 plant species to test for plant phylogenetic signals within rhizosphere fungal communities.We used the ribulose-bisphosphate carboxylase (rbcL) and maturase K (matK) genes, as published by De Vere et al. (2012).One accession was selected for each gene for each species (Table S1), genes were aligned using MUSCLE (Edgar, 2004) in MEGA X (Kumar et al., 2018), using default parameters and gene alignments were concatenated.A maximum likelihood and general time reversible model was used to construct the phylogeny with a discrete gamma distribution.The rate variation model allowed for some sites to be evolutionarily invariable (GTR+G+I).A well resolved phylogeny (Fig. 3, see later) resulted that was selected on the basis of the highest log likelihood score.Pairwise distances between species were calculated by estimating the number of substitutions per site between sequences using a maximum composite likelihood model with the same parameters as above.The tree and distances were used for downstream phylogenetic analysis.

Statistical analysis
All analyses were conducted in R v.3.5.1 (R Core Team, 2018).Data were log transformed when necessary before analysis to improve the distribution of residuals within our models.We constructed a Bray-Curtis dissimilarity matrix of fungal ASV relative abundances to represent fungal community dissimilarity in our analyses below.

Trait plasticity
We assessed whether our environmental treatments (control, shaded, water limited) successfully induced root trait plasticity.We built linear mixed effects models using the 'lmer' function within the LME4 package (Bates et al., 2015) with a trait as the response variable and total plant biomass and the environmental conditions as explanatory variables, including the interaction.Species and experimental block were random effects.If the treatment or the interaction term were significant, we inferred that the trait had a 'plastic' response to the environmental conditions, not just differing due to plant biomass or ontogeny (de Vries et al., 2016).

Plant species identity
To assess the effects of plant species identity, the environmental treatments and any interaction between them in determining rhizosphere fungal community dissimilarity, we used permutational multivariate analysis of variance (PERMANOVA) with 1000 permutations, using the 'adonis' function within the VEGAN package (Dixon, 2003).We used the 'strata' argument to restrict the permutations within experimental blocks.We tested whether the effect of plant species identity on fungal community composition was mediated via plant functional groups.We used PERMANOVA to test how much variation in fungal community composition was explained by plant species identity after accounting for the effects of plant functional groups in a sequential model.Data were checked by factor for homogeneity of multivariate dispersions and all data satisfied this requirement.

Plant traits
We generated a nonmetric multidimensional scaling (NMDS) ordination of fungal communities based on the Bray-Curtis dissimilarity matrix.We used the 'ordisurf' function within the VE- GAN package (Dixon, 2003) to fit trait data to the fungal community ordination using generalised additive models (GAMs).These models fitted the trait data as a smooth response surface over the fungal community ordination accounting for both NMDS axes.We used GAMs for two reasons: first, they allow for both a linear or nonlinear fit of the traits to fungal community composition; and second, they enabled us to map the raw trait data over the ordination of fungal community composition.This greatly increases the interpretability of the resultant traitfungal correlations by contrast with methods such as Mantel tests, which first transform trait data into distance matrices.However, this analysis does not handle interaction terms, so to test for any interactions between traits and environmental treatments on fungal community dissimilarity we used PERMANOVAs, as above.
Where we identified significant relationships between traits and fungal community dissimilarity, we further tested whether relationships between plant traits and fungal communities were driven by specific trophic guilds of fungi.We used FUNGuild (Nguyen et al., 2016) to annotate ASVs as either pathotrophs or saprotrophs, at the lowest taxonomic group available, with an upper limit of family.We considered sequences classified within the Glomeromycota as AMF (Sch€ ußler et al., 2001;Lekberg et al., 2018).All remaining ASVs were left unassigned.Any ASVs that had multiple trophic assignments were not included.We calculated the relative abundance of pathotrophs, saprotrophs and AMF as the proportion of reads identified as belonging to a trophic guild relative to the total number of reads identified in the sample.Trophic guild richness was calculated as the number ASVs assigned to a given guild.We used linear mixed effects models with experimental block as a random effect to assess whether the relative abundance or richness of pathotrophs, saprotrophs and AMF was driven by plant traits, environmental treatments and the interaction between them.We used the 'r.squaredGLMM' function in the MUMIN package (Bart on, 2014) to derive the conditional and marginal r 2 values.
Finally, to identify specific ASVs that were driving trait-fungal correlations, we built random forest regression models using the 'RANDOMFOREST' (Breiman & Cutler, 2012) package.Model performance was tested by randomly permuting the trait data 1000 times and comparing the permuted and observed fit using the 'RFUTILITIES' package (Murphy et al., 2010).We cross validated the models using leave-one-out cross-validation, implemented in the 'CARET' package (Kuhn, 2019).We considered an ASV to be important if it was assigned an increase in mean squared error of greater than 2% percent (%MSE).Larger %MSE values indicated increasing importance of an ASV to the model as it means that there is a greater increase in error of the model should the abundance of that ASV be randomly permuted across the dataset.For each ASV identified as important we built a linear model of the abundance of that ASV, when present in the samples, and regressed it against the variable in question.We used this model coefficient sign and value as the indicator of the directionality and strength of the relationship between the ASV and measured variable.

Edaphic conditions
We tested whether the relationships between plant traits and soil fungi we found were a function of edaphic conditions, as soil fungal communities are known to vary with soil variables (V etrovsk y et al., 2019), particularly nutrient availability (Lauber et al., 2008).We conducted a PCA on scaled (z-score) soil variables (plant available N, water extractable nitrate and water extractable total N, total soil C and N, dissolved organic C and microbial biomass C).We extracted the loadings of the axes from this PCA and used them as continuous variables.Using GAMs within the 'ordisurf' function in the VEGAN package (Dixon, 2003), as with the trait data, we assessed the fit of the soil PCs to the NMDS ordination of fungal community dissimilarities.We tested whether any significant relationships between traits and fungal guilds were a function of the soil principal components that were significantly related to soil fungal community dissimilarity.

Plant phylogenetic relatedness
To test for correlations between plant phylogenetic relatedness and fungal community dissimilarity we used Mantel tests with 1000 permutations using the 'mantel' function within the VEGAN package (Dixon, 2003).To test whether fungal trophic guild abundances were structured by plant phylogeny we checked for phylogenetic signal (M€ unkem€ uller et al., 2012) using the 'phyloSignal' function within the 'PHYLOSIGNAL' package (Keck et al., 2016) to calculate Blomberg's K value (Blomberg et al., 2003) and Pagel's k value (Pagel, 1999).

Plant trait variation and plasticity
The species used in this study represented a range of trait syndromes across the root economics space (Fig. 1; Table S1).The growth of plants under either shaded or water limited conditions resulted in significant plant trait plasticity when compared with the control (Fig. 2; Table S2).The direction and magnitude of changes in these trait values compared with the control varied, with shading inducing stronger plastic trait responses than the water limited treatment (Fig. 2).We observed significant phylogenetic signal within the trait data, although plant trait plasticity was not structured by plant phylogenetic relatedness (Table S3).

Fungal communities
In total, 7528 fungal ASVs were identified in our analysis, which were primarily composed of Ascomycota, Basidiomycota and Mortierellomycota, with an average relative abundance of 37%, 30% and 20%, respectively, across all samples (Fig. S1).Fungal community composition was comparable with that found in similar studies on soils from the same field site (Leff et al., 2018).In total, 14.6% of ASVs were annotated with a functional guild with 361, 98 and 643 ASVs identified as saprotrophs, pathotrophs or AMF, respectively, with an average relative abundance of 4.06%, 0.96% and 4.70%.Despite inducing significant trait plasticity, the imposed environmental treatments explained little (1.2%) variation in rhizosphere fungal community dissimilarity (Table 1; Fig. S2).By contrast, plant species identity explained c. 12% of variation in rhizosphere fungal community dissimilarity and captured variation above and beyond that explained by plant functional groups (Table 1).Given the low amount of variation in fungal community dissimilarity explained by the environmental treatments, and the absence of interactive effects (Table 1) with species identity and plant traits (aside from leaf nitrogen, although in no further analyses was leaf nitrogen identified as correlating with fungal community structure) (Table S4), we pooled data across the environmental treatments for the phylogenetic and trait analyses.

Plant phylogenetic relatedness and fungal communities
We detected a positive correlation between fungal community dissimilarity and phylogenetic distance between plant species (Mantel's r = 0.383, P < 0.001; Table S5).The correlation between plant phylogenetic relatedness and fungal community dissimilarity was stronger across plant functional groups than within them (Table S5).The proportion of fungal pathotrophs in the rhizosphere was weakly correlated with plant phylogeny (k = 0.331, P = 0.060; K = 0.269, P = 0.082).This was due to a greater relative abundance of pathotrophs in the rhizosphere of grasses compared with both forbs and legumes (Fig. 3).Strong phylogenetic signals were also found for the relative abundance of rhizosphere AMF (k = 0.836, P < 0.001; K = 0.571, P < 0.001), which were attributed to a lower relative abundance of AMF in the rhizosphere of grasses and higher relative abundance of AMF in the rhizosphere of legumes compared to forbs (Fig. 3).

Plant traits and fungal communities
Root diameter, root N, specific root length, specific root area and leaf C were predictive of rhizosphere fungal community composition, represented by the NMDS ordination of these communities (Fig. 4; Table 2).Root traits were stronger determinants of rhizosphere community structure, with leaf C explaining a low amount of deviance within our models (Fig. 4; Table 2).
We sought to identify which fungal guilds were driving the observed relationships between plant traits and rhizosphere fungal community composition.Root traits, but not leaf traits, were determinants of the relative abundance and richness of AMF and fungal pathotrophs and saprotrophs (Fig. 5; Table S6).We observed a positive relationship between root diameter (r 2 = 0.339, P < 0.001) and root N (r 2 = 0.211, P < 0.001) with AMF relative abundance, while specific root length was negatively related to AMF relative abundance (r 2 = 0.111, P < 0.001).Increased saprotroph abundance was explained by increased root diameter (r 2 = 0.020, P = 0.007), increased root N (r 2 = 0.015, P = 0.015) and decreased specific root length (r 2 = 0.011, P = 0.040).Root N (r 2 = 0.025, P = 0.003) was negatively related to fungal pathotroph abundance.We did not observe any relationships between the proportion of AMF, pathotrophs or saprotrophs and specific root area and leaf C. Increased AMF richness was explained by increases in root diameter (r 2 = 0.239, P < 0.001) and root N (r 2 = 0.087, P < 0.001), and decreased specific root length (r 2 = 0.084, P < 0.001).Less variation in AMF richness was explained by traits compared with AMF relative abundance, however the opposite was true for pathotrophs.Root diameter (r 2 = 0.031, P = 0.001) and root N (r 2 = 0.079, P < 0.001) were both negatively, and specific root length (r 2 = 0.019, P = 0.016) positively, related to pathotroph richness (Fig. 5).There were no interactions detected between the environmental treatments and the trait-fungal relationships described (Table S7), although pathotroph relative abundance was increased in the water limited treatment and the proportion and richness of AMF were reduced in the shaded treatment (Table S8).
We identified which ASVs were driving the observed trait-fungal relationships (referred to as 'important ASVs'see Materials and Methods).Random forest regression models were built for the traits that best predicted rhizosphere fungal communities, namely root diameter, root N, specific root length, specific root area and leaf C. We identified strong, significant relationships between root diameter, root N, specific root length, specific root area, and specific fungal ASVs (Tables 3, S9; Datasets S1-S4), but not leaf C (Table S9).Here, c. 50% of important ASVs were identified to family level or below, and only 47% could be assigned to trophic guilds (Table S10), making functional assessments of these results challenging.
The important ASVs within the root diameter models were predominately Glomeromycota.Within the respective models, important Glomeromycota ASVs were always positively correlated with root diameter and root N, and negatively correlated with specific root length (Table 3).Patterns amongst other phyla were less clear, with ASVs belonging to the Ascomycota being both positively and negatively correlated with all four traits.Despite representing c. 30% and 20% of ASVs respectively, the Basidiomycota and Mortierellomycota were poorly represented amongst the top 20 most important ASVs driving the trait-fungal models (Table 3).
The ASVs identified as important were not ubiquitous across all samples.Across the four traitsroot diameter, root N, specific root length, specific root areaimportant ASVs were present on average in 140 out of 305 samples.We found that the majority of the top 20 most important ASVs was different for each trait (Table S11).Only nine ASVs were found to be important across all four traits, but these ASVs were not well taxonomically or functionally characterised (Table S12).All of these ASVs exhibited the same pattern: if they were positively correlated with root diameter and root N, they were negatively correlated with specific root length and specific root area and vice versa.

Soil variables as predictors of fungal community dissimilarity
The PCA of soil variables revealed two main PC axes (Fig. S3) explaining c. 60% of variation in the soil data.Soil pH and moisture loaded positively and plant available N, water extractable nitrate and water extractable total N, loaded negatively on PC1.
Total soil C and N, dissolved organic C and microbial biomass C loaded positively on PC2.PC1 was influenced by both plant species identity and environmental treatment, but not plant functional group (Table S13).PC2 was unaffected by plant species identity, functional group or environmental treatment (Table S13).As measured with GAMs, PC1 (deviance explained = 6.39%,P = 0.004), but not PC2, was associated with fungal community dissimilarity.The effect of root diameter on the proportion and richness of AMF, and of specific root length on the proportion of AMF, was affected by PC1 (Fig. S4;   S14, S15), although by comparison of marginal r 2 values, the interaction explained little extra (< 3.5%) variation in the models (Tables S6, S14, S15).

Discussion
Our results showed that plant phylogeny and species identity were important determinants of rhizosphere fungal community composition, which is consistent with previous studies (Barber an et al., 2015;Schroeder et al., 2019).Our data also indicated that plant traits, especially root traits, are strong determinants of rhizosphere fungal community composition.We found that rhizosphere fungal community composition was driven by root diameter, root N, specific root length, specific root area and leaf C.However, we did not find leaf C to be predictive of the relative abundance or richness of any fungal trophic guilds, suggesting that root traits, rather than the leaf traits we measured, are the primary determinants of rhizosphere fungal communities among grassland plants species.Although our environmental manipulations induced significant plant trait plasticity, we found that the effect of plant traits on fungal community composition was independent of the environmental conditions.We observed that fungal communities were driven by traits that had both plastic (root diameter, specific root length/area) and nonplastic (root N) responses to environmental treatments.Importantly, this suggests that intraspecific variation in trait values was not sufficient to overwhelm interspecific relationships between plant traits and rhizosphere fungal communities.Our findings provide evidence that intraspecific trait plasticity in response to benign environmental change has negligible impacts on the rhizosphere fungal communities of grassland plants, and that interspecific trait variation is of greater importance.The robustness of these relationships to trait plasticity induced by environmental stress remains to be evaluated.
We found that the relative abundance of AMF was strongly driven by root traits, and we propose that this underpins all of the trait-fungal correlations we detected.We found that AMF abundance was positively related with root diameter, and negatively related with specific root length.Increases in root diameter (and therefore decreases in specific root length) offer increased opportunity for mycorrhizal colonisation due to the associated increase in root cortical tissue (Reinhardt & Miller, 1990) and the positive correlation between root diameter and mycorrhizal colonisation rate is well documented (Kong et al., 2014;Ma et al., 2018;McCormack & Iversen, 2019).Given that AMF hyphae are rich in N, compared with plant roots (Hodge & Fitter, 2010), it is likely that the increased relative abundance of AMF with increasing root diameter results in increased root N, possibly causing the positive association between root N and AMF relative abundance that we observed.These results suggested that the 'collaboration gradient' within the root economics space (Bergmann et al., 2020) is a significant driver of rhizosphere fungal communities at the local level and further highlights the importance of constructing higher diameter roots to accommodate AMF partnerships (Eissenstat et al., 2015;Li et al., 2017).This likely explains why we did not find the leaf traits we measured to be of importance for rhizosphere fungal communities, as the 'collaboration gradient' of root resource acquisition is independent of 'fast'-vs-'slow' economics as captured by leaf traits in this study primarily by specific leaf area and leaf N (Wright et al., 2004;Bergmann et al., 2020).We further suggest that the positive correlation between root diameter and AMF relative abundance can explain the observed correlations between plant traits and fungal pathotroph and saprotroph relative abundance and richness.The richness of fungal pathotrophs was negatively correlated with root diameter, and both the relative abundance and richness of pathotrophs were negatively correlated with root N. AMF are known to increase plant resistance to fungal pathogens (Cameron et al., 2013;Berdeni et al., 2018).Thus, an increased abundance of AMF within the rhizosphere could lead to a reduced abundance of fungal pathotrophs.These results are supported by our phylogenetic analysis, which revealed that increased pathotroph and reduced AMF relative abundance was associated with the functional group of grasses that typically have high specific root length and low root diameter.This finding is also in line with Semchenko et al. (2018), who found root systems of high specific root length were associated with an increased diversity of fungal pathogens.Furthermore, we found that both root diameter and root N were positively correlated with the relative abundance and richness of fungal saprotrophs, respectively.While AMF have no saprotrophic ability (Frey, 2019), they do promote faster plant litter decomposition rates (Gui et al., 2017), although the precise mechanisms involved are not yet known (Frey, 2019).Thus, the increased abundance of AMF due to increased root diameter, may explain the positive correlations between root diameter and root N with fungal saprotrophs.The negative correlation between specific root length with both AMF and saprotroph relative abundance we found further supports this hypothesis.
Overall, these results suggested that a plant's strategy to partner with AMF plays a fundamental role in determining rhizosphere fungal communities (Bergmann et al., 2020).It is of note that we used universal fungal primers to characterise the entire fungal community, which inevitably leads to a low relative abundance of AMF taxa (Lekberg et al., 2018).Given the importance of AMF to trait-fungal relationships, methods using AMF-specific primers should be used to establish whether the patterns reported here are consistent across AMF clades.Furthermore, as studies are increasingly using molecular data to infer fungal functional data (Semchenko et al., 2018;Che et al., 2019;Phillips et al., 2019), future studies are needed to assess how well molecular-derived relative abundances and richness estimates relate to functional data, such as percentage root colonisation by AMF.
Our finding that rhizosphere fungal communities can be explained by plant traits contrasts with previous studies of temperate grassland plants where relationships have not been detected (Leff et al., 2018).This may be due to differences in the timescale of studies, in that Leff et al. (2018) examined trait-microbial relationships formed over several growing seasons, whereas the plants in this study were comparatively young.Fungal communities are known to vary both temporally (Hannula et al., 2019;Kivlin & Hawkes, 2020) and with ontogeny (Gosney et al., 2014).Therefore, further study is required to monitor how seasonal, temporal and ontogenic variation affects the traitfungal relationships detected here.However, methodological differences between studies may also have contributed.Here, we used GAMs to fit the trait data to the ordination of our fungal communities.Importantly, these allow for a nonlinear fit of the traits to the fungal community data (see trait splines in Fig. 2).While we can pull out important and significant linear interactions between plant traits and specific components of the fungal communities (Fig. 3), these interactions are complex and vary in both direction and magnitude.The additive effect of this is that across the complex and multivariate nature of fungal community structure, the influence of a trait over a given community may not be linear.Allowing for this within our attempts to model trait-microbial interactions may allow us to uncover key plantsoil relationships.
Using random forest modelling we found that important Glomeromycota ASVs were strongly and positively correlated with root diameter and root N, and negatively correlated with specific root length.These results reinforce the findings of our correlative analyses in revealing a significant link between these root traits and AMF.However, across the traits modelled we could only assign 47% of the ASVs identified as important to a trophic guild.This means 53% of the ASVs identified as important for trait-fungal relationships were not represented in our correlative analyses.For example, ASV_1028 was found to have the highest %MSE within the root N model and was identified to the order Sordariomycetes (Phylum: Ascomycota).Sordariomycetes consists of several thousand species, which may be pathotrophic, saprotrophic or endophytic (Zhang et al., 2006).Consequently, this ASV remains functionally uncharacterisable and not represented within models of trait-fungal trophic guild relationships.Attempts to understand plant-fungal relationships are clearly limited by the scarcity of functional data available for soil fungal communities.Furthermore, most of the ASVs identified as important within our random forest analyses were not ubiquitous and were present on average in c. 45% of samples.In support of the findings of Leff et al. (2018), this suggests that less generalist, more specialist, taxa are driving relationships between plants and rhizosphere fungi.However, it remains unknown whether these more specialist taxa are performing plant hostspecific functions, displaying strong plant host preferences, or simply stochastically selected from a functionally redundant population of fungi originating from the diverse initial field soil community.With the increasing availability of plant trait data through global databases (Iversen et al., 2017;Guerrero-Ramirez et al., 2020;Kattge et al., 2020) combined with the ever-growing literature on the global diversity and biogeography of soil fungi (Tedersoo et al., 2014;Egidi et al., 2019;Delgado-Baquerizo et al., 2020), we move closer to the development of predictive frameworks of how plant phylogenetic and functional diversity shape global plant-fungal interactions.While much variation in rhizosphere fungal community composition remained unexplained by plant traits, our results demonstrate the potential for root traits in particular to be incorporated into a predictive trait-based framework of rhizosphere fungal communities.Moreover, given the ecological importance of soil fungi in terrestrial ecosystems, our findings provide new insights into the role of plant trait-fungal community relationships as key drivers of ecosystem functioning and plant community dynamics.Table S2 Results from linear mixed effects models assessing environmental treatment effects on plant trait plasticity.
Table S3 Phylogenetic structuring of plant functional traits and the calculated trait plasticity index.
Table S4 Results from PERMANOVA showing the effect of traits, the environmental treatments and their interactions on the whole fungal community.
Table S5 Mantel tests testing the correlation between plant phylogenetic relatedness and rhizosphere fungal community composition.
Table S6 Correlations between traits and the trophic guilds identified within the fungal community.
Table S7 Results from models explaining the proportion and richness of fungal trophic guilds by plant traits and the environmental conditions the plants were grown under.
Table S8 The effect of the environmental treatments upon the identified functional guilds of fungi.
Table S9 Summary of the performance of random forest modelling and the model validation procedure.
Table S10 Summary of the taxonomic and functional depth to which the ASVs identified as important within our random forest modelling were identified.
Table S11 Matrix showing the number of shared ASVs between traits with significant correlations to fungal community structure as determined by random forest regression models.
Table S12 Summary of the nine ASVs that were found to be important across the four significant random forest analyses.
Table S13 Drivers of variation in principal components derived from soil chemical and biological variables.
Table S14 Fungal guild abundance and richness explained by plant traits and first principal component of soil variables.
Table S15 Fungal guild abundance and richness explained by plant traits and second principal component of soil variables.
Please note: Wiley Blackwell are not responsible for the content or functionality of any Supporting Information supplied by the authors.Any queries (other than missing material) should be directed to the New Phytologist Central Office.

Fig. 2
Fig.2Plasticity of aboveground belowground plant traits.Log ratio of average trait values per species, per treatment relative to the control (above, shaded, below, water limited).Values > 0 indicate higher trait values in the treatment.Values < 0 indicate higher trait values in the control.Green bars means traits are exhibiting a 'true' plastic response as determined by linear mixed models(de Vries et al., 2016).Blue bars means traits are not exhibiting a plastic response, instead variation may be the result of ontogeny or changes in plant biomass.A value of zero means no change between treatment and control.Biomass are not coloured as used to determine plasticity but included for comparison between treatments.Box size represents the interquartile range.Whiskers extend to the largest or smallest value no more than 1.5 times the interquartile range.Data outside this range are represented by points.

Fig. 3
Fig.3The relative abundance of pathotrophs, saprotrophs and arbuscular mycorrhizal fungi (AMF) as structure by plant phylogeny.Upper plot shows species level data with Pagel's k and Blomberg's K measures of phylogenetic signal above the respective bar chart.Lower plot shows average relative abundance of pathotrophs, saprotrophs and AMF per plant functional group with SE presented.Lettering denotes significance.ns, not significant.Colours denote plant functional groups.

Fig. 4
Fig. 4 Results from the generalised additive models fitting measured trait data across the nonmetric multidimensional scaling (NMDS) ordination of measured fungal communities.(a) Root diameter, (b) Root nitrogen, (c) Specific root length, (d) Specific root area, (e) Leaf carbon, (f) NMDS ordination of fungal communities, without trait overlay.Points are coloured by plant functional group.Green splines show the fit of the trait data from high values (dark green) to low values (light green) over the ordination.The trait overlay indicates that communities, as represented by the points on the NMDS, are associated with higher or lower trait values in line with the coloured trait gradient.Note, the gradient splines would be parallel if the relationship between the trait and the community was linear.Nonlinear relationships between traits and fungal community composition are represented by curved splines.'de' shows the deviance explained by the respective model.
Geange et al. (2017)12)ironmental treatments of varying light intensity and water availability to induce trait variation and plasticity, without inducing plant stress.Semchenko et al. (2012)detected increased plant growth under 50%Geange et al. (2017)found that reductions in water addition of c. 25% are sufficient to induce plant trait plasticity.Therefore, the three treatments were as follows: control; 60% WHC with no shade cover, shaded; 60% WHC with 55% shade cover, and water limited; 45% WHC with no shade.Plants were watered every other day by weight, to the desired soil water content.The shade covers were measured to reduce light intensity by c. 52%.On average, our water limited treatments received c. 23% less water than the control treatment.Five replicates were established per species, per treatment, totalling 315 pots.

Table 1
Results from PERMANOVA showing the relative importance of the environmental treatments, plant species identity and plant functional groups on fungal community composition.

Table 2
The role of plant traits in determining rhizosphere fungal community composition as measured with generalised additive models (GAMs).

Table 3
Correlations between plant traits and the top 20 ASVs identified as 'important' for the respective random forest models regressing plant traits against fungal ASVs.consistent across traits for comparison.'% MSE increase' assesses importance of a given ASV to the model, higher values indicate higher importance.'Coefficient' is the result of a linear model regressing the relative abundance of the ASV, when it was present in the sample, against the trait, to indicate the direction of the relationship.'Samples present' indicates the number of samples a given ASV was present in, with a potential maximum of 305.Species selection with average measured trait values and GenBank accessions used to create the phylogeny.