Trait‐based responses to land use and canopy dynamics modify long‐term diversity changes in forest understories

Aim: Land use is the foremost cause of global biodiversity decline, but species do not respond equally to land-use practices. Instead, it is suggested that responses vary with species traits, but long-term data on the trait-mediated effects of land-use on communities is scarce. Here we study how forest understory communities have been affected by two land-use practices during 4-5 decades, and whether changes in plant diversity are related to changes in functional composition.
Location: Finland.
Time period: 1968-2019.
Major taxa studied: Vascular plants.
Methods: We resurveyed 245 vegetation plots in boreal herb-rich forest understories, and used hierarchical Bayesian linear models to relate changes in diversity, species composition, average plant size, and leaf economic traits to reindeer abundance, forest management intensity, and changes in climate, canopy cover and composition.
Results: Intensively managed forests decreased in species richness and had increased turnover, but management did not affect functional composition. Increased reindeer densities corresponded with increased leaf dry matter content, evenness and diversity, and decreased height and specific leaf area. Successional development in the canopy was associated with increased specific leaf area and decreased leaf dry matter content and height in the understory over the study period. Effects of reindeer abundance and canopy density on diversity were partially mediated by vegetation height, which had a negative relationship with evenness. Observed changes in climate had no discernible effect on any variable.
Main conclusions: Functional traits are useful in connecting vegetation changes to the mechanisms that drive them, and provide unique information compared to turnover and diversity metrics. These trait-dependent selection effects could inform whether species benefit or suffer from land use changes and explain observed biodiversity changes under global change.


| INTRODUC TI ON
There is an urgent need to develop better tools for monitoring biodiversity change, as the most commonly used metrics, species richness and compositional turnover, do not capture its full extent (Blowes et al., 2019;Hillebrand et al., 2018). One way to move towards understanding the processes behind biodiversity change is to establish links between the functional composition of communities and potential selective forces, such as land-use change. Changes in average functional traits can be directly informative of changing selection pressures on communities and resulting ecosystem consequences (Lavorel & Garnier, 2002;Violle et al., 2007). That is, if we knew how the functional composition of communities responded to anthropogenic pressures, we could use traits to monitor the effects of human actions on natural systems, and perhaps predict the resulting changes in ecosystem properties.
The above-ground traits of vascular plant species and communities have been found to vary on at least two important and independent axes: plant size and the leaf economics spectrum (LES).
Size-related traits determine competitive hierarchies in relation to light, while leaf economic traits describe, for example, leaf construction costs and maximum photosynthetic output Díaz et al., 2016). Plant size is often described with vegetative height, and the LES with traits such as specific leaf area (SLA) and leaf dry matter content (LDMC) (van der Plas et al., 2020).
When external drivers such as land-use change affect the optimal allocation of resources to vertical growth and to leaves, they are also expected to exert directional selection pressure on the composition of plant communities, which should manifest as changes in average functional traits.
Climate has traditionally been seen as the dominant driver of boreal biodiversity, while the impacts of other drivers such as large herbivores have been underestimated; yet they play a prominent role in community assembly (Pausas & Bond, 2018). Large herbivores have strong impacts on vegetation that include favouring traits that decrease susceptibility to grazing, such as small size and unpalatable leaves (e.g. high LDMC and low SLA, Díaz et al., 2007). Nowadays wild megafauna has been partly replaced by human livestock (Hempson et al., 2017), and humans regulate herbivore densities. In many parts of Eurasia and North America reindeer/caribou (Rangifer tarandus) is herded as semi-domesticated free-ranging livestock. Reindeer grazing is known to cause declines in total plant biomass and in favoured forage species such as large herbs and deciduous shrubs in boreal and tundra environments (Bråthen & Oksanen, 2001;Olofsson et al., 2010;Sundqvist et al., 2019). Evidence of the impact of large herbivores on plant traits exists from other systems (Cingolani et al., 2007;Díaz et al., 2007), but is very limited from the boreal zone and especially from its forested systems.
One third of global forested area belongs to the boreal zone, and a majority of it is used for industrial wood production (Gauthier et al., 2015). Management effects on ecosystems can lag decades behind changes in management practices (Muurinen et al., 2019), underlining the importance of long-term resurvey data for detecting biodiversity changes. Even though boreal forests are estimated to be highly sensitive to modern climate change (IPCC, 2018), the species compositions of both canopies and understories have shown relatively small long-term responses to climate change compared to land use (Danneyrolles et al., 2019;Tonteri et al., 2016). A possible reason for the weak climate change responses is that closed forest canopies can buffer microclimatic conditions in the understorey against macroclimatic warming (De Frenne et al., 2019). Trends in canopy cover can be affected by land-use histories, as disturbed tree-layers densen naturally when left alone, while new disturbances, such as logging, can increase temperature and light availability in the understorey (Tonteri et al., 2016). Changes in light availability should be reflected in SLA, which increases in response to shade (Dahlgren et al., 2006). A recent resurvey from European temperate forests shows increasing canopy shade to be the primary cause of changes in functional composition of herb-layers in these systems (Depauw et al., 2020), but the generality of this finding in other biomes remains untested.
The monitoring of average functional traits alongside diversity is also important because diversity is not independent of functional composition. It has long been known that high productivity can lead to decreased diversity (Grime, 1973). The unimodal productivitydiversity relationship is in part a mechanistic consequence of productivity increasing average species size, since communities composed of larger species have room for fewer plant individuals and therefore exhibit decreased diversity (Oksanen, 1996). Similar hypotheses regarding size-diversity relationships have been put forward and tested for other taxa as well (Siemann et al., 1996). Thus, in systems where human actions modify the selective landscape in relation to species size, anthropogenic disturbances are expected to also affect local diversity. For example, as discussed above, increased grazing pressure by livestock should proportionately favour smaller species and thus increase diversity (Jia et al., 2018). species benefit and which suffer from land-use changes and explain observed biodiversity changes under global change.

