Evaluation of the Community Land Model-Simulated Specific Leaf Area with Observations over China: Impacts on Modeled Gross Primary Productivity

Abstract

Specific leaf area (SLA) is a key leaf functional trait associated with the ability to acquire light. Substantial variations in SLA have not been well described in the community land model (CLM) and similar terrestrial biosphere models. How these SLA variations influence the simulation of gross primary productivity (GPP) remains unclear. Here, we evaluated the mismatch in SLA between the CLM4.5 and observed data collected from China and quantified the impacts of SLA variation calculated from both observations and the default values across seven terrestrial biosphere models on modeled GPP using CLM4.5. The results showed that CLM4.5 tended to overestimate SLA values at the top and gradient of the canopy. The higher default SLA values could cause an underestimation of the modeled GPP by 5–161 g C m⁻² yr⁻¹ (1%–7%) for temperate needleleaf evergreen tree (NET), temperate broadleaf deciduous tree (BDT), and C3 grass and an overestimation by 50 g C m⁻² yr⁻¹ (2%) for temperate broadleaf evergreen tree (BET). Moreover, the observed SLA variation among species ranged from 21% to 59% for 14 plant functional types (PFTs), which was similar to the variation in default SLA values across models (9%–60%). These SLA variations would lead to greater changes in modeled GPP by 7%–19% for temperate NET and temperate BET than temperate BDT and C3 grass (2%–9%). Our study suggested that the interspecific variation in SLA and its responses to environmental factors should be involved in terrestrial biosphere models; otherwise, it would cause substantial bias in the prediction of ecosystem productivity.

1. Introduction

Plant functional traits are core attributes that reflect important ecological strategies in plants [1,2]. As a leaf functional trait, specific leaf area (SLA) is an important indicator of resource trade-off strategies under environmental changes, which determines the ability of a plant to capture light [3,4,5]. Specific leaf area refers to the leaf area per unit of leaf biomass, which is of great significance for vegetation photosynthesis and ecosystem productivity [6,7].

Substantial variation in SLA exists among plant functional types (PFTs) (e.g., needleleaf evergreen tree and broadleaf deciduous tree) and species for a given PFT. Numerous surveys have shown that SLA varies with different plant functional types (PFTs), following a general pattern of herb > shrub > tree, broadleaved > coniferous, and deciduous > evergreen [8,9,10]. Moreover, plant SLA among different species also shows significant variability. Based on the TRY Global Plant Traits Database, the variation coefficients (CV) of SLA data within the same PFT reach up to 23%–78% [9]. Such a large interspecific variation in SLA for a given PFT has also been reported in China and other regions [11]. For example, Liu et al. investigated plant SLA for 76 natural communities in China, the CV of different forest types ranged from 40% to 60%, and those of herb and shrub were 66% and 57%, respectively [12]. In the study by Reich et al., the CV of SLA was 34.3% for the broadleaf deciduous tree and 44.7% for the needleleaf evergreen tree based on field data across sites in the USA [13]. However, such a leaf trait diversity has not been well described in state-of-the-art terrestrial biosphere models.

In most terrestrial biosphere models, vegetation is grouped into different plant functional types with specific plant functional traits. This simple scheme is an effective way to characterize the functional diversity of ecosystems on the regional and global scales [14,15]. As for the expression of SLA, each plant functional type has a corresponding mean SLA value that can be calculated from observations. Few models consider the vertical variability of SLA within the canopy as in the community land model (CLM) [16,17]. Nevertheless, the default trait value settings for each PFT in different terrestrial biosphere models vary with each other and are different from the observations made at multiple spatial scales [9,18]. The uncertainties in plant trait parameters would considerably impact the prediction of vegetation productivity [19,20]. For example, the uncertainties in leaf longevity resulted in a greater than 30% change in vegetation biomass for temperate broadleaf evergreen tree in the Lund–Potsdam–Jena model [21]. However, how the SLA variation within each plant functional type influences the terrestrial gross primary productivity (GPP) simulation still remains unclear.

