Plant belowground organs perform essential functions, including water and nutrient uptake, anchorage, vegetative reproduction and recruitment of mutualistic soil microbiota. Recently, multivariate analyses showed that root traits of species can largely be linked to a ‘conservation’ and a ‘collaboration’ gradient. Here, we tested whether this species-level bidimensional belowground trait space also exists at the community level in grasslands. Furthermore, we tested whether the position of grassland communities in belowground trait space relates to environmental variables.

For a total of 313 species, we collected data on eight belowground traits in greenhouse and common garden experiments and supplemented it with data on bud-bank size and specific leaf area from databases. We calculated community weighted means (CWMs) of these ten traits for 150 temperate grassland plots to investigate belowground plant-trait dimensionality and its variation along ten soil and land-use parameters.

Using PCA, we found that about 55% of variance in CWMs was explained by two main dimensions, corresponding to a mycorrhizal ‘collaboration’ and a resource ‘conservation’ gradient. Frequently overlooked traits such as rooting depth, bud-bank size and root branching intensity were largely integrated in this trait space. The two plant-strategy gradients were partially dependent on each other, with communities that do ‘outsourcing’ of resource uptake to mycorrhizal fungi along the collaboration gradient also being more ‘slow’ along the conservation gradient. (i.e. high root tissue density and high root weight ratio). ‘Outsourcing’ communities were also more often deep-rooting and associated with soil parameters, such as low moisture and sand content, high topsoil pH, high C:N and low δ15N. ‘Slow’ communities had large bud banks and were associated with low land-use intensity, high topsoil pH, and low nitrate but high ammonium concentration in the soil. Surprisingly, we did not find an association of phosphorus availability with the mycorrhizal ‘collaboration’ gradient.

In conclusion, the ‘collaboration’ and ‘conservation’ gradients previously identified among species scale up to the community level in grasslands, encompass more traits than previously described, and vary with the environment.

Data on grassland vegetation composition

