Invasive legume affects species and functional composition of mountain meadow plant communities

Plant invasions are among the key drivers of global biodiversity and ecosystem change. They often cause reductions in native species richness and overall biodiversity. Nitrogen-ﬁxing plants are prob-lematic as they affect soil nutrient availability and outcompete species of nutrient-poor sites. Here we assessed the impacts of the legume Lupinus polyphyllus on species and functional diversity of mountain meadow communities in the UNESCO Biosphere Reserve Rho¨n. We compared species diversity (rich-ness, evenness and effective species number), functional diversity (functional richness, evenness, divergence and dispersion) and similarity of plots in three characteristic vegetation types ( Nardus grassland, mesic and wet mountain hay meadows) between different lupine cover classes. We calculated community weighted means (CWMs) of single plant traits and plotted them against lupine cover classes. The invasion of L. polyphyllus homogenizes vegetation composition since the similarity among plots of the different vegetation types increased with increasing lupine cover. It signiﬁcantly affected species diversity in terms of richness and effective species number and the functional divergence of the vegetation. The trait set of species occurring together with lupine was shifted towards more competitive trait values. We demonstrate strongly negative impacts of L. polyphyllus on different mountain meadow vegetation types since L. polyphyllus , fosters the growth of competitive species and leads to overall more productive plant communities.


Introduction
Plant invasions are considered one of the major drivers of ecosystem modification and biodiversity change at the global scale (Davis et al. 2011;Duraiappah et al. 2005;Hejda and Pyšek 2009;Keller et al. 2011). The immediate and long-term consequences of invasions may often be detrimental for the invaded ecosystems, as plant invasions have led to substantial declines in biodiversity and ecosystems functioning. Invasive plant species may also affect native plant communities through altered habitat structure and related changes in light conditions (Otte and Maul 2005) or effects on water availability (Drenovsky et al. 2012) and nutrient supply (Ehrenfeld 2010). As a result, invasive species may alter the species composition of native communities (Hejda 2013;Otte and Maul 2005;Thiele et al. 2010). They are generally considered to reduce species richness or even erode whole gene pools by the extinction of endemics (Vilà et al. 2010). Consequently, dominance of invasive plants may lead to biotic homogenization of the resident communities, which manifests as an increase in genetic, taxonomic or functional similarity (Tordoni et al. 2019). However, Sax and Gaines (2003) pointed out that non-native species might also increase species richness at regional scales.
Species diversity, measured as species richness and species evenness, is related to productivity and population dynamics and, thus, ecosystems functions and services (Mace et al. 2012). Measuring changes in these indices therefore constitutes a first approach to assess impacts of invasive species on native communities. Functional traits reflect species interactions with the biotic and abiotic environment (Brym et al. 2018) and variation in traits is therefore strongly related to species composition (Tordoni et al. 2019). As invasive plant species may increase community similarity between invaded sites, they may also change the functional trait composition of invaded communities and thereby affect ecosystem functioning. A high functional diversity as ''the range of the functional traits of the organisms in a given ecosystem'' (Tilman 2001) may be beneficial for resistance to invasion (Mason et al. 2005). Functional diversity indices, calculated from plant traits, are thus tools that adequately reflect the different components of functional diversity and are, therefore, useful to detect impacts of plant invasions.
One of the most common alien invasive species in Europe is Lupinus polyphyllus Lindl. Introduced as an ornamental and for soil melioration, it has spread in many countries and has become dominant in various ecosystems. It has become invasive in several European countries like Norway, Lithuania, Latvia and Germany, but also on other continents e.g. in New Zealand (Fremstadt 2010). Further, L. polyphyllus is among the 15 most common plant invaders and listed on the blacklist of invasive species in Germany (Nehring et al. 2013). As a legume, it may increase the nitrogen availability of soils (Holdaway and Sparrow 2006), which will deteriorate conditions for species adapted to nutrient poor sites. It may furthermore outcompete smaller species for light due to its tall stature (Otte and Maul 2005) and thus reduce species richness (e.g. Ramula and Pihlaja 2012).
In the study area, the Rhön UNESCO Biosphere reserve, L. polyphyllus was introduced in spruce forests in order to meliorate soil conditions. A change of the mowing regime due to the foundation of the Biosphere reserve in 1991 was followed by a rapid spread of L. polyphyllus. Today, it dominates the landscape, building stands with covers up to 90% on many meadows. Especially affected by the invasion are Nardus grasslands (Habitats Directive 92/43/EEC, habitat type 6230: species-rich Nardus grasslands) and mesic and wet mountain hay meadows (habitat type 6520: mountain hay meadows; Klinger et al. 2019;Otte and Maul 2005;Volz 2003), which constitute the typical vegetation types of the Biosphere reserve.
Assessing the impacts of plant invasions on the community level, but also on species and functional diversity of resident communities, is crucial in order to coordinate restoration efforts and utilize resources as efficient as possible (Hejda 2013;Vilà et al. 2010). Particularly in sensitive areas, which host many rare and endangered plant species and are hence vulnerable to biodiversity loss, information on the impacts of invasive plant species are necessary. While previous studies have mainly focused on the impacts of L. polyphyllus on species diversity and community composition (Hejda 2013;Otte and Maul 2005), the simultaneous effects of L. polyphyllus on species and functional diversity are still poorly understood. To our knowledge, few studies addressed potential habitatspecific responses to the invasive L. polyphyllus (Thiele et al. 2010).
Therefore, we present a comprehensive study testing the impacts of L. polyphyllus on community composition, species diversity and functional diversity of grasslands. Since the impact of invasive plants may vary among invaded plant communitites, and different communities may be driven by different environmental constraints, we explored the effects of L. polyphyllus on the vegetation separately for Nardus grasslands, mesic and wet mountain hay meadows.
We ask the following research questions: a. Does L. polyphyllus affect the species composition of mountain meadows?
We hypothesize that L. polyphyllus will increase the similarity across plant communities of the different meadow types.
b. Does L. polyphyllus affect the species diversity of the mountain meadows?
We hypothesize that L. polyphyllus increases species diversity (expressed by species richness, evenness, and effective species number) at low cover of L. polyphyllus by creating new niches, while all three variables will decrease at high cover.
c. Does L. polyphyllus affect multivariate functional diversity as well as community-weighted means of single traits of mountain meadows?
We hypothesize that with increasing lupine cover, sites become more productive, which will be reflected by the respective plant traits (SLA, seed releasing height, LDMC, seed number, flowering duration, life persistence). Functional diversity as well as the proportion of certain plant traits in the vegetation types will decrease with increasing cover of L. polyphyllus.