The objectives of this study were to investigate (1) the mismatch between simulated SLA in the CLM (version 4.5) and observed data for 1056 species in China and the differences in default SLA values across seven terrestrial biosphere models and (2) how the variation in SLA values influences GPP modeling for different ecosystems. Here, we first compared the default SLA values in the CLM4.5 model with published plant SLA observation data collected in China from 2005 to 2022. The variation in observed SLA data for each PFT and the difference in default parameter values of SLA in seven terrestrial biosphere models (i.e., BEPS, JULES, Hybrid, BIOME-BGC, SiBCASA, CLM4.5, and IBIS) were then quantified. Taking three forest ecosystems and one grassland ecosystem as an example, we finally quantified the uncertainties in modeled GPP caused by the mismatch in SLA between the mean observed values and original default values in CLM4.5, and the SLA variations among species and across different models for a given PFT.

2. Materials and Methods

2.1. Data

We collected 2632 records of observed specific leaf area (SLA) data for 1056 species from published papers during 2005–2020 by searching for keywords (i.e., leaf traits, specific leaf area, leaf mass per area, China) in the Web of Science and the China National Knowledge Infrastructure. SLA values estimated from models were excluded. These SLA data were distributed in different climates zones (Figure 1).

According to the plant functional types setting in the CLM4.5 model and observations in China, we statistically analyzed the SLA data for 15 plant functional types (PFTs), including temperate and boreal needleleaf evergreen tree (temperate and boreal NET), boreal needleleaf deciduous tree (boreal NDT), tropical and temperate broadleaf evergreen tree (tropical and temperate BET), tropical and temperate and boreal broadleaf deciduous tree (tropical and temperate and boreal BDT), temperate broadleaf evergreen shrub (temperate BES), temperate and boreal broadleaf deciduous shrub (temperate and boreal BDS), C3 grass, C4 grass, and rainfed crop (Table S1). Here, we combined C3 arctic grass and C3 grass with the same default SLA value in the CLM4.5 model into one plant functional type. We used the coefficient of variation (CV) to quantify the variation in SLA among species within the same plant functional type.

The default SLA parameter values for plant functional types were collected from seven state-of-the-art terrestrial biosphere models, including Boreal Ecosystems Productivity Simulator (BEPS), the Joint UK Land Environment Simulator (JULES), Hybrid model, Biome Biogeochemical model (BIOME-BGC), Simple Biospher and Carnegie-Ames-Stanford Approach model (SiBCASA), Community Land Model (CLM4.5), and Integrated Biosphere Simulator (IBIS) [18,22,23,24,25,26]. Considering the differences in PFT assignment among the seven models, we just compared the differences in default SLA values among the models for eight plant functional types, i.e., NET, NDT, BDT, BET, BES, BDS, grass, and crop). For the models with subgroups for each of the eight PFTs, we used the average default SLA values of all subgroups.

Half-hourly climate data (i.e., downwelling short-wave radiation (W m⁻²), downwelling long-wave radiation (W m⁻²), air temperature (K), precipitation (mm s⁻¹), surface pressure (Pa), relative humidity (%), and wind speed (m s⁻¹)) at four ChinaFLUX sites (http://www.nesdc.org.cn/ accessed on 15 June 2022) were collected to drive the CLM4.5 model. The four sites included Qianyanzhou subtropical coniferous forest (QYZ), Changbaishan temperate broadleaved Korean pine mixed forest (CBS), Dinghushan subtropical evergreen coniferous forest (DHS), and Haibei alpine shrub-meadow (HBG). The CBS temperate broadleaf deciduous forest is an old growth forest in northeastern China. The HBG is a typical grassland site in the Qinghai–Tibet Plateau. QYZ and DHS are two evergreen forests in southern China. A brief description of the sites’ characteristics is listed in

Table 1. Site information of four ChinaFLUX sites.

Site Name	QYZ	CBS	DHS	HBG
latitude (E)	26.74	42.40	23.17	37.67
longitude (N)	115.06	128.10	112.57	101.33
elevation (m)	102	738	300	3327
plant functional type	temperate NET	temperate BDT	temperate BET	C3 grass
simulated years	2003–2008	2003–2008	2003–2008	2003–2008

.

Missing climate data were filled in using observations from meteorological stations at the same site in the Chinese Ecosystem Research Network (CERN) (http://ww.cnern.org.cn/ accessed on 25 June 2022). A seven-day running mean diurnal cycle was used under the absence of station data. Other input data (e.g., plant functional type, soil depth, and texture) were also collected at these sites. More details were given by Zhang et al. [27].

2.2. Model

CLM4.5 is a state-of-the-art land surface process model in the Community Earth System Model (CESM1.2), which couples terrestrial biogeophysical processes, biogeochemical processes, and hydrological processes. CLM models have been evaluated widely for carbon fluxes and pools, evapotranspiration, leaf area index, land water storage, and soil moisture at different temporal and spatial scales [28,29,30,31,32,33,34]. The CLM4.5 model showed a good performance when simulating GPP at the above four sites [27].

In CLM4.5, photosynthesis in C3 and C4 plants was simulated based on the models of Farquhar et al. and Collatz et al. [35,36], respectively. Leaf net photosynthesis (${\mathrm{A}}_{\mathrm{n}}$) is defined as follows:

$${\mathrm{A}}_{\mathrm{n}}=\min \left({\mathrm{A}}_{\mathrm{c}},{\mathrm{A}}_{\mathrm{j}},{\mathrm{A}}_{\mathrm{e}}\right)-{\mathrm{R}}_{\mathrm{d}}$$

(1)

where ${\mathrm{A}}_{\mathrm{c}}$ is the RuBP carboxylase (Rubisco) limited rate of carboxylation ($\mathsf{\mu }\mathrm{mol} {\mathrm{CO}}_2{\mathrm{m}}^{-2} {\mathrm{s}}^{-1}$), ${\mathrm{A}}_{\mathrm{j}}$ is the maximum rate of carboxylation allowed by the capacity to regenerate RuBP ($\mathsf{\mu }\mathrm{mol} {\mathrm{CO}}_2 {\mathrm{m}}^{-2} {\mathrm{s}}^{-1}$), ${\mathrm{A}}_{\mathrm{e}}$ is the product-limited rate of carboxylation for C3 plants and the PEP carboxylase-limited rate of carboxylation for C4 plants, and ${\mathrm{R}}_{\mathrm{d}}$ is the leaf dark respiration rate ($\mathsf{\mu }\mathrm{mol} {\mathrm{CO}}_2 {\mathrm{m}}^{-2} {\mathrm{s}}^{-1}$).

Photosynthesis is calculated separately for sunlit and shaded leaves to scale up carbon flux from the leaf to canopy levels.

$${\mathrm{A}}_{\mathrm{c}}={\mathrm{A}}_{\mathrm{sun}}{\mathrm{LAI}}_{\mathrm{sun}}+{\mathrm{A}}_{\mathrm{shade}}{\mathrm{LAI}}_{\mathrm{shade}}$$

(2)

where ${\mathrm{A}}_{\mathrm{sun}}$ and ${\mathrm{A}}_{\mathrm{shade}}$ are the photosynthesis for sunlit and shaded leaves (CO₂ m⁻² s⁻¹), ${\mathrm{LAI}}_{\mathrm{shade}}$ is the sunlit and shaded leaf area indices, and ${\mathrm{A}}_{\mathrm{c}}$ is expressed as follows:

$${\mathrm{A}}_{\mathrm{c}}=\left\{\begin{array}{c}\begin{array}{c}\frac{{\mathrm{V}}_{\mathrm{cmax}}\left({\mathrm{c}}_{\mathrm{i}}-{\mathsf{\Gamma }}_{*}\right)}{{\mathrm{c}}_{\mathrm{i}}+{\mathrm{K}}_{\mathrm{c}}\left(1+{\mathrm{o}}_{\mathrm{i}}/{\mathrm{K}}_{\mathrm{o}}\right)} \mathrm{for} \mathrm{C}3 \mathrm{plants}\\ {\mathrm{V}}_{\mathrm{cmax}} \mathrm{for} \mathrm{C}4 \mathrm{plants}\end{array}\end{array}\right\}$$

(3)

where ${\mathrm{V}}_{\mathrm{cmax}}$ is the maximum rate of carboxylation (${\mathsf{\mu }\mathrm{mol} \mathrm{m}}^{-2 }{\mathrm{s}}^{-1}$); ${\mathrm{c}}_{\mathrm{i}}$ is the internal leaf CO₂ partial pressure (Pa); ${\mathrm{o}}_{\mathrm{i}}$ is the O₂ partial pressure (Pa); ${\mathrm{K}}_{\mathrm{c}}$ and ${\mathrm{K}}_{\mathrm{o}}$ are the Michaelis–Menten constants (Pa) for CO₂ and O₂, respectively; and ${\mathsf{\Gamma }}_{*}$ is the CO₂ compensation point (Pa).

${\mathrm{V}}_{\mathrm{cmax}}$ depends on temperature and is calculated from the value at 25 °C (${\mathrm{V}}_{\mathrm{cmax}25}$). ${\mathrm{V}}_{\mathrm{cmax}25}$ varies with foliage nitrogen concentration and specific leaf area.

$${\mathrm{V}}_{\mathrm{cmax}25}={\mathrm{N}}_{\mathrm{a}}{\mathrm{F}}_{\mathrm{LNR}}{\mathrm{F}}_{\mathrm{NR}}{\mathrm{a}}_{\mathrm{R}25}$$

(4)

where ${\mathrm{N}}_{\mathrm{a}}$ is the area-based leaf nitrogen concentration (g N m⁻² leaf area), ${\mathrm{F}}_{\mathrm{LNR}}$ is the fraction of leaf nitrogen in Rubisco (g N in Rubisco g⁻¹ N), ${\mathrm{F}}_{\mathrm{NR}}$ is the mass ratio of total Rubisco molecular mass to nitrogen in Rubisco (g Rubisco g⁻¹ N in Rubisco), and ${\mathrm{a}}_{\mathrm{R}25}$ is the specific activity of Rubisco ($\mathsf{\mu }\mathrm{mol} {\mathrm{CO}}_2 {\mathrm{g}}^{-1} {\mathrm{Rubiscos}}^{-1}$). ${\mathrm{N}}_{\mathrm{a}}$ is calculated from the mass-based leaf nitrogen concentration and specific leaf area:

$${\mathrm{N}}_{\mathrm{a}}=\frac{1}{{\mathrm{CN}}_{\mathrm{L}}{\mathrm{SLA}}_0}$$

(5)

where ${\mathrm{CN}}_{\mathrm{L}}$ is the leaf carbon-to-nitrogen ratio (g C g⁻¹ N) and ${\mathrm{SLA}}_0$ is specific leaf area at the canopy top (m² leaf area g⁻¹ C).

The vertical variability in SLA within the canopy was simulated by a linear function of SLA and the canopy depth as follows:

$$\mathrm{SLA}\left(\mathrm{x}\right)={\mathrm{SLA}}_0+\mathrm{mx}$$

(6)

where ${\mathrm{SLA}}_0$ is a fixed value of SLA at the top of the canopy (m² g⁻¹), m is a linear slope coefficient, and x is the canopy depth expressed as overlying leaf area index (m² overlying one-sided leaf area m² ground area). ${\mathrm{SLA}}_0$ and m are both assumed to vary with plant functional type. More details of the CLM4.5 model can be found in Oleson et al. and Thornton et al. [25,37]. Although the leaf nitrogen content and V_cmax25 are updated from static values in CLM4.5 to variates simulated by the leaf utilization of nitrogen for the assimilation (LUNA V1.0) model in the new version of the CLM model (CLM5.0), the quantification of SLA and related main processes in simulating photosynthesis remain the same as in CLM4.5 [38,39].

2.3. Analysis of the Impact of SLA Variation on Gross Primary Productivity

We conducted six model experiments to examine the impacts of SLA variation on gross primary productivity at the QYZ (temperate NET), CBS (temperate BDT), DHS (temperate BET), and HBG (C3 grassland) sites. Different SLA values were used in the six experiments. The default SLA values in the CLM4.5 model (p₀) were used in experiment S1. In experiment S2, we used the mean value of the observed SLA data collected in this study for each PFT (p_obs). In experiments S3 and S4, we added the observed SLA variation to the mean SLA value used in experiment S2 (i.e., p_obs − SD_obs and p_obs + SD_obs, respectively). In experiments S5 and S6, the variation in default SLA values across the models was added to the mean default SLA values (p_mod) of the seven terrestrial biosphere models (i.e., p_mod − SD_mod and p_mod + SD_mod, respectively).

For each model experiment, we first ran the CLM4.5 model to obtain the initial values of the state variables under equilibria using a two-stage spin-up approach. The first stage of spin-up followed the accelerated decomposition for 600 years; then, a normal decomposition was implemented for 200 years, with a repeating cycle of 6 years (2003–2008) with dynamic meteorological forcing and constant land use types, CO₂ concentration (284 ppmv), and N deposition (0.5 g N m⁻² yr⁻¹) at the pre-industrial level. A transient run was operated from 1850 to 2008 after reaching equilibrium. The CO2 concentration and N deposition data for the four sites were downloaded from the global dataset described in Thornton et al. [37].

The relative changes in gross primary productivity (GPP), RuBP-limited photosynthesis rate (Ac), and leaf area index (LAI) in model experiments S3 and S4 were compared with the results in experiment S2 (Equation (7)). The relative changes in GPP, Ac, and LAI in model experiments S4 and S5 were compared with the results in experiment S6 (Equation (8)).

$$\mathrm{R}1=\mathrm{MAX}(\frac{\left|{\mathrm{M}}_{{\mathrm{S}}_3}-{\mathrm{M}}_{{\mathrm{S}}_2}\right|}{{\mathrm{M}}_{{\mathrm{S}}_2}},\frac{\left|{\mathrm{M}}_{{\mathrm{S}}_4}-{\mathrm{M}}_{{\mathrm{S}}_2}\right|}{{\mathrm{M}}_{{\mathrm{S}}_2}})\times 100\%$$

(7)

$$\mathrm{R}2=\mathrm{MAX}(\frac{\left|{\mathrm{M}}_{{\mathrm{S}}_5}-{\mathrm{M}}_{{\mathrm{S}}_1}\right|}{{\mathrm{M}}_{{\mathrm{S}}_1}},\frac{\left|{\mathrm{M}}_{{\mathrm{S}}_6}-{\mathrm{M}}_{{\mathrm{S}}_1}\right|}{{\mathrm{M}}_{{\mathrm{S}}_1}})\times 100\%$$

(8)

where R1 and R2 are the relative changes in modeled GPP, Ac, and LAI caused by SLA variation from observations and models and ${\mathrm{M}}_{{\mathrm{S}}_{\mathrm{i}}}$ is the model outputs of GPP, Ac, and LAI in model experiments S1, S2, S3, S4, S5, and S6.

3. Results

3.1. Comparison of SLA between the CLM Model and Observations over China

We compared the SLA of mean observed values in China with default values in the CLM4.5 for 14 plant functional types (PFTs), as shown in Figure 2. The default SLA values were shown as broadleaf deciduous tree (BDT) > needleleaf deciduous tree (NDT) > broadleaf evergreen tree (BET) > needleleaf evergreen tree (NET), and broadleaf deciduous shrub (BDS) > broadleaf evergreen shrub (BES), which were consistent with the observations. However, the CLM4.5 model overestimated SLA values by 0.009 m²/g on average for most PFTs, except for tropical BET, temperate BET, and temperate BES. The positive bias of the default SLA values was highest (0.013 m²/g) in tropical BDT and lowest (0.003 m²/g) in temperate NET. By contrast, the CLM4.5 model underestimated the SLA values by 0.002 m²/g in tropical BET and temperate BET and by 0.005 m²/g in temperate BES.

We also evaluated the vertical gradients in simulated SLA varying with light from the bottom to top of the canopy for temperate BDT, temperate NET, and temperate BET in CLM4.5 (Figure 3). Compared with the observations, the CLM4.5 model overvalued canopy gradients of SLA in temperate BDT (Figure 3a) and temperate NET (Figure 3b), mainly because of the overestimation of SLA values at the top canopy. Specifically, the decrease of 0.016 m²/g in simulated SLA in temperate BDT was larger than the mean reduction of 0.014 m²/g observed for Phellodendron amurense Rupr., Fraxinus mandschurica Rupr., and Juglans mandshurica Maxim. when light intensity increased from 15% to 100%. For temperate NET, the simulated SLA decreased by 0.005 m²/g as light intensity increased from 15% to 100%, which was the same as the observed change in Picea asperata Mast., but larger than the slight reduction of 0.002 m²/g for Pinus koraiensis Siebold & Zucc. In addition, the simulated vertical gradients in SLA in temperate BET were within the range of SLA values observed for the three species (Figure 3c). However, the response of simulated SLA value was underestimated by 0.006 m²/g compared with the observed variation for Ficus tinctoria G.Forst. and Serianthes nelsonii Merr. when light intensity varied from 0 to 100%. Moreover, the CLM4.5 model ignored the interspecific diversity of plants in the response of SLA to light gradients, especially for the temperate BDT (Figure 3a) and BET (Figure 3c).

3.2. Interspecific Variation in Observed Plant SLA within Plant Functional Types

The SLA among different species within a given plant functional type in China showed large variability. The observed plant SLA varied from 0.0002 m²/g to 0.0997 m²/g, with a mean variation coefficient (CV) of 42% across different PFTs (

Table 2. The variation in specific leaf area (SLA) observations among different plant functional types (PFTs) in China. Abbreviations: NET, needleleaf evergreen tree; NDT, needleleaf deciduous tree; BET, broadleaf evergreen tree; BDT, broadleaf deciduous tree; BES, broadleaf evergreen shrub; BDS, broadleaf deciduous shrub.

PFT	Maximum SLA Value (m2/g)	Minimum SLA Value (m2/g)	Coefficient of Variation (%)	Number of Samples (n)
temperate NET	0.002	0.016	42.86	69
boreal NET	0.003	0.008	40.0	10
boreal NDT	0.005	0.024	40.0	11
tropical BET	0.004	0.027	28.57	103
temperate BET	0.0002	0.037	42.86	277
tropical BDT	0.007	0.032	41.18	9
temperate BDT	0.002	0.081	50.0	501
boreal BDT	0.005	0.028	25.0	4
temperate BES	0.004	0.048	47.06	221
temperate BDS	0.002	0.010	59.09	454
boreal BDS	0.004	0.025	21.05	6
C3 grass	0.002	0.082	57.14	822
C4 grass	0.003	0.052	47.37	126
rainfed crop	0.002	0.042	52.63	19

). The interspecific variations in SLA were relatively small for tropical BET, boreal BDT, and boreal BDS, with the CV values varying from 21% to 29%, but were large for the remaining PFTs, with the CV values ranging from 40% to 59%. Among the different tree and shrub types, temperate BDT and temperate BDS had the largest CV of SLA, respectively. Moreover, temperate plants with a great interspecific diversity showed higher SLA variability than tropical and boreal plants in China. Specifically, the CV of SLA for temperate BDT was two times that for boreal BDT.

3.3. Variation in the Parameter Values of SLA among Different Terrestrial Biosphere Models

Figure 4 displayed the variation in default SLA values among seven terrestrial biosphere models within eight plant functional types, compared with the observed SLA data. Almost all default SLA value settings in the models for each PFT were in the range of the observed SLA values in China but had large differences within the same PFT. The CV of the default SLA values across models within one PFT varied from 8.7% for crop to 60.0% for NET. Particularly, the SLA value of grass in the BIOME-BGC model was higher than that in the other models by 0.024 m²/g on average, although SLA variation in grass across the seven models was relatively low. The SLA values for the PFTs of BES and NDT in the SiBCASA model were high, which were greater than that of the others by 0.016 m²/g and 0.011 m²/g, respectively. By contrast, the SLA values in the JULES model were generally lower than that in the other models by 0.011 m²/g on average for all PFTs. The variation in SLA value assigned among models for a given PFT was lower than that from the observations, except for NET, BET, and BES. Specifically, the CV of SLA values among models for NET was higher than the CV of observed values by 17.1%.

The mean default values across models for a given PFT was greater than the observed mean by 0.004 m²/g on average (Figure 4). The difference was largest for crop, reaching 0.009 m²/g, followed by NDT, with 0.0087 m²/g. The mean default SLA values across different models were approximately equal with the observed mean for BES. In addition, the default SLA values in few models were close to the mean observed values. For example, the SLA value of BET in the CLM4.5 model was 0.012 m²/g and the average observation was 0.014 m²/g. The difference in SLA values between the hybrid model and the observed mean was only 0.001 m²/g.

3.4. Impacts of Variation in SLA on Modeled Gross Primary Productivity

We first estimated the differences in CLM4.5-modeled gross primary productivity (GPP) at the QYZ, CBS, DHS, and HBG sites between model experiment S2 using the observed SLA values and model experiment S1 with the default SLA values. Figure 5 presented the changes in annual GPP, mean photosynthesis rate (Ac), and mean leaf area index (LAI) at four sites in the experiment S1 compared with those in experiments S2. The mismatch in SLA between default values in CLM4.5 and observations had a larger influence on GPP simulation for temperate NET (QYZ) and temperate BET (CBS) than that for temperate BDT (DHS) and grass (HBG). The overestimation of the default SLA value by 0.003 m²/g for temperate NET could result in a lower annual GPP estimation at QYZ by 161 g C m⁻² yr⁻¹. This weakened productivity was mainly caused by a decrease in the photosynthesis rate in spite of an increase in LAI. Similarly, the higher default SLA value for temperate BDT by 0.008 m²/g, induced a lower modeled GPP at CBS by 69 g C m⁻² yr⁻¹. On the contrary, a 0.002 m²/g underestimation of the SLA value for temperate BET in the CLM4.5 caused a decrease in LAI and a slight increase in Ac, which led to a higher annual GPP by 51 g C m⁻² yr⁻¹. In addition, a large, overvalued SLA in CLM4.5 for C3 grass had a small impact of 5 g C m⁻² yr⁻¹ on simulated GPP and Ac at HBG, although there was an increase in LAI by 0.3 m²/m².

We then quantified the impacts of SLA variation calculated from the observations and models on the modeled GPP at these sites.

Table 3. The relative changes in GPP, Ac, and LAI modeled with different SLA values. R1 is the relative changes in modeled GPP, Ac, and LAI caused by the observed SLA variation based on model experiments S2, S3, and S4 (Equation (7)). R2 is the relative changes in modeled GPP, Ac, and LAI caused by the SLA variation across models based on model experiments S1, S5, and S6 (Equation (8)). Abbreviations: PFT, plant functional type; NET, needleleaf evergreen tree; BDT, broadleaf deciduous tree; BET, broadleaf evergreen tree.

Site	PFT	Relative Change	CV of SLA (%)	GPP (%)	Ac (%)	LAI (%)
QYZ	temperate NET	R1	42.9	7.0	43.9	14.1
		R2	60.0	18.5	31.0	61.7
CBS	temperate BDT	R1	50.0	6.3	30.9	16.9
		R2	29.4	8.8	24.5	43.8
DHS	temperate BET	R1	42.9	8.0	37.3	14.7
		R2	42.8	14.1	24.6	63.7
HBG	C3 grass	R1	57.1	3.3	57.9	48.4
		R2	34.3	1.6	3.3	46.2

presents the CV of SLA quantified by observed data (Table 2) and default values across seven terrestrial biosphere models (Figure 4), and the corresponding relative changes in modeled GPP, Ac, and LAI with the involvement of these variation in SLA calculated by Equations (7) and (8), respectively. The results showed a larger influence on modeled GPP at QYZ and DHS than that at the CBS and HBG sites due to the variation of the SLA values. The SLA values in temperate NET and temperate BET both changed by 43% at the mean observed level, which caused the modeled GPP at QYZ and DHS to change by 7% and 8%, respectively. The observed SLA variations in temperate BDT and C3 grass had little impact on the GPP simulation at CBT and HBG, although Ac and LAI changed greatly. In particular, the SLA value of C3 grass changed by 57% at the mean observed level and could lead to Ac and LAI changes by 48% and 58%, while GPP at the HBG site only changed by 3%. Moreover, a large SLA variation among models in temperate NET and temperate BET could result in greater impacts on the modeled GPP than the effects from variations in the observation. For example, SLA variation among the models with a CV of 60% in temperate NET caused a GPP change of 19%, resulting from changes in Ac with 31% and LAI with 62% at QYZ.

4. Discussion

Although current terrestrial biosphere models generally consider the differences in SLA among plant functional types (PFTs), there still remain mismatches in SLA values between model and observations. Our results suggested a remarkable mismatch of SLA between the CLM4.5 model and observations collected in China, especially for tropical broadleaf deciduous tree overestimated by 0.013 m²/g (Figure 2). It is necessary for us to revise the SLA values when we use the CLM4.5 model to simulate terrestrial carbon cycle dynamics at the regional scale. Our results also showed that the observed SLA values for China plant were higher than the global average for all PFTs (Figure 6), which is supported by a recent work showing that SLA in Asia was higher than that in Europe and North America by about 0.004 m²/g and 0.001 m²/g for a given leaf dry matter content of 0.25 g/g [44]. The higher values of SLA in Asia and China might be caused by intense resource competition stress, where the communities are dominated by more acquisitive species with high SLA [45]. Therefore, simply referring to other models or global traits data to set SLA values in a regional gross primary productivity simulation study is unreliable. Furthermore, since the observed SLA data for some plant functional types (e.g., temperate BDT, temperate BET, and C3 grass) in China and around the world do not follow normal distributions [9,46], the default model values and mean observed values remain uncertain. We recommend using the probability density distribution of observed trait data within a given PFT rather than setting trait parameter values based on the mean values.

The SLA among species within a given plant functional type has substantial variability in China and other regions due to differences in the genotype and environmental changes [47,48], especially climate and soil factors, e.g., temperature, light, precipitation, and soil nutrient [49,50]. As trait variation was revealed by previous studies as inducing large effects on the simulations of ecosystem productivity (e.g., GPP, NPP, and NEP) and biomass (e.g., vegetation biomass and litter pool) [21,23,51], our results also suggested that SLA observation variation could bring great uncertainty to the GPP simulation, especially for temperate NET and temperate BET. As shown in Table 3, the SLA value varying by 43% at the mean observation level could induce a change in the simulated annual GPP by 161 g C m⁻² yr⁻¹ (7.0%) for temperate NET. It is necessary to describe trait variation within each plant functional type to reduce model uncertainties.

The responses of plant SLA to environmental changes, as shown in previous studies, have also not been well quantified in terrestrial biosphere models. The SLA vertical variability under the change of light in canopy involved in the CLM4.5 model has not been considered in many other terrestrial biosphere models (Figure 2). Moreover, the variation in SLA with development stages was only simulated in some crop models, such as the WOFOST (World Food Studies) model [52], as shown in Figure 7a. These variations in SLA at different growth stages for both crops and trees (Figure 7b) should be added to terrestrial biosphere models in the future. In addition, SLA also varies with other factors, such as water stress (Figure 8a) and soil nitrogen content (Figure 8b), which need to be paid more attention in terrestrial biosphere models.

Within a given plant functional type, plants grouped by morphology and structure may better describe the variability in traits [66,67]. In the Functionally Assembled Terrestrial Ecosystem Simulator (FATES) in CLM5.0, the plant population in each patch is divided first into plant functional types, and each plant type is presented as numerous height classes, according to tree height, diameter, canopy layer, and other variables [68]. Trait–climate relationships have also been analyzed in few modelling studies. For instance, SLA, V_cmax25, and J_cmax25 within the PFTs were re-parameterized yearly depending on the local climatic conditions (e.g., MAT, MAP, and soil moisture) in the JSBACH model [69]. Yang et al. used trait–climate relationships to predict the spatial patterns of LMA, N_mass, and LAI and then simulated vegetation distributions and vegetation responses to climate changes in China using a Gaussian mixture model trained with these trait data [70]. However, the underlying mechanisms behind these trait–climate relationships require more long-term observations, to be able to simulate vegetation responses to future climate change.

The collected SLA data in this study mainly distributed in temperate and subtropical forest and grassland ecosystems, SLA observations in boreal ecosystems, woody plants in the Tibetan Plateau region, and crops were relatively few, and need to be supplemented by further research. Different protocols for measuring SLA (e.g., all leaves vs. the topmost leaf, with vs. without petioles) have been used in published studies [15], which may cause bias in data statistics and analysis. The methods, time, and positions of sampling should be standardized in the future to enhance the representativeness of the plant trait database, especially at the region scale. This paper only quantified the impacts of SLA variation on GPP simulation with the CLM4.5 model for four PFTs; its effects in other PFTs and terrestrial biosphere models need to be further investigated.

5. Conclusions

In this study, we evaluated the CLM4.5-simulated SLA against the observed data collected from China and examined the impacts of SLA variation on GPP simulation using the CLM4.5 model. The results showed that CLM4.5 overestimated the default SLA values at the top of canopy for 11 PFTs and the canopy gradient of SLA for temperate BDT and temperate NET. The higher default SLA values in temperate NET, temperate BDT, and C3 grass caused an underestimation in the modeled GPP at the QYZ, CBS, and HBG sites compared with the results from the mean SLA observations. Substantial SLA variation could cause great changes in modeled GPP, especially for temperate NET and temperate BET. Our study suggested that the interspecific variation in SLA within a given PFT and its responses to environmental changes should be considered in terrestrial biosphere models to reduce the uncertainties in GPP and LAI estimations. More efforts are needed to make full use of the plant trait database to understand the underlying mechanisms of trait variation and to promote model development so as to enhance the prediction ability of ecosystem responses to future climate changes.