Evaluating alternative hypotheses behind biodiversity and multifunctionality relationships in the forests of Northeastern China

Qingmin Yue, Minhui Ho, Yn Geng, Xuerui Wng, Klus von Gdow, Chunyu Zhng,Xiuhi Zho, Lushung Go,＊

a Research Center of Forest Management Engineering of State Forestry and Grassland Administration, Beijing Forestry University, Beijing, 100083, China

b Faculty of Forestry and Forest Ecology, Georg-August-Universitat, G¨ottingen, D-37077, Germany

c Faculty of AgriSciences, Stellenbosch University, Matieland, 7600, South Africa

Keywords:Biodiversity and ecosystem multifunctionality Biomass Functional traits Mass ratio effect Niche complementarity effect Vegetation quantity effect

ABSTRACT

1. Background

Forested landscapes are important in providing a variety of ecosystem functions and services, including carbon sequestration, timber production, and the provision of habitats for plants and animals (Pan et al.,2011;Mori et al.,2017;Gadow et al.,2021).A fundamental goal in forest ecological research is thus to deepen our understanding of how biodiversity affects ecosystem functioning and processes. This is particularly relevant in the light of climate change and biodiversity loss (Sax and Gaines, 2003; Bennett et al., 2015). Recent studies have shown that species-rich communities are more effective in providing certain ecosystem services, such as biomass productivity and ecological resilience than communities with few species (Ruiz-Benito et al., 2014;Chiang et al., 2016; Ratcliffe et al., 2016; Adair et al., 2018; de Avila et al., 2018). The biodiversity effect should be more important when multiple ecosystem functions are considered simultaneously (i.e.,ecosystem multifunctionality). Because species differ in their contributions to different ecosystem functions, high biodiversity may promote ecosystem multifunctionality through a “portfolio” or “synergy” effect(Gamfeldt et al.,2008;Lefcheck et al.,2015;Gamfeldt and Roger,2017).However, such evidence is mostly based on small-scale studies (Meyer et al., 2018; Schuldt et al., 2018; Yuan et al., 2020). Effects on the biodiversity and ecosystem multifunctionality (BEMF) relationships at greater scales, including large landscapes or even entire countries, are largely unknown.Since many policy decisions on land management tend to focus on the landscape, province, and country scales (Ratcliffe et al.,2017;Adair et al.,2018;van der Plas et al.,2018;Ouyang et al.,2019),a better understanding of the mechanisms behind the BEMF relationships at the macroscale is crucial for conserving and enhancing multifunctional forest ecosystems.