Study area
The study was carried out in the mountainous region of the UNESCO Rhön Biosphere Reserve in Central Germany, which was founded in 1991 and comprises an area of approx. 2400 km 2 . Tertiary sands and clays shaped the geology by building an elevated shelf, which is covered by basaltic rocks (Klausing 1988). Soils over basaltic bedrocks are usually well supplied with cations. However, in the Rhön region, high precipitation and traditional land use resulted in low nutrient availability and low pH values (Puffe and Zerr 1988). The mean annual precipitation at the highest elevation in the area, Mt. Wasserkuppe (950 m a.s.l.), amounts to 1,135 mm (mean of 1981-2010; DWD 2019) and the mean annual temperature reaches no more than 5.5°C (mean of 1981-2010; DWD 2019).
Traditional land-use such as regular mowing and pastoral sheep-herding in combination with low fertilizer input shaped the landscape and formed extended semi-natural grasslands with high conservation value (Otte and Maul 2005). The Rhön Biosphere Reserve comprises 8,900 ha of low-intensively used species-rich grasslands and thus plays an important role in nature conservation (Grebe 1995).

Data analysis 1. Community composition
To obtain the main floristic gradients, we performed a non-metric multidimensional scaling (NMDS) ordination. We chose NMDS as a robust, distance-based method that accurately displays the vegetation data. NMDS was calculated based on Bray-Curtis distances as dissimilarity measure with 20 random starts and three dimensions using the metaMDS function as implemented in the R vegan package 2.5-5 (Oksanen et al. 2019). We included L. polyphyllus in the ordination (for an NMDS without L. polyphyllus, see Fig. 7 in ''Appendix''). The vegetation types and the lupine cover classes were used as grouping variables in the NMDS analysis. To avoid noise in the dataset we excluded species with less than three occurrences. We included weighted Ellenberg L (light), N (nutrients), R (reaction) and M (moisture) indicator values (including L. polyphyllus) as well as species diversity indices in the analysis in order to further evaluate the effects L. polyphyllus has on the vegetation composition.
To quantify the effects of L. polyphyllus on the community composition in the NMDS, we calculated the average distance of all plots in one cover class to the respective cover class centroid and compared these distances among the cover classes. Furthermore, we compared the distances of lupine cover class centroids with each other. In order to determine whether increasing cover of L. polyphyllus homogenizes the community composition, we estimated the similarity based on quantitative Sörensen dissimilarity (similarity = 1 -Sörensen dissimilarity) between all vegetation types within the cover classes. In order to account for environmental effects on similarities, we analyzed environmental alongside lupine effects on the similarity by means of a random intercepts linear mixed effect model without interaction term using the lme4 package (Bates et al. 2015). We calculated multivariate Gower distances between environmental variables (topographic wetness index and topographic position index derived from a digital elevation model and pH measured in 2015) of each plot using the FD package (Laliberté and Legendre 2010). L. polyphyllus cover classes and environmental distances served as fixed effects, whereas plot, i.e. vegetation type combinations, constituted the random effect. Statistical significance was obtained via bootstrapped p values based on 500 bootstrap samples from likelihood ratio tests between the full model with effects in question against reduced models without the effects in question.
Additionally, we calculated permutation-based analysis of variance (PerMANOVA, adonis function, Oksanen et al. 2019) in order to test for differences in species composition between lupine cover classes. Additionally, we used pH, topographic wetness index and position index as covariates to account for potential environmental effects. Subsequently, pairwise differences in the species composition of the different lupine cover classes in each of the three vegetation types were tested using the pairwise.adonis function (Martinez Arbizu 2020). We ran PerMANO-VAs and following pairwise comparisons with 1000 permutations and Bray-Curtis distances as dissimilarity measure. We obtained adjusted p values using Bonferroni corrections for the pairwise comparisons.