The plant-community data used as a baseline for Central European agricultural grassland vegetation originate from the ‘Biodiversity Exploratories’ project (Fischer et al., 2010). In each of three regions of Germany, the Schwäbische-Alb (south-western Germany), Hainich-Dün (central Germany) and Schorfheide-Chorin (north-eastern Germany), 50 grasslands covering a wide range of land-use intensities were selected. From 2008 to 2019, the vegetation composition of a 4 m × 4 m plot in each of the 150 grasslands was assessed annually in May/June by identifying all vascular plant species and visually estimating their aboveground cover. To align the species names between the vegetation and trait datasets, we standardized the species names according to the accepted names in The Plant List (www.theplantlist.org, accessed 15 June 2019, using the Taxonstand R package (Cayuela, La Granzow-de Cerda, Albuquerque, & Golicher, 2012). In total, 319 vascular plant species have been identified in the 150 grassland plots.

Plant species traits

We obtained mean species values for eight traits from four pot experiments that we performed, and for two further traits from already existing databases. For 291 of the 319 grassland species, we were able to obtain seeds from commercial seed suppliers or botanical gardens. We then performed four pot experiments to measure species traits. Taraxacum spp. are abundant in the grassland plots, though, due to their complex taxonomy, rarely identified at the species level. We here used trait values of Taraxacum campylodes for Taraxacum spp. The trait values are part of a previously published dataset (Lachaise, Bergmann, Rillig, & van Kleunen, 2021) and an unpublished dataset (Bergmann et al. unpublished data), and comprehensive descriptions of the experiments are provided in Appendix S1. In brief, we did one greenhouse experiment in which we grew 2659 individual plants, representing 216 species, for four weeks after which we weighed the roots and analysed scanned images of the roots with WinRHIZO 2017a software (Regent Instruments Inc., Canada) to determine root tissue density, specific root length, fine root diameter, root weight ratio and root branching intensity (Lachaise et al., 2021). Because these traits were measured on young root systems, most of the roots could be considered fine roots with principally a resource uptake function rather than a transport or storage function. We did a second greenhouse experiment using 2007 plants, representing 196 species, to determine the N content of fine roots (Fine roots %N) using isotope-ratio mass spectrometry. In a third greenhouse pot experiment, we determined mycorrhizal colonization rates for 225 plants, representing 75 species that are among the most common ones in the grasslands plots (mean cover of 65%, Appendix S3). Six weeks after inoculation with spores of Rhizophagus irregularis (Bergmann et al. unpublished data), roots were harvested and washed, and the percentage of mycorrhizal colonization was determined using the line-intersect method (McGonigle, Miller, Evans, Fairchild, & Swan, 1990). In a fourth experiment, we grew 752 plants, representing 183 species, outdoors in growth-tubes to determine the depth above and below which plants have 50% of their root biomass (Rooting depth 50%, see Appendix S1 or Schenk & Jackson, 2002b for the calculation method) for about 16 weeks. In addition, to have an estimate of the belowground regeneration potential, we extracted bud-bank size, including stem and root-derived buds occurring belowground or at the soil surface, from the CLO-PLA database (Klimešová, Danihelka, Chrtek, de Bello, & Herben, 2017) for 313 of the 319 species. Finally, to also have a reliable indicator of the plant communities’ acquisitive side of the plant economics spectrum (Allan et al., 2015; Busch et al., 2019), we extracted specific leaf area, the one and only aboveground trait in our analyses, for 279 of the 319 species from the LEDA database (Kleyer et al., 2008).

 Environmental variables of grassland plots

To relate the different dimensions of variation in trait CWMs to the abiotic environment, we used ten environmental variables related to land-use intensity and soil conditions. The goal was to capture a relatively independent set of descriptors likely to drive the belowground functioning of plants. A detailed description of each variable can be found in Appendix S2. We used the land-use-intensity index (Blüthgen et al., 2012), which aggregates information on the intensity of mowing, fertilization and grazing, and is a major driver of ecosystem properties (Allan et al., 2015). We used a variety of physicochemical indicators related to soil fertility of the topsoil (0-20 cm). Soil-moisture content and sand content were measured to capture soil water availability and texture, respectively. Soil pH was chosen, as it affects the availability of essential plant nutrients such as P in soils. We used soil extractable NO₃, extractable NH₄ and δ15N as indicators of soil nitrogen availability and related processes (Robinson, 2001; Kleinebecker et al., 2014), and the C:N ratio as a coarse indicator of stoichiometry and organic matter decomposability (Schachtschabel, Blume, Brümmer, Hartge, & Schwertmann, 1998). We further made use of resin-bag-adsorbed P and the N:P ratio to capture phosphorus availability in soil (Güsewell, 2004). Because soil volume is a central element in soil fertility and root-system distribution, we used data on bulk densities to convert per-mass nutrient concentrations to per-volume concentrations (Appendix S2). Few of the grassland-site descriptors were measured for each of the years for which we had vegetation-composition data (i.e. for the period 2008-2019). However, we tried to maximize the coverage for this period by using all available census dates for these variables (see Appendix S2 for years covered) and averaging the values per plot.

Statistical analyses

All the statistics were done using R v 4.0.1 (R Core Team, 2020).

Community weighted trait means

To characterize the plant communities of each of the 150 grassland plots based on values of functional traits of their species, we calculated community weighted means (CWMs) as

_CWMTrait= _j=1^Sp_jTrait_j

Here p_j is the relative aboveground cover of species j in the community, Trait_j is the trait value of species j, and S is the number of species in the community with available trait data. Because some plots had patches of bare soil in some of the annual vegetation surveys, and because for some species trait data were missing, we normalized plant cover to cumulate to 100% for all species with available trait data in each plot before calculating the CWMs. As trait data for most of the abundant grassland species was available, this analysis includes about 90% of the total plant cover in most plots, for most traits (Appendix S3). The only exception is mycorrhizal colonization, which is only available for 78 species, but, even for that trait, the average cover of species included is 65% (range 32 - 87%, Appendix S3).

Principal components of CWMs variation

As the CWMs of several traits were correlated (Appendix S9), we performed principal component analyses (PCA) to reduce the dimensionality of the data. To assess how robust the resulting dimensions are to the inclusion of additional information, we performed four separate PCAs. Each of these PCAs included all nine belowground traits, but they differed in that we also included or excluded _CWMSpecific leaf area, as one of the major traits associated with the aboveground ‘fast’ side of the plant economics spectrum, and that we included or excluded plant-functional-type information i.e. the percentages cover of grasses (Poales), N-fixing forbs (Fabaceae) and non-N-fixing forbs. So, one PCA included CWMs of belowground traits only (“Belowground PCA”), one additionally included _CWMSpecific leaf area (“Above-Belowground PCA”), one additionally included the proportions of Poales, Fabaceae and non-N-fixing forbs, and one included all. To increase the separation of the variable loadings (the trait CWMs) on the two first axes, we performed an ‘oblimin’ rotation on these axes for the Belowground PCA and the Above-Belowground PCA. To complement the information provided on taxonomic or phylogenetic influence on community trait values, we also looked at the ten most dominant species or taxa in the trait space formed by PC1 and PC2 and the indicator species or taxa that associated with each quadrant of the two-dimensional space formed by PC1 and PC2 (Appendix S14). CWMs are mainly determined by the values of the abundant species in a plot, which may differ in some of their average trait values from less abundant species (Lachaise et al., 2021). As measure of abundance, we used the aboveground cover of species which only provides a two-dimensional estimate of abundance (i.e. area instead of volume). Moreover, it has recently been argued that the relative aboveground cover of a species might deviate from its relative belowground cover (Ottaviani et al., 2020). Therefore, to assess how robust our analyses are with regard to weighting the species trait values, we also did our four PCAs using Community Arithmetic Means (CArMs), where the trait values are not weighted by the species aboveground cover in the community (Appendix S16). Furthermore, to compare the relationships observed at the community level and at the species level, we also did the Above-Belowground PCA using trait means of the species instead of CWMs (Appendix S15). For each PCA, _CWMRoot tissue density was log₁₀ transformed and for each trait or proportion of plant functional type, data was standardized by subtracting the mean and dividing by the standard deviation to conform to the multinormality requirements.

Associations of the principal components of CWMs with environmental variables

To test for associations between the principal components of CWMs of the grassland plots and the environmental variables, we performed multiple regressions. The PC1 and PC2 scores from each of the four PCAs on CWMs of the functional traits were used as response variables, and the environmental variables were used as predictors. Soil C:N, N:P, sand content, NH₄, NO₃, and δ15N were log-transformed before analysis to get a more regular (less clumped) distribution of the predictor values. To account for the fact that the grassland plots are located in three different regions of Germany, we also included region as a predictor in the models. However, to avoid overfitting of the models, we did not include interactions between regions and other predictors. For model reduction, backward stepwise model selection based on AIC was performed using the function step(). This procedure selects a parsimonious set of predictors while minimizing the variance inflation factor (max VIF = 3.6 for Above-Belowground PCA). Because the two first axes (PC1 and PC2) of the four PCAs produced similar scores for the CWMs of the grassland plots (all pairwise correlations of the PC1s were >0.98 and those of the PC2s were >0.67), we present the results of the analysis of the “Above-Belowground PCA” in the main text (based on the PC axes of Fig. 2; see Fig. 3), and the results for the other three PCAs in Appendix S7. We did the same for the PC3 to PC6 scores from the Above-Belowground PCA (Appendix S12), and for each of the ten _CWMTraits (Appendix S13). We further tested if the proportion of the three plant functional types (Poales, Fabaceae, non-N-fixing forbs), as related to the trait dimensions, responded to environmental variables in a similar way, and ran the same models with the proportion of plant functional types as the response variables (Appendix S11).

Please contact tom.lachaise@gmail.com for details not mentioned in methods,

otherwise please refer to Soil conditions drive belowground trait space in temperate agricultural grasslands | bioRxiv

for details about Experimental design and environmental variables.