Understanding the role of phenology and summer physiology in controlling net ecosystem production: a multiscale comparison of satellite, PhenoCam and eddy covariance data

Understanding the temporal and spatial variability (SV) of net ecosystem productivity (NEP) is critical for coupling ecosystem carbon (C) cycle and climate system. Previous studies have shown responses of NEP to changes of plant phenology, but impacts of summer physiological status on annual NEP and how this may vary across different ecosystems and spatial scales were largely unknown. Combining large regional satellite derived indices (MODIS), 676 site-year local data (FLUXDATA) covering seven vegetation types, and 57 site-year regional data (PhenoCam), we found that phenological metrics and summer physiological indicators were significantly correlated with their respective gross primary production-based estimates. The interannual variability of NEP was mostly explained by summer physiology than phenology for most ecosystems, while phenology showed a better performance than summer physiology only for grassland sites. In comparison, inconsistent results were derived from three scales when explaining the SV of NEP. Summer physiology showed more potential in interpreting the SV of NEP at large regional scale, but both physiology and the length of growing season exhibited similar performances at local scale. Observations from regional scale were not able to explain NEP, given that the green chromatic coordinate signal cannot track photosynthesis in summer with a high canopy closure. The finding highlights the important role of summer physiology in controlling the C accumulation of terrestrial ecosystems and understanding the responses of summer physiology to environmental drivers is of great significance for improving the simulation accuracy of C sequestration under the global climate change.