Species diversity
To test for changes in species diversity among lupine cover classes, we calculated species richness, effective species number based on Shannon diversity (Jost 2006), and evenness (Hill 1973). After visually inspecting diagnostic plots (Zuur et al. 2010), differences in effective species number and evenness were tested through one-way ANOVAs separately for each vegetation type and Post-hoc Tukey tests. In the case of species richness, we computed generalized linear models (GLMs) for Poisson distributions for the single vegetation types and pairwise post-hoc Tukey tests.

Functional diversity and traits
We preselected 14 traits that we considered being responsive to changing lupine cover (life persistence, leaf dry matter content, seed number, flowering onset, end of flowering, reproduction type, leaf persistence, flowering duration, strategy type, life form, canopy height, releasing height, seed mass, specific leaf area). We obtained trait data from the open source trait databases LEDA and BIOLFLOR (Kleyer et al. 2008;Klotz et al. 2002). Several studies have shown that using trait data from databases provides valuable estimates and meaningful results (Bernhardt-Römermann et al. 2011;Bernhardt-Römermann et al. 2008;Busch et al. 2018;Hattermann et al. 2019). In order to gain a reduced set of traits, we performed a double canonical correspondence analysis (double CCA) following the approach of Kleyer et al. (2012), using lupine cover as environmental gradient. By visual inspection of the resulting ordination plots, we identified six traits that best reflected functional responses to increasing lupine cover: (1) seed releasing height as a proxy for plant height, (2) LDMC (leaf dry matter content), (3) seed number, (4) SLA (specific leaf area), (5) flowering duration and (6) life persistence. Missing trait values were either obtained from literature or the mean value of the whole genus was taken. This was necessary in 23 (5 times canopy height, 10 times LDMC, 6 times seed releasing height, 2 times seed number) of 1020 cases. Four species (Tephroseris helenitis, Ranunculus polyanthemos, Ranunculus nemorosus, Crepis mollis) had to be omitted because of lacking data. In order to achieve normal distribution of the data, seed number was log-transformed.
We calculated the diversity indices functional richness, i.e. the amount of niche space filled by the species in the community, functional evenness, i.e. the degree to which the biomass of a community is distributed in niche space, functional divergence, i.e. the degree of niche differentiation (Mason et al. 2005) as well as functional dispersion, i.e. the spread of the species in the trait space (Laliberté and Legendre 2010). Calculations were done using the dbFD function of the R package FD 1.0-12 (Laliberté et al. 2014) and included the most important dimensions of functional diversity (Mason et al. 2005). In order to test for differences in functional diversity indices between the lupine cover classes for each vegetation type, we performed ANOVAs and post-hoc Tukey tests after visual inspection of diagnostic plots (Zuur et al. 2010).
In addition, we examined if single traits were influenced by lupine cover through calculating the community weighted mean (CWM) for each trait, i.e. the averaged trait values, weighted by the relative abundance of each species (Garnier et al. 2015). In order to detect the direction of change within the single traits we plotted the CWMs against the lupine cover classes and added weighted linear least squares regression curves. We calculated linear models to assess the strength and significance of the relation between CWMs and L. polyphyllus cover.
Since we were furthermore interested in the proportions of functional groups on plots with different lupine cover we estimated the cover fractions of legumes (Fabaceae), separated into L. polyphyllus and other legumes, grasses (Poaceae, Juncaceae, Cyperaceae), and herbs (remaining species including some dwarf shrubs). For each vegetation type, we summarized the cover of all species belonging to each functional group within each vegetation type and cover class and used this as 100% cover, so that the results are represented as relative cover values of each functional group within each vegetation type.
All statistical analyses were carried out in the R statistical environment version 4.0.0 (R core team 2020).

