Effects of Forest Canopy Vertical Stratification on the Estimation of Gross Primary Production by Remote Sensing

Gross primary production (GPP) in forests is the most important carbon flux in terrestrial ecosystems. Forest ecosystems with high leaf area index (LAI) values have diverse species or complex forest structures with vertical stratifications that influence the carbon–water–energy cycles. In this study, we used three light use efficiency (LUE) GPP models and site-level experiment data to analyze the effects of the vertical stratification of dense forest vegetation on the estimates of remotely sensed GPP during the growing season of two forest sites in East Asia: Dinghushan (DHS) and Tomakomai (TMK). The results showed that different controlling environmental factors of the vertical layers, such as temperature and vapor pressure deficit (VPD), produce different responses for the same LUE value in the different sub-ecosystems (defined as the tree, shrub, and grass layers), which influences the GPP estimation. Air temperature and VPD play important roles in the effects of vertical stratification on the GPP estimates in dense forests, which led to differences in GPP uncertainties from−50% to 30% because of the distinct temperature responses in TMK. The unequal vertical LAI distributions in the different sub-ecosystems led to GPP variations of 1–2 gC/m2/day with uncertainties of approximately −30% to 20% because sub-ecosystems have unique absorbed fractions of photosynthetically active radiation (APAR) and LUE. A comparison with the flux tower-based GPP data indicated that the GPP estimations from the LUE and APAR values from separate vertical layers exhibited better model performance than those calculated using the single-layer method, with 10% less bias in DHS and more than 70% less bias in TMK. The precision of the estimated GPP in regions with thick understory vegetation could be effectively improved by considering the vertical variations in environmental parameters and the LAI values of different sub-ecosystems as separate factors when calculating the GPP of different components. Our results provide useful insight that can be used to improve the accuracy of remote sensing GPP estimations by considering vertical stratification parameters along with the LAI of sub-ecosystems in dense forests.


Introduction
As a key component of the global carbon cycle, the terrestrial gross primary production (GPP) in forests plays an important role in the global carbon, water, and energy cycles [1][2][3].GPP is the amount of illumination energy from the sun that is transformed into biomass energy, and is dependent on the illumination, heat, and water in terrestrial ecosystems [4,5].Precisely estimating GPP is essential for evaluating terrestrial carbon cycles, which also influence climate change [6].Approximately 40% of GPP is directly or indirectly used by humans [7].Vegetation is the key component in the ability of humans to adjust the climate and mitigate increasing atmospheric CO 2 [8].Currently, several GPP remote sensing products, such as MPI-BGC [9], MODIS [10], JULES [11,12], BESS [13], and BEPS [14,15], which estimate GPP over large regions, exhibit different estimation accuracies for different vegetation types.
Forest ecosystems cover more than 30% of the land surface and represent 80% of the aboveground biomass on Earth [16,17].In primeval forest ecosystems, the aboveground parts of the forest system can be divided into two components, the overstory and understory.The overstory component is the tree layer that intercepts direct solar radiation.The understory is defined as the vegetation structures under the upper layer of the canopy, and it is mainly composed of species of vegetation that prefer shaded conditions.In tropical and subtropical dense forests, the understory component can be divided into shrub and grass layers, but in temperate forests, this component consists of mostly grasses and herbs [18,19].
Both overstory and understory species are important contributors to the carbon uptake in terrestrial ecosystems [20,21].Some studies have shown that understory vegetation contributes 32% of the total GPP in tropical forests [22,23], but this contribution may be 25% in temperate forests [24].A recent study showed that the shaded GPP has increased by 1.1% over the past three decades [25].This variability is due to variations in the light, temperature, and humidity conditions in the understory and overstory.Shaded leaves also have unique characteristics, with understory and overstory vegetation exhibiting different responses to the environment and various photosynthetic traits [26][27][28][29][30].However, in regions with high biodiversity and a high leaf area index (LAI) values in the forest zone, such as evergreen broadleaf forests (EBF) and artificial forests, the estimation of GPP exhibits high variance [31,32].Furthermore, GPP estimations by remote sensing methods combined with canopy photosynthesis models do not consider the understory incident radiation, and thus may underestimate the total GPP [33][34][35].
In dense forest ecosystems, canopy species show multiple leaf functional traits that are important for GPP estimation.Leaf traits vary greatly among overstory and understory species.When scaling up different leaf traits to whole ecosystems [36,37], the structural parameters of ecosystems are some of the most important controlling factors [38,39].However, species in the same vertical stratification, such as the tree and grass layers in temperate forests and the trees, shrubs, and grasses in tropical rainforests, have similar living environments, and thus similar leaf traits [40][41][42][43][44][45][46][47][48].
Remote sensing is a useful approach for estimating ecosystem GPP at larger spatial scales [8,49].In forest sites, many studies have demonstrated that micrometeorological parameters such as temperature, humidity, photosynthetically active radiation (PAR), and longwave incident radiation exhibit vertical gradients [50][51][52][53][54].Moreover, the physiological parameters of vegetation such as chlorophyll contents also exhibit vertical gradients [42,55,56].However, it is difficult to capture these signals from deep in the canopy.Many scientists have used models and methods to retrieve vertical environmental and physiological parameters such as temperature, humidity, PAR, the fraction of absorbed PAR (FPAR), the LAI, and the light use efficiency (LUE) [57][58][59][60][61][62].Wu et al. and Kosugi et al. [63,64] also determined the rules of the leaf-scale maximum carboxylation rates (Vcmax25) at different heights.All of these models or methods were developed in specific regions; therefore, it may be difficult to determine the parameters on a global scale.In this article, we introduce three GPP models based on LUE models that are based on the Monteith [65] methods.Monteith's function is calculated as: where LUE is the daily LUE, and FPAR is the fraction of the absorbed PAR estimated from the LAI and zenith angle.The daily LUE is calculated by different factors that are used to determine LUEmax, and the numerical exercises have specific input parameters.Different models have various parameterizations of environmental stress for LUE; for example, the vegetation photosynthesis model (VPM; Xiao et al. [66]) assumes that temperature is of low importance to LUE at optimal growing temperatures, but it is of high importance at high (50 • C) and low temperatures (0 • C).The MOD17 (MOD, Running et al. [4,67]) model was designed to use remote sensing data for GPP estimation, and this model utilizes a higher stress factor (less than 10% of the LUE for optimal growth) at low daily minimal temperatures (−8 • C) and a low stress factor (100% of the LUE for optimal growth) at high daily minimal temperatures (9 • C).The model developed by Wu et al. [63] assumed that high temperatures were important for LUE (more detailed information can be found in the Supplementary Material).
However, the models mentioned above did not consider the effect of vegetation stratification on the forest canopy.The goal of this study was to investigate the effect of vertical stratification on remotely sensed GPP estimations in selected dense forest.We raised the following two questions: (a) how do environmental drivers vary with height in multilayer dense forests; and (b) does adding vertical stratification parameters to remotely sensed GPP algorithms improve the performance of the estimations in forest canopies?To address these questions, we used eddy covariance flux data from two sites, including one EBF and one artificial Japanese larch forest, in three remotely sensed LUE models to evaluate the responses of GPP to different factors during the growing season.These sites have meteorological data measurements at multiple heights combined with vertical stratification information, such as temperature, vapor pressure deficit (VPD), and PAR.

Study Sites
In this research, we selected Dinghushan (DHS) and Tomakomai (TMK) as the research areas.Both of these areas have dense forests and a conspicuous vertically stratified vegetation structure during the growing season from April to October.DHS, which is located in a monsoon subtropical forest climate zone in the southeast part of China, has a 15 m-high broadleaved forest canopy, while TMK, which is located in the northeast part of Japan, has a 15 m-high deciduous needle leaf (DNF) canopy.Vertical meteorological and canopy eddy covariance (EC) data were measured at both sites.A certain amount of understory LAI (>0.9 m 2 /m 2 ) was found in these two sites during the growing season.Two selected sites' detailed information is showed in Figure 1 and Table 1.

Tomakomai (TMK)
The TMK study area is in the Tomakomai National Forest, which is approximately 15 km northwest of Tomakomai, Hokkaido, in northern Japan.The area of larch forest is approximately 100 ha.From 2001-2003, the average temperature was 6.5 • C, and the site was warm in summer and snowy in winter.The mean annual rainfall was 1055 mm.The growing season extends from April to October, while there is little vegetation cover in other months.
The flux tower was located in the center of the larch forest (42 • 44 13.1" N, 141 • 31 7.1" E, 115 masl), and the surrounding trees were approximately 45 years old at the time of this study.The height of the vegetation canopy was 18 m in 2003.The dominant species near the tower were classified as either overstory or understory species.The overstory species included Japanese larch, birch, spruce, and Japanese elm, and the understory species were ferns and other herbs.The LAI values near the flux tower in the overstory and understory layers were 5.4 m2 /m 2 and 3 m 2 /m 2 , respectively, in August 2005 [72].
The flux tower contained eight levels of meteorology instrumentation for measuring different variables such as temperature and humidity above and within the canopy, which were measured with a platinum resistance thermometer and a capacitive hygrometer (HMP-45D, Inc., Logan, UT, USA) at 2 m, 5 m, 8 m, 14 m, 18 m, 22 m, 27 m, and 40 m.The PPFD was measured at 2 m, 5 m, 14 m, 18 m, and 40 m by a quantum sensor (LI-190S, LI-COR Inc., Lincoln, NE, USA) in the tower and on the forest floor.At 27 m, the CO 2 flux between this ecosystem and the atmosphere was measured by open and closed-path EC instruments (LI-7500 and LI-6262, LI-COR Inc., Lincoln, NE, USA).The tower was destroyed by a typhoon in 2004, and monitoring therefore ceased.All of these data were recorded at 30-min intervals.LAI was measured every month from July to December 2003.During the study period from 1 July 2005 to 31 December 2005, in DHS, the meteorological instruments on the sixth and seventh levels malfunctioned, so all of the data from 31 m and 36 m were missing for this period.Therefore, we did not compare the meteorology data from 1 July 2005 to 4 December 2005.For the other data during the study period at the two sites, the temperature and VPD data were gap-filled by the data processor, and low-quality and missing data were removed [68,69,71].

PPFD Data Preprocessing
The PPFD data had some missing entries during the study period that were gap-filled by the moving average method at 30-min intervals, and the PPFD data for different vertical layers in the DHS site were provided [71].To convert the mean PPFD data from µmol/m 2 /s to W/m 2 /s or MJ/day for the PAR input for each layer, we used a linear regression between shortwave radiation at the top of the tower and the PPFD of the top layer.These two factors were very highly correlated (y = 1.4907x − 0.2526, x is the PPFD, unit in µmol/m 2 /s; y is the shortwave radiation, unit in W/m 2 /s; r 2 = 0.9951, root mean square error (RMSE) = 16.37 W/m 2 , bias (Bias is defined as the predicted result minus the observational data.)<2.5 W/m 2 ).In TMK, the PPFD was also calculated using this method.

GPP Data Processing
EC instruments directly measure the carbon exchange between terrestrial ecosystems and the atmosphere, but the data must be transformed to determine the GPP.The net ecosystem exchange (NEE) in these two sites included exchange items (Fc) and storage items (Fs), so these values would be more accurate at these two sites under high LAI.
In this article, we converted the FLUXNET2015 flux GPP data at DHS and the gap-filled Asia Flux NEE at TMK to GPP using the flux analysis tool [74].After obtaining the quality-controlled NEE, we partitioned the NEE into GPP and RE.For unreliable data, we used the harmonic analysis of time series (HANTS) algorithm [75] to gap-fill the GPP data at 30-min intervals.The 30-min GPP data were excluded when the solar radiation in the PPFD was less than 10 µmol/m 2 /s.All of the data, including the carbon flux data, micrometeorological data, and LAI, were from the same source to ensure accurate representation of the carbon exchange near the flux tower.

Methods
In this section, we used three LUE-based GPP models and parameterization methods, which applied unequal temperature and water scalars to LUE.These weights can be applied to remote sensing methods and site-level GPP modeling.
Diverse species inhabit the overstory and understory layers of the forest canopy, which have unique ecophysiological parameters.The controlling environmental factors, such as temperature, humidity, and light radiation, have vertical gradients corresponding to different heights in the canopy.However, for certain vegetation groups, such as trees, shrubs, and grasses, and their living heights in the forest canopy, the diverse ecophysiological properties of individual species are much more likely to be affected by vertical gradients than the whole forest ecosystem [43,44].To simplify the photosynthetic processes of multiple species in this study, we defined the concept of sub-ecosystems, which divided the entire terrestrial ecosystem into two or three parts.Examples include the overstory and understory sub-ecosystems in TMK and the arbor, shrub, and grass sub-ecosystems in DHS.Different sub-ecosystem LAI values can be used as the criteria to separate ecosystems into two or three individual types, where GPP would be the aggregation of three different types of photosynthetic processes, because leaves are the only photosynthetic organs in a plant species.

Theoretical Background of the GPP Model
In this research, we used three LUE-based GPP models, namely, VPM [66], MOD17 [10,67], and Wu [63], to model the GPP according to the following formula: In Equation (2), PAR is the photosynthetically active radiation, and FPAR is the fraction of the absorbed PAR, which is calculated as: where k is the extinction coefficient, and LAI is the leaf area index in the selected area at a specific canopy height.LUEmax is the maximum LUE, which is scaled by temperature and the VPD response under diverse situations.Therefore, LUE can be expressed as: Detailed information on the three GPP models can be found in the Supplementary Material.The models have different parameterizations for f(T) and g(VPD).

GPP Estimation Considering the Vertical Structure
In situ meteorology data for the height of each sub-ecosystem are available for the DHS and TMK sites.Thus, the GPP of every sub-ecosystem can be described by the GPP model with the meteorology data and model parameters from each height.In this research, we assumed that each sub-ecosystem would have a different LUEi, which is defined as the LUE of the sub-ecosystem.LUEi is affected by only certain environmental factors, such as temperature and the VPD at the corresponding height.LUEi was calculated using Equation (4), and the LUEmax of the sub-ecosystems (such as those of shrubs and grasses), temperature, and the VPD stress functions are provided in the Supplementary Material.
At the ecosystem scale, using temperature and the VPD scalar produces the LUEmax, which characterizes carbon sequestration from the LUE.Moreover, different vegetation types, such as the tree, shrub, and grass layers, have various LUEmax values.Table 2 presents the various LUEmax inputs for the GPP estimation results from four different articles.The results of Zhou and Madani stated that overstory layers, such as the overstory tree layer, have an average LUEmax of 0.28 gC/m 2 /MJ, which is higher than the LUE of the understory layer (grass).
Table 2. Tree, shrub, and grass maximum light use efficiency (LUEmax) parameterization values (units: gC/m 2 /day) (EBF and DNF were defined as the LUEmax values of the tree layer in the GPP simulations).

Tree (DHS) Overstory (TMK) Shrub Grass
Zhou [32] 1.120 1.180 0.930 1.520 MODIS [67] 1.268 1.051 0.841 0.860 Yuan [31] 1.680 1.640 0.660 1.520 Madani [76] 0.980 1.171 0.631 1.294 To analyze the effect of vertical structure on GPP estimation, we calculated GPP considering vertical layers using the following function: where GPPm is defined as the modeled sum of GPP in two (overstory and understory) or three (tree, shrub, and grass) sub-ecosystems.The ecosystems were divided into overstory and understory groups at TMK and into tree, shrub, and grass groups at DHS using in situ LAI measurements.The meteorology data can provide in situ air temperature and VPD stress values for each sub-ecosystem.Temperature sensors at different heights measured the air temperature, while the VPD was calculated from the relationship between VPD and both temperature and humidity [77,78].APARi is the absorbed PAR for the LUE of a specified height, and was calculated from FPAR and PAR.At the DHS site, the whole ecosystem can be divided into three sub-ecosystems, namely, the tree, shrub, and grass layers; thus, the GPP was calculated for the EBF, shrub, and grass layers (sub-ecosystems) using the in situ temperature, VPD, PAR, and long-term LAI data for each layer, and specific LUEmax input values for each sub-ecosystem.At the TMK site, we divided the ecosystem into DNF and grass sub-ecosystems, and used the GPP input parameters for each layer, as used in the approach for DHS: GPPs represents the GPP of a specific ecosystem that was estimated using the big-leaf model, which assumes that the entire ecosystem is a single leaf.At the DHS site, the LUE was calculated as an EBF ecosystem with an EBF-based LUEmax (Table 3), and the temperature and VPD stress were given by the mean canopy temperature and VPD.At TMK, the DNF-based LUEmax (Table 3) was used, and the temperature and VPD stress were given by the mean canopy temperature and VPD.
The GPP bias (GPPdif in Equation ( 7)) is the multilayer GPP result (GPPm) minus the single-layer GPP result (GPPs), which represents the difference between the two models.Therefore, GPPdif is the difference between considering and not considering vertical stratification.TMK and DHS have overstory and understory vegetation, so using representative GPPs for the whole canopy would lead to estimation uncertainties.The modeled GPP values produced by a variety of models were given their own terms.For example, if a GPP result was produced without considering the vertical structure in a VPM, this value was represented as GPP_s_VPM.If a GPP result was produced based on a MODIS model that considered the vertical structure, it was represented as GPP_m_MOD.

Modeling Vertical PAR with a Radiative Transfer Model
In this study, we also compared the measured PAR and modeled PAR values to simulate PAR in the vertical direction.According to the light transfer function in the vegetation canopy [58], the vertical PAR at different heights impacts the LAI through light absorption, reflectance, and transmittance.In the vertical direction, incident PAR for direct light follows the Beer-Lambert law: Remote Sens. 2018, 10, 1329 where PAR 0 is the incident PAR at the top of the canopy, z is the length of light transmission through the canopy, and the τ b is the transmitted ratio, which is calculated as: where K b is the extinction coefficient, which is affected by the zenith angle φ [79] Therefore, PAR(z) is the understory incident PAR, which can be expressed as follows: To test the proportional impact of the LAI in different sub-ecosystems on the GPP estimation in dense forests, we used the discrete anisotropic radiative transfer (DART [80,81]) model to simulate APAR based on clear sky conditions in the tree, shrub, and grass layers under different vertical LAI distributions.

Statistical Indicators
At the same time, the maximal air temperature in the canopy was usually found at the top of the canopy, while the minimal value was often found at the ground layer.We evaluated the largest variance in the vertical temperature and VPD profiles using two parameters: Tdif and VPDdif.Tdif is defined as the air temperature at the top of the canopy minus the temperature at the grass layer.VPDdif is the VPD of the top layer of the canopy minus the VPD of the grass layer.

Analysis of Variance in Forest Canopy Vertical Stratification
In forest canopies, vertical temperature and VPD stratification vary with forest type and season.Figure 2A,B show the differences in the daily average vertical temperature from day of year (DOY) 210 to 365 at the DHS site (EBF) in 2005, and from DOY 182 to 365 in 2003 at the TMK site (DNF).
Variances in air temperature occur in forest canopies during the growing season.In most cases, the air temperature in the tree layer is higher than that in the grass layer, with an average temperature difference of approximately 0.4 • C at the DHS site.In the DNF vegetation type at the TMK site, the mean temperature difference was approximately 0 • C before DOY 240.From DOY 240 to 280, the temperature of the tree layer was 0.5 • C higher than that of the grass layer.However, during the non-growing season after DOY 290, the temperature differences of the forest layer ranged from 0.5 • C to 1.5 • C, which means that the temperatures in the tree layer were higher than those in the grass layer.
Similarly, the VPD of the tree layer in the forest canopy was higher than that of the grass layer.VPD varies with air temperature.During the summer to autumn period (DOY 180 to 300), the daily average VPD of the tree layer was approximately 120 Pa higher than that of the grass layer at the DHS site, and approximately 70 Pa higher than that at the TMK site.With the decrease in temperature during the winter period (DOY 300 to 365), the VPD of the tree layer was 50 Pa higher than that of the grass layer at the DHS site, and 25 Pa higher than that at the TMK site, which was not as great as the difference during the summer period.Figure 2B,D shows that Tdif and VPDdif have a robust linear relationship in months during the growing season.

Analysis of Variance in Forest Canopy Vertical Stratification
In forest canopies, vertical temperature and VPD stratification vary with forest type and season.Figure 2A,B show the differences in the daily average vertical temperature from day of year (DOY) 210 to 365 at the DHS site (EBF) in 2005, and from DOY 182 to 365 in 2003 at the TMK site (DNF).In a forest canopy, overstory incident PAR is much higher than that in the understory layer because the incident energy in the understory is affected by sky conditions and the overstory forest structure parameters.At the TMK site, the understory incident PAR was approximately 10% of the overstory incident PAR during the growing season when the LAI of the overstory was high (LAI = 5).However, when the overstory LAI was low during the winter period (DOY 300 to 365), the understory incident PAR was 30% of the overstory incident energy (Figure 3A).However, in EBFs such as those at the DHS site, the overstory LAI did not intensely change throughout the year.Figure 3B shows that the incident PAR of the shrub layer at the DHS site was 40% that of the tree layer, and the incident PAR of the grass layer was 20% that of the tree layer.After DOY 250, the incident PAR of the overstory layer was steadier than that in the summer period (before DOY 250) because there were fewer cloudy days; thus, the incident PAR of the grass layer was approximately 15% that of the overstory layer.
The vertical forest structure parameters also influence the amount of understory incident PAR.Table 3 indicates that the total amount of incident PAR in the understory layer was underestimated by 2% to 3% at these two sites.
such as those at the DHS site, the overstory LAI did not intensely change throughout the year.Figure 3B shows that the incident PAR of the shrub layer at the DHS site was 40% that of the tree layer, and the incident PAR of the grass layer was 20% that of the tree layer.After DOY 250, the incident PAR of the overstory layer was steadier than that in the summer period (before DOY 250) because there were fewer cloudy days; thus, the incident PAR of the grass layer was approximately 15% that of the overstory layer.The vertical forest structure parameters also influence the amount of understory incident PAR.Table 3 indicates that the total amount of incident PAR in the understory layer was underestimated by 2% to 3% at these two sites.

Influence of Vertical Stratification on GPP Estimation
Two main factors determine the estimation of GPP in LUE models: (1) the influences of temperature, VPD, and LUEmax on LUE, and (2) the effects of PAR and LAI on the absorbed shortwave energy in different components of the whole ecosystem.In this study, the major difference between considering and not considering the vertical structure during the GPP estimation is the difference between the overstory and understory LUE and APAR.Vertical variations in temperature, VPD, and the LUEmax of sub-ecosystems led to differences in LUE at various forest canopy heights.Some studies have shown that the vegetation in the understory layer such as grass has a lower LUEmax.In addition, the vertical differences in temperature and VPD led to variations during different periods of the year.Figure 4 shows the modeled LUE, which combined the vertical variations in temperature and VPD in the sub-ecosystems.At the DHS site, high layers in the canopy, such as the tree and shrub layers, had LUEs that were 0.3 gC/MJ to 0.6 gC/MJ lower than those of understory layers, such as grass layers, throughout the study period.Prior to DOY 320, the LUE of the shrub layer was 0.1 gC/MJ lower than that of the tree layer.After DOY 330, the LUEs of the shrub layer and the tree layer were almost the same.Similarly, the understory LUE at TMK was approximately 0.4 gC/MJ higher than that of the overstory layer throughout the growing season (DOY 210 to 300).

Influence of Vertical Stratification on GPP Estimation
Two main factors determine the estimation of GPP in LUE models: (1) the influences of temperature, VPD, and LUEmax on LUE, and (2) the effects of PAR and LAI on the absorbed shortwave energy in different components of the whole ecosystem.In this study, the major difference between considering and not considering the vertical structure during the GPP estimation is the difference between the overstory and understory LUE and APAR.Vertical variations in temperature, VPD, and the LUEmax of sub-ecosystems led to differences in LUE at various forest canopy heights.Some studies have shown that the vegetation in the understory layer such as grass has a lower LUEmax.In addition, the vertical differences in temperature and VPD led to variations during different periods of the year.Figure 4 shows the modeled LUE, which combined the vertical variations in temperature and VPD in the subecosystems.At the DHS site, high layers in the canopy, such as the tree and shrub layers, had LUEs that were 0.3 gC/MJ to 0.6 gC/MJ lower than those of understory layers, such as grass layers, throughout the study period.Prior to DOY 320, the LUE of the shrub layer was 0.1 gC/MJ lower than that of the tree layer.After DOY 330, the LUEs of the shrub layer and the tree layer were almost the same.Similarly, the understory LUE at TMK was approximately 0.4 gC/MJ higher than that of the overstory layer throughout the growing season (DOY 210 to 300).The mean canopy air temperature influenced the estimation of GPP.The trend in GPPdif in Figure 5 shows that with decreasing canopy mean temperature, the estimation based on GPPm becomes much lower than that based on GPPs, and with increasing mean temperature, the estimation based on GPPm becomes larger than that based on GPPs.The estimates of the three models are similar: GPPdif reaches 0 at approximately 16 °C.The optimal growing temperatures of different The mean canopy air temperature influenced the estimation of GPP.The trend in GPPdif in Figure 5 shows that with decreasing canopy mean temperature, the estimation based on GPPm becomes much lower than that based on GPPs, and with increasing mean temperature, the estimation based on GPPm becomes larger than that based on GPPs.The estimates of the three models are similar: GPPdif reaches 0 at approximately 16 • C. The optimal growing temperatures of different species in one sub-ecosystem, such as the tree layer, would not exhibit much variance.When the mean canopy temperature reaches approximately 16 • C, the grass ecosystem reaches its optimal growing temperature, and the LUE of grass would be 10% higher than that at 10 • C or lower.At the same time, the LAI of the understory increases with increasing temperature, and the grass understory has a higher photosynthetic ability.This situation causes the value of GPPm to be higher than that of GPPs.The canopy vertical VPD stratification has a pattern that is similar to that of the canopy vertical temperature stratification: there are variations in photosynthesis in the vertical layers that lead to GPP estimation bias.Figure 6 shows that the GPPdif values are completely different at the DHS and TMK sites.The trend of the change in GPPdif in Figure 6B shows that the estimation from GPPm becomes much lower than that from GPPs as the mean canopy VPD decreases, and with an increase in the mean VPD, the estimation from GPPm becomes larger than that from GPPs.At the DHS site, the GPPdif varied nonsignificantly from approximately −5% to 5%.GPPdif declines with increasing canopy mean VPD, which is opposite the pattern observed at the TMK site, where GPPdif increases with VPD, because the VPD increases from 0 to a higher value at the DHS site.When the VPD is greater than 0, the photosynthetic capacities of different vegetation components, including the tree, shrub, and grass layers decrease, and the total GPP declines.GPPdif changes from −50% to 30% when the mean canopy VPD increases from approximately 0 Pa to 500 Pa at the TMK site.The canopy vertical VPD stratification has a pattern that is similar to that of the canopy vertical temperature stratification: there are variations in photosynthesis in the vertical layers that lead to GPP estimation bias.Figure 6 shows that the GPPdif values are completely different at the DHS and TMK sites.The trend of the change in GPPdif in Figure 6B shows that the estimation from GPPm becomes much lower than that from GPPs as the mean canopy VPD decreases, and with an increase in the mean VPD, the estimation from GPPm becomes larger than that from GPPs.At the DHS site, the GPPdif varied nonsignificantly from approximately −5% to 5%.GPPdif declines with increasing canopy mean VPD, which is opposite the pattern observed at the TMK site, where GPPdif increases with VPD, because the VPD increases from 0 to a higher value at the DHS site.When the VPD is greater than 0, the photosynthetic capacities of different vegetation components, including the tree, shrub, and grass layers decrease, and the total GPP declines.GPPdif changes from −50% to 30% when the mean canopy VPD increases from approximately 0 Pa to 500 Pa at the TMK site.
the GPPdif varied nonsignificantly from approximately −5% to 5%.GPPdif declines with increasing canopy mean VPD, which is opposite the pattern observed at the TMK site, where GPPdif increases with VPD, because the VPD increases from 0 to a higher value at the DHS site.When the VPD is greater than 0, the photosynthetic capacities of different vegetation components, including the tree, shrub, and grass layers decrease, and the total GPP declines.GPPdif changes from −50% to 30% when the mean canopy VPD increases from approximately 0 Pa to 500 Pa at the TMK site.

Vertical Stratification of APAR Influenced GPP Estimation by the Multilayer Model
Forest structure parameters such as LAI influenced the absorption of radiation by the forest canopy, so the total APAR of the forest canopy can be calculated using LAI and incident PAR.In this section, we analyze LAI and incident PAR as influenced by the APAR in different sub-ecosystems and the effect of separate APARs in the sub-ecosystems on GPP modeling.
In a natural environment, the variety of vertical LAI distributions leads to different absorbed PAR distributions in the tree, shrub, and grass layers, which have diverse photosynthetic capacities that lead to GPP modeling biases.Figure 7 shows the total GPP in different LAI distributions.We made two assumptions.The first assumption is that the LUE on a daily scale in different layers is the mean value, which is 1.05 gC/m 2 /MJ in the tree layer, 0.83 gC/m 2 /MJ in the shrub layer, and 1.47 gC/m 2 /MJ in the grass layer.These values are the average values of the results of three models considering optimal growing environments.The APAR for each layer is influenced by only the LAI and the total incident PAR in the top layer.The second assumption is that the LUEmax in the shrub layer is lower than that in the tree and grass layers [32].Similar to the results at the DHS site, when the LAI ratio in the tree, shrub, and grass layers is 2.5:0.6:0.9 (at DOY 249), GPPm is 5% higher than GPPs.At the TMK site during the summer, the LAI ratio in the tree and grass layers is 2:1, and GPPm is 15% higher than GPPs.As in the shrub and grass layers, the understory LAI would thus be the most important factor for GPPdif, and its ratio to total LAI leads to a range of −40% to 25% for the GPP estimation uncertainties.
Incident PAR in the forest canopy also determined the GPP modeling difference whether or not vertical stratification was considered.Figure 8A,C show that considering the effect of vertical stratification on GPP leads to results that are similar to the GPP values that are measured by the flux towers because understory species such as shrubs and grasses have high photosynthetic capacities despite receiving little incident PAR [82].When the incident PAR of the top of the canopy is lower than 175 W/m 2 at these two sites, the value of the flux tower-based GPP is 1-2 gC/m 2 lower than GPP, and is closer to GPPm than GPPs.When the incident PAR is higher than 175 W/m 2 , GPPm is 0.5-2 gC/m 2 higher than the flux tower-based GPP, which has a lower bias than GPPs.
mean value, which is 1.05 gC/m 2 /MJ in the tree layer, 0.83 gC/m 2 /MJ in the shrub layer, and 1.47 gC/m 2 /MJ in the grass layer.These values are the average values of the results of three models considering optimal growing environments.The APAR for each layer is influenced by only the LAI and the total incident PAR in the top layer.The second assumption is that the LUEmax in the shrub layer is lower than that in the tree and grass layers [32].Similar to the results at the DHS site, when the LAI ratio in the tree, shrub, and grass layers is 2.5:0.6:0.9 (at DOY 249), GPPm is 5% higher than GPPs.At the TMK site during the summer, the LAI ratio in the tree and grass layers is 2:1, and GPPm is 15% higher than GPPs.As in the shrub and grass layers, the understory LAI would thus be the most important factor for GPPdif, and its ratio to total LAI leads to a range of −40% to 25% for the GPP estimation uncertainties.Ternary contour for the differences between considering and not considering the vertical structure during the GPP estimation with LAI separated into tree, shrub, and grass layers.This graph illustrates that when the total LAI of the canopy is four, there is little temperature and humidity stress on the LUE.We assumed that the LUE of the grass layer was 1.28 gC/m 2 /day/MJ, that of the shrub layer was 0.85 gC/m 2 /day/MJ, and that of the tree layer was 1.02 gC/m 2 /day/MJ.
Incident PAR in the forest canopy also determined the GPP modeling difference whether or not vertical stratification was considered.Figure 8A,C show that considering the effect of vertical stratification on GPP leads to results that are similar to the GPP values that are measured by the flux towers because understory species such as shrubs and grasses have high photosynthetic capacities despite receiving little incident PAR [82].When the incident PAR of the top of the canopy is lower Figure 7. Ternary contour for the differences between considering and not considering the vertical structure during the GPP estimation with LAI separated into tree, shrub, and grass layers.This graph illustrates that when the total LAI of the canopy is four, there is little temperature and humidity stress on the LUE.We assumed that the LUE of the grass layer was 1.28 gC/m 2 /day/MJ, that of the shrub layer was 0.85 gC/m 2 /day/MJ, and that of the tree layer was 1.02 gC/m 2 /day/MJ.
Remote Sens. 2018, 10, x FOR PEER REVIEW 14 of 25 than 175 W/m 2 at these two sites, the value of the flux tower-based GPP is 1-2 gC/m 2 lower than GPP, and is closer to GPPm than GPPs.When the incident PAR is higher than 175 W/m 2 , GPPm is 0.5-2 gC/m 2 higher than the flux tower-based GPP, which has a lower bias than GPPs.Moreover, at high temperatures at the DHS site, the estimation of GPPm is lower than the estimation of GPPs, and the difference between GPPm and GPPs increases with increasing incident PAR (Figure 8B,D).At the TMK site, the difference between GPPm and GPPs increases with increasing PAR, and reaches a maximum of approximately 70% of GPPs, which is much higher than that at the DHS site.When the incident PAR increases, the PAR reaching the understory also increases.Understory layers have higher LUE values, and GPPm is capable of considering the higher photosynthesis from understory components, which explains why GPPm is higher than GPPs.At the TMK site, GPPm is lower than GPPs when the temperature is low, which is the opposite pattern of that at the DHS site.The high temperatures at the DHS site and the low temperatures at the TMK site represent non-optimal growing temperature conditions, under which the biases between GPPm and Moreover, at high temperatures at the DHS site, the estimation of GPPm is lower than the estimation of GPPs, and the difference between GPPm and GPPs increases with increasing incident PAR (Figure 8B,D).At the TMK site, the difference between GPPm and GPPs increases with increasing PAR, and reaches a maximum of approximately 70% of GPPs, which is much higher than that at the DHS site.
When the incident PAR increases, the PAR reaching the understory also increases.Understory layers have higher LUE values, and GPPm is capable of considering the higher photosynthesis from understory components, which explains why GPPm is higher than GPPs.At the TMK site, GPPm is lower than GPPs when the temperature is low, which is the opposite pattern of that at the DHS site.The high temperatures at the DHS site and the low temperatures at the TMK site represent non-optimal growing temperature conditions, under which the biases between GPPm and flux tower-measured GPP are always lower than those between GPPs and flux tower-measured GPP with different degrees of variation at different sites.When the temperature enters the optimal growing temperature range, the estimation of GPPm is higher than that of GPPs, which considers the high photosynthesis capacity of the understory.

Performance of GPP Estimation Results Based on Flux Tower-Based GPP
Estimating GPP with stratification information in the multilayer GPP model results in higher estimation precision than with the single-layer method (Table 4).Table 4 illustrates that the GPPm estimates at the DHS and TMK sites had higher accuracies than the GPPs estimations, and generally had a lower RMSE (GPPm is 0.2 gC/m 2 /day lower than GPPs).At the DHS site, the modeled GPP that considered vertical stratification showed a similar correlation when vertical layers were considered, but with the higher LAI of the understory, the R 2 value was 2% to 5% higher with GPPm than with GPPs at the TMK site.The low mean bias in the multilayer model (DHS_m and TMK_m) indicates that the multilayer model performs better than the single-layer model.These selected GPP estimation models, including VPM, Wu's model, and MOD, showed high correlations with the flux tower-based GPP.Wu's model and the VPM showed higher R 2 values than the MOD model at DHS, because they consider scalar temperature at high temperatures.The MOD17 model showed lower RMSE and higher R 2 values than both Wu's model and the VPM at the TMK site, indicating that the MOD model exhibits better performance at that site.The estimation results for the two sites were combined, which indicated that Wu's model exhibited the best model performance with the highest mean R 2 (R 2 = 0.685) and the lowest mean RMSE at these two sites (RMSE = 1.66 gC/m 2 /day) when the vertical stratification information was included in the GPP estimation.Figure 9 indicates that Wu's model can track the variations in growing season GPP when the vertical stratification information in the canopy is incorporated.
site, indicating that the MOD model exhibits better performance at that site.The estimation results for the two sites were combined, which indicated that Wu's model exhibited the best model performance with the highest mean R 2 (R 2 = 0.685) and the lowest mean RMSE at these two sites (RMSE = 1.66 gC/m 2 /day) when the vertical stratification information was included in the GPP estimation.Figure 9 indicates that Wu's model can track the variations in growing season GPP when the vertical stratification information in the canopy is incorporated.Here, we used three different types of LUEmax parametric approaches to model GPP using Wu's method.GPPmW is the multilayer result using two or three different LUEmax inputs for GPP estimation.GPPsW is the single-layer result using an EBF-based LUEmax.GPPaW is the average LUEmax input in GPPmW, using the average LUEmax for the different components.
Simplifying the vertical stratification parameters in the GPP model will also yield high estimation accuracies.To simplify the parameters when considering vertical stratification during GPP estimation, the adjusted method in Table 4 used the same temperature and VPD for every layer as inputs, the LAI values in different sub-ecosystems were separated from that for the total ecosystem, and LUEmax was used for different sub-ecosystems.The results showed that among all the methods, the adjusted method had the smallest mean bias at the TMK and DHS sites.Additionally, the R 2 value differed depending on whether or not multilayer models were considered.This method requires the canopy mean temperature and the VPD and LAI proportions in the understory and overstory, the LUEmax set for different sub-ecosystems, and the incident PAR at the top of the canopy.The data used here are easier to obtain than those required by the multilayer GPP estimation method, and the performance results were also better than those when not considering vertical stratification.
Modeled GPP and flux tower-based GPP values still differ to some degree, as illustrated in Figure 10.On a daily scale, the multilayer GPP model has a lower mean bias than the single-layer GPP model.Figure 10b demonstrates that the multilayer method performs well during the growing season (DOY 180 to 280) when the understory LAI is large.On some days when the shortwave incident radiation is high, GPPm was overestimated because the understory PAR was overestimated, but GPPm showed little variance when the understory LAI was low (Figure 10a,b DOY 290 to 360).Most of the results estimated by GPPs were lower than the GPP measured by the flux towers because the understory LUEmax was underestimated.Therefore, the GPP estimation method that considered vertical structure would be useful when the understory LAI is high (>1 m 2 /m 2 ).Generally, adding the understory vegetation GPP to the total GPP reduces the underestimation when the understory LAI Figure 9. Differences in GPP estimations in different seasons.Subfigure (A) showed the seasonal GPP variance in DHS, while subfigure (B) showed the condition in TMK.Here, we used three different types of LUEmax parametric approaches to model GPP using Wu's method.GPPmW is the multilayer result using two or three different LUEmax inputs for GPP estimation.GPPsW is the single-layer result using an EBF-based LUEmax.GPPaW is the average LUEmax input in GPPmW, using the average LUEmax for the different components.
Simplifying the vertical stratification parameters in the GPP model will also yield high estimation accuracies.To simplify the parameters when considering vertical stratification during GPP estimation, the adjusted method in Table 4 used the same temperature and VPD for every layer as inputs, the LAI values in different sub-ecosystems were separated from that for the total ecosystem, and LUEmax was used for different sub-ecosystems.The results showed that among all the methods, the adjusted method had the smallest mean bias at the TMK and DHS sites.Additionally, the R 2 value differed depending on whether or not multilayer models were considered.This method requires the canopy mean temperature and the VPD and LAI proportions in the understory and overstory, the LUEmax set for different sub-ecosystems, and the incident PAR at the top of the canopy.The data used here are easier to obtain than those required by the multilayer GPP estimation method, and the performance results were also better than those when not considering vertical stratification.
Modeled GPP and flux tower-based GPP values still differ to some degree, as illustrated in Figure 10.On a daily scale, the multilayer GPP model has a lower mean bias than the single-layer GPP model.Figure 10b demonstrates that the multilayer method performs well during the growing season (DOY 180 to 280) when the understory LAI is large.On some days when the shortwave incident radiation is high, GPPm was overestimated because the understory PAR was overestimated, but GPPm showed little variance when the understory LAI was low (Figure 10a,b DOY 290 to 360).Most of the results estimated by GPPs were lower than the GPP measured by the flux towers because the understory LUEmax was underestimated.Therefore, the GPP estimation method that considered vertical structure would be useful when the understory LAI is high (>1 m 2 /m 2 ).Generally, adding the understory vegetation GPP to the total GPP reduces the underestimation when the understory LAI is high.

Discussion
In this research, we examined the addition of vertical stratification parameters such as temperature, VPD, PAR, LAI, and LUEmax to support remote sensing GPP estimation at the site level.Temperature and VPD influenced the LUEmax at various canopy heights because the vertical meteorology environments have greater variations at forest sites than those in other vegetation types, such as shrub and grass sites.PAR and LAI were used to determine the APAR of different subecosystems, which would lead to a great deal of uncertainty in energy allocation.
Vertical LUE is influenced by the vertical distributions of temperature, VPD, and LUEmax, but there are uncertainties when modeling the LUE in different sub-ecosystems.LUEmax values differ in sub-ecosystems, as explained in Section 4.2.1 of this article.At the TMK site, the GPPdif changes from −50% to 30% when the air temperature at the canopy height increases from approximately 7 °C

Discussion
In this research, we examined the addition of vertical stratification parameters such as temperature, VPD, PAR, LAI, and LUEmax to support remote sensing GPP estimation at the site level.Temperature and VPD influenced the LUEmax at various canopy heights because the vertical meteorology environments have greater variations at forest sites than those in other vegetation types, such as shrub and grass sites.PAR and LAI were used to determine the APAR of different sub-ecosystems, which would lead to a great deal of uncertainty in energy allocation.
Vertical LUE is influenced by the vertical distributions of temperature, VPD, and LUEmax, but there are uncertainties when modeling the LUE in different sub-ecosystems.LUEmax values differ in sub-ecosystems, as explained in Section 4.2.1 of this article.At the TMK site, the GPPdif changes from −50% to 30% when the air temperature at the canopy height increases from approximately 7 • C to 20 • C. The variation is more significant at the TMK site than that at the DHS site, and one reason for this difference is that the temperature differences among layers are more obvious at the TMK site than at the DHS site, as shown in Figure 2. The effect of temperature stress on LUE differs for tree, shrub, and grass layers under the same mean canopy weather conditions.As the optimal temperatures of the tree, shrub, and grass layers differ by 2-5 • C at the same air temperatures in the forest canopy, the temperature stresses of the layers on LUEmax differ, varying by 5-20% around the optimal temperature [78,83].However, when the temperature is above 23 • C, GPPdif changes only slightly because the canopy and the understory layers experience high temperature stress, and the LUE would decrease rapidly at this temperature range.On the other hand, the variation in GPPdif did not significantly change with the increase in VPD, as Figure 6B shows that GPPdif at the DHS site decreases only slightly when the VPD increased.Nepstad et al. [84] showed that the water contents of plants are more related to soil moisture than to the water content of the air in seasonally moist forests.The plants in the overstory and understory layers share the soil water content.Thus, the VPD does not have a significant effect on the GPP estimation in dense forests.All of these effects result in different LUE values at distinct canopy heights, as shown by previous research [85].
Uncertainties in the vertical LAI distribution lead to uncertainties in APAR in the sub-ecosystems.LAI, which is a key parameter for separating an ecosystem into different photosynthetic components, determines the contribution of different species to total ecosystem photosynthesis, the light transfer process in the canopy, and the micrometeorological environment inside the forest [58].As many studies have shown that the distribution of PAR through the canopy structure is highly affected by LAD, LAI, and the characteristics of different species in forest canopies with high heterogeneity, the PAR distribution must consider the vertical structure [39,58,82].The overstory tree layer determines the intercept ratio of the incident PAR, so if the LAI in the overstory is high, a large amount of PAR would be absorbed during the first turn of intercepted light.The vertical LAI distribution in different vegetation covers follows its own rules (Figure 11), and therefore, remote sensing-based understory and overstory LAI estimations are important for separating the APAR of sub-ecosystems and adding unique LUE values when mapping the global GPP [86][87][88][89].For all of the cases in a natural environment, when the total LAI is four, the highest GPP appears when all of the leaves are grass species.As shown in Figure 7, the LAI of the grass layer would not be higher than 2 m 2 /m 2 , and the highest GPP would thus appear when the LAI of both the tree and grass layers is two.The highest GPP in a dense natural forest would have a ratio of 0.5:0:0.5 for the tree, shrub, and grass layers.Under such circumstances, GPPdif reaches its highest value, and increases in the LUE of both the tree and grass layers would produce a greater GPP.The greatest GPP difference is 25% higher than that estimated via the single-layer method, and if the canopy has a higher proportion of shrubs, the GPP difference is 20% lower.When the proportion of LAI is high in the shrub layer, the total GPP is lower than that estimated by the single-layer method because the shrub LUE is low.As discussed in Section 4.2.2, a high LAI in the understory produces a high GPP because the understory component absorbs diffuse light and provides a better growing environment than that of the overstory.Figure 10 shows that different ecosystems have different vertical LAI ratios.In forest sites such as EBF, DBF, ENF, and DNF, the tree LAI ratio can be 60% or more, and the understory LAI can be as a great as 20% during the growing season.With different LAI distributions in the forest ecosystem, the total APAR allocation of the canopy in a sub-ecosystem would lead to uncertainties in the GPP estimation.[71,90]; DNF, deciduous needleleaf forest [72]; ENF, evergreen needleleaf forest [91,92]; DBF, deciduous broadleaf forest [26]; SAV, savanna [93]).
Wu's model exhibited the best performance of the three models.This model better considers the effect of temperature on GPP in overstory and understory layers.Wu's model and the VPM consider the temperature stresses when the temperature is high.Therefore, these models exhibit better variance when the canopy temperature is high in both the multilayer model and single-layer model.In the VPM, phenology was proven to be an important indicator.However, at the site level, the VPM cannot obtain accurate spatially matched enhanced vegetation index (EVI) values, and the LAI of the canopy changes only slightly during the growing season; thus, we omitted the effect of phenology.In future research on multilayer GPP models, the phenology variance between understory and overstory layers would be a potential way to improve the GPP estimation accuracy in the VPM and Wu's model.
The continuous vertical acquisition of LAI is also challenging.On the one hand, many researchers have shown that remote sensing methods underestimate understory LAI, indicating that other GPP components are also underestimated in dense forests [39,94].On the other hand, most remote sensing-based LAI products currently do not separate the overstory from the understory in forest sites.Without the LAI of different layers, the GPP of distinct vertical layers cannot be modeled; thus, understory GPP is often underestimated from the high understory LUE.
Future versions of remote sensing-based GPP products should add forest vertical stratification information to GPP models.For example, the use of remote sensing-based understory and overstory LAI [90][91][92][93] to separate sub-ecosystems and the addition of prior LUE when combining the specific LUEmax from sub-ecosystem-based temperature and VPD response functions can improve GPP estimations compared to the performance of approaches that do not consider the vertical structure ("adjusted" method in Table 4).The 'adjusted' method in Table 4 showed lower systemic biases and a RMSE that was similar to that of GPPm, which indicated that this method has an estimation accuracy that is similar to that of the multilayer method.However, this method requires only the two main data sources as prior information, the LUEmax in height of the canopy, the canopy mean air temperature, VPD, and the LAI values of the sub-ecosystems, which are simple to obtain via remote sensing.This method can potentially be applied in remote sensing-based multilayer GPP models.
Wu's model exhibited the best performance of the three models.This model better considers the effect of temperature on GPP in overstory and understory layers.Wu's model and the VPM consider the temperature stresses when the temperature is high.Therefore, these models exhibit better variance when the canopy temperature is high in both the multilayer model and single-layer model.In the VPM, phenology was proven to be an important indicator.However, at the site level, the VPM cannot obtain accurate spatially matched enhanced vegetation index (EVI) values, and the LAI of the canopy changes only slightly during the growing season; thus, we omitted the effect of phenology.In future research on multilayer GPP models, the phenology variance between understory and overstory layers would be a potential way to improve the GPP estimation accuracy in the VPM and Wu's model.
The continuous vertical acquisition of LAI is also challenging.On the one hand, many researchers have shown that remote sensing methods underestimate understory LAI, indicating that other GPP components are also underestimated in dense forests [39,94].On the other hand, most remote sensing-based LAI products currently do not separate the overstory from the understory in forest sites.Without the LAI of different layers, the GPP of distinct vertical layers cannot be modeled; thus, understory GPP is often underestimated from the high understory LUE.
Future versions of remote sensing-based GPP products should add forest vertical stratification information to GPP models.For example, the use of remote sensing-based understory and overstory LAI [90][91][92][93] to separate sub-ecosystems and the addition of prior LUE when combining the specific LUEmax from sub-ecosystem-based temperature and VPD response functions can improve GPP estimations compared to the performance of approaches that do not consider the vertical structure ("adjusted" method in Table 4).The 'adjusted' method in Table 4 showed lower systemic biases and a RMSE that was similar to that of GPPm, which indicated that this method has an estimation accuracy that is similar to that of the multilayer method.However, this method requires only the two main data sources as prior information, the LUEmax in height of the canopy, the canopy mean air temperature, VPD, and the LAI values of the sub-ecosystems, which are simple to obtain via remote sensing.This method can potentially be applied in remote sensing-based multilayer GPP models.Thus, a high GPP estimation accuracy could be obtained in ecosystems with high understory LAI values by considering the vertical structure (multilayer model).
The accuracy of GPP estimates in a complicated, vertically stratified forest canopy where the understory LAI is high (e.g., in a tropical rainforest or old-growth forest) can be improved by using multilayer methods.Additionally, Wu's model showed that a potential model could estimate the GPP with multilayer data inputs, but further global validation of the GPP model parameters is needed.Some articles have added vertical stratification data to GPP estimations in cropland [95,96], but unlike the homogeneous landscape in croplands, a vertical structure would lead to variations in the physiological parameters of vegetation in forest sites.Thus, there are still some limitations to this study.

Conclusions
In this article, we used three different LUE-based GPP models to analyze the effect of vertical vegetation stratification on GPP calculations.Additionally, in situ temperature, VPD, and PAR data were used to analyze the environmental factors that control GPP modeling.The results showed the following.(1) Vertical temperature and VPD stratification at the selected sites varied with the forest type, time of day, and season, but the difference was not significant.Moreover, incident PAR followed the Beer-Lambert law; there was a significant difference between the overstory and understory layers, and this difference was affected by the distributions of LAI in different sub-ecosystems.
(2) The stratifications of temperature, VPD, PAR, LAI, and LUEmax have an interactive influence on the improvement of GPPm.Due to the importance of PAR in the estimation of GPP, the stratification of PAR becomes obvious when the incident PAR increases and the difference between GPPm and GPPs exceeds 50% of GPPs.GPPdif is approximately 30% when temperature stratification is considered, and this value is 20% when LAI stratification is considered.Multilayer LUEmax inputs are important for accurate GPP estimations, and different inputs may lead to contrary results concerning the difference between GPPm and GPPs.Therefore, quantitatively corrected understory LUEmax reduces the bias in GPP estimations.(3) The accuracy of the GPP estimation can be improved by using multilayer models, but the improvement is related to the forest type, the stratification of the canopy, and the growing conditions.Furthermore, the improvement is more obvious when the multilayer model is applied to canopies with obvious vertical stratification under optimal growing conditions.In this article, GPPm was the most accurate at the TMK site from July to September.Since a multilayer GPP estimation model considers the high photosynthetic ability and LAI of grass, GPPm is always larger than GPPs.However, when the conditions are not optimal for growth-for example, if the temperature becomes very low or high-the GPPm will be smaller than the GPPs.
In future studies on GPP remote sensing products, vertical stratification should be considered in the models.Adding overstory and understory information such as LAI and LUEmax for diverse sub-ecosystems can improve the estimation accuracy of the GPP product.

Figure 1 .
Figure 1.Locations of the two study sites in East Asia.

Figure 1 .
Figure 1.Locations of the two study sites in East Asia.

Figure 2 .
Figure 2. Seasonal patterns in the vertical variance in temperature and vapor pressure deficit (VPD) at the DHS and TMK sites.In the subplots (A,C), the blue bars (VPDdif in Pa) are the biases between

Figure 2 .
Figure 2. Seasonal patterns in the vertical variance in temperature and vapor pressure deficit (VPD) at the DHS and TMK sites.In the subplots (A,C), the blue bars (VPDdif in Pa) are the biases between the tree layer VPD and the grass layer VPD.The red line (Tdif in • C) is the bias between the air temperatures of the top and grass layers.The green dotted line is the 0 • C line, indicating no difference between the temperatures of the tree and grass layers.Subplots (B,D) show the relationships between Tdif and VPDdif at the two sites.The dot colors are clustered by month.

Figure 3 .
Figure 3. Bar chart for the daily average incident PAR in sub-ecosystems at the selected sites; each bar represents the daily mean incident PAR in each layer (A) TMK; (B) DHS.

Figure 3 .
Figure 3. Bar chart for the daily average incident PAR in sub-ecosystems at the selected sites; each bar represents the daily mean incident PAR in each layer (A) TMK; (B) DHS.

4. 2 . 1 .
LUE in Different Sub-Ecosystems Is Influenced by the Vertical Stratification of Temperature, VPD, and LUEmax

4. 2 . 1 .
LUE in Different Sub-Ecosystems Is Influenced by the Vertical Stratification of Temperature, VPD, and LUEmax

Figure 4 .
Figure 4. Modeled LUE values in the different sub-ecosystems.

Figure 4 .
Figure 4. Modeled LUE values in the different sub-ecosystems.

Figure 5 .
Figure 5. Differences in gross primary production (GPP) under different mean canopy temperatures.The GPP difference is calculated as (GPPm-GPPs)/GPPs.DHS_VPM, DHS_W, and DHS_MOD are the values calculated at the DHS site and by VPM, Wu's method, and the MOD17 model, respectively.

Figure 5 .
Figure 5. Differences in gross primary production (GPP) under different mean canopy temperatures.The GPP difference is calculated as (GPPm-GPPs)/GPPs.DHS_VPM, DHS_W, and DHS_MOD are the values calculated at the DHS site and by VPM, Wu's method, and the MOD17 model, respectively.

Figure 6 .
Figure 6.Effects of both considering and not considering the vertical differences in VPD on GPP.The black line is the fit line for the GPP difference scatterplot at different sites.Subfigures (A,B) are the GPP difference comparisons based on different canopy mean VPDs at TMK and DHS, respectively.

Figure 6 .
Figure 6.Effects of both considering and not considering the vertical differences in VPD on GPP.The black line is the fit line for the GPP difference scatterplot at different sites.Subfigures (A,B) are the GPP difference comparisons based on different canopy mean VPDs at TMK and DHS, respectively.

Figure 7 .
Figure 7. Ternary contour for the differences between considering and not considering the vertical structure during the GPP estimation with LAI separated into tree, shrub, and grass layers.This graph illustrates that when the total LAI of the canopy is four, there is little temperature and humidity stress on the LUE.We assumed that the LUE of the grass layer was 1.28 gC/m 2 /day/MJ, that of the shrub layer was 0.85 gC/m 2 /day/MJ, and that of the tree layer was 1.02 gC/m 2 /day/MJ.

Figure 8 .
Figure 8. GPP differences considering multilayer and single-layer estimations of incident PAR.Subfigure (A,C) show the GPP simulation results with flux tower measurements of GPP.The green points represent multilayer results (GPPmMOD), and the red points represent single-layer results (GPPsMOD).The black points represent the flux tower-based GPP results (GPPactual).Subfigure (B,D) show GPP differences under various incident PAR.

Figure 8 .
Figure 8. GPP differences considering multilayer and single-layer estimations of incident PAR.Subfigure (A,C) show the GPP simulation results with flux tower measurements of GPP.The green points represent multilayer results (GPPmMOD), and the red points represent single-layer results (GPPsMOD).The black points represent the flux tower-based GPP results (GPPactual).Subfigure (B,D) show GPP differences under various incident PAR.

Figure 9 .
Figure 9. Differences in GPP estimations in different seasons.Subfigure (A) showed the seasonal GPP variance in DHS, while subfigure (B) showed the condition in TMK.Here, we used three different types of LUEmax parametric approaches to model GPP using Wu's method.GPPmW is the multilayer result using two or three different LUEmax inputs for GPP estimation.GPPsW is the single-layer result using an EBF-based LUEmax.GPPaW is the average LUEmax input in GPPmW, using the average LUEmax for the different components.

Figure 10 .
Figure 10.Differences in GPP estimation using different models with carbon flux tower measurements.Subfigure (a) showed the condition in DHS, and subfigure (b) showed it in TMK.VPMm results from using the multilayer VPM GPP minus the flux tower-based GPP.VPMs results from using the single-layer VPM GPP minus the flux tower-based GPP.The Y-axis on the right is the daily mean shortwave incident radiation.

Figure 10 .
Figure 10.Differences in GPP estimation using different models with carbon flux tower measurements.Subfigure (a) showed the condition in DHS, and subfigure (b) showed it in TMK.VPMm results from using the multilayer VPM GPP minus the flux tower-based GPP.VPMs results from using the single-layer VPM GPP minus the flux tower-based GPP.The Y-axis on the right is the daily mean shortwave incident radiation.

Table 3 .
Understory photosynthetically active radiation (PAR) simulation: parameter value inputs and results.

Table 4 .
Comparison of flux tower-based GPP and different model output results when considering and not considering the vertical layer.Here, sitename_s is the single-layer GPP, sitename_m is the multilayer GPPm, and sitename_adjust uses the same temperature and VPD inputs (root mean square error (RMSE) and mean bias unit is gC/m 2 /day).