Introduction
Terrestrial ecosystems sequestered nearly 30% of anthropogenic carbon dioxide emissions (Ballantyne Vegetation phenology describes the recurring seasonal changes of plant development status in response to environmental changes (Lieth 1974), and includes key events such as sprouting, leaf development, flowering, and defoliation (Ivits et al 2012). Previous studies emphasized the impact of spring and autumn phenology on C uptake, and found the SV of NEP can be explained by growing season length (GSL) but sensitivities vary among ecosystems (Churkina et al 2005, Wu et al 2012. Nevertheless, the relationship between GSL and NEP is not always significant and predictable (Richardson et al 2013). Due to the uncertainty of the impact of climate-driven shifts in phenology on NEP, taking summer physiological status into consideration, when drought always happens, is consequently of great significance to understand the process of C sequestration. Summer physiological status (e.g. leaf pigmentation, leaf mass per unit area, water and nutrient content, etc (Ma et al 2011), which could be indicated by vegetation index or production (Xia et al 2015, Zhou et al 2016, reflect the greatest potential of canopy photosynthesis, thereby affecting carbon accumulation (Hilker et al 2011, Keenan et al 2014. Previous studies showed that the maximum gross primary production (GPP) (GPP max ) has the potential to estimate the IAV of GPP (Xia et al 2015), and found the better performance of GPP max than phenological metric (i.e. the start and the end of growing season; SOS and EOS) in forest ecosystems (Zhou et al 2016, Xu et al 2019. However, the comparison of photosynthetic phenology and physiology in controlling NEP has not been explored across vegetation types, and for different scales.
In general, the dynamic of canopy development can be quantitatively tracked at three spatial scales, including large regional satellite-based observations (Zhang et  Thus, using multiscale data are a critical need for a more profound and full-scale understanding of phenology and physiology in explaining the variability of NEP. Using 67 FLUXNET sites (676 site-year in total) flux data across 7 vegetation types, 16 PhenoCam sites (57 site-year in total), as well as MODIS imageries, our specific objectives were (1) to compare the phenological metrics (the start, peak, end, and length of growing season; SOS, POS, EOS and GSL) and physiological metrics (the maximum of growing season, max) obtained from large region satellite platforms, regional near-surface and local eddy covariance, (2) to investigate the potential of phenological and physiological metrics in interpreting IAV and SV of NEP from multi-spatial scales.

Study sites
We selected 67 sites from the FLUXNET community, including 15 deciduous broadleaf forests (DBF), 18 evergreen needleleaf forests (ENF), 7 mixed forests (MF), 9 croplands (CRO), 12 grasslands (GRA), 2 woody savannas (WSA) and 4 wetlands (WET), considering the data quantity and representativeness (figure 1). All sites had at least 3 years of concurrent high-quality C mass fluxes and MODIS data. Also, we used PhenoCam data for flux data that were also available, and only Type I PhenoCam sites were considered to guarantee the high data quality . Sixteen PhenoCam sites were selected finally and details of these sites are shown in supplementary  Table 2). For each study site, all VIs were extracted from 3 × 3 MODIS pixels (1.5 km × 1.5 km) centered on flux tower, which is confirmed to be reliable considering both footprints and land cover (Chen et al 2011).

PhenoCam data
To obtain the near-surface canopy development, a subset of PhenoCam Dataset V2.0 was used in this study (https://daac.ornl.gov/cgi-bin/ dsviewer.pl?ds_id=1674). Green chromatic coordinate (GCC) was used to tracking seasonal change of vegetation in the region-of-interest (Richardson 2019). To reduce unwanted variability in GCC time series (equation (1)), the 90th percentile of GCC values obtained in a 3-day moving window was defined to the central day (Sonnentag et al 2012).
where DN R , DN G and DN B represent the red, green and blue digital number (DN), respectively.

Flux data
Flux data of all study sites were obtained from the FLUXNET2015 Dataset (https://fluxnet.fluxdata.org/ data/fluxnet2015-dataset/). Several improvements were applied to data quality control protocols and the data processing pipeline (Pastorello et al 2014). NEP, which represents the total amount of organic Figure 1. The spatial distribution of study sites in this study. The circle size represents the available data years when MODIS and FLUXNET data were included simultaneously.

Phenological and physiological metrics of multiple scales
The comparison of multiscale data time series of seven plant function types at seven representative FLUXNET sites were shown in Supplementary figure 1. The phenological and physiological metrics of the three scales were listed in table 1. The specific extraction methods of all phenological and physiological metrics at three scales were provided in Supplement, and the schematic diagram of extraction methods is shown in figure 2. EVI, OSAVI, GCC , gross primary productivity, respectively. SOS, POS, EOS, GSL and max represent the start of growing season, the peak of growing season, the end of growing season, the length of growing season and the maximum of growing season.

Analyzing the role of phenological and physiological metrics in controlling the variability of NEP
The yearly anomalies between NEP and both phenological and physiological variables at three scales were calculated to evaluate their potential in explaining the IAV of NEP across sites (Richardson et al 2010).
We further used the analysis of variance technique (Duncan's multiple range test, p = 0.05) to analyze whether the reduction in residual sum of squares between these correlations are statistically significant or not. We also compared the ability of phenological and physiological variables in explaining the SV of NEP across all sites. The mean values and standard errors of all variables were calculated at each site. Then, we used a simple linear regression to analyze the relationship between NEP and these variables across all sites in three scales. The correlation was evaluated by the Pearson correlation coefficient (R) and a p-value threshold of 0.05 (R ≥ 0.6, 0.6 > R ≥ 0.4 and R < 0.4 are defined as 'high ' , 'medium' and 'low' correlation, respectively). In addition to the overall correlation, we also calculated the correlation separately for seven vegetation types, considering the C sequestration characteristics of different vegetation types.

Results
3.1. The multiscale comparison of physiological and phenological metrics 3.1.1. Phenology and physiology at large regional (satellite) and local (GPP) scales. All three satellite-based physiology and phenology were significantly correlated with GPP-based estimates, and phenological metrics (SOS, POS, EOS and GSL) showed relatively stronger relationships with GPP-based estimates than physiological metric (VI max ) (figure 3). Among four phenological metrics, satellite-based VIs exhibited relatively weak ability in estimating GPP-based GSL than that of SOS, POS and EOS for all sites (figure 3(d)). The satellite-based VIs had overall comparable performance in predicting GPP-based SOS, POS and EOS for all sites (figures 3(a), (b) and (c)). Furthermore, VIs differed in their ability in tracking GPP-based phenological metrics. Three VIs showed similar performances in predicting GPPbased SOS (figure 3(a)). EVI-derived POS was more tightly correlated with GPP-based POS than NDVI and OSAVI, with high R of 0.68 (p < 0.001) for all sites (figure 3(b)). However, NDVI and OSAVI performed better in estimating GPP-based EOS and GSL than EVI (figure 3(c), (d)).
The potential of satellite-based VIs differed substantially in tracking GPP-based phenology across vegetation types (Supplementary Table 3). Among seven vegetation types, satellite-based phenology was more tightly linked to GPP-based phenology at DBF, GRA, WSA and WET sites. The overall lowest correlation appeared at ENF sites. In addition, satellite-based VIs showed limited potential in estimating EOS and GSL at CRO sites (R = 0.29, p < 0.001 and NS), and POS at MF sites (NS).
Compared with phenology, three satellite-based VIs showed relatively limited ability in indicating GPP-based physiology (figure 3(e)). These potentials also varied among vegetation types that better results were found in grassland ecosystems (GRA and WSA) than in forest (DBF, ENF and MF) and cropland ecosystems, while no significant correlation was observed at WET sites (Supplementary Table 3).

Phenology and physiology at regional (PhenoCam) and local (GPP) scales
GCC-derived phenology was more tightly link to GPP-based estimates than physiology, and all GCC-derived phenological metrics (SOS, POS, EOS and GSL) were significantly correlated with GPPbased estimates (figure 4). Among them, GCCderived EOS and SOS performed relatively better than POS and GSL, with high R of 0.87 (p < 0.001) and 0.79 (p < 0.001), respectively (figure 4(a), (c)). However, no significant correlation was observed between GCC max and GPP max at all PhenoCam sites (figure 4(e)).

Explaining the IAV of NEP using multi-scales physiological and phenological metrics
At large regional scale, all satellite derived VIs yielded roughly similar results (figure 5, Supplementary Table 4). We found that three VI max showed more potential in explaining the IAV of NEP than phenology (R = 0.23, p < 0.001), and SOS and GSL were more tightly correlated with the IAV of NEP than POS and EOS (figure 5). The significant superiority of VI max also illustrated by the Duncan's multiple range test (figure 7). Among different vegetation types, physiology also performed much better than phenology at DBF, ENF, CRO, WSA sites. However, phenological metrics, especially GSL, showed higher ability at GRA and MF site. Furthermore, the spring and summer phenology were better than autumn phenology in interpreting the IAV of NEP at GRA sites. Interestingly, the GSL of MF was negative correlated with the IAV of NEP, which was contrary to that of other vegetation types.
At local scale, physiology also exhibited more potential in interpreting the IAV of NEP than phenological metrics, with a medium R of 0.45 (p < 0.001) (figure 6). The better performance of physiology also was illustrated by inter-comparing statistical test (figure 7). The significant GPP-derived SOS and GSL showed relatively stronger relationships with the IAV of NEP than POS and EOS. Among different vegetation types, GPP max exhibited a higher correlation than any other GPP-derived phenological metrics in interpreting the IAV of NEP for DBF, ENF, MF, CRO and WET ecosystems (Supplementary Table 5). However, phenology had a better performance than physiology only for GRA sites, and SOS, POS and GSL were significantly correlated with the IAV of NEP. In addition, although weaker than summer physiology, phenology still made a nonnegligible contribution to regulate the IAV of NEP at CRO sites, and SOS, POS and GSL were also powerful indicators of NEP than EOS.

Explaining the SV of NEP using multi-scales physiological and phenological metrics
In contrast with the dominant role of summer physiology in controlling the IAV of NEP, phenology seemed to play an equally important role in affecting the SV of NEP. At large regional scale, physiology deriving from all VIs exhibited slightly higher potential than phenology ( figure 8(e)). The GSL also played an important role in interpreting the SV of NEP, with a slightly lower R of 0.33 (p < 0.01) for all three VIs ( figure 8(c)). At regional scale, GCC-derived phenology performed better than physiology in controlling SV of NEP, and GSL was the only metric that was able to provide marginally significant correlation with the SV of NEP (R = 0.47, p = 0.06) (figure 9). At local scale, both phenology and physiology exhibited similar relationships with the SV of NEP with equaled R of 0.40 (p < 0.001) (figure 10(d), (e)). However, no significantly correlation was observed at POS.

The match and mismatch of phenological and physiological metrics across multiple scales
We found that phenological metrics derived from large regional and regional scales (satellite and PhenoCam) tended to be more consistent with respective to local scale (GPP) than that for physiology (figures 3 and 4).
The phenological metrics (SOS and EOS), usually occur in the periods of relatively low vegetation biomass, when VIs and GCC are particularly well-suited to tracking the canopy greenness and may provide an explanation for the better consistency of phenology across the scales. However, differences also existed due to the mismatch of the rate of leaf expansion and the extent of photosynthetic ability ( The discrepancies of physiological metricsat three scales can be explained by the limitation of VIs. The VIs have the ability to indicate the photosynthetic capacity to a certain extent, owing to the responsiveness to canopy structural variations and biophysical properties (Rondeaux et al 1996). However, satellitebased VIs tends to be less sensitive when canopy grows to high closure (D'Odorico et al 2015). EVI was expected to be more effectiveness than NDVI in high canopy closure, which was supported by the higher correlation between GPP max and EVI max than NDVI max in our study (figure 3(e1), (e2)). Moreover, OSAVI max was more tightly link to GPP max than EVI max and NDVI max (figure 3(e3)), implying its underestimated ability in tracking photosynthesis at high vegetation biomass. In addition, no significant correlation was observed between GCC max and GPP max (figure 4(e)), because the peak of GCC was nonlinearly correlated with the peak of the total chlorophyll concentration (Yang et al 2014, Wingate et al 2015.

The potential of phenology and physiology in controlling IAV of NEP
At both large regional and local scales, the physiology played a more important role than phenology in controlling the IAV of NEP (figures 5 and 6). However, such conclusion varied among vegetation types (Supplementary Table 4, 5). Consistent results have been found at both two scales that physiology was more important for explaining the IAV of NEP at deciduous and needleleaf forest, cropland and wood savannas ecosystems, while phenology showed a greater potential in grassland ecosystem. The inconsistency also existed that the local scale physiology was significantly correlated to the IAV of NEP at mixed forest and wetland sites, but large regional phenology dominated the IAV of NEP in mixed forest ecosystem. The consistency and inconsistency indicate the different strategies of C sequestration among different vegetation types (

The role of physiology in regulating IAV of NEP among vegetation types.
Both VI max and GPP max are highly correlated with the IAV of NEP in deciduous forest (Supplementary Table 4, 5), implying that the summer photosynthesis is more crucial in controlling annual C sequestration than other vegetation types. Physiology (both VI max and GPP max ) showed the relatively lowest correlation with the IAV of NEP in needleleaf forest. We suggest that the summer photosynthetic rate of needleleaf forest is relatively less important in regulating C uptake than other ecosystems. It can be explained by its physiological characteristic with the lower rates of photosynthesis and longer GSL as the survival strategies in adapting to the environmental changes (Piao et al 2007, Richardson et al 2010, Wu and Chen 2013. For cropland, the dominate role of physiology was also observed by a previous study that the IAV of the peak NDVI rate was the main cause of NEP in a maize cropland ecosystem (Du et al 2019), which emphasizes the importance of summer growth in crop development.

The role of phenology in regulating IAV of NEP among vegetation types.
Compared with physiology, phenology exhibited relatively limited ability in controlling the IAV of NEP at both two scales, while high potential of phenology was found at grassland and cropland site. This was also supported by other study that annual NEP of GRA and CRO was more sensitive to GSL (7.9 g C m −2 d −1 ), while the NEP of forest exhibited less sensitive (5.8 g C m −2 d −1 for DBF and 3.4 g C m −2 d −1 for ENF) (Churkina et al 2005). Phenology was better than physiology in determining the IAV of NEP in GRA ecosystem, given similar positive correlations for GSL (R = 0.27, p < 0.01 for NDVI and R = 0.26, p < 0.01 for GPP) (Supplementary Tables 4 and 5). The biophysical characteristics of GRA may explain the higher sensitivities to phenology. As a water restricted ecosystem, the variability of GSL determined by precipitation events, which thereby appreciably changes annual carbon sequestration (Poulter et al 2014, Ahlstrom et al 2015, Tang et al 2015, Zhang et al 2020. Additionally, our results showed that SOS rather than EOS was significantly correlated with the IAV of NEP at these two ecosystems, which suggests that the spring phenology had a more important impact on annual C sequestration than autumn phenology. In addition, it has reported that an advanced POS is associated with an increase of GPP and vegetation greenness throughout the growing season (Xu et al 2016, Park et al 2019, Wang and Wu 2019, especially in temperate GRA (Yang et al 2019), yet the effect of POS in controlling the IAV of NEP has not been documented. We found the crucial role of POS for IAV of NEP in grassland and cropland ecosystems. The possible biophysical mechanism is that an early POS implies the favorable spring meteorological conditions, and subject to more solar irradiance and longer day length, which would enhance vegetation photosynthesis, and thus increasing C assimilation (Duveneck and Thompson 2017, Gonsamo et al 2018).

The potential of physiology and phenology in controlling SV of NEP
We found that both physiology and phenology were significant correlated with the SV of NEP at both large regional and local scale, however, physiology seemed to show relatively limited potential at regional scale (figure 9). Considering that GCC is a simple calculation of visible bands, and therefore it lacks an inherent connection with photosynthesis (Serbin et al 2015, Liu et al 2020, especially in summer with a high canopy closure (Brown et al 2017). Previous studies have shown that the chlorophyll concentration is still rising when 'summer greendown' appears on GCC (Elmore et al 2012, Yang et al 2014, Toomey et al 2015. Keenan et al (2014) also found that canopy greenness rises very quickly to the maximum before measured physiological and morphological maturity. Therefore, it may provide the explanation that no significant correlation was observed between NEP and GCC max in our study (figure 9(e)).
As a more biophysical mechanism supported index, satellite-based vegetation index showed completely different performances. The VI max exhibited relatively higher potential than phenology in explaining the SV of NEP ( figure 8(e)). Furthermore, EVI max and OSAVI max performed better than NDVI max , which is in line with the view that NDVI trends to be saturated than EVI and OSAVI in summer (Huete et al 2002, Jin andEklundh 2014). The phenological metrics, including SOS, POS and GSL, were also important in explaining the spatial difference of NEP (figures 8(a), (b) and (d)). Among them, GSL showed the highest correlation, which is also confirmed by other studies (Wu et al 2012, Fu et al 2017. The negative correlation between POS and NEP implies that an early POS usually corresponds to a higher NEP, which agrees with that the earlier occurrence of peak plant activity will cause increased plant productivity in extratropical ecosystems (Gonsamo et al 2018). The generally consistent results were obtained from the local scale that both GSL and GPP max have comparable strong ability in regulating the SV of NEP (figures 10(d), (e)).

Conclusion
Using satellite observation, 676 site-year eddy covariance measurements from 67 sites across 7 vegetation type, as well as 57 site-years PhenoCam data from 16 sites in North America and Europe, we explored the role of phenology and physiology in controlling the variability of carbon sequestration from multiple scales. The relatively consistent results were drawn from both large regional and local scales. We found that physiology exhibited greater potential than phenology in controlling the IAV of NEP in most ecosystems, while phenology, especially spring and summer phenology, performed better in grassland ecosystem. Additionally, physiology exhibited better performance in interpreting the SV of NEP at large regional scale, and phenology also played an important role at the local scale. The inconsistent of regional scale results can be explained that GCC lacks biochemical and biophysical mechanism to characteristic process of photosynthesis, especially in the summer with high canopy closure. Our study highlights that the underestimated summer physiology is critical in determining the variability of carbon sequestration. Understanding the response of summer physiology to climate drivers is of great significance for simulating terrestrial ecosystem carbon cycle under the global climate change, which can be the focus of future research.