Biodiversity-ecosystem functioning relationships in a long-term non-weeded ﬁ eld experiment

. Many grassland biodiversity experiments show a positive relationship between biodiver- sity and ecosystem functioning, however, in most of these experiments plant communities are established by sowing and natural colonization is prevented by selective weeding of non-sown species. During ecosystem restoration, for example on abandoned fields, plant communities start on bare soil, and diversity is often manipulated in a single sowing event. How such initial plant diversity manipula- tions influence plant biodiversity development and ecosystem functioning is not well understood. We examined how relationships between taxonomic and functional diversity, biomass production and sta- bility develop over 16 yr in non-weeded plots sown with 15 species, four species, or that were not sown. We found that sown plant communities become functionally similar to unsown, naturally colo- nized plant communities. However, initial sowing treatments had long-lasting effects on species composition and taxonomic diversity. We found only few relationships between biomass production, or stability in biomass production, and functional or taxonomic diversity, and the ones we observed were negative. In addition, the cover of dominant plant species was positively related to biomass production and stability. We conclude that effects of introducing plant species at the start of secondary succession can persist for a long time, and that in secondary succession communities with natural plant species dynamics diversity – functioning relationships can be weak or negative. Moreover, our findings indicate that in systems where natural colonization of species is allowed effects of plant dominance may under- lie diversity – functioning relationships.


