Animal-Mediated Ecosystem Process Rates in Forests and Grasslands are Affected by Climatic Conditions and Land-Use Intensity

Decomposition, vegetation regeneration, and biological control are essential ecosystem functions, and animals are involved in the underlying processes, such as dung removal, seed removal, herbivory, and predation. Despite evidence for declines of animal diversity and abundance due to climate change and land-use intensification, we poorly understand how animal-mediated processes respond to these global change drivers. We experimentally measured rates of four ecosystem processes in 134 grassland and 149 forest plots in Germany and tested their response to climatic conditions and land-use intensity, that is, grazing, mowing, and fertilization in grasslands and the proportion of harvested wood, non-natural trees, and deadwood origin in forests. For both climate and land use, we distinguished between short-term effects during the survey period and medium-term effects during the preceding years. Forests had significantly higher process rates than grasslands. In grasslands, the climatic effects on the process rates were similar or stronger than land-use effects, except for predation; land-use intensity negatively affected several process rates. In forests, the land-use effects were more pronounced than the climatic effects on all processes except for predation. The proportion of non-natural trees had the greatest impact on the process rates in forests. The proportion of harvested wood had negative effects, whereas the proportion of anthropogenic deadwood had positive effects on some processes. The effects of climatic conditions and land-use intensity on process rates mirror climatic and habitat effects on animal abundance, activity, and resource quality. Our study demonstrates that land-use changes and interventions affecting climatic conditions will have substantial impacts on animal-mediated ecosystem processes. Electronic supplementary material The online version of this article (10.1007/s10021-020-00530-7) contains supplementary material, which is available to authorized users.


Supplementary material for Ambarlı et al. "Animal-mediated ecosystem process rates in forests and grasslands are affected by climatic conditions and land-use intensity"
Authors: Didem Ambarlı*, Nadja K. Simons, Katja Wehner, Wiebke Kämper, Martin M. Gossner,Thomas Nauss,Felix Neff,Sebastian Seibold,Wolfgang Weisser,Nico Blüthgen *Corresponding author;phone: +49.8161.71.4861 Table of contents S1: Relation of seed removal rate with short-term precipitation S2: Details of herbivory measurements S3: Details of processing explanatory data S4: GLMM results in detail S5: Comparison of the mean effect sizes of short-term vs. medium-term land-use variables in grasslands S6: Process rates in detail S7: Effects of the vegetation structure on process rates S8: Linear correlations between the process rates S9: Within-plot variation of the process rates SUPPLEMENTARY MATERIAL S1 Relation of seed removal rate with short-term precipitation Figure S1 Scatterplot of mean seed removal rate per plot (original data) versus total precipitation in the survey periods in grasslands and forests. Complete removal after 30 mm precipitation implies removal due to wash-away. Therefore data from plots with higher than 30 mm total precipitation were excluded from the analyses.

A. Quantifying herbivory in grasslands
In grasslands, herbivory assessments were based on biomass samples that were taken from a temporarily fenced area where early vertebrate grazing and mowing were excluded. A metal frame of 10 cm × 45 cm × 2 cm was placed randomly at two sampling locations in the fenced area. The vegetation was cut at the edge of the frame above 2 cm. In the lab, the proportion of grasses and stems or leaves of herbaceous plants were estimated visually to account for differences in herbivory between plant functional groups.
Then an overall subset of 100 leaves was randomly picked from the two groups (herbs and grasses) according to their proportion in the biomass. Leaf damage was estimated as damaged surface area of each leaf in mm² by comparing the damaged areas in the leaves with circular and square templates ranging in size from 1 to 500 mm² (see Gossner and others, 2014). For each leaf, the total damaged area from different damage types (sucking, scraping, mining, chewing etc.) were summed per leave to get one value for damaged area per leaf. We reported regularly shaped, whitish points, which were often clearly visible against the light, as sucking damage caused by sap-sucking insects. Though this method may fail to include all sucking damage especially by small auchenorrhynchans, it captures the majority and is therefore used commonly in insect herbivory studies (Loranger and others, 2014). Leaf damage at the edges was estimated by reconstructing the removed area based on leaf shape. The remaining leaf area was measured with a leaf-area-meter (LI-COR area meter (LI-3100C, Lincoln (NE) USA)). Herbivory rates were calculated as the total damaged area divided by the total leaf area, i.e. damaged and remaining. We accounted for the damage types also measured by leaf-area-meter so included in the remaining area, i.e. sucking, scraping, mining and galling, by subtracting the total area of those damages from the measured area to find the remaining area.