K E Y W O R D S
community composition, forest management, large herbivores, plant size-diversity relationship, reindeer, species diversity, temporal change, vegetation resurvey Because the processes that drive changes in plant communities may be spatially structured, changes in community properties may be as well. For example, diversity changes might be positive in some areas and negative in others, resulting in zero average change (Vellend et al., 2013). This needs to be taken into account when interpreting average changes, lest we fall victim to the ecological fallacy of assuming that all members of a population behave in the same way as the population on average (Robinson, 1950).
Here we study the temporal and spatial variation in diversity and community composition by analysing the long-term effects of two main land-use types (reindeer husbandry and forestry) and changes in canopies (canopy SLA and cover) on boreal herb-rich forest understories. By comparing diversity-, dissimilarity-, and functional traitbased approaches using hierarchical Bayesian regression modelling, we seek answers to four questions: (a) have species diversity, functional composition (measured as average height, SLA and LDMC), and species composition changed in boreal herb-rich forests over 40-50 years; (b) are there spatial differences in these changes; (c) are the changes related to changes in semi-domesticated reindeer abundance, forest management intensity, or changes in canopy SLA and cover; and finally (d) are the changes in functional composition correlated with changes in diversity across both time and space, hinting at processes that link traits to diversity? 2 | MATERIAL S AND ME THODS

| Vegetation data and resampling
We conducted a vegetation resurvey in boreal herb-rich forests in northern Finland (Figure 1). The resurveyed herb-rich forest sites were originally surveyed in 1968-1975 by Eero Kaakinen (Kaakinen, 1971(Kaakinen, , 1974. The purpose of the original surveys was to identify phytosociological species associations in mature herb-rich forests. For this reason, the surveys were conducted in relatively undisturbed habitats. These forests harbour higher species richness than the surrounding less fertile and homogenous boreal forests (Maliniemi et al., 2019), are an important habitat for many threatened species (Kouki et al., 2018), and are thus important to monitor. They occur typically as small patches and cover only a fraction of the forested area in northern Scandinavia reflecting the scattered distribution of calcareous bedrock and soil. The field layer is dominated by herbs such as Geranium sylvaticum and Filipendula ulmaria, ferns such as Gymnocarpium dryopteris and Athyrium filixfemina, and graminoids such as Milium effusum and Elymus caninus.
Occasionally, they also host species typical of less fertile boreal forests, such as the dwarf-shrubs Vaccinium myrtillus and Vaccinium vitis-idaea. The three most abundant tree species Picea abies, Alnus F I G U R E 1 Map of the study area. The circles depict the biogeographical district (BGD) division used in this study. The line in the middle is the southern edge of the reindeer herding area. The images on the right depict typical vegetation in reindeer-grazed (upper) and ungrazed (lower) herb-rich forests. Basemap CC-BY National Land Survey of Finland incana and Betula pubescens formed more than 80% of tree cover in our study plots, during both the original survey and the resurvey.
Using detailed notes on the locations of the original sites, we relocated and resurveyed a majority of the original 336 sites in 2013-2019. Most surveys were done in teams of ≥ 2 investigators. We omitted sites that lacked detailed enough descriptions of location, and some sites that were located very far from other sites. The final number of plots, after applying the exclusion criteria described below, was 245. Relocation error is inherent in vegetation resurveys (Kapfer et al., 2017), but was minimized as much as possible by the relocation help of the original surveyor, the patchy distribution of the studied vegetation type, and the fact that herb-rich forests are easily told apart from the boreal forest matrix. This study thus applies all of the measures listed by Verheyen et al. (2018) for minimizing observer and relocation error. During both the original survey and the resurvey, species composition was estimated from 5 m × 5 m quadrats for the field layer and from 10 m × 10 m quadrats for shrub and canopy layers. The smaller quadrats were nested in the larger quadrats. The abundances of field and shrub-layer species were estimated as absolute percentage covers, and the covers of tree-layer species were estimated as relative covers, that is, they always summed to 100%.
In addition, total canopy cover was estimated visually on a threepoint ordinal scale (0-30, 31-70, 71-100%). In this study, we focus only on the vascular plants of the field and tree layers. Species in the shrub layer (woody species generally over 0.5 m but under 2 m tall) were excluded from analyses because they were measured at a different plot size compared to the field layer. We deemed this justified because the cover of such species was low (c. 2% on average). For species list and trait values, see Table S1.1 in Supporting Information Appendix S1.

| Environmental changes
We estimated forest management intensity during the resurvey based on the protocol used for the assessment of threatened habitat types in Finland (Kouki et al., 2018, p. 180). This protocol assigns sites on an ordinal scale of habitat quality ranging from 0 to 4, describing absence of forest management impacts. Here, we inverted the scale to describe management intensity. We also omitted one site that was destroyed by deforestation and land conversion to arrive at a management impact variable ranging from 1 to 4. The management intensity assessment criteria include the structure of the canopy and vegetation layers, disturbance of the soil, and the amount of deadwood, and are described in Supporting Information Appendix S2, Table S2.2. These criteria are applied in a one-out-allout manner, meaning that if any of the criteria for a more intensive management class were fulfilled, the site was assigned in that class even if it would have been a less intensively managed site based on all other criteria. The criteria are independent of the other variables used in this study, meaning that species composition and diversity were not used as measures of management intensity. Signs of forest management were present in c. 70% of the resurveyed plots, whereas the other c. 30% were deemed unmanaged (Supporting Information Figure S2.1). Our observations in the field (including old stumps, plantations of c. 30-year-old trees, and old ditches) confirm that forest management was not concentrated near the second resampling, but was spread over a longer time period.
Reindeer are the major large grazer in our study area. In Finland, semi-domesticated reindeer occur only inside a designated reindeer herding area. Data on reindeer abundances in the different reindeer herding districts were provided by the Natural Resource Institute Finland. Reindeer densities in the herding area had increased c. 40% during the study period, from 1 to 1.4 reindeer/km² (Supporting Information Figure S2.2). We calculated 10-year average reindeer densities prior to each sampling time to be used as a predictor of vegetation changes. The population sizes of other large herbivores in the study area, as inferred from hunting statistics, are negligible (Natural Resources Institute Finland, 2020, Supporting Information Figure S2.1), with the exception of the moose (Alces alces). Although the moose population in Finland has grown due to changes in forestry practices and hunting, the growth has been rather uniform across the country (Nygrén, 2009). Moreover, herbaceous field layer vegetation constitutes less than 5% of moose summer diets (Hjeljord et al., 1990). Thus, while moose density changes might have affected forest regeneration and shrub layer composition, the direct effects of increased moose densities should be rather small and uniform across the study area. Contrarily, fertile herb-rich forests are a key foraging habitat for reindeer during summer, and we expect that changes in reindeer abundances over time can strongly affect trait composition and diversity of plant communities.
We included climate change in our analyses to test its potential effect on vegetation. Ten-year average minimum and maximum summer temperatures in our study plots have increased by 0.6 and 0.8 ℃, respectively, in 45 years (1970-2015, Supporting Information Figure   S2.4), and changes experienced by the field layer will have been buffered by the canopy layer (De Frenne et al., 2019). To avoid choosing between highly correlated maximum and minimum temperatures and implying that these interpolated variables carry more information than they do, we chose to explain vegetation responses to climate using 10-year average vapour pressure deficit (VPD) during the summer months from the Terra Climate dataset (Abatzoglou et al., 2018).
VPD calculated in this way has a correlation of 0.9 and 1.0 with average minimum and maximum summer temperatures, respectively, from the same dataset, and has the added benefit of being defined in terms

| Community properties
All data processing and analyses were done in R (version 4.0.5, R Core Team, 2019).
As species diversity measures, we calculated species richness, Shannon diversity, and species evenness. By Shannon diversity we mean the so-called true diversity of order 1, which is a unit conversion of the Shannon diversity index from entropy to species (Jost, 2006). As species evenness, we used Pielou's J, or the ratio of Shannon entropy to log-transformed species richness, which is a measure of relative evenness ranging from 0 to 1 (Jost, 2010).
We calculated community weighted means (CWMs) of logtransformed trait values for one size-structural trait and two leaf economic traits: vegetative height, specific leaf area (SLA) and leaf dry-matter content (LDMC). Trait values were log-transformed because generally, traits follow a log-normal distribution (Bjorkman et al., 2018). Species absolute covers were transformed to relative covers and used as weights. Trait values were retrieved from the TRY database (Kattge et al., 2020(Kattge et al., , version 5, 2011 and the LEDA database (Kleyer et al., 2008), and supplemented with our own measurements for common species that were not found in the databases.
In four plots, some of the traits were available for less than 80% of total cover. These plots were excluded from analyses. In the remaining plots, trait values were available for > 99.5% of total cover for all traits. References to the original trait datasets are listed in Supporting Information Appendix S3.
We also calculated a CWM for the SLA of the canopy layer. SLA was not available for the non-native Larix sibirica (present on one site), for which we used the SLA of Larix decidua instead. We used the R package vegan (version 2.5-6, Oksanen et al., 2019) to calculate temporal turnover in community composition using the version of Jaccard distance that takes into account species abundances. Plot-scale turnover was logit-transformed before analyses to satisfy the assumption of homoskedasticity. We also calculated temporal turnover for the entire dataset.

| Modelling
We used those sites for which we had the full set of both response variables (species richness, species diversity, species evenness, height, LDMC and SLA) and predictor variables (canopy cover and canopy SLA, reindeer densities, management intensity and climate).
The final dataset thus consisted of 245 resampled sites. These sites were then assigned to five biogeographical districts (BGDs, Figure 1) based on the Finnish biogeographical province division. Assignment to BGDs was done to account for spatially structured confounding variables in modelling. The original provinces were slightly modified for this study by merging districts with very few plots with neighbouring districts, and by moving plots between districts so as to make each district either completely within or outside the reindeer herding area.
We used hierarchical Bayesian regression models to analyse plot-level changes in species diversity and composition. Two kinds of models were deployed based on the nature of the response variable.
For species richness, evenness, diversity, and the three functional trait composition variables, the response variables could be modelled as a function of time, whereas temporal turnover could only be modelled with spatially varying effects, because turnover already represents the difference between two timepoints. We built four models with increasing complexity to answer our stated questions.  The structure for the expectation function of different models is presented in Table 1. Note that responses (with the exception of turnover) were fed into the model separately for each time period, that is, changes in diversity and functional composition were not pre-calculated. This was done to allow pooling of information across time and space, and to describe the data-generating process as accurately as feasible.
All models assumed a Gaussian error distribution. In the models in Table 1 Bayesian hierarchical models were built using the R package for Bayesian modelling rethinking (version 2.11, McElreath, 2020), an interface to the probabilistic programming language Stan (Carpenter et al., 2017). In all models, coefficients for continuous variables were given a Normal(0,0.5) prior distribution, coefficients for categorical variables were given a Normal(0,1) prior distribution, all standard deviations were given an Exponential(1) prior distribution, and temporal correlations in the BGD-level random effect were given an LKJ (2) prior. There were no divergent transitions in any of the Markov chains in any of our models, which would have indicated unreliable estimates of the parameter posterior distributions (Carpenter et al., 2017). In addition, all parameters converged to a stationary posterior distribution (split chain R-hat < 1.01).
To test whether vegetation height is related to evenness in time (modelled above) as well as in space, we modelled evenness as a function of height using a generalized additive model (Wood, 2011). The effect of height was entered as a thin plate spline with basis dimension 20. We also tested whether the height-evenness relationship  (1) to north (5) was different between sampling times (whether there was a timeheight interaction), but a comparison of Akaike information criterion (AIC) values indicated this was not the case.

| Question a: Overall vegetation changes
Examining average values across all biogeographical districts we detected no changes in species richness, but observed increases in plot-level Shannon diversity and evenness during the sampling interval. Average LDMC increased, but there was no detectable change in average species SLA or height ( Figure 2). Average plot-scale turnover between the sampling periods was about 71% (Figure 2). For the full dataset (i.e. the relative abundances of all species across the study area), temporal Jaccard dissimilarity was 37%, which, when divided over the average sampling interval (44 years) using the equation for exponential decay (long-term turnover = 1 − (1 − annual turnover) 44 ), corresponded to an annual turnover of 2.2%.

| Question b: Geographical variation in vegetation changes
Vegetation changes across the study area were not representative of changes in all the BGDs (Figure 2). Increases in species diversity and evenness had a clear geographical trend, with no changes in the south and clear increases in the north. Height changes also had a weak latitudinal gradient, with uncertain increases in height associated with more southern areas, and vice versa for the north. Patterns in species richness and SLA change and turnover were more idiosyncratic with respect to geography. Changes in LDMC did not differ geographically, but increased uniformly in all the BGDs.

| Question c: Causes of vegetation changes
There was no discernible effect of observed climate change on any of the studied variables (Figure 3a).
Increasing reindeer densities were associated with increasing species diversity, evenness and LDMC, and decreased SLA and height (Figure 3b). In addition, plots within the reindeer herding area had more positive trends in diversity and evenness and more negative trends in height compared to those outside the reindeer herding area.
Forest management negatively affected species richness, and increased temporal turnover compared to unmanaged forests ( Figure 3c). The responses were reinforced by increasing management intensity. Temporal trends in trait composition were only weakly related to forest management. SLA showed a tendency to increase and LDMC a tendency to decrease in the most intensively managed sites compared to natural sites. Unmanaged forests had an uncertain positive richness trend, with a mean change of +0.3 standardized units and an 80% probability of the trend being positive.
Changes in canopy SLA correlated negatively with field-layer SLA and temporal turnover (Figure 3d). Canopy SLA changes were almost completely explained by changes in the dominance of the evergreen tree Picea abies (Supporting Information Figure S2.5), with P. abies dominance corresponding to decreasing canopy SLA.
Canopy cover changes also affected these communities, with increasing and decreasing canopy covers changing communities in different directions (Figure 3e). Changes in canopy cover were positively correlated with changes in diversity, evenness and SLA, and negatively with changes in LDMC and height.

| Question d: Effect of height on diversity change
Height had a nonlinear but mostly negative spatial relationship with evenness ( Figure 4a). The temporal relationship between height and evenness was also negative (Figure 4b). Including height as a predictor somewhat decreased the inferred effects of reindeer and canopy cover changes on evenness (Supporting Information Figure S2.7).

| D ISCUSS I ON
We found that land-use, but not climate change, alters the functional composition of forest understorey plant communities and explains part of the variation in temporal biodiversity trends both at regional

| Vegetation trends vary in space
Our analysis revealed that many community characteristics were stable over the 4-5 decade study interval across the study area, but displayed trends when vegetation changes were conditioned on geographical location. Diversity, evenness, and vegetation height had different trends in different parts of the latitudinal gradient.
Apparent stability at larger scales can thus be the sum of opposing trends in different regions. We report that across the study area, species richness was on average stable as in a previous meta-analysis (Vellend et al., 2013), even though species turnover was on average very high (c. 70%). Our results are thus consistent with the hypothesis that local species richness is highly constrained by regional processes such as compensatory colonization-extinction dynamics, meaning that lost species are readily replaced by new species from the metacommunity (Supp & Ernest, 2014). Average turnover within plots was high, most likely attributable to high species richness (Supporting Information Table S2.3), a long sampling interval, and inherent relocation error. Nonetheless, temporal turnover for the full dataset was of the same magnitude as average turnover in terrestrial ecosystems in a recent meta-analysis (2% per year, Blowes et al., 2019), showing that the rate of compositional change in our data is representative of previously reported biodiversity changes.

| Land-use and canopy dynamics change the composition of field-layer vegetation
Our results provide clear evidence that reindeer control forest community composition by diminishing the relative abundance of tall plants with high SLA, thereby increasing species diversity, expanding the earlier findings from tundra and mountain forests (Kaarlejärvi et al., 2017;Sundqvist et al., 2019) to boreal forests. These findings are consistent with theoretical expectations, since SLA has F I G U R E 3 Standardized direct effects of climate change (a), reindeer (b), forest management intensity (c), canopy specific leaf area (SLA) increase (d), and canopy cover change (e) on community-level properties and turnover in boreal herb-rich forest understories, as inferred from Model 3. Reindeer effects are displayed as the direct effect of reindeer density change and the combined effects of density change and being situated inside the reindeer herding area (average difference in trend compared to a plot outside the area). Management intensity effects are contrasts against unmanaged forests (intensity class 1). The points are medians, the narrow and thick lines are 2/3 and 20/21 credible intervals, respectively, corresponding to 2/1 and 20/1 odds for the true parameter value to lie within the given interval. Response abbreviations as in Figure 2. Traits were log-transformed and turnover logit-transformed before analyses been shown to correlate positively with nutrient concentrations in leaves (i.e. forage quality, Díaz et al., 2016), while large plants have been shown to be more susceptible to mammalian herbivory . Furthermore, reindeer presence in itself seemed to affect vegetation changes even when animal densities remained the same, since changes in vegetation height and diversity were more pronounced in BGDs that were inside the reindeer herding area (BGDs 3, 4 and 5, Figure 3). However, these effects could also be related to other confounding factors that change along the latitudinal gradient.
Surprisingly, forest management intensity did not have clear effects on functional composition. Instead, it decreased species richness and increased turnover. Forest management effects on turnover and diversity were mostly direct, that is, not mediated by management effects on canopy cover and composition, as evidenced by the small difference between management effects in models that included versus excluded canopy variables. The finding of higher turnover in managed forests agrees with previous studies from boreal forests (Brice et al., 2019;Kaarlejärvi et al., 2021) and other biomes (Barlow et al., 2016;Lake et al., 2000) reporting accelerated biodiversity change following human disturbance. Taken together, our results suggest that species richness changes in vascular plant communities in response to forest management do not seem to be driven by a filtering of species based on plant height or leaf traits. However, the highest management intensity class (class 4) that included recent clearcuts had more negative SLA and more positive LDMC trends compared to other managed forests, which supports previous findings of harvesting having a signal on the light-interception strategy of the understorey community (Tonteri et al., 2016), but suggests that this effect is transient and undetectable in later forest development phases.
Canopy development also affected functional composition.
Canopy cover had increased on average (see Supporting Information Figure S2.6), leading to increased shading and decreased abundance of tall, light-demanding species with high SLA and low LDMC. SLA is a succession trait, which increases in response to shade (Dahlgren et al., 2006), whereas LDMC has probably decreased because of its negative covariance with SLA . Canopy composition, expressed as canopy SLA, also influenced the functional composition of the understorey. We would have expected changes in canopy SLA to correlate positively with nutrient availability and thus also with fast leaf economic traits (i.e. high SLA, low LDMC). Instead, we found the opposite. In our data, canopy SLA was almost perfectly negatively correlated with the relative abundance of P. abies (see Supporting Information Figure S2.5), which is the only abundant evergreen tree in these systems and frequently forms monodominant stands. This lack of diversity at one end of the trait spectrum limits the usefulness of functional traits in inferring the mechanisms by which canopy composition affects understories, since average trait metrics become perfectly correlated with all qualities of P. abies. We hypothesize that in our data, the true mechanism behind the connec-

| Changes in diversity are partially mediated by vegetation height
Changes in species diversity were mostly caused by changes in evenness, and those in turn were directly related to changes in vegetation height and its environmental drivers, namely reindeer herbivory and changes in canopy cover. Furthermore, the relationship between evenness and plant height was consistent across both space and time, providing further evidence that changes in functional composition can be a direct mechanism behind diversity changes. Reindeer grazing and canopy cover effects are most likely results of the mechanistic link between size and diversity: tall species are being eaten or suppressed by lack of light, which reduces dominance, which in turn increases diversity (Grime, 1973;Oksanen, 1996). In these forests, both canopy cover and reindeer densities have increased on average. According to our results, both these changes favour shorter plants and result in increased evenness and diversity, thus partially explaining the observed large-scale changes in diversity. These results agree with Depauw et al. (2020), who also found that increased shading is responsible for increased evenness in the herb-layers of temperate European forests.
The effects of reindeer and canopy variables on evenness were not fully mediated by vegetation size (Supporting Information Figure   S2.7), because reindeer still had effects on diversity after adjusting for changes in plant height (Model 4). One explanation may be that our study does not account for intraspecific changes in height, which could happen in response to altered light availability or grazing pressure (Jessen et al., 2020). Thus, the effects of reindeer herbivory and canopy cover on evenness via vegetation height are probably even stronger than detected here.

| CON CLUS IONS
Herbivory and canopy properties affect the diversity and composition of understorey plant communities in boreal herb-rich forests over several decades ( Figure 5). Our results illustrate that land-use F I G U R E 5 A conceptual model of vegetation change drivers in these boreal forests. Red lines correspond to positive effects, and blue lines to negative effects. Positive trends in land-use and canopy variables are marked with red plus signs. While based on our modelling results, the arrows are our subjective interpretation of their implications. LDMC = leaf dry matter content; SLA = specific leaf area and canopy changes impact diversity indirectly by affecting average plant height. By favouring smaller species and limiting the size of plant individuals, increased shade and grazing pressure increase diversity. Therefore, plant size has a mechanistic connection with diversity. In addition, high light availability and reindeer densities favour slower leaf economic traits. Our results thus confirm that functional traits are both indicators and mediators of land-use change effects on plant communities, and underline the importance of including functional metrics in biodiversity monitoring programmes.

ACK N OWLED G M ENTS
The

DATA AVA I L A B I L I T Y S TAT E M E N T
The data and scripts used to produce these results have been deposited in Zenodo (Happonen et al., 2021).