Forindividualecosystemfunctions,twohypotheseshavebeenproposed toexplainbiodiversityeffects.Thenichecomplementarityhypothesisstates that niche differences among species, such as interspecific differences in resource use,should lead to more efficient acquisition of limiting resources and therefore promoting ecosystem functioning(Tilman et al.,1997).The mass ratio hypothesis postulates that the effects of functional traits of individual species are related to their relative abundance and that the most dominant values of plant functional traits will be the main determining factors of ecosystem functioning (Grime, 1998). Measures of multifunctionality are usually based on multiple individual functions(Gamfeldt et al.,2008;Byrnes et al.,2014).It is therefore reasonable to speculate that these effects influencing individual functions should also underlie relationships between BEMF. Both niche complementarity and mass ratio effects can be quantified using functional traits based on two complementaryapproaches,functionaldiversityandfunctionalidentity(Lohbecketal.,2015;Chiang et al.,2016;Hao et al.,2020).Functional diversity represents the ranges,values,and distributions of traits withinmultidimensionalniche spaces(Lalibert'e and Legendre,2010),and thus provides a proxy for niche complementarity.Functional identity(also known as functional composition)is usually expressed in terms of community weighted mean traits,that is, the average of species traits weighted by their relative abundance or biomass(Violle et al.,2007).Accordingly,functional identity is often used to evaluate mass ratio effect.

Recent studies have demonstrated that forest ecosystem functioning is affected not only by vegetation quality (i.e., functional diversity and functional identity) but also by vegetation quantity (i.e., stand density,represented by stand biomass or basal area per unit area;Finegan et al.,2015;Lohbeck et al.,2015).The vegetation quantity hypothesis proposes that higher stand biomass means greater canopy volume, light absorption,and water interception,thus fostering forest ecosystem functioning(Lohbeck et al., 2015; Yuan et al., 2018; Hao et al., 2020; Yue et al.,2020). Specifically, these positive stand biomass effects on ecosystem functioning are more evident in forests at the early development stage following disturbance,where the sharp increase in biomass may override more subtle changes in biodiversity(Lohbeck et al.,2015).However,as successional development proceeds, greater stand biomass is associated with more severe competition, which may suppress plant growth and forest productivity (Finegan et al., 2015). Therefore, the vegetation quantity hypothesis may also be valid in explaining the BEMF relationships. However, this has never been confirmed in forest ecosystems,although it has potential implications for practical management.

The forests of northeastern China represent 37%of the country's total forest land area(Wu et al.,2019)and account for more than 27.5%of the total carbon stocks of forests in China(Peng et al.,2009;Dai et al.,2018;Huang, 2019; Luo et al., 2020). These forests also constitute one of the most biodiverse temperate forest regions in the world due to the physiographical heterogeneity and well-preserved from a lack of Quaternary glaciation (Axelrod et al., 1996; Qian and Ricklefs, 2000; Qian et al.,2005).However,during the past century,excessive logging and absence of professional management resulted in substantial decreases in the natural forest area and quality.More than 70%of the natural forests were converted into secondary forests(Hao et al.,2000;Yu et al.,2011;Bryan et al., 2018). Towards the end of the 20th century, these forests were protected for ecological restoration,especially since the launch of China's Natural Forest Protection Project in 1998 (Zhang and Liang, 2014).Nevertheless, secondary forests are still the major vegetation type in northeastern China. Forest restoration and reforestation of degraded forests have enjoyed a high priority of the local government.However,a major practical challenge is to preserve forest biodiversity while at the same time enhancing ecosystem multifunctionality. Therefore, understanding the mechanisms that regulate BEMF relationships is of increasing practical relevance.

Accordingly,the objective of this study is to examine how ecosystem multifunctionality is affected by niche complementarity,mass ratio,and vegetation quantity effects. A structural equation modeling approach is used to explore the direction and strength of each effect (Fig. 1). We expect that forests with higher functional diversity exhibit greater ecosystem multifunctionality due to interspecific facilitation and resource partitioning (H1, niche complementarity hypothesis). We also expect that forests dominated by fast-growing acquisitive traits display greater ecosystem multifunctionality,because acquisitive traits will lead to higher photosynthetic efficiency and carbon sequestration rate (H2,mass-ratio hypothesis). Furthermore, we hypothesize that higher stand density (biomass per unit area) would promote multifunctionality.Because in secondary forest, higher density is usually associated with larger canopy and root volume, and thus with greater resource acquisition potential(H3,vegetation quantity hypothesis).All three effects may be crucial for maintaining ecosystem multifunctionality. Their relative contribution may covary with forest development stage and environmental condition. Therefore, we considered environmental conditions and stand age as covariates in the analyses.

2. Methods

2.1. Study area and sampling design

This study was conducted in a large forested landscape throughout northeastern China,stretching across the provinces of Heilongjiang,Inner Mongolia, Jilin, and Liaoning. The total area spanned approximately 1,000 km in an east-west direction(119.80°–134.01°E),and 1,500 km in a north-south direction (39.71°–53.36°N), accounting for 37% of the country's total forest land area (Wu et al., 2019). The regional climate varies from a moderate temperate zone to a cool temperate zone, with a mean annual temperature range from-5.8 to 9.9°C.The annual precipitation ranges from 430 to 1,134 mm, with a predominant rainy season from June to September.The study area includes all the major forest types in temperate northeastern Asia. The dominant broad-leaved species are

Quercus mongolica,Betula platyphylla,Tilia amurensis,Betula davurica,Ulmus japonica, and Acer mono; the main coniferous species are Pinus koraiensis,Abies nephrolepis,and Larix gmelinii(see full species list in Table S1).

Fig. 1. Conceptual model linking stand biomass, age, functional diversity,functional identity, environmental factors, and ecosystem multifunctionality.The solid arrows represent the hypothesized causal relationships. Specifically,the three alternative hypothetical pathways represent the three key hypotheses tested in this study: niche complementarity (H1), mass ratio (H2), and vegetation quantity (H3). Biomass: stand biomass; Age: stand age; FD: functional diversity; FI: functional identity; MF: ecosystem multifunctionality.

Fig. 2. Locations of study plots and elevation ranges.

A systematic grid of 412 circular plots, each with an area of 0.1 ha,was established in 2017 and 2018 (Fig. 2). These plots are distributed along and near eight mountain ranges, including Changbai, Greater Khingan, Hadaling, Laoyeling, Lesser Khingan, Longgang, Wanda, and Zhangguangcailing, with elevations ranging between 79 and 1,255 m.The distances between neighboring plots vary between 20 and 60 km.The geographical data (longitude, latitude, and altitude) were recorded for each plot (Zhang et al., 2020). All live trees greater than 5-cm diameter outside bark at breast height (1.30 m, DBH) were measured and mapped.A total of 32,307 individuals of 67 species were recorded.The 15%tallest trees in each plot were selected for sampling stem cores.Stem cores were extracted using an increment borer (5 mm, Suunto,Finland) from the north side of trees at a height of 1.3 m (Wu et al.,2019). The ring counts of these samples were averaged to estimate the stand age, which reflects the development stage of a secondary forest after disturbance(Ouyang et al.,2019).

2.2. Measurement of functional diversity and identity

Seven traits that have been suggested as having a great effect on individual tree growth and survival, and thus on ecosystem functioning,were measured according to standard protocols (Perez-Harguindeguy et al.,2013).Included are six leaf traits(leaf area,specific leaf area,leaf dry matter content, leaf carbon content, leaf nitrogen content, and leaf carbon-to-nitrogen ratio), and one stem trait (wood density). More detailed descriptions of the data collection processes for the different traits are provided in the supplementary materials.

Rao's quadratic entropy (RaoQ) was used to quantify the functional diversity. RaoQ describes the variation in pairwise trait dissimilarities among all individuals within a community(niche complementarity;Rao,1982).Community weighted mean traits(CWMs)were used to quantify the functional identity (also known as the functional composition).CWMs were calculated as the average of trait values weighted by the biomass of each species in a community(mass ratio;Violle et al.,2007).Since plant traits were correlated,a principal component analysis(PCA)was performed on the CWMs using the method presented by Ruiz-Benito et al.(2017),with the first two axes representing the functional identity.The first axis (CWM1), which explained 55% of the total variation, was positively related to the CWMs of the leaf area,specific leaf area,and leaf nitrogen content,but negatively correlated with leaf dry matter content,leaf carbon content, and leaf carbon-to-nitrogen ratio (Fig. S1). The second axis (CWM2), which explained 19% of the total variation, was positively related to CWMs of the leaf area, but negatively related to wood density (Fig. S1). Therefore, CWM1 and CWM2 represent the gradients of functional identities from species with slow growth and conservative resource use strategies (negative values) to species with rapid growth and acquisitive resource use strategies(positive values).

2.3. Environmental factors

To consider the effects of environmental factors on the BEMF relationships,key climate and soil properties were used as predictors.Four climate factors, including mean annual temperature, mean annual precipitation, temperature seasonality (i.e., the coefficient of temperature variation within 12 months), and precipitation seasonality (i.e., the coefficient of precipitation variation within 12 months), were extracted from the WorldClim site (www.worldclim.org) with a resolution of 1 km2. Two soil factors,including soil depth and soil pH,were measured.Within each circular plot,one soil profile was excavated in the center of the plots and the soil depths were measured.Five samples of the topsoil(depths between 0 and 20 cm) were taken. The samples were evenly spaced along a straight line extending from one plot perimeter to the opposite perimeter.In a standard 0.1 ha circular plot with radius 17.85 m,the distance between soil samples is thus 8.93 m.A cylindrical metal corer was used to extract the soil samples; soil pH was measured following the procedures recommended by the Soil Science Society of China(Bao,1999).The mean pH value of the five samples from each plot was used to represent the soil pH.

To eliminate possible correlations while retaining effective environmental information, another PCA was conducted for climate and soil properties. Again, the first two axes were retained for further examination. The first axis (ENV1), which explained 47% of the total variation in the dataset, was positively related to mean annual temperature and mean annual precipitation, but negatively related to temperature and precipitation seasonality (Fig. S2). The second axis(ENV2) explained 16% of the total variation and was positively related to soil depth, but negatively related to soil pH. Thus, ENV1 and ENV2 represented the gradients of the environment from harsher conditions with limited resources (negative values) to fertile conditions with sufficient resources (positive values).

2.4. Stand biomass

Total tree biomass(including stem,branch,root,and foliage biomass)was calculated based on the available regional-specific allometric equations.The basic model was developed by Wang(2006),Dong et al.(2014,2015),and Luo et al.(2020),and includes 80%of all the species found in the study area:

where ^TB represents individual tree biomass;B(cm)is DBH;H(m)is tree height; a, b, and c are estimated coefficients; the subscripts s, t, r, and f refer to stem,branch biomass,root biomass,and foliage.For the species with unknown coefficients (20% of the species), a generic allometric equation was chosen. The stand biomass of each plot was calculated as the sum of the individual tree biomass values and standardized to a 1-ha value.

2.5. Ecosystem multifunctionality

Seven ecosystem functions were used in the analysis: (1) forest productivity; (2) shrub density; (3) herb coverage; (4) herb richness; (5)litter biomass; (6) litter thickness; and (7) soil carbon stock. The supplementary materials provide detailed descriptions of the data collection process for different functions.

During the analysis, based on the individual functions, two complementary approaches were adopted to calculate ecosystem multifunctionality for each plot(following Gamfeldt et al.,2008;Byrnes et al.,2014). The first method involves an averaging approach, in which ecosystem multifunctionality was quantified as the average value of the seven standardized functions (hereafter referred to as MFave). The second is a threshold approach,in which ecosystem multifunctionality was quantified as the number of functions that exceeded 10%, 20%, 30%,40%,50%,60%,70%,80%,and 90%of the maximum of each function(hereafter MF10-90).The maximum value was calculated as 97.5%of the observed values for each function.

2.6. Structural equation models

Structural equation models were used to test the effects of biodiversity on ecosystem multifunctionality.A structural equation model(SEM)is a probabilistic model containing or specifying multiple causal pathways. The model has several distinctive characteristics. First, SEM attempts to satisfy the criteria for drawing causal inferences and the results can be used to verify certain hypotheses. Second, SEM permits specific endogenous variables to be functions of other endogenous variables,thus containing possible indirect effects (Lamb et al., 2014). Therefore, SEM has the potential to disentangle the complex direct and indirect causal relationships between biodiversity and ecosystem functions, in which a series of dependent and independent variables may be correlated (van der Sande et al.,2017;Hao et al., 2018). Our conceptual model(Fig.1)was used to evaluate the direct effects of the environmental factors and stand age on ecosystem multifunctionality,as well as the indirect effects via interactive functional diversity (RaoQ), functional identity (CWMs)and forest density(biomass).

The explanatory and response variables have different dimensions and cannot be compared directly, prior to the SEM analyses. All the variables were therefore standardized to mean zero and unit standard deviation to eliminate the dimension difference and to make the ranges of all variables comparable on a similar scale.For each of the ecosystem multifunctionality indices (i.e., MFave and MF10-90; Table 1), separate SEMs were fitted via fixed model structures using the R package“l(fā)avaan”(Rosseel, 2012). In addition, to further understand the BEMF relationships, we fitted separate SEMs for seven individual functions. The goodness of fit for the SEMs was evaluated using Bentler's comparative fit index(CFI)and the standardized root mean square residual(SRMR).The indirect effects were calculated by multiplying the path coefficients for the effects of certain predictors on the median predictor with the path coefficients for the effects of the median predictor on ecosystem multifunctionality. The total effects of each predictor on ecosystem multifunctionality were determined as the sum of the direct and indirect effects. All the analyses were implemented in R (Version 4.0.3; R Development Core Team,2019).

3. Results

The SEMs were used to investigate the BEMF relationships following our conceptual model.The conceptual model was found to provide good fits to our data(on average CFI=1.000;SRMR=0.001).In terms of the averaging approach, the SEM accounted for 10% of the standardized variation in the averaged ecosystem multifunctionality (i.e., MFave;Fig.3).Our results indicate that stand biomass had a direct positive effect(r = 0.21; p ＜0.05) on MFave. Stand age had no direct effect but an indirect positive effect via stand biomass(r=0.09;p ＜0.05)and CWM1(r=0.04;p ＜0.05)on MFave(Fig.3 and Table 2).For the biodiversity properties,CWM1(positively related to leaf area,specific leaf area,and leaf nitrogen content; but negatively related to leaf dry matter content,leaf carbon content, and leaf carbon-to-nitrogen ratio) had a direct negative effect on MFave(r=-0.18;p ＜0.05),while CWM2(positively related to leaf area, but negatively related to wood density) and RaoQ had no significant relationships with MFave (p ＞0.05, Fig. 3 and Table 2). ENV1 (positively related to mean annual temperature and precipitation, but negatively to temperature and precipitation seasonality)had a direct positive effect(r=0.13;p ＜0.05)on MFave and mixed indirect effects through biomass, CWM1, and RaoQ. While ENV2 (positively related to soil depth but negatively to soil pH) had no significant effect on MFave(Fig.3 and Table 2).

In terms of the threshold approach, our models accounted for 9.0%(1%–25%) of the standardized variation in the threshold-based ecosystem multifunctionality (i.e., MF10-90). The SEMs showed that stand biomass had consistently direct positive effects on ecosystem multifunctionality except for MF80 and MF90, while stand age had indirect positive effects caused by greater stand biomass(Table S2).Stand biomass and age jointly explained up to 36.4% (18.6%–49.9%) of the accounted variance (Fig. 4). For the biodiversity properties, CWM1 had direct negative effects, particularly at higher ecosystem multifunctionality thresholds (MF30, MF60-90); while CWM2 had direct negative effects at lower ecosystem multifunctionality thresholds(Table S2; MF10-20). RaoQ had no significant effects on ecosystem multifunctionality, with the exception that a direct positive effect was observed when MF10 was considered(Table S2).Biodiversity properties jointly explained on average 40.0% (23.2%–50.7%) of the accounted variance (Fig. 4). For the environmental factors, ENV1 had significant positive total effects in the lower ecosystem multifunctionality thresholds(MF10-40). ENV2 had a positive effect at lower ecosystem multifunctionality threshold(MF10)and negative effects at higher ecosystem multifunctionality thresholds (MF80-90; Table S2). Environmental factors explained on average 23.6% (8.9%–43.2%) of the accounted variance.

To better explain the relationships between BEMF, a series of SEMs that used the same model structure are shown in Fig. 3, but each time with an individual ecosystem function being examined.It was found that the relative importance of different factors may change in various ecosystem functions,as detailed in Fig.5 and Table S3.In summary,the ecosystem functions related to stand productivity, herb coverage, and herb richness were dominated by stand biomass and age.The ecosystem functions related to litter biomass and thickness were mainly attributed to biodiversity properties(Fig.5 and Table S3).The ecosystem functions related to shrub density and soil carbon stock were driven by environmental factors(Fig.5 and Table S3).

Table 1 Description of the explanatory and response variables.

Fig. 3. Structural equation model of stand biomass,stand age, biodiversity, and environmental factors as predictors of averaged ecosystem multifunctionality(MFave). Solid arrows represent significant paths (p＜0.05), and dashed arrows represent nonsignificant paths(p ≥0.05).The path coefficients are reported as standardized effect sizes. The values on the arrows represent the standardized path coefficients. The values of R2 represent the proportions of the response variations explained by the observed variables.CWM1 and CWM2:the first two principal component analysis axes of the community weighted mean traits; RaoQ:Rao's quadratic entropy;ENV1 and ENV2:the first two principal component analysis axes of the environmental factors.

4. Discussion

A clear understanding of how multiple ecological mechanisms interact to shape BEMF relationships is essential for managing ecosystem functions at the macroscale while facing possible changes in global landuse and climate conditions. Therefore, this study used SEMs to test therelative importance of three alternative hypotheses (niche complementarity, mass ratio, and vegetation quantity effects) that are known to affect the BEMF relationships.Our results provide evidence in support of the mass ratio hypothesis and the vegetation quantity hypothesis.

Table 2 Direct, indirect, and total standardized effects of predictors on the averaged ecosystem multifunctionality based on the structural equation model.

4.1. Functional identity is more important than functional diversity

Functional diversity and identity have been widely used in ecological research to describe two complementary aspects of biodiversity(Tobner et al., 2016; Adair et al., 2018). Functional diversity measures the distribution and dispersion of traits in the multidimensional niche space,while functional identity describes the dominance of traits in the niche space. Therefore, functional diversity and identity are key variables related to the niche complementarity and mass ratio hypotheses,respectively (Ratcliffe et al., 2016; Hao et al., 2020). If the niche complementarity hypothesis is valid,functional diversity should increase ecosystem multifunctionality. This is because in functionally diverse forest ecosystems,greater trait dispersion allows species to stably coexist through niche partitioning and efficient utilization of resources(Tilman et al., 1997). However, in this study, we found that functional diversity generally had no significant effect on ecosystem multifunctionality,which implies that the niche complementarity effect may not be the main driver of ecosystem multifunctionality in the forests in northeastern China. These findings differ from earlier studies that found strong support for the niche complementarity hypothesis for productivity and other individual ecosystem functions(Hao et al.,2018,2020).Such divergent results would be due to the trade-off among ecosystem functions, as previous studies have demonstrated that a species’ capacity to support some functions at high levels will compromise its ability to support others(Felipe-Lucia et al., 2018; Le Bagousse-Pinguet et al., 2019; Yuan et al.,2021).On the other hand,if functional identity is found to play a crucial role in ecosystem multifunctionality, it may override the effect of functional diversity (Le Bagousse-Pinguet et al., 2019). In this analysis, two PCA axes (CWM1 and CWM2) of seven CWM traits were extracted to represent the functional identity.We found that CWM1,which explained 55% of the total variation in functional identity, had consistently negative effects on ecosystem multifunctionality. CWM1 was positively related to the functional identity of the leaf area, specific leaf area, and leaf nitrogen content,but negatively correlated with the leaf dry matter content,leaf carbon content,and leaf carbon-to-nitrogen ratio(Fig.S1).Thus, CWM1 reflects the gradient of forest composition from climax species with slow-growing conservative traits to pioneer species with fast-growing acquisitive traits(van der Sande et al.,2017,2018;de Avila et al., 2018; Fotis et al.,2018). Compared with previous results that acquisitive traits increase productivity and other ecosystem functions in temperate or tropical forests(Lohbeck et al.,2015;de Avila et al.,2018;Fotis et al., 2018), our findings show that acquisitive traits have consistently negative effects on ecosystem multifunctionality. However, this unexpected result is reasonable. Most of the forests in our study are secondary forests(average age=40.1 years),in which both pioneer and climax species coexist. The forests dominated by pioneer species with acquisitive traits have lower ecosystem multifunctionality, because pioneer species allocate more energy and resources to growth, thereby neglecting other ecosystem functions. In contrast, forests dominated by climax species with slow-growing conservative traits may have more robust ecosystem multifunctionality (Prado-Junior et al., 2016; Yuan et al., 2020). Therefore, from a forest management perspective, the results suggest that species with conservative traits should be favored to maintain higher ecosystem multifunctionality values. In summary,higher functional identity of acquisitive traits tended to decrease ecosystem multifunctionality,while functional diversity had not significant effects on ecosystem multifunctionality.

Fig. 4. Relative importance of stand biomass, stand age, biodiversity, and environmental factors as predictors for multiple threshold-based ecosystem multifunctionality. The relative importance of the predictors is expressed as the explained percentages of the variance and is calculated based on the absolute value of their standardized regression coefficients.Biomass: stand biomass; Age: stand age; CWM1 and CWM2: the first two principal component analysis axes of the community weighted mean traits; RaoQ:Rao's quadratic entropy;ENV1 and ENV2:the first two principal component analysis axes of the environmental factors; 10 to 90: 10%–90% based thresholds of the ecosystem multifunctionality.

Fig. 5. Relative importance of stand biomass, stand age, biodiversity, and environmental factors as predictors of individual ecosystem functions. The relative importance of the predictors is expressed as the explained percentages of the variance and is calculated based on the absolute value of their standardized regression coefficients. Biomass: stand biomass;Age: stand age; CWM1 and CWM2: the first two principal component analysis axes of the community weighted mean traits;RaoQ:Rao's quadratic entropy;ENV1 and ENV2: the first two principal component analysis axes of the environmental factors; PRO: forest productivity; SD: shrub density; HC: herb coverage; HR: herb richness; LB: litter biomass; LT:litter thickness; SC: soil carbon stock.

4.2. Stand biomass is a key factor linking biodiversity and ecosystem multifunctionality

Previous studies have shown that stand density(biomass or basal area per unit area) is a potential predictor of certain ecosystem functions(Lohbeck et al.,2015;Yuan et al.,2018,2019;Hao et al.,2020;Yue et al.,2020). Lohbeck et al. (2015) found that stand biomass had a positive effect on productivity in tropical secondary forests and explained this result using the “vegetation quantity hypothesis”. Corral Rivas et al.(2016)presented strong evidence of the same effects in Mexican forests.

In this study, we assess vegetation quantity effects on ecosystem multifunctionality,not only on ecosystem productivity.Our results show that increasing stand biomass strongly increased ecosystem multifunctionality,confirming the vegetation quantity effect in this region that is dominated by secondary forests.We speculate that,at the early stage of forest development, stands with higher biomass are characterized by greater crown and root volumes(Lohbeck et al.,2015;Yuan et al.,2018).Greater crown and root volumes make fuller use of available space and resources, thereby promoting higher levels of productivity (Hao et al.,2020; Yue et al., 2020). In addition, an accumulation of biomass may promote litter production and ameliorate soil conditions(Lohbeck et al.,2015), which in turn provide favorable abiotic conditions (higher relative humidity and soil nutrients) for the understory vegetation. These processes jointly contribute to improved ecosystem multifunctionality.Our results, based on individual ecosystem functions, show that herb coverage and richness increased with forest development, and thus provide additional support for this view(Fig.5 and Table S3).

Our results indicate that stand density(i.e.,vegetation quantity)may have equal or even stronger total effects than functional diversity or functional identity (i.e., vegetation quality) on ecosystem multifunctionality (Fig. 4). Previous studies have demonstrated that the magnitude and direction of density effects,including stand biomass and basal area, may vary with the forest development stage (Lohbeck et al.,2015; Yuan et al., 2018; Hao et al., 2020; Yue et al., 2020). During the early stage of forest development, stand biomass can increase forest productivity and other ecosystem functions as a result of efficient resource (light, water, and nutrients) acquisition (Lohbeck et al., 2015;Hao et al., 2020), while during the late stage of forest development,biomass may decrease forest productivity as a result of competition and negative density effects(Finegan et al.,2015).Most of the forests used in our study are secondary forests (average age = 40.1 years); therefore,positive relationships between stand biomass and ecosystem multifunctionality were detected. Our results indicate that these forests of northeastern China have the potential to deliver more ecosystem functioning and services at this stage of forest development.In addition,our results show that stand age had no direct effect on ecosystem multifunctionality,but indirect effects by increasing stand biomass.In a forest management context, it would thus be advisable to maintain a balance between the minimum stand biomass that ensures a certain level of ecosystem multifunctionality and the maximum biomass where competitive effects outweigh facilitation effects(Bauhus et al.,2017;Hao et al.,2020).

Our results based on the threshold approach show that the magnitude of the effect of biomass density tends to be smaller,while the magnitude of the effect of biodiversity properties tends to be greater when higher thresholds were considered. This also has implications for the management of secondary forests. If the aim is to sustain all functions at their maximum values,greater attention should be given to vegetation quality instead of vegetation quantity.

4.3. Environmental conditions

Within the context of climate change, it is important to understand the direct and indirect effects of environmental factors on ecosystem multifunctionality. Accordingly, climatic and soil conditions were included in this study as predictors. Two PCA axes represented the environmental gradients from harsher conditions with limited resources to fertile conditions with sufficient resources. We found that the first principal component of environmental conditions (ENV1, representing climate conditions)contributed to greater ecosystem multifunctionality,specifically higher temperature and increasing precipitation. These results are not surprising,and agree with previous studies that found that temperature and precipitation are the most crucial abiotic variables for tree productivity at the regional scale (Ratcliffe et al., 2016). However,there is evidence that soil conditions are strongly correlated with the availability of light, water, and nutrients (Balvanera et al., 2006).Different from studies at the local scale (Hao et al., 2018; Yuan et al.,2018),we found that the second principal component of environmental conditions (ENV2, representing soil conditions) did not influence ecosystem multifunctionality at the macroscale.

5. Conclusion

The objective of this study was to assess the key mechanisms affecting the relationships between biodiversity and ecosystem multifunctionality.One of the most significant findings to emerge from this study is that stand density (live biomass per ha) has a strong positive effect on ecosystem multifunctionality.Functional diversity generally had no significant effect on ecosystem multifunctionality, in contrast, functional identity had significant direct and indirect effects on ecosystem multifunctionality. Our results provide strong evidence for mass ratio and vegetation quantity effects, but little or no evidence of niche complementarity effects. The study thus illustrates the importance of accounting for stand biomass when analyzing forest diversity-multifunctionality relationships.

Funding

This research is supported by the Program of National Natural Science Foundation of China (No. 31971650), the Key Project of National Key Research and Development Plan (No. 2017YFC0504005), and the National Natural Science Foundation of China(No.31800362).

Availability of data and material

The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.

Authors’contributions

Q. Yue and M. Hao designed the conceptual idea. C. Zhang and X.Zhao designed the fieldwork and owned field data.Q.Yue analyzed the data and wrote the manuscript with input from M.Hao,L.Gao,Y.Geng,X.Wang,C.Zhang,and K.V.Gadow through multiple rounds of revision.All authors have approved the final article and have no conflict of interest.

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Declaration of competing interest

The authors declare that they have no competing interests.

Acknowledgements

We would like to thank all the students who took part in the forest field survey.

Appendix A. Supplementary data

Supplementary data to this article can be found online at https://do i.org/10.1016/j.fecs.2022.100027.