INTRODUCTION
Positive relationships between plant diversity and ecosystem functioning have been well established in many biodiversity experiments (e.g., Hooper et al. 2005, van Ruijven and Berendse 2005, Tilman et al. 2006, Ives and Carpenter 2007, Lefcheck et al. 2015, Venail et al. 2015, Weisser et al. 2017) and often become stronger over time (van Ruijven and Berendse 2009, Reich et al. 2012, Ravenek et al. 2014, Weisser et al. 2017). These biodiversity experiments allow for changes in the relative abundances of plant species, which may lead to shifts in plant community evenness or functional diversity (Weisser et al. 2017). However, the natural process of species colonization, establishment, and succession is prevented by selective weeding of non-sown species, resulting in increased importance of selection and complementarity effects, and a possibly reduced influence of competitive exclusion (Jiang et al. 2009). Consequently, diversity-functioning relationships found in diversity experiments might deviate from those that can be found in nature (Grace et al. 2016) but see Duffy et al. (2017). Yet, in natural systems, the occurrence of low-or high-diversity plant communities is not random, but will be determined by environmental filters (Thompson et al. 2005). Therefore, understanding how diversity-ecosystem-functioning relationships develop in the presence of natural species dynamics could benefit from studying non-weeded, synthetic, plant communities with different levels of initially sown diversity that have been allowed to further develop under the same environmental conditions.
Compared to the many studies in weeded biodiversity experiments, relatively few studies have addressed how diversity and ecosystem functioning may develop in synthetic, non-weeded, plant communities (Pfisterer et al. 2004, Fukami et al. 2005, Bezemer and van der Putten 2007, Doherty et al. 2011. Such synthetic, non-weeded communities may be established as an ecosystem restoration measure, for example, to develop grasslands communities on abandoned fields (van der Putten et al. 2000). In these communities, initial plant diversity is controlled in a single sowing event, and from then on both plant diversity and productivity may respond to plant community developments. How diversity-functioning relationships develop under conditions where plant diversity is allowed to respond to changes in plant diversity and productivity is not well understood.
Previous work showed that in synthetic, non-weeded communities plant trait composition converges, while taxonomic composition remains divergent (Fukami et al. 2005), but it was not studied how these patterns link to ecosystem functioning. Bezemer and van der Putten (2007) found long-term positive diversity-functioning relationships in sown communities, while unsown communities also reached high diversity, but with low biomass production and stability. In other non-weeded biodiversity-functioning experiments, diversity-functioning relationships disappeared over time, due to high rates of colonization in species-poor communities (Pfisterer et al. 2004) and increased dominance of productive species over time (Doherty et al. 2011). Although each of these studies provides valuable insights into diversity-functioning relationships in non-weeded communities, it has remained unresolved how changes in taxonomic and functional diversity, as well as in composition and cover of dominant species under conditions of unrestricted natural colonization of non-sown species relate to the amount and stability of biomass production.
The aims of the present study were to quantify temporal changes in taxonomic and functional biodiversity-ecosystem-functioning relationships in experimental plant communities that are undergoing natural development by unrestricted natural colonization of unsown species. Our first hypothesis was that biodiversity-ecosystem-functioning relationships may be positive in general (Tilman et al. 2014), but that they are stronger for functional diversity than for taxonomic diversity, as a result of functional complementarity (Hooper et al. 2005, Schumacher and). Our second hypothesis is that relationships between biodiversity and ecosystem functioning will become weaker when time since sowing increases. First, strong competitors may outcompete subordinate plant species in productive communities, thereby reducing biodiversity but not productivity (Smith and Knapp 2003). Second, and in contrast to weeded experiments, diversity differences between the treatments may become lower over time due to the colonization of new species in low diversity treatments (Pfisterer et al. 2004). Third, diversity-functioning relationships may disappear over time as a result of convergence of functional diversity and functional composition of plant communities (Fukami et al. 2005). In addition to testing our hypotheses, we explored whether ecosystem functioning was related to the cover of abundant plant species, because plant identity may be a strong predictor of plant community productivity regardless of biodiversity (Grime 1998, Lep s et al. 2001, Smith and Knapp 2003, Doherty et al. 2011.
In a 16-yr old plant diversity experiment, we compared high-diversity sown (HDS), low-diversity sown (LDS) and unsown, naturally colonized (NC) grassland plant communities, which were allowed to develop without any other human interference except mowing and aboveground biomass removal once a year (van der Putten et al. 2000, Fukami et al. 2005, Bezemer and van der Putten 2007. First we investigated temporal patterns in taxonomic and functional diversity and composition, as well as biomass production and stability in HDS, LDS, and NC. Then, we related community diversity patterns, and the cover of individual plant species to temporal patterns in biomass production and stability of biomass production. We discuss how sowing-associated priority effects, competitive exclusion, and effects of individual plant species may explain our findings.

Plant species composition
In 1996, we initiated a long-term grassland biodiversity field experiment (CLUE: changing land usage, enhancement of biodiversity and ecosystem development) at Mossel, the Netherlands (52°04 0 N, 5°45 0 E). The experiment was set up on a former maize field and consisted of five randomized blocks that each contained three 10 9 10 m plots (van der Putten et al. 2000). After preparing the soil with a rotavator, within each block one plot was sown with 15 mid-successional plant species (high-diversity sown; HDS), one with four mid-successional species (low-diversity sown; LDS) while the third plot was not sown (natural colonization; NC). There were five blocks. Seeds were obtained from Cruydt Hoeck (Nijeberkoop, Netherlands), a small company that collects seeds locally. Details about the sowing treatments and species combinations are presented in Appendix S1: Table S1. After sowing, plant communities were allowed to develop naturally and were not weeded. At the end of each growing season (late September/early October) each plot was mown and all aboveground biomass was removed.
In each plot, we established 12 permanent quadrats of 1 9 1 m. From 1997 to 2012, plant community composition was measured yearly in July, at peak standing biomass, by visually estimating the percent cover of each of the plant species in each quadrat. Percent cover was estimated in six different classes (1%, 2%, 5%, 15%, 25%, and 50%) until 2001 and from 2002 onward the actual estimate of the percent cover was recorded as the real cover value. Adjacent to each permanent quadrat, aboveground biomass was clipped annually at the end of the growing season (August) in 0.25 9 0.25 m subplots. Plant biomass in the subplots was dried at 70°C until constant weight.

Plant traits
For each of the plant species we obtained the values of five key plant traits that are related to plant biomass production (Diaz et al. 2004): leaf dry matter content (LDMC), specific leaf area (SLA), plant height, leaf size, and leaf mass. These traits were obtained from the LEDA traitbase (Kleyer et al. 2008). Using LEDA-based trait data allowed us to include trait values for the majority of the plant species in our experiment, including the ones that had disappeared from the plant community over time. We are aware that trait values can be context dependent and that traits measured in the field may link more closely to local plant community processes and ecosystem functioning (Cordlandwehr et al. 2013). Therefore, in order to compare how our field-based measurements related to LEDA-obtained trait values we determined four trait values (SLA, plant height, leaf size, and leaf mass) from five dominant plant species in the field in September 2013 according to Cornelissen et al. (2003) and P erez- Harguindeguy et al. (2013). For each species, we collected a maximum of 10 adult plants in each plot (see Appendix S1: Table S1 for the total number of individuals used per species per trait value). For each individual plant, we measured plant height (m) and collected two young, fully-grown leaves for which we determined leaf area (mm 2 ), leaf mass (mg), and specific leaf area (mm 2 /mg). We measured leaf area on fresh leaves within 4 h after collection by scanning the leaves. Leaf mass was determined after drying at 60°C for 48 h. More than one-half of the field-obtained measures differed significantly from the LEDA-derived values (Appendix S2: Table S1). However, high field trait values often were related to high LEDA trait values and significant differences were often due to single August 2018BIODIVERSITY-ECOSYSTEM FUNCTIONING 1837 LEDA values or very low variation in the field. None of the field-obtained measures differed among the three sowing treatments (Appendix S2: Table S2).

Data analysis
Plant species and functional diversity were calculated per plot using the average cover for each species over the 12 permanent quadrats of 1 m² each. We omitted tree species for calculating functional diversity indices, because all trees were in the sapling stage and functional traits were not representative. Trait values were available for 83 out of 116 herbaceous plant species that were recorded in the plots, and these 83 species made up 94% of the total plant cover on average. Species for which no trait values could be collected were omitted when calculating trait diversity and trait composition. Species diversity was calculated based on all species recorded and expressed as plant species richness (the number of plant species) and Shannon's diversity index (H 0 ; H 0 = À∑ p i ln(p i ), where p i is the relative cover of plant species i; Magurran 2003). We calculated functional diversity as functional richness and functional divergence (Petchey and Gaston 2006, Villeger et al. 2008, Mouchet et al. 2010 using the R script FDind . Functional richness represents the multivariate trait space that is occupied by a community without taking into account plant cover (Villeger et al. 2008). For a single-trait approach this means that functional richness is the difference between the maximum and minimum trait value in a community (Villeger et al. 2008). Functional divergence is the deviance of species from the mean distance to the center of gravity in the multivariate trait space weighted by the relative cover of the species (Villeger et al. 2008). Low functional divergence indicates that abundant species have trait values that are close to the community-weighted trait average, while high functional divergence indicates that the most abundant species have more extreme trait values (Villeger et al. 2008). Before calculating functional richness and divergence plant trait data were standardized using a Z-transformation to obtain a mean of zero and a standard deviation of one because they were measured in different units.
Stability of biomass production was calculated as the mean biomass production in each plot divided by the standard deviation (Tilman et al. 2006). Stability needs to be calculated across multiple years, or cohorts of years. We therefore divided the 16-yr period in four cohorts of 4 yr each (1997-2000, 2001-2004, 2005-2008, 2009-2012). Using four cohorts of 4 yr allowed us to include the full 16 yr of experimental data and examine temporal changes in stability, while periods were long enough to exclude overriding effects of extreme climate events (Isbell et al. 2015). We also calculated stability using different lengths of the cohort time periods resulting in overall qualitatively similar temporal stability patterns (Appendix S3).
We tested the effects of treatment and sampling year on plant species and functional diversity, and on biomass production and stability of biomass production with a general linear mixed model. Treatment (HDS, LDS, NC) was used as fixed factor and year as a predictor variable, and we used a random factor to account for repeated measures. Residuals were tested for normality with a Shapiro-Wilk test and data for homogeneity of variances with a Levene's test. The relationships between species richness, Shannon's diversity, functional richness, and functional divergence, representing plant community diversity, with biomass production and stability were determined using a general linear mixed model with biomass or stability as a response variable and the diversity indices and year (or time period in the case of stability) as predictor variables. We included year as a factor to allow for non-linear temporal changes in the diversity-productivity relationship. Block was included as a nested random effect and we accounted for the repeated sampling design using a compound symmetric correlation structure of the residuals for each plot (Zuur et al. 2009). Moreover, to determine the relationship between diversity and biomass production or stability independent of the effect of sowing treatments, we performed a similar general linear mixed model that included sowing treatment as a fixed factor, in addition to the diversity indices and year. Finally, we used multiple regression analyses to test how the cover of individual species related to biomass production and stability of biomass production. We focused on dominant species of which the cover was higher than the mean cover of all species.
We used multivariate redundancy analysis (RDA) to test how species and trait composition varied over time for the three different sowing treatments (van den Brink andter Braak 1999, Lep s andSmilauer 2003). Plot identity was used as covariate. We plotted the canonical regression coefficients from the RDA analyses against time to graphically present changes in species and trait composition over time for each of the treatments. This is similar to principle response curves (van den Brink and ter Braak 1999) but without one of the treatments being used a reference (Pierik et al. 2011). We tested the significance of treatment 9 year interactions on species and trait composition using a Monte Carlo permutation test with 999 permutations. Permutations were restricted using a split-plot design, where permutations were allowed to occur across field plots within years, but not across years, i.e., field plots were used as split plots and years as whole plots. Split plots were permuted freely within whole plots, and whole plots were not permuted. In addition, we used NMDS analyses to visualize how the plant community changes over time.
General linear mixed models, multiple regression analyses and NMDS analyses were carried out in R version 3.4 (R Development Team 2013). RDA analyses were performed using CANOCO 5 (ter Braak and Smilauer, 2012).

Temporal patterns in diversity and ecosystem functioning
The effects of the three sowing treatments on Shannon diversity and functional divergence changed over time (significant interaction between treatment and period) and a similar trend was observed for taxonomic richness (Table 1; Fig. 1). Taxonomic richness and Shannon diversity were highest in NC plots and initially lowest in LDS plots, but richness and diversity increased over time in LDS plots (Fig. 1a, c). Functional divergence in HDS plots was initially lowest, but increased over time toward the levels in the NC and LDS plots (Fig. 1d). This means that abundant plant species in the HDS plots initially had functional traits close to the community average, but that traits of the abundant species diverged more from the community average over time in HDS plots. Functional richness differed between periods, but was not affected by sowing treatments (Table 1, Fig. 1b).
The sowing treatments had a long-lasting effect on plant species composition (first RDA axis: F = 61.1, P = 0.005, explained variation = 24.5%; second RDA axis: F = 35.8, P < 0.001, explained variation = 12.1%) and the species composition between the treatments remained different over time (Fig. 2a, Appendix S4: Fig. S1). NMDS analysis showed that the cover of species like Anthoxanthum odoratum and Taraxacum officinale was higher in NC treatments, while the cover of Poa pratensis or Plantago lanceolata was higher in HDS treatments (Appendix S4: Fig. S1). Still, in all treatments, the species composition developed in the same direction, i.e., cover of Senecio jacobaea (syn Jacobaea vulgaris) and Cirsium arvense decreased over time, while cover of Festuca rubra and Holcus lanatus increased (Appendix S4: Fig. S1). In contrast, plant trait composition converged between the three treatments (first RDA axis, F = 48.5, P < 0.001, explained variation = 20.5%; second axis, F = 61.1, P = 0.096, explained variation = 19.5%). In particular, SLA, leaf area, and leaf mass of the plant community decreased over time in NC plots toward levels in HDS and LDS plots (Fig. 2b).  , 1997-2000; 2001-2004; 2005-2008; 2009-2012). The table shows F values (with P values in parentheses). Values in boldface type represent significant differences at P < 0.05. Generally, biomass production was highest in HDS plots and lowest in NC plots. Biomass production in sown plots (HDS and LDS) declined over time, whereas biomass production remained constant in non-sown plots (NC). The LDS and NC treatments converged over time (Table 1, Fig. 3a). Stability was initially highest in HDS plots and initially higher in sown than in unsown plots. However, stability of sown plots declined over time and the three treatments converged (Table 1, Fig. 3b). Results were qualitatively similar when the data were analyzed across other period lengths, instead of 4-yr periods (Appendix S3: Table S1).
The cover of bare soil in each plot was affected by interaction between sowing treatment and year (Appendix S5: Fig. S1), indicating that the difference between sowing treatments was not significant in all years. In early years, the amount of bare soil was relatively similar between sowing treatments while, in later years, NC had the highest percentage of bare soil and HDS plots the lowest (F 20,128 = 1.82, P = 0.025) (Appendix S5: Fig. S1). Over time the percentage of bare soil generally declined for all treatments (F 10,128 = 14.51, P < 0.001).

Diversity-ecosystem functioning
There was a negative relationship between taxonomic richness and biomass production (F 1,204 = 24.93, P < 0.001, Fig. 4). Biomass production was not affected by an interaction between year and taxonomic richness (F 15,204 = 1.06, P = 0.393) indicating that over time this relationship did not change significantly. Both Shannon diversity (F 1,204 = 1.79, P = 0.182) and functional richness (F 1,204 = 0.31, P = 0.580) did not affect biomass production. Also the interaction between Shannon diversity and year (F 15,204 = 0.90, P = 0.564), and between functional richness and year (F 15,204 = 1.46, P = 0.124) did not affect biomass production (Fig. 4). The interaction between functional divergence and year was significant (F 15,204 = 1.86, P = 0.029), indicating that the relationship between functional divergence and biomass production changed over time. In early years, the relationship between functional divergence and biomass production was generally negative, while it became neutral to slightly positive in late years (Fig. 4). When the NC plots were removed from the analyses (i.e., when only the sown plots were compared) the relationship between taxonomic richness and biomass, between Shannon diversity and biomass, and between functional divergence and biomass were negative, and there were no significant interactions between the diversity indices and year (data not shown). When sowing treatment was included in the model, the diversity-ecosystem functioning relationships generally disappeared (data not shown), indicating that the diversity-functioning relationships were largely dependent on the effects of the sowing treatments on biomass production.  Vol. 99,No. 8 Stability of biomass production was negatively affected by increasing taxonomic diversity (F 1,48 = 2.28; P = 0.012) and by increasing functional divergence (F 1,48 = 25.12, P < 0.001), but not by functional richness (F 1,48 = 1.49, P = 0.2328) (Fig. 5). Also, there were no interactions between taxonomic richness, taxonomic diversity, functional diversity, functional richness, and time period. When NC plots were removed from the analyses (i.e., when only sown plots were included) stability was significantly affected by the interaction between taxonomic richness and time period (F 3,28 = 3.78; P = 0.021) and stability still decreased with higher functional divergence (F 1,28 = 15.30; P < 0.001). When sowing treatment was included in the model, the interactive effect of species richness and year on stability of biomass production was significant, as well as the interactions between functional divergence and stability (data not shown), indicating that diversity effects on stability were not only driven by sowing treatments.

Cover-ecosystem-functioning
We found relationships between the cover of dominant plant species and biomass production and stability (Table 2). For example, biomass production increased with increasing cover of the three most dominant species, i.e., the forbs Tanacetum vulgare and Senecio jacobaea and the legume Lotus corniculatus. In contrast, the cover of other highly abundant species, such as Festuca rubra was not correlated with biomass production. Also, there were dominant species with a negative relationship between cover and biomass production, such as Taraxacum officinale (Table 2). We found that stability of biomass production increased with increasing cover of Lotus corniculatus, while it decreased with increasing cover of Holcus lanatus and Taraxacum officinale (Table 2). When sowing treatment was included in the analyses, many of the relationships between plant cover and ecosystem functioning disappeared, indicating that the relationships were strongly dependent on differences in plant cover or biomass production between sowing treatments (data not shown). However, some of the relationships were maintained. For example, biomass production increased with increasing cover of Cirsium arvense (t = 6.15, P < 0.001) and Senecio jacobaea (t = 3.83, P < 0.001) and stability increased with increasing cover of Lotus corniculatus (t = 2.99, P = 0.005).

August 2018
BIODIVERSITY-ECOSYSTEM FUNCTIONING 1841 vegetation is allowed. Moreover, in contrast to many other biodiversity studies, both under experimental (see, e.g., Tilman et al. 2014, Weisser et al. 2017) and natural conditions (Grace et al. 2016, Duffy et al. 2017, in our experiment we observed very few or negative relationships between diversity and biomass production. Instead, biomass production in our plant communities was better explained by the cover of dominant plant species (Lep s et al. 2001) than by the diversity of the plant community. Although, the LEDA trait values represented the plant traits in our experiment relatively well, deviances from field-obtained measures in our site indicate the relationship between functional diversity and ecosystem functioning might have been stronger if we would have had field-measured trait values available for the whole experimental period (Cordlandwehr et al. 2013).

Temporal patterns in diversity and ecosystem functioning
Species and functional diversity were generally higher in the natural colonization treatment (NC), than in the high diversity sown (HDS) and low diversity sown (LDS) treatments. Natural colonization is often suppressed in sown plant communities, resulting in more homogeneous communities that are less species rich than those originating from colonization without sowing (van der Putten et al. 2000, Pywell et al. 2002, Lep s et al. 2007). Particularly, plant communities sown with highly diverse seed mixtures are often less impacted by colonizers than communities sown with low-diversity mixtures resulting in increased predictability of functioning and stability of high diversity grassland communities (van der Putten et al. 2000, Pywell et al. 2002, Roscher et al. 2013). Our study supports these findings, and may be partly explained by higher openness and therefore higher potential for plant colonization in unsown plots in some experimental years (Appendix S5: Fig. S1). In addition, we also provide evidence that, in the long term, functional diversity and composition converge (see also Fukami et al. 2005, Petermann et al. 2010, Roscher et al. 2014, as well as that biomass production and the stability of biomass production in plots sown with a high-diversity seed mixture decline toward the level of unsown plots. These results indicate that, in the long run, sown non-weeded communities become functionally similar to unsown naturally colonized plant communities, even though their taxonomic composition can remain different. Our findings are complementary to classical biodiversity-functioning experiments, which often do not allow natural colonization of new plant species into the sown communities (but see Petermann et al. 2010).

Biodiversity-ecosystem-functioning relationship
In contrast to the positive diversity-ecosystem-functioning relationship often reported in grasslands biodiversity FIG. 5. The relationship between stability of biomass production and species richness (first row), Shannon diversity (second row), functional richness (third row), and functional divergence (fourth row) for the each of the experimental periods (1997-2000, 2001-2004, 2005-2008, 2009-2012). Each dot in each panel represents a plot, resulting in 3 treatments 9 5 replicate plots = 15 dots. Red dots represent HDS (high-diversity sown) treatments, green dots represent LDS (low-diversity sown) treatments, and blue dots represent NC (natural colonization) treatments.

1842
G. F. VEEN ET AL. Ecology,Vol. 99,No. 8 experiments (e.g., Lehman and Tilman 2000, van Ruijven and Berendse 2005, Tilman et al. 2006, Hector et al. 2010, Weisser et al. 2017 or natural conditions (Duffy et al. 2017, Gross et al. 2017, we found very few relationships between diversity and ecosystem functioning and the ones we found were generally negative, for both taxonomic and functional diversity. Our experiment differs fundamentally from weeded biodiversity experiments, as weeded experiments are targeted to study species extinction (Weisser et al. 2017), while our experiment allows colonization and hence the reappearance of "extinct" species. Therefore, our study resembles the successional development that may occur during the restoration of plant communities on abandoned, exarable fields. In addition, in our experiment, all treatments developed a relatively high diversity compared to weeded biodiversity experiments, because of the establishment of species in low-diversity and unsown communities (Doherty et al. 2011). Yet, weeded experiments typically observe the largest diversity effects between the monocultures and mixtures (Weisser et al. 2017). These differences in the experimental approach taken in our study could explain the different relationships between diversity and ecosystem functioning found in our experiment and those reported in weeded experiments. Although, with the few relationships found it is hard to detect a temporal pattern, our findings suggest that in line with our second hypothesis and previous findings the relationships between functional diversity and species richness and ecosystem functioning became weaker over time when natural colonization of species was allowed (Petermann et al. 2010, Doherty et al. 2011. However, previous studies investigating biodiversity-functioning relationships over time, showed that the original positive relationships with diversity often weakened after weeding was abandoned (Petermann et al. 2010, Doherty et al. 2011), while we found a negative relationship.
The diversity-functioning relationships in our experiment were largely driven by long-lasting sowing treatment effects on plant diversity. High-diversity sown (HDS) communities were harder to invade by new plant species than unsown (NC) communities van der Putten 2007, Lep s et al. 2007). As a result, sown species had strong priority effects by influencing the establishment of immigrants, thereby driving long-term divergence in taxonomic diversity and community composition (Fukami et al. 2005). These results indicate that priority effects can have a strong impact on plant community diversity and composition (Egler 1954), and can shape diversity-functioning relationships (Kardol et al. 2013, von Gillhaussen et al. 2014). This will be particularly important in situations where species colonize new sites, such as in the case of biological invasions, ecosystem restoration, or re-establishment of communities after environmental disturbances. Therefore, it will be crucial to enhance understanding on how species arrival sequences may drive diversity-functioning relationships.
The negative relationships observed in our study could be due to an increased role for competitive exclusion in the absence of weeding (Jiang et al. 2009). This is because, in highly productive communities, such as in our study, dominant species can outcompete subordinates, thereby reducing diversity, but not productivity (Grime 1998, Smith and Knapp 2003, Thompson et al. 2005, Fraser et al. 2015, leading to low-diversity but highly productive communities. This pattern is, however, not universal and the shape of diversityfunctioning relationships in natural systems may depend on ecosystem productivity or changes in environmental conditions (Flombaum and Sala 2008, Hautier et al. 2014, Grace et al. 2016. Therefore, in order to understand the shape of diversity-functioning relationships in communities where natural species dynamics is allowed, we need integrative consideration of multiple mechanisms that drive both diversity and ecosystem functioning (Grace et al. 2016).
The negative functional diversity-functioning relationship in the beginning of the experiment may be the result from species selection in sown communities in combination with high biomass production of these species. In most sown plant communities selected species have a broad range of plant functional traits resulting in strong complementarity effects and positive diversity-functioning relationships (Tilman et al. 1997, Jiang et al. 2009, Flynn et al. 2011, Polley et al. 2013. In our experiment, however, the species were selected for natural co-occurrence in later successional grasslands that are typical for the experimental area (van der Putten et al. 2000). Hence the selection and sowing of these species may inherently have resulted in environmental filtering (Keddy 1992, Weiher et al. 1998, which could lead to low initial functional diversity in the sown communities. Over time, the negative relationship disappeared, probably as a result of functional convergence, indicating that community assembly was driven by deterministic rules (Fukami et al. 2005). It is also important to note that, in our experiment, the originally Notes: The t and P values are from the regression analysis. Cover is the mean cover of plant species (%) in the experiment. Values in boldface type are significant at P < 0.05.

August 2018
BIODIVERSITY-ECOSYSTEM FUNCTIONING 1843 manipulated biodiversity turned into a dependent variable itself as we allowed free community development in all plots. This may have resulted in plant species filling up all niches until the maximum functional diversity is reached, leading to disappearance of the causal relationship between functional diversity and ecosystem functioning.

Cover-functioning relationships
We found that the cover of individual species strongly correlated with ecosystem functioning (Grime 1998, Lep s et al. 2001, Smith and Knapp 2003, Doherty et al. 2011. Stronger relationships between the cover of certain species and ecosystem functioning than between diversity and ecosystem functioning may be partly explained by larger variation in species cover than in the plant diversity. Still, productive species, such as Tanacetum vulgare (Grime 1998, Smith and Knapp 2003, Thompson et al. 2005, or species with unique functions, such as the nitrogen-fixer Lotus corniculatus (Thompson et al. 2005) may drive biomass production and stability, regardless of diversity (Smith and Knapp 2003). Strong positive correlations between abundances of less dominant species with ecosystem functioning could be due to preference of these species for certain vegetation types or biotic or abiotic conditions, while negative correlations may result from increased abundance of subordinates in unproductive plots. Our results suggest that specific characteristics of dominant species can be as important for ecosystem functioning as diversity per se (Mokany et al. 2008, Finegan et al. 2015. Moreover, our findings that several community properties can be related to ecosystem functioning, show that a complex network of underlying ecological processes jointly drives productivity and diversity (Hautier et al. 2014, Grace et al. 2016.

CONCLUSION
We show that in an experimental design including both non-sown and sown plant communities, sown plant communities become functionally similar to unsown naturally colonized plant communities, while sowing still had long-lasting effects on species composition and diversity-functioning relationships. This suggests that priority effects can play an important role in shaping diversity-functioning relationships in natural communities, for example within the context of exotic invasions, ecosystem restoration, or recolonization after natural disturbances. In contrast with previous biodiversity experiments (van Ruijven and Berendse 2005, Tilman et al. 2014, Weisser et al. 2017, we found very few and negative diversity-ecosystem-functioning relationships. This indicates that, in communities where plant diversity is manipulated by single-sowing events, such as during ecosystem restoration, biodiversity-ecosystem-functioning relationships may differ from those found in weeded experiments. In the absence of weeding, biodiversity itself might turn into a response variable and obscure causal relationships between diversity and ecosystem functioning. To further enhance insights in biodiversity-ecosystem functioning relationships, we propose that including natural colonization and priority effects in plant diversity-related community studies is necessary.