Community composition
The NMDS ordination plot (Fig. 1a) depicted floristic variation among and within the three main vegetation types. Across vegetation types, the centroid of lupine cover class 1 was closer to cover class 2 than to class 3 and 4, while the latter were situated very close together (Table 1).
In relation to the respective class centroid, mean distance of all single plots within each lupine cover class was largest in class 1 (0.789) and smallest in class 4 (0.297). One-way ANOVA showed significant differences in distances to centroids (F 3,80 = 21.1; p \ 0.001) among cover classes. Pairwise significant differences (Tukey test) between classes 1 and 2 (p \ 0.01), 3 and 1 (p \ 0.001), as well as 4 and 1 (p \ 0.001) and 4 and 2 (p \ 0.01) showed that the vegetation composition within the plots becomes more similar with higher lupine cover. In addition, the three vegetation types showed strong overlaps for plots with higher lupine cover (Fig. 1a). In Nardus grasslands, high lupine cover plots were less scattered than those with low cover (Fig. 1b). Richness and lupine cover correlated strongest with axis one (Richness: R 2 = 24%, lupine cover: R 2 = 61%). Effective species number and evenness on the other hand showed the strongest correlation with axis two (Effective species number: R 2 = 24%, Evenness: R 2 = 26%). In mesic hay meadows, high lupine cover plots were situated very closely together whereas low lupine cover plots were rather scattered (Fig. 1c). Effective species number and evenness in mesic hay meadows had an R 2 of 19% and 36%, respectively, and correlated with axis one, as did lupine cover and the N value (R 2 = 68%, R 2 = 38%, respectively). Mesic hay meadows had higher Ellenberg R and L values in plots with lower lupine cover. They correlated strongest with axis two (R 2 = 50%, and 17% respectively). Species richness was also correlated with axis two (R 2 = 19%). Similar to Nardus grasslands and mesic hay meadows also wet hay meadows showed a clear pattern, with plots with high lupine cover being located close together (Fig. 1d). Significant variables in wet hay meadows were Ellenberg indicator values N, R, L and M as well as lupine cover, all of which correlate with axis one (N: R 2 = 23%, R: R 2 = 61%, L: R 2 = 38%, M: R 2 = 64%, lupine cover: R 2 = 31%).
Mean similarity among vegetation types increased with increasing lupine cover, roughly doubling from plots without lupine to plots with [ 70% cover (Fig. 2). Similarity was significantly affected by lupine cover (bootstrapped p \ 0.01), while the environmental distance among the plots had no effect (p = 0.93).
Community weighted means of single traits showed different responses among vegetation types (Fig. 5). In Nardus grassland communities, LDMC (R 2 = 0.52) In all three vegetation types, the proportion of grasses was considerably smaller (0.24-0.26) in the highest L. polyphyllus cover class than in the lowest cover class (0.34-0.50), whereas the fractions of herbs declined to a minor extent. Other legumes were only affected and reduced in mesic and wet hay meadows (Fig. 6).

Discussion
Our results suggest that the invasive L. polyphyllus has significant effects on species composition of mountain meadow communities. When L. polyphyllus becomes dominant, many species of the native communities drop out (e.g. Thiele et al. 2010). The similarity among different vegetation types increased, which was reflected in the ordination and through direct comparison of centroid distances. This homogenization of native vegetation communities after plant invasion has been reported in several studies (La Sorte et al. 2014). This is most probably because only a specific set of other plants are able to co-exist with dominant non-native species (Hejda 2013), because they have a certain suite of functional traits (see below) matching with the conditions created by the dominant species (mass ratio hypothesis, Grime 1998). In contrast to Ramula and Pihlaja (2012), who could not show an effect of L. polyphyllus on community composition, using pooled relevés from meadows, forests, road verges and wastelands, we found a homogenization of the vegetation with increasing cover of L. polyphyllus among vegetation types. Moreover, our results showed that the species composition differed significantly within the vegetation types as lupine cover increased.
In Nardus grassland and wet hay meadow, the Ellenberg M value was lower in sites with a dense lupine cover. As the species itself has an Ellenberg moisture indicator value of 5 (Ellenberg et al. 1992) and hence an ecological optimum on mesic sites, L. polyphyllus apparently prefers sites with lower moisture in Nardus grassland and wet hay meadows. In mesic and wet hay meadows Ellenberg L value was higher were lupine cover was lower. Due to its tall growth, it fosters shade-adapted species that have a rather low light indicator value (Otte and Maul 2005;Thiele et al. 2010). Accordingly, species with high light value occur only in low lupine cover sites, where illumination is relatively high. Surprisingly, also the N value was higher in low lupine cover stands, while one would expect a higher N value in dense lupine stands, since L. polyphyllus increases plant available nitrogen (Hiltbrunner et al. 2014). The supposed discrepancy between expected and actual N value may result from the calculation of weighted Ellenberg values since L. polyphyllus has not been assigned an Ellenberg N value. Other species that co-occur with lupine and show high covers in mesic and wet mountain hay With df, degrees of freedom; F, F-statistic, R 2 and p. adj. = adjusted p value. Asterisks mark the significance level (* p \ 0.05; **p \ 0.01; ***p \ 0.001) meadows have rather low or no N value, such as Deschampsia cespitosa, Anthoxantum odoratum, Alopecurus pratensis, Cirsium palustre or Plantago lanceolata. While our study was not designed to examine the soil conditions of lupine stands, these results may suggest that the change in vegetation composition is rather related to changing light and not soil conditions. Many studies have shown that the character of the invaded community is decisive for the magnitude of the impact of invasive species (e.g. Hejda 2013; Mason et al. 2009). As pointed out by Thiele et al. (2010) especially L. polypyllus strongly interacts with the respective habitat type, which is in line with our findings. Species diversity was mainly affected in Nardus grasslands and wet hay meadows, where the introduction of L. polyphyllus increased the overall number of species. Similar results were found in studies conducted in Arizona (US), Australia and Uruguay, where exotic species increased local plant species diversity (Sax and Gaines, 2003 and see Davis, 2009). However, species richness and effective species numbers tended to decrease in the highest cover class, probably because small species will be shaded and outcompeted in the long term (Hejda 2013;Thiele et al. 2011). Mesic hay meadows were least affected by L. polyphyllus, as no significant impacts on species diversity indices were found. Species that co-occur here with lupine are often invasive in other regions of the world, such as Anthoxantum odoratum, Trisetum flavescens, Festuca rubra and Holcus lanatus (Gross et al. 2010;Hejda 2013). Thus, these species have high competitive abilities and can probably better cope with high lupine cover (Hejda 2013).
Multivariate functional diversity was not much affected by the lupine invasion. Functional divergence in mesic hay meadows and wet hay meadows showed increasing values in high lupine cover classes. This pattern of functional divergence may indicate a higher degree of niche differentiation in stands with high lupine cover (Funk et al. 2017) due to the increasing competition for light through lupine. Functional b Fig. 3 Box whisker plots of the species diversity indices, with dots representing outliers. Whiskers represent minimum and maximum of the data except for the outliers, the box upper and lower quartiles and the black dot the median. Equal letters indicate homogenous groups based on Tukey post hoc tests dispersion showed a similar pattern in mesic hay meadows which may indicate a change of the spread of species in trait space (Laliberté and Legendre 2010) that may be caused by L. polyphyllus changing the location of the center of this trait space due to its high abundance. Concerning the CWMs of single traits, all three vegetation types were significantly influenced by L. polyphyllus. With high cover of lupine, they harbored taller species with a lower LDMC. Large plant height implies high competition capacity (Bernhardt-Römermann et al. 2008), highlighting the increased competition capacity of co-occurring species in dense lupine stands. Species with low LDMC are associated with rather productive environments (Kleyer et al. 2008). In Nardus grasslands, lupine cover also significantly increased SLA, which may indicate higher growth rates, higher leaf N concentrations (Westoby et al. 2002), and a higher productivity (Duru et al. 2014). Moreover species with a high SLA quickly invest nutrients in high-quality biomass and growth (Westoby et al. 2002). Taken together, increasing seed releasing height, decreasing LDMC in all high lupine cover vegetation types and increasing SLA in Nardus grasslands present evidence for higher productivity of lupine-dominated sites. In accordance with Tordoni et al. (2019), our study suggests that a homogenization of the vegetation communities as well as significant shifts of single functional traits among the vegetation types occurred due to the L. polyphyllus invasion; yet the multivariate functional space itself was not affected. We suggest that since L. polyphyllus originates from mountain meadows (Starfinger and Kowarik 2003) it might, generally, have a similar set of functional traits as the Our results showed that the proportion of functional plant groups changed substantially and consistently in all vegetation types with increasing L. polyphyllus cover. In unfertilized one-cut meadows the fraction of legumes is typically 10% (Voigtländer and Jacob 1987), while grasses and herbs have both a fraction of 45%. In the present study, the total proportion of legumes rose considerably with increasing L. polyphyllus cover, mainly to the disadvantage of grasses, which strongly declined. All three vegetation types hosted more than 30% of legumes and are thereby outside the range of unfertilized farmed grassland (Voigtländer and Jacob 1987). The changes of functional groups induced by L. polyphyllus again mirrors its shaping and homogenizing of the vegetation types.

Conclusions
Increasing cover of L. polyphyllus changed community composition and led to a homogenization of all three vegetation types. However, plots with low cover of lupine showed higher local species diversity. The primary effect of the L. polyphyllus invasion on functional diversity was a significant shift of certain species' traits to a more competitive suite of traits. Compared to mountain grasslands without lupine, the proportion of grasses dropped significantly in plots dominated by L. polyphyllus, whereas legumes increased. Since the study comprises vegetation types of high conservation value, hosting many rare and endangered species, measures to control the invasion and reduce the L. polyphyllus number are urgently needed. Availability of data and materials Data will be provided on request.
Code availability Code can be provided on request.

Compliance with ethical standards
Conflict of interest We declare that we have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The length of the arrows represents the relationship between ordination and gradient with a significance level of p B 0.05