Pollen-based reconstruction reveals the impact of the onset of agriculture on plant functional trait composition.

The onset of agriculture improved the capacity of ecosystems to produce food, but inadvertently altered other vital ecosystem functions. Plant traits play a central role in determining ecosystem properties, therefore we investigated how the onset of agriculture in Europe changed plant trait composition using 78 pollen records. Using a novel Bayesian approach for reconstructing plant trait composition from pollen records, we provide a robust method that can account for trait variability within pollen types. We estimate an overall four-fold decrease in plant size through agriculture and associated decreases in leaf and seed size. We show an increase in niche space towards the resource-acquisitive end of the leaf economic spectrum. Decreases in leaf phosphorus might have been caused by nutrient depletion through grazing and burning. Our results show that agriculture, from its start, has likely been gradually impacting biogeochemical cycles through altered vegetation composition.

INTRODUCTION

Recent studies have demonstrated that humans began to alter the world's ecosystems earlier and to a greater extent than previously thought (Ellis et al., 2021; Mottl et al., 2021; Stephens et al., 2019). Agriculture is considered to be the dominant activity through which humans have induced change in ecosystems. In northwestern Europe, forest clearance, burning and land management for agriculture increased openness, heterogeneity and biodiversity of landscapes from 5500 (Central Swiss Alps) and 6000 years ago (British Isles) (Colombaroli et al., 2013; Woodbridge et al., 2020). Agriculture also alters ecological processes. For instance, agriculture may increase erosion and loss of soil organic carbon with such impacts noted as early as 7500 years ago (Hoffmann et al., 2013; Houben, 2008; Notebaert et al., 2009). Plants are central to many ecosystem processes and their variable forms and functions mean that they influence these processes in different ways, including decomposition (Cornwell et al., 2008), carbon storage (Deyn et al., 2008) and erosion control (Zuazo & Pleguezuelo, 2009). Yet little is known about how agriculture changes plant functional composition on timescales of centuries to millennia.

The onset of agriculture may alter plant functional composition in multiple ways. First, the intentional removal of woody species to create agricultural land may have changed plant functional composition. In pre‐agricultural Europe, closed forest most likely dominated fertile upland areas, although open vegetation, maintained by large herbivores and fire, was common on floodplains and infertile soils (Svenning, 2002). Disturbances by agriculture create opportunities for fast‐growing species which disperse rapidly. Especially In the early stages of agriculture, fallowing would have created a heterogeneous agricultural landscape, with areas in use and in various stages of succession (Mazoyer & Roudart, 2007). Second, agriculture can change plant functional composition through the deliberate introduction of species and associated spread of weeds. Crop plants have been selected over millennia for traits that make them more suitable for human consumption and production, such as higher net photosynthesis rates, higher leaf nitrogen content and higher seed production (Roucou et al., 2018).Weed species typically have high growth rates and produce a high number of seeds (Navas, 2012). Finally, agricultural practices can alter plant community functional composition by changing biotic and abiotic factors (Laliberte & Tylianakis, 2012). Even Europe's earliest farmers implemented intensive land management practices, such as livestock manuring and water management (Bogaard et al., 2013).

Whilst land use appears to have strongly controlled northwestern European vegetation composition since the mid‐Holocene, climate was a main driver of vegetation change throughout the Holocene as well (Marquer et al., 2017). 12,000 years ago, mean annual temperature on the Northern Hemisphere was about 2°C lower than the temperature in the 19th century (Kaufman et al., 2020). The increase in temperature at the beginning of the Holocene led to the expansion of mixed deciduous forest, while boreal forest and tundra and steppe vegetation shifted northwards (Roberts, 1998). Maximum forest cover in Europe was reached between 8000 and 6000 BP (i.e. years before 1950) (Roberts et al., 2018). The replacement of tundra and steppe vegetation by forest consequently changed plant functional composition by increasing plant height and other size‐related traits such as leaf size and seed mass (Moles et al., 2009). Further changes in size‐related traits may have occurred through succession within the woodland primarily driven by light competition (Douma et al., 2012). Furthermore, as temperature and water availability are main constraints on the primary productivity of terrestrial ecosystems, rising temperatures in the early Holocene could also have facilitated more competitive species that typically have large leaves and high leaf nutrient concentration (Wright et al., 2017).

The large number of palynological records produced over the last century are invaluable archives of temporal vegetation dynamics. Pollen records have been used to study local or regional history and compiled to reconstruct continental scale vegetation patterns using landscape reconstruction algorithms (see for instance Roberts et al., 2018). Many pollen records are now available in open repositories, most prominently via the Neotoma database (Williams et al., 2018). This allows meta‐analyses which can consider site‐specific factors like bioclimatic zones and timings of different land‐use activities. Classifying pollen according to their functional characteristics or ecological attributes allows comparisons and generalisations across taxa and biomes. In palaeoecology, the classification of pollen types into plant functional types is commonplace, the simplest division being made between arboreal and non‐arboreal pollen (Birks, 2020). More refined plant functional types have been used in vegetation reconstruction by biomisation (Davis et al., 2015) and climate reconstructions using the modern analogue technique (Mauri et al. 2015). However, because plant functional types are relatively coarse categorisations, the range in trait values can be larger within types than between them (Wright et al., 2005). Furthermore, functional types may not capture key plant characteristics that govern relationships between environmental and ecosystem change (Funk et al., 2017; Thomas et al., 2019). In ecology and biogeography, the use of plant traits has become popular in the past two decades, aided by the upsurge in the availability of data on plant traits (Kattge et al., 2020). Functional traits can provide a better understanding of community assembly as well as ecosystem functioning (Hevia et al., 2017). Here, we specifically use the term ‘functional trait’ for traits that influence an organism's fecundity, growth, development or survival and that have been previously defined as representing key dimensions of variation both at the species and the community level (Bruelheide et al., 2018; Díaz et al., 2016; Volaire et al., 2020).

Recently, the trait‐based approach has been applied in several palaeoecological studies to investigate, among other things, the effects of climate and fire on community trait composition (Blaus et al., 2020; Brussel et al., 2018; Carvalho et al., 2019; Reitalu et al., 2015; van der Sande et al., 2019). However, reliable reconstructions of past trait composition are hindered by the general low taxonomic resolution of the palynological record, which can rarely be resolved to the species level (Birks, 2020). As trait data are usually collected at the species level, ascribing trait values to taxa in the pollen record is not straightforward. Trait variation within species, let alone genera or families, can comprise a substantial part of the total trait variation in communities (Siefert et al., 2015). Because of the low taxonomic resolution of pollen data, previous attempts to reconstruct functional diversity from pollen data using mean trait values to represent a whole taxon might have introduced a large amount of uncertainty (Blaus et al., 2020; Carvalho et al., 2019; van der Sande et al., 2019). Here, we offer a new approach to reconstruct functional composition from pollen records. We use Bayesian modelling to consider the full trait distribution within each pollen taxon and account for the consequences of low pollen taxonomic resolution in community‐level trait reconstruction (Bjorkman et al., 2018; Funk et al., 2017).

By reconstructing plant trait composition in multiple pollen records from a variety of different northwestern European locations, we aim to determine trends in plant trait composition within landscapes undergoing the transition to agriculture. We also assess the relative importance of climate versus agriculture in driving these trends. More specifically, we aim to (1) determine whether functional composition exhibits spatial and/or temporal patterns throughout the Holocene in northwestern Europe, (2) examine whether agriculture triggered a change in functional composition and (3) assess the additional role of climate in driving functional composition. We analyse 78 northwestern European pollen records with agricultural histories, spanning the time period from 10,000 BP to the present. Our results provide a first assessment of the changes that the early establishment of agriculture induced in plant trait composition at macroecological scales, which is fundamental to understanding ecological processes of community assembly under human modification.

METHODS

Selection of study sites

To identify relevant studies, we conducted a structured search in Web of Science using search terms related to agriculture, the time period and the spatial coverage (full list of search terms is presented in Appendix S1 in Supporting Information). Two hundred and ten publications were retrieved after filtering for relevance. To be included, agriculture had to be identified as one of the drivers of vegetation change in the article, the studies also had to cover more than 500 years and include a dated pollen record. Both arable and pastoral agriculture were considered in the search. We only included studies from which we could retrieve raw pollen data, which left us with 78 sites. The Neotoma R package was used to retrieve 74 pollen records from the Neotoma database, in particular its constituent databases the European Pollen Database and the Alpine Pollen Database (Fyfe et al., 2009; Goring et al., 2015; Williams et al., 2018). Pollen records for an additional four sites were supplied by a study author (Davies, 1999). We recorded site characteristics and the start and type of agriculture as defined by the authors of the original studies. An overview of the sites is presented in Appendix S2.

Preparation of the pollen data

To overcome inconsistencies in pollen nomenclature that exist between laboratories, the pollen nomenclature was harmonised based on the nomenclature of Mottl et al. (2021) and the European Pollen Database (Flantua et al., 2021; Giesecke et al., 2019).

Pollen counts do not precisely represent species abundance because of taxon‐specific differences in pollen productivity and dispersal (Dawson et al., 2016; Seppä, 2013). To correct for this, pollen counts were divided by pollen productivity estimates, which are correction factors obtained from comparison between modern vegetation and pollen surface sediment samples (Bunting et al., 2013). Only taxa with pollen productivity estimates (PPEs) were selected for further analyses (Wieczorek & Herzschuh, 2020, tables 3 and 4) (Appendix S3). On average, this resulted in the exclusion of 2.7% of the total pollen count.

To standardise the radiocarbon dating calibration method across cores, we constructed new chronologies for all records using the Bchron R package and the chronological information published in Neotoma (Parnell & Haslett, 2021).

Attributing species to pollen taxa

Trait data are usually collected at the species level, thus to infer trait values at the pollen taxonomic level, we needed to make assumptions about the species belonging to each pollen taxon in our dataset (Figure 1). For this purpose, we made a conversion table for pollen taxon to species from the current distribution of terrestrial species in each pollen taxon in our entire study area. Species distribution data for every pollen taxon in our study area were downloaded from the Global Biodiversity Information Facility database the rgbif R package using the pollen taxonomic name and study area's country names as search criteria (Chamberlain et al., 2021; GBIF, 2019). Possible duplicates and incorrectly geo‐referenced observations were cleared from the species distribution data using the coordinatecleaner R package, which flags occurrence records to databases of common spatial errors in biological collection data, such as country and province centroids or the location of biodiversity research institutions (Zizka et al., 2019). To improve matching of species names to the trait data, the taxonstand R package was used to standardise the species nomenclature according to ‘The Plant List’ (Cayuela et al., 2012; The Plant List, 2013). To filter out (semi)‐aquatic species, we used Ellenberg indicator values for moisture which were downloaded from the TRY database (Ellenberg & Leuschner, 2010; Kattge et al., 2020). Species for which the Ellenberg values were unknown (68%) were retained in the dataset. The resulting conversion table was used for every site in the dataset, and thus included all species found in the geographic range of this study.

FIGURE 1

Diagram representing the structure of analysis. 1) Species were assigned to pollen types using the current geographical distribution of the pollen taxon in the GBIF database. 2) Trait distribution at the taxon level was estimated using gap‐filled data trait data from the TRY database. 3) Pollen percentages were corrected by pollen productivity estimates. 4) Trait distribution at the community level was estimated by weighting the taxon estimates by the corrected pollen percentages 5) Community trait values were modelled using generalised additive models.

Selection of traits

Ten traits were selected for the analysis: plant height, specific leaf area (SLA), leaf area, leaf nitrogen, leaf phosphorus, leaf carbon, leaf dry matter content (LDMC), seed length, seed mass and seed count. These traits represent key axes of trait variation in plants, as well as main plant strategies that we expect to be affected by the onset of agriculture (Díaz et al., 2016; Grime, 1988; Pierce et al., 2017). Plant height was included because the size of plants and their parts reflect a key dimension of variation in plants (Díaz et al., 2016). The trade‐off between SLA and leaf nutrients on the one hand and LDMC on the other is well established and relates to the resource acquisition strategy of species, referred to as the leaf economic spectrum (Wright et al., 2004). Finally, we included seed traits to represent ruderal species which might benefit from increased disturbance by agriculture. These species are characterised by large leaves with high nutrient content, combined with small and numerous seeds (Grime, 1998; Pierce et al., 2017; Westoby et al., 2002).

We used trait data from the TRY database (Kattge et al., 2020), which is currently the largest global database of plant functional traits. This trait data were previously gap‐filled using Hierarchical Bayesian Probabilistic Matrix Factorisation (HBPMF), an approach specifically developed for plant functional trait imputation (Schrodt et al., 2015). All analyses were done using log‐transformed trait values to approach normally distributed data.

Calculating community‐weighted mean values

We used a novel Bayesian modelling approach for reconstructing plant trait composition from pollen records. Using the gapfilled trait data of the species that were assigned to the pollen taxa, we first modelled the trait distributions at the pollen taxonomic level. Traits are inherently correlated because organisms have to balance their allocation to survival, growth and reproduction with the availability of resources (Díaz et al., 2016). To allow correlation between traits, we modelled the trait distribution of pollen taxa using a Multivariate Normal distribution for the likelihood, so that where i represents a trait observation from pollen taxon j from the gap‐filled trait data, and k the number of dimensions of the multivariate distribution, corresponding to ten traits in this analysis. is a vector of mean trait values of length k. is a covariance matrix of k by k. Vague priors were used for the taxon mean and standard deviation estimation, so the posterior is strongly informed by the data. The trait data were Z‐score standardised before modelling.

We then modelled the trait distribution of a community (community weighted mean, CWM) by weighting the taxon trait distributions of each pollen taxon by the corrected pollen percentages at site s and moment in time y, so that where n is the total number of taxa and π j the abundance of pollen taxon at a given site s and moment in time t. We carried through the mean (CWM ) and the standard deviation (CWM ) of the modelled distribution in the subsequent analyses to account for inherent uncertainties in estimations of the mean weighted trait value.

All Bayesian models were run in JAGS (v4.3.0), a program for Bayesian analysis using Markov Chain Monte Carlo simulation (MCMC). The runjags R package provided an interface for running JAGS in R (Denwood, 2016). The performance of the MCMC random walk, that is the chain, to sufficiently represent the posterior was reviewed by checking effective sample size and chain convergence using the coda R package (Plummer, 2004; Plummer et al., 2006). The JAGS code for all models can be found in Appendix S4. We also performed the analysis using a univariate Beta for LDMC and log‐Normal likelihood distribution for the other traits (Appendix S5).

Principal component analysis

To explore variation in the community trait values of all 10 traits, a principal component analysis (PCA) was performed on the mean estimates of CWM (CWM ) using the prcomp function in R. To investigate spatial and temporal differences in multivariate trait composition, we plotted the sample scores of the first and second principal component for six clusters that were based on the proximity of sites to each other. We ran generalised additive models (GAMs) with a smooth for time and random effect of site for the first two principal components and the six clusters using the mgcv package. To test the sensitivity of the analysis to the inclusion of plant height, as a trait that is obviously affected by agriculture, we recalculated CWM without plant height, and performed the PCA on these values (Appendix S6).

Trait composition change over time

To examine trends in individual trait values over time, we used Bayesian GAMs. In this analysis, both the estimates of CWM and CWM were used. The JAGS code for the GAM was generated using the jagam function of the mgcv R package and adjusted to include the estimates of CWM and a random effect of site (α ) (Wood, 2016, 2017). where f(time) is the smoothing function for time before present. To show the dispersion of the reconstructed trait values and give some insight in to species compositional changes, we plotted the reconstructed community trait values as points together with the modelled trends and coloured them by the proportion of trees. JAGS code for the GAMs can be found in Appendix S4. Model fit was evaluated by simulating datasets from the posterior as done in Simpson (2018) and Wood (2016). Robustness of the trends to the choice of sites was tested by running the GAMs, while leaving one site out at the time. Effect of uncertainty in the radiocarbon dating was evaluated by running the GAMs using 50 randomly sample age distributions from the posterior of the Bchron age models. These simulations can be found in Appendix S7.

Effects of agriculture and climate on trait composition

We assessed relationships between the arrival of agriculture and changes in CWM by fitting another GAM with a smoothing function for the time since the arrival of agriculture (f[agriculture]). Where possible, the start of agriculture at each site was retrieved from the associated publication and was thus based on the original author's expert knowledge of the pollen record, other palaeoecological records and archaeological findings, depending on the study. When this information was not clearly defined in the publication, the EUROEVOL and the ArchaeoGLOBE databases were used to distinguish the before and after agriculture time period in the records. The EUROEVOL dataset consists of radiocarbon‐dated archaeological findings from Neolithic Europe and the ArchaeoGLOBE holds consensus data on the extent of agriculture in that region (Manning et al., 2016; Stephens et al., 2019). Eighteen sites were removed from the dataset for this analysis as they already had agriculture present before the start of the pollen record. To account for climatic change through time, we included a smoothing function for temperature in the model (f[temperature]). Mean annual temperature simulations were obtained from CHELSA‐TraCE21k(Karger et al., 2020) (Appendix S8). CHELSA‐TraCE21k is a downscaling algorithm of the TraCE21k palaeoclimate simulations that can generate global temperature estimates at a temporal resolution of 100 years and 30 arcsec spatial resolution for the last 21,000 years. Site was included as a random effect (α ). To show the dispersion of the reconstructed trait values and give some insight in to species compositional changes, we plotted the reconstructed community trait values as points together with the modelled trends and coloured them by the proportion of crop species. R version 4.0.5 was used throughout this study (R Core Team, 2021).

RESULTS

Site characteristics and pollen taxa

The start of agriculture at the sites ranged between 7000 BP and 750 BP. Arable farming was identified in the majority of sites (69/78). Tree cultivation was performed in 15 sites, mainly in the south of the study area. In four of the 78 sites, primarily in upland areas, pastoralism was the sole type of farming identified. Other important drivers of vegetation change as identified in the original publications were climate (40/78), fire (anthropogenic and natural; 34/78) and woodland clearance for iron smelting (14/78).

Forty‐seven pollen taxa were included in the analysis and a pollen‐type‐species conversion table was made for all of them. In total, the trait data for 2357 species were included in the analysis. An overview of the number of species included in a pollen taxon and the number of trait observations can be found in Appendix S3.

Spatial and temporal differences in multiple traits

The first two components of the PCA account for 75.8% of the total variation CWM (Figure 2, Table 1). The first principal component (47.5%) mainly represents the variation in plant size and their organs, with the strongest contribution of leaf area, seed size and plant height. The second principal component (28.3%) expresses variation in traits of the leaf economic spectrum, where communities of high specific leaf area and high leaf nutrients are found on one end of the spectrum and communities of high leaf carbon and LDMC at the opposite end. Seed count is the strongest contributor to the third component (10%) (Appendix S6). The strongest positive correlation was found between size‐related traits, such as seed mass and seed length (r = 0.96, 95% credibility interval; [0.957, 0.961]) and seed mass and leaf area (r = 0.88 [0.87, 0.89]). Correlation was also high between leaf traits, such as leaf phosphorus and leaf area (r = 0.75 [0.74, 0.76]). Leaf carbon and specific leaf area were negatively correlated (r = −0.71 [−0.72, −0.7]) (Appendix S9).

FIGURE 2

Principal component analysis of community level trait values. Every coloured point represents the reconstructed community trait value. The colour indicates whether the observation is from before or after agriculture. Table 1 Presents the loadings of the traits on the principal components.

TABLE 1

Loadings of traits on the principal components

Trait	PC1	PC2	PC3
Leaf area	−0.449	0.017	0.034
Seed mass	−0.419	−0.198	0.037
Seed length	−0.416	−0.17	0.001
Plant height	−0.375	−0.288	−0.103
Leaf phosphorus	−0.365	0.262	0.173
Leaf nitrogen	−0.317	0.272	0.064
Leaf dry matter content	−0.13	−0.463	0.311
Specific leaf area	−0.098	0.52	0.107
Seed count	0.155	−0.047	0.914
Leaf carbon	0.157	−0.468	−0.097
Proportion of variance explained (%)	47.5	28.3	10

The lowest scores on the first principal component are found in the early Holocene (10000–6000 BP). Trends in the first principal component are similar and increasing across the study area, indicating a general decrease in community whole plant size. The trend in the second principal component changes from about 4000 BP in the most Southern sites and 2000 BP in the rest of the study area, meaning communities become higher in SLA and lower in leaf carbon (Figure 3).

FIGURE 3

Trends in the first and second principal component over time for six clusters of proximate sites (see Table 1 and Figure 2). The first principal component comprises the variation in plant, leaf and seed size, while the second principal component expresses variation in traits of the leaf economic spectrum.

Changes in individual traits over time

Community plant height, leaf area, seed size increased in the first half of the Holocene and decreased from about 7500 BP (Figure 4). Community plant height was on average highest at 7700 BP with 4.3 m [3.5, 5.3], and lowest at the present with 0.76 m [0.6, 1.0]. Leaf phosphorus and LDMC decrease significantly over time. LDMC decreases from 0.33 g/g [0.31, 0.34] to 0.30 g/g [0.29, 0.31]. SLA, leaf carbon, leaf nitrogen and seed number remain stable over time.

FIGURE 4

Smooth lines represent the generalised additive model (GAM) of the relationship of CWM of 10 traits with time (calibrated years before present), for all 78 sites. 95% credibility interval is indicated by the shaded area. The reconstructed trait values are plotted and coloured by the corrected proportion of tree pollen at that time point.

Plant height and seed size decrease around the onset of agriculture and onwards (Figure 5), whereas climate had little effect on these traits in the time period investigated in this study (Figure 6). Seed mass changes from 23.8 mg [14.0, 40.4] 5000 years before the arrival of agriculture to 9.5 mg [6.7, 13.4] at the start of agriculture, and to 1.8 mg [1.2, 2.6] 5000 years after the arrival of agriculture. Leaf traits were influenced by climate, with increases in SLA, leaf area and leaf nutrients corresponding with increasing temperature. Leaf area and leaf phosphorus decrease with agriculture, leaf phosphorus from 1.69 mg/g [1.65, 1.73] at the start of agriculture to 1.55 mg/g [1.51, 1.60] at the present. No relationship between the arrival of agriculture and SLA, leaf carbon, leaf nitrogen and LDMC could be detected.

FIGURE 5

Partial component plots of the relationship between the 10 traits and years since arrival of agriculture. Eighteen sites were removed in this analysis, because agriculture was already present before the start of the records. 95% credibility interval is indicated by the shaded area. The reconstructed trait values are plotted and coloured by the corrected proportion of crop pollen (log transformed) at that time point.

FIGURE 6

Partial component plots of the relationship between the 10 traits and temperature. Eighteen sites were removed in this analysis, because agriculture was already present before the start of the records. 95% credibility interval is indicated by the shaded area.

DISCUSSION

By combining functional trait and pollen data, we have identified changes in functional composition of northwestern European sites after the onset of agriculture. These results are especially important for the understanding human impacts on carbon and nitrogen cycles. We find a general decrease in whole plant size alongside a shift in trait space related to the leaf economic spectrum after the onset of agriculture. The changes in the leaf economic spectrum emerge from a gain of communities with a resource‐acquisitive strategy and a loss of communities with a resource‐conservative strategy (Figure 2). Both agriculture and climate play a role in explaining the changes in plant community trait composition, but the contribution of these factors varies with the different traits (Figures 4, 5 and 6). Plant height and seed size appear to be especially influenced by agriculture and little by climate, whereas leaf traits appear to be influenced by both agriculture and climate.

Decrease in whole plant size is triggered by the onset of agriculture

The trend in the first principal component demonstrates a common decrease in the size of whole plants, seeds and leaves across the study area (Table 1, Figures 2 and 3). This decrease is likely solely attributable to the onset of agriculture as we found little effect of climate on seed size or plant height (Figure 6). A decrease in plant height is the most intuitive change due to the opening of forest for agriculture, shown in many previous studies (Rey et al., 2019; Roberts et al., 2018). Here, we find an overall four‐fold reduction in plant height in our study area since the start of agriculture (Figure 5). Plant height is the most important predictor of vegetation carbon storage (Conti & Díaz, 2013), thus quantification of the effect of the onset of agriculture on plant height is central to the understanding of human impact on the carbon cycle. Although decreasing plant height with agriculture is not terribly surprising, the fact that we can identify and quantify this trend demonstrates that pollen‐based functional reconstructions can provide reliable insights, even though some pollen types vary widely in plant height (Appendix S10).

We hypothesised that opening up of ecosystems for new agricultural land by fire or deforestation, provides an advantage to species that are fast growing and rapidly dispersed. These species are characterised by low seed size and high seed count to facilitate effective dispersal (Pierce et al., 2017). Although our results show a reduction in seed size with the arrival of agriculture, seed count remains stable over time (Figure 5). We thus find only limited support for this hypothesis. The decrease in seed mass and length with the arrival of agriculture is most likely also due to removal of forests, as these traits are highly correlated (Díaz et al., 2016).

Communities became more resource‐acquisitive

We show a trend towards the acquisitive end of the leaf economic spectrum, especially in the last 2000 years (Figures 2 and 3). This is likely both driven by agriculture and climate (Figures 5 and 6). With increasing temperature at the beginning of the Holocene, limits on primary productivity were released, giving a competitive advantage to species with traits that support rapid growth and efficient use of available resources (Wright et al., 2017).

We expected that due to the onset of agriculture, communities would become more resource‐acquisitive because of the deliberate manipulation of soil fertility and the characteristics of crop and early successional species (Roucou et al., 2018). The gain in resource‐acquisitive species could have consequences for ecosystem functioning (Lavorel & Grigulis, 2012; Reich, 2014). Acquisitive communities have highly digestible litter and a high availability of nutrients leading to rapid decomposition rates and low litter accumulation (Cornwell et al., 2008). Although some communities after agriculture occupy trait space outside the range found before agriculture, there is a significant overall decrease in leaf phosphorus and leaf area (Figure 5). We think that this seemingly contradictory result is due to differences in the trajectories of agriculture among locations. Crop species as well as other agricultural indicator species (e.g. Rumex, Plantago) are high in leaf phosphorus content (Appendix S10). However, the depletion of nutrients by grazing, burning and removal of trees favoured the establishment of heathland, especially in Scotland and Ireland, which is low in leaf phosphorus and leaf area (Hjelle et al., 2010; Trondman et al., 2015; Webb, 1998).

Uncertainties and future directions

Our Bayesian modelling approach for functional composition reconstruction allowed us to account for the uncertainty introduced by the low taxonomic resolution of the pollen records. The method can furthermore be expanded to other traits and can be used for the quantification of functional diversity, which is likely to be an important contributor to ecosystem functioning alongside the dominance of trait values (Garnier et al., 2015). Palaeoecological studies provide a rich array of proxies on past ecosystem functioning, including erosion control and nutrient retention (Bennike et al., 2021; Jeffers et al., 2015). Combining these proxies and functional composition reconstructions could open a new avenue for research into the relationships between plant composition and ecosystem functioning, and the legacies of human presence. While we were able to incorporate the uncertainty caused by the low taxonomic resolution of pollen data, a limitation of this study is that we have not been able to evaluate the effect of trait dispersion within pollen taxa on the community‐level trait reconstruction. Further work will be needed to examine this effect. Better ways of making conversion tables to assign species to pollen taxa could reduce the amount of uncertainty in the reconstructions. For this more information on past species distributions is necessary, for instance from macrofossil records (Birks & Birks, 2000), ancient DNA from lake sediments (Parducci et al., 2017) and species distribution modelling (Svenning et al., 2011). It is also important to reflect that this landscape scale analysis may have masked finer‐scale trends and that we have not been able to study the effect of differences in agricultural method and intensity between sites. Lack of independent data regarding the variable agricultural trajectories at sites made us unable to study these factors. Furthermore, here we included mean annual temperature as the only climate proxy, while other climatic factors as well as soil differences could also be important in explaining plant trait composition (Joswig et al., 2022). However, the approach we developed for reconstruction of plant functional composition can easily be applied in more nuanced analyses. These kinds of analyses will be significantly aided when more pollen records and other palaeoecological proxies are published in public databases.

The palaeoecological record provides a unique long‐term perspective into the legacies of human impact and can, by applying concepts from ecology and biogeography, provide us with valuable insights on the maintenance of ecosystem functioning in a domesticated world. This analysis is the first empirical demonstration of common changes in the plant functional trait composition of landscapes in the transition to agriculture. Our results show that early agricultural might had significant impacts on biogeochemical cycles through modifying vegetation composition, and we show that these impacts can be demonstrated early on.

AUTHORSHIP

AV designed the study supervised by MJS, SM and FS. AV conducted the analysis with help of FS. AD provided data. AV wrote the manuscript and all authors contributed substantially.

PEER REVIEW

The peer review history for this article is available at https://publons.com/publon/10.1111/ele.14063

Data S1

Click here for additional data file.