B. Quantifying herbivory
In forests, herbivory rate was estimated on the 10 (in ALB and HAI regions) and 8 (in SCH region) most abundant plant species of all forest layers in terms of percent cover. This also included geophyte species sampled in spring, which were excluded from the subsequent analyses to enable use of the climatic data in the short-term timescale. We selected those species based on plant survey data from former years (Boch and others, 2013, D.Prati, pers. comm). The leaf material was collected along the outer border of the plots in summer (from mid-June (ALB) to mid-August (SCH)). Leaves from the tree layer were collected using elongated lopping shears (max. height: 7 m) to cut branches from the shaded canopy of at least three trees.
Where shears were not long enough, branches were collected using shooting and/or climbing. When individuals of the tree species were present in the shrub layer, leaves were collected from that layer as well in an amount approximately representing the distribution of tree and shrub layer within the plot.
When a species could not be found on the plot margin or when sampling was not possible out of other reasons, this species was replaced by the next abundant species on the species list. All plant material was collected in plastic bags containing a moist cloth to prevent the leaves from drying out and transferred to a fridge until further processing. In the lab, leaf material was further processed. A total of 12 to 200 randomly picked leaves per species per plot were processed (depending on leaf size and composition, TableS1). For each single leaf, the area was measured using a LI-COR area meter (LI-3000C, Lincoln (NE) USA). For coniferous trees (Abies, Picea, Pinus sp.), sets of 10 leaves were measured together and average values were taken for a single leaf area for better accuracy. For grass species with needle-shaped leaves (Deschampsia flexuosa), no reliable area measurements were possible using the leaf area meter.
Thus, leaf length was measured using a ruler and the approximate area was determined by multiplying the length by an estimated leaf width of 1 mm. As for grassland herbivory, the total damaged area from different damage types (sucking, scraping, mining, chewing, etc.) was determined for each leaf using circular and square templates ranging in size from 1 to 500 mm² (see Gossner and others, 2014). Damage types were summed and herbivory rate per leaf was calculated as the proportion of herbivory-affected leaf area compared to the corrected total area of the leaf (i.e., the sum of measured leaf area and chewing damage). From those leaf-level herbivory rates, we calculated mean herbivory rates per plant species and plot. To determine community-level herbivory rates per plot, we used plant cover values obtained from plot core areas by vegetation surveys, the same data as used for the selection of the most abundant plant species. Plant-species herbivory rates were weighted by their relative contribution to the total cover of assessed plants to calculate community-level herbivory rates per plot.  When data were missing for specific dates of a plot, we substituted it with the mean value of the specific habitat of the region for that date. Other parameters that might affect the processes rates, specifically mean wind velocity and sunshine duration, were not used because of the high number of missing data in the dataset. We did not use the number of frost days in summer periods, as in each summer, there was only one plot with frost days. Ordination of the first two axis of PCA on medium-term climatic variables used in the dung removal, seed removal and predation analyses. PCA was conducted with "factoextra" package (Kassambara and Mundt, 2017) and plotted with ggbiplot package in R. Tmax: Maximum temperature, CoolD: number of cool days in that period, Prec: Precipitation. S2016 and S2017 indicate variables for spring and summer seasons of 2016 and 2017, respectively. September_Tmin and Wnt_Tmin indicate minimum temperatures in September 2016 and across winter between 2016 and 2017, respectively. Letters "c" and "m" at the end of each variable indicate cumulative or mean of that variable for that period, respectively. Plots of each region are colored and indicated in the legend as A: Schwäbische Alb in red, H: Hainich-Dün in green, and S: Schorfheide-Chorin in blue.

Figure S4
Ordination of the first two axis of PCA on medium-term climatic variables used in the herbivory analyses. It was conducted with "factoextra" package (Kassambara and Mundt, 2017) and plotted with ggbiplot package (Vu, 2011) in R. Tmax: Maximum temperature, CoolD: number of cool days in that period, Prec: Precipitation. S2016 indicates variables for the spring and summer of 2016. September_Tmin and Wnt_Tmin indicate minimum temperatures in September 2016 and across winter between 2016 and 2017, respectively. Letters "c" and "m" at the end of each variable indicate cumulative or mean of that variable for that period, respectively. Plots of each region are colored and indicated in the legend as A: Schwäbische Alb in red, H: Hainich-Dün in green, and S: Schorfheide-Chorin in blue.    E. Histograms of land-use variables

Figure S11
Histograms of the land-use data used in the dung removal, seed removal and predation analyses. The x-axes of the first row include variables at short-term grazing (livestock unit*days/ha), mowing (number of cuts) and fertilization (kg nitrogen/ha), respectively. The second row has the same variables for the medium-term timescale, standardized relative to its mean within that year, and then calculated the average of two preceding years (2015 and 2016).

Figure S12
Correlations among grassland explanatory variables used for the dung removal, seed removal and predation analyses. Numbers indicate correlation coefficients based on Pearson's test. Blue color indicates a positive correlation whereas red color indicates a negative correlation. Positive correlations which are smaller than 0.25 appear white. ST: short-term variables corresponding to two days of survey period. MT: medium-term. MT variables are represented with the first two axes of the PCA. Short-term precipitation is the sum of daily precipitation of 2 days. Tmean: mean of daily temperature, Tmax: mean of daily maximum temperature.

Figure S13
Correlations among forest explanatory variables used for dung removal, seed removal and predation analyses. Numbers indicate correlation coefficients based on Pearson test. Blue color indicates a positive correlation whereas red color indicates a negative correlation. Positive correlations smaller than 0.35 appear white. ST: short-term variables corresponding to two days of survey period. MT: mediumterm. MT variables are represented with the first two axes of the PCA. Short-term precipitation is the sum of daily precipitation of 2 days. Tmean: mean temperature, Tmax: maximum temperature.

Figure S14
Correlations among grassland explanatory variables considered for the herbivory analyses. Numbers indicate correlation coefficients based on Pearson test. Blue color indicates a positive correlation whereas red color indicates a negative correlation. Correlations that are very close to zero appear white. ST: short-term variables corresponding to two days of survey period. MT: medium-term. MT variables are represented with the first two axes of the PCA. Precipitation is the sum of daily precipitation of the short-term period. Tmean: mean temperature, Tmax: maximum temperature. cum: cumulative.

Figure S15
Correlations among forest explanatory variables considered for the herbivory analyses. Numbers indicate correlation coefficients based on Pearson test. Blue color indicates a positive correlation whereas red color indicates a negative correlation. ST: short-term variables corresponding to two days of survey period. MT: medium-term. MT variables are represented with the first two axes of the PCA. Precipitation is the sum of daily precipitation of the short-term period. Tmean: mean temperature, Tmax: maximum temperature, cum: sum of daily values over the period.

SUPPLEMENTARY MATERIAL S4 GLMM results in detail
Significance codes indicated in bold are as follows: ("." : p≤0.1, *: p≤0.05,**: p≤0.01,***: p≤ 0.001.   Figure S16 Effects of the interaction between medium-term climatic variables on seed removal in grasslands. Both variables were standardized to mean = 0 and standard deviation = 1.  Figure S17 Effects of interaction between the short-term climatic variables on seed removal in forests. Both variables were standardized to mean = 0 and standard deviation = 1.       Figure S18 Effects of the interaction between the short-term climatic variables on predation by arthropods in forests. Both variables were standardized to mean = 0 and standard deviation = 1.

D. Herbivory
Models were conducted at the plot level. Due to collinearity, short-term precipitation and medium-term fertilization were not included in the models. ST: short-term, MT: medium-term.

SUPPLEMENTARY MATERIAL S5 Comparison of the mean effect sizes of short-term vs. medium-term land-use variables in grasslands
This section aims to compare the effect sizes of land-use variables in the short-term vs. medium-term timescales. As those variables were correlated within and among each variable set, we did not analyze them together but we obtained the mean effect sizes from separate models using either short-term or medium-term variables together with the climatic conditions, following the main model structure indicated in the Methods. As we focused on the short-term effects in the main findings, we used the results presented in the Results and the Supplementary material S6 for the short-term effects of dung removal, seed removal, predation by birds or rodents. In case that those results included the medium-term effects such as the case for some predation and herbivory tests, we also run a model only with the shortterm effects and then only with the medium-term effects of the corresponding variables. The mediumterm variables included grazing and fertilization intensities but not of mowing due to collinearity. The short-term variables included all three land-use components except for herbivory as the subplots that herbivory measurements took place were not mown. We then calculated the mean effect size of land-use variable sets, the absolute of each effect size weighted inversely by its standard error. The results show that the mean effect sizes of short-term land-use variables were higher than that of the medium-term variables for all processes except for the seed removal and the predation by the arthropods (Fig. S19).

Figure S19
Mean effect sizes of short-term (blue bars) vs. medium-term (red bars) land-use components in grasslands. The effect sizes were obtained from the GLMM models and weighted inversely by their standard errors.                Figure S20 Boxplots showing plot-level mean rates of major processes in forests and grasslands. G: Grassland, F: Forest.

Figure S21
Boxplots of plot-level predation rates by each predator group in two habitats. Predation by each group was calculated as the proportion of caterpillars with at least one predation mark by each group. G: Grassland, F: Forest.

Figure S22
Boxplots of plot-level herbivory rates on two plant groups in grasslands.

SUPPLEMENTARY MATERIAL S7: Effects of vegetation structure on process rates
We recorded vegetation structure parameters at the subplot level: vegetation height was measured for the majority (90%) of the ground cover in cm. Vegetation cover of ground layer and canopy density were recorded as percentages. Tables S41 & 42 show the results of GLMM analyses for the effect of vegetation on the process rates at the subplot level. We used all vegetation parameters and region (representing environmental differences between the regions) as fixed effects. We included plot as a random effect (1 | Plot). All numerical explanatory variables were standardized. Significance (Sig.) codes: ("." : p≤0.1, *: p≤0.05,**: p≤0.01,***: p≤ 0.001. No vegetation structure data are available from herbivory surveys.

Table S46
Correlations between within-plot level variations of process rates (CV within in %) and land-use intensities. Correlation coefficients of Pearson correlation tests were given and significant results are indicated in bold with p values in parentheses. For land-use intensity in forests, a compound measure, ForMI (Kahl and Bauhus, 2014), was used which combines harvest intensity, proportion of non-natural trees in the tree composition, and anthropogenic deadwood. For grasslands, a compound measure, LUI (Blüthgen and others, 2012), combining the intensity of livestock grazing, the number of cuts and fertilization was adopted and the dataset was produced based on the short-and medium-terms of the study. Herbivory in forests was measured at the plot level so no CV within was available.