Response of Growing Season Gross Primary Production to El Niño in Different Phases of the Pacific Decadal Oscillation over Eastern China Based on Bayesian Model Averaging

Gross primary production (GPP) is the total amount of carbon dioxide "fixed" by land plants per unit time through the photosynthetic reduction of CO₂ into organic compounds (Gough, 2011). The capacity of terrestrial GPP to partly offset anthropogenic CO₂ emissions is especially important for China, where the emissions of fossil fuel CO₂ into the atmosphere are high (Gregg et al., 2008). The subtropical forests in eastern China can take in a considerable amount of CO₂ (Yu et al., 2014), and contribute remarkably to the total carbon sequestration of Chinese ecosystems (Zhang et al., 2019a). The climate in eastern China is dominated by the East Asian Summer Monsoon (EASM), which is one of the most active components in the global climate system and can significantly affect the movement of local rain belts (Ding et al., 2008). Given the ecological importance of the monsoon climate in China and the high sensitivity of the carbon uptake of ecosystems in EASM region to rainfall (Dan et al., 2015; Zhang et al., 2019a), it is important to develop a deeper understanding of the response of GPP to the variability in rainfall during the EASM.

Different modes of sea surface temperature (SST) can explain most of the inter-annual and inter-decadal variations of rainfall in China (Yang and Lau, 2004). Two modes of SST that reflect strong climatic teleconnections with the EASM rainfall are El Niño and the Pacific Decadal Oscillation (PDO) (Ma, 2007). An El Niño event is characterized by anomalous warming in the central and eastern equatorial Pacific Ocean and occurs, on average, every two to seven years. By contrast, the PDO (Mantua et al., 1997; Qian and Zhou, 2014) is marked by an inter-decadal cycle, and is mainly present in the North Pacific (poleward of 20°N). During the warm phase of the PDO, SST tends to be anomalously cool in the central North Pacific, coincident with anomalously warm SST in the eastern North Pacific, whereas the spatial pattern of SST is reversed during the cool phase of the PDO.

Because the inter-annual variability in the atmospheric concentrations of CO₂ induced by El Niño is mainly attributed to tropical terrestrial ecosystems (Yang and Wang, 2000; Cox et al., 2013; Wang et al., 2013; Liu et al., 2017), previous studies have primarily focused on how terrestrial carbon fluxes in the tropics responding to changes in the Earth’s climate (especially changes in temperature and rainfall) caused by El Niño (Gu and Adler, 2011; Bastos et al., 2013; Wang et al., 2016; Fang et al., 2017). The most common way to investigate the decadal variability in the terrestrial carbon budget caused by the PDO is to use the long-term teleconnection between the PDO indices and specific carbon cycle variables (Wharton and Falk, 2016; Zhang et al., 2018). However, this method is unable to reflect the impact of the variability of local chaotic weather events (Ito, 2011), especially when the responses of terrestrial carbon cycle to different El Niño events are contrasting (Liu et al., 2017; Wang et al., 2018). Hence, the mixed effects of El Niño and the PDO on the carbon cycle in the extratropical monsoon regions of China are still poorly understood.

Even though the spatiotemporal patterns of rainfall anomalies over East Asia associated with individual El Niño episodes are variable (Wang et al., 2017), Feng et al. (2014) still found relatively stable relationships between the EASM rainfall and El Niño under different phases of the PDO. When the PDO is in a cool phase, with a high-pressure area in the subtropical western Pacific clearly experiencing two northward shifts, positive rainfall anomalies move from southern to northern China. By contrast, when the PDO is in a warm phase, a greater amount of rainfall is observed in central China, whereas less rainfall is observed in southern and northern China. How GPP responds to these specific rainfall patterns is still unclear.

The purpose of this study was to determine the relationship between the GPP and the changes in rainfall associated with El Niño under different phases of the PDO in eastern China. In particular, we attempted to identify the regions driving this relationship and to demonstrate the response of GPP to El Niño and the PDO. The results may have implications for how terrestrial carbon dynamics responding to variabilities in the Earth’s climate.

The rainfall patterns described in the fourth paragraph of the Introduction were revealed and explained by data spanning from 1957 to 2010. Given the data consistency, this study also focused on the same period and selected the same El Niño years (see Table 1) as Feng et al. (2014). Compared with other approaches, the Terrestrial Biosphere Model (TMB) is a preferred method to estimate carbon fluxes over such a long period because of the availability of climate forcing data and process-based simulation. The monthly GPP products of 15 TMBs in the Multi-scale synthesis and Terrestrial Model Intercomparison Project (MsTMIP) (Huntzinger et al., 2013; Wei et al., 2014) were applied here. The MsTMIP version 1.0 dataset includes four global (0.5° × 0.5° resolution) sensitivity simulations (SG1, SG2, SG3 and BG1) which can systematically assess the impacts of different forcing factors. The model-driven data in MsTMIP version 1.0 are climate variables (Climate), land use/land cover changes (LUCC), the atmospheric CO₂ concentration (CO₂) and nitrogen deposition (N). The time-varying drivers in the four global sensitivity simulations are Climate (SG1), Climate + LUCC (SG2), Climate + LUCC + CO₂ (SG3), and Climate + LUCC + CO₂ + N (BG1). BG1 only includes a limited number of models as a result of the lack of a nitrogen cycle in other MsTMIP members. SG3 and BG1 are presumed to be closer to reality than SG1 and SG2 and were therefore used in this study. The CLASS-CTEM-N model was excluded because its GPP estimates for China are much lower than is reasonable (Shao et al., 2016).

The available climatology datasets from the Climate Research Unit (CRU) (Harris et al., 2014) and the National Centers for Environmental Prediction (NCEP)/National Center for Atmospheric Research (NCAR) (Kalnay et al., 1996) cannot fully meet the spatial and temporal requirements of MsTMIP. Therefore, a new dataset named CRU-NCEP was produced by combining the strengths of the datasets from the CRU and NCEP/NCAR as driver data for MsTMIP. The CRU–NCEP dataset not only maintains a good match to the observed rainfall in China (Zhao and Fu, 2006), but also corrects the known biases in temperature and downward shortwave radiation in the NCEP/NCAR reanalysis products (Wei et al., 2014). Therefore, the globally gridded (0.5° × 0.5°) and monthly CRU-NCEP dataset was chosen to determine the relationships between GPP and the particular meteorological factors (rainfall, temperature, and downward shortwave radiation).

Eight typical ChinaFLUX sites (see Appendix A.1) were selected to validate modeled GPP. After comparing with the observed GPP from ChinaFLUX, it was found that there were large uncertainties in the MsTMIP GPP (see Appendix A.2). Selecting only one model may lead to statistical biases and an underestimation of the uncertainties (Neuman, 2003; Raftery et al., 2005). Bayesian model averaging (BMA) (Hoeting et al., 1999) can effectively solve this problem. After comprehensively considering the conditional distribution of each ensemble member, BMA takes the maximum of the probability density function (PDF) as an individual model’s weighted coefficient following a training process. Thus, the final weighted average results reflect each model’s optimum contribution to the quantities of interest (Wasserman, 2000). We used the MATLAB toolbox provided by Vrugt (2016) to implement BMA. The core mathematical expressions are as follows. If $y$ is the quantity of interest on the basis of the training data y^T, and ${{M}} = ({{{M}}_1}, \cdot \cdot \cdot,{{{M}}_K})$ is the set of all models with mean $\overline {{M}} $ and variance ${\sigma ^2}$, then ${\rm{Pr}}$ is expressed as the PDF, the weighted coefficient of ${{{M}}_k}$ is defined as ${\beta _k}$, and the PDF of BMA at a given time $t$ is:

The optimum of the likelihood function is:

After implementation of training with data ${y^T}$, ${\beta _k}$ from the solutions of Eq. (2) should satisfy the conditions:

The plant function types (PFTs) were simply merged into three classes (trees, grasses, and crops) when the BMA was trained during 2003–08 (see Appendix A.3). Before the weighted coefficients for different models were used at the regional scale, they needed to be further modified depending on the proportion of land cover types in each grid cell. If the BMA weighted coefficients trained in the ChinaFLUX sites for ${{{M}}_k}$ are ${\beta _{{\rm{Trees}},k}}$, ${\beta _{{\rm{Grasses}},k}}$, and ${\beta _{{\rm{Crops}},k}}$, and the proportion of trees, grasses, and crops in each grid cell are ${{\rm{frac}}_{{\rm{Trees}}}}$, ${{\rm{frac}}_{{\rm{Grasses}}}}$, and ${{\rm{frac}}_{{\rm{Crops}}}}$, respectively, then the modified weighted coefficient $\beta _k^ * $ for ${{{M}}_k}$ in the target grid cell is defined as:

where ${w_k} \;=\; {\beta _{{\rm{Trees}},k}}{\rm{fra}}{{\rm{c}}_{{\rm{Trees}}}}\;\;{\rm{ + }}\;\;{\beta _{{\rm{Crops}},k}}{\rm{fra}}{{\rm{c}}_{{\rm{Crops}}}}\;\;{\rm{ + }}\;\;{\beta _{{\rm{Grasses}},k}} $$ {\rm{fra}}{{\rm{c}}_{{\rm{Grasses}}}}$. The BMA GPP in the target grid cell was given as:

To highlight the impact of the EASM on GPP, May and September were defined as the start and end months of the growing season corresponding to the onset and withdrawal months of the EASM. When analyzing the response characteristics of GPP to El Niño over China and the corresponding mechanism, GPP and the meteorological data were not only limited to the growing season in each year, but also were detrended to minimize the possible influence of long-term trends. Trends were calculated by simple linear regression.

The MsTMIP and BMA GPP data in the grid cells containing ChinaFLUX sites from 2009 to 2010 were chosen to evaluate the rationality of using BMA. The efficiency criteria used in this section are E_f indices (see Appendix A.4). The E_f indices for different models are listed in Table 2. BMA performed better than the MsTMIP models in the qualified grid cells corresponding to Changbaishan, Qianyanzhou, Yucheng, and Inner Mongolia, thereby proving its strength for the typical PFTs at these sites. Even though the E_f indices calculated in the grid cells for Dinghushan and Haibei were not statistically significant, the relatively small absolute values still indicated that the magnitude of the BMA GPP was acceptable, especially in the subtropical region where Dinghushan is located. The evergreen broad leaf trees in Xishuangbanna and the grasses in Dangxiong are rainforests and alpine meadows, respectively, and the negative E_f indexes suggested that the manner in which BMA was used in this study was unable to reduce the uncertainties from these two specific PFTs.

Figure 1 shows the average spatial anomalies of the BMA GPP and the CRU–NCEP rainfall for the selected El Niño years with different phases of the PDO. The GPP and rainfall with large anomalies were mainly distributed in northern China (32°–38°N, 111°–122°E) and the Yangtze River valley (28°–32°N, 111°–122°E). The similar spatial distribution of the GPP and the rainfall anomalies shown in Fig. 1 indicate that rainfall is the connection between the El Niño and GPP in eastern China (28°–38°N, 111°–122°E).

Figure 1.  Spatial anomalies of BMA GPP during the growing season (May–September) averaged from El Niño years under (a) the cool phase of the PDO and (b) the warm phase of the PDO in China. (c) and (d) are the same as (a) and (b) respectively, but for rainfall. The anomalies are relative to the mean of the time period 1957–2010. The dashed lines denote negative anomalies and the colored areas denote positive anomalies. The target zones (28°–32°N, 111°–122°E; 32°–38°N, 111°–122°E) are shown by the rectangles. GPP and rainfall were detrended in advance.

The spatial anomalies of the GPP averaged from multiple El Niño years (Fig. 1a, Fig. 1b) leads to the following hypotheses: (1) the zone in which GPP was significantly affected by El Niño was in eastern China (28°–38°N, 111°–122°E); (2) when the El Niño years were in the cool phase of the PDO, GPP was higher in northern China (32°–38°N, 111°–122°E) and lower in the Yangtze River valley (28°–32°N, 111°–122°E); and (3) when the El Niño years were in the warm phase of the PDO, the spatial anomaly patterns of the GPP in these two regions were reversed.

To check whether these hypotheses are valid, the response of GPP to each El Niño year was analyzed. Figure 2 shows the spatial anomalies of GPP in eight typical El Niño years under the cool phase of the PDO. Except for 1968 and 1994, the graphs for other years in Fig. 2 all support these hypotheses and the most positive anomalies in northern China (32°–38°N, 111°–122°E) were statistically significant.

Figure 2.  Spatial anomalies of BMA GPP during the growing season (May–September) for eight selected El Niño years under the cool phase of the PDO in China. The anomalies are relative to the mean of the time period 1957–2010. The blank areas denote negative anomalies and the colored areas denote positive anomalies. The target zones (28°–32°N, 111°–122°E; 32°–38°N, 111°–122°E) are shown by the rectangles. σ is the standard deviation of the BMA GPP anomalies, where ≥1.5σ means that the positive anomalies are significant. GPP was detrended in advance.

Figure 3 shows the spatial anomalies of GPP in eight typical El Niño years under the warm phase of the PDO. Excluding 1957 and 1976, the positive anomalies of GPP for the remaining El Niño years were mainly distributed in the Yangtze River valley (28°–32°N, 111°–122°E), which further verified the hypotheses from the averaged El Niño years. Less significant positive GPP anomalies were found in the Yangtze River valley. This difference suggests that plants in northern China are more sensitive to El Niño than those in the Yangtze River valley.

Figure 3.  Spatial anomalies of BMA GPP during the growing season (May–September) for eight selected El Niño years under the warm phase of the PDO in China. The anomalies are relative to the mean of the time period 1957–2010. The blank areas denote negative anomalies and the colored areas denote positive anomalies. The target zones (28°–32°N, 111°–122°E; 32°–38°N, 111°–122°E) are shown by the rectangles. σ is the standard deviation of the BMA GPP anomalies, where ≥1.5σ means that the positive anomalies are significant. GPP were detrended in advance.

The evolution of the zonal-time anomaly (May–September) of rainfall and GPP averaged from 111° to 122°E for typical El Niño years in the cool and warm phases of the PDO are shown in Fig. 4 and Fig. 5. Even when there are some deflections, it is still clear that the large anomalies in the rainfall and GPP are synchronized both spatially and temporally. This confirms that the close relationship between GPP and El Niño in eastern China is related to the north–south movement of the rain belt.

Figure 4.  Latitude–time cross-sections of the (a) rainfall (mm yr⁻¹) and (b) BMA GPP (gC m⁻² yr⁻¹) anomalies during the growing season (May–September) for El Niño years under the cool phase of the PDO in eastern China (28°–38°N, 111°–122°E). The dashed lines denote negative anomalies and the colored areas denote positive anomalies. GPP and rainfall were detrended in advance.

Figure 5.  Latitude–time cross-sections of the (a) rainfall (mm yr⁻¹) and (b) BMA GPP (gC m⁻² yr⁻¹) anomalies during the growing season (May–September) for El Niño years under the warm phase of the PDO in eastern China (28°–38°N, 111°–122°E). The dashed lines denote negative anomalies and the colored areas denote positive anomalies. GPP and rainfall were detrended in advance.

The inter-annual relationships between meteorological factors and vegetation activity can be changed in different decades (Piao et al., 2014). The EASM weakened from the late 1970s but has gradually recovered since the early 1990s (Liu et al., 2012). Two abrupt climate change points are recognized in 1978 and 1992 (Ding et al., 2008). To determine whether the relationships between GPP and the particular meteorological factors (rainfall, temperature, and the downward shortwave radiation) are consistent under different EASM intensities, they were compared in different time periods (1957–78, 1979–92, and 1993–2010).

To confirm the role of key meteorological factors for GPP, partial correlations between GPP and rainfall, downward shortwave radiation, and temperature in China are given in Fig. 6. This analysis shows that rainfall dominated the inter-annual changes of the GPP in eastern China from 1957 to 2010. By contrast, the roles of temperature and downward shortwave radiation in the same region were rather limited, although the effect of temperature was clearer in the time period 1957–78. The bar diagrams in Fig. 7 show the proportions of partial correlations between GPP and the particular meteorological factors in northern China and the Yangtze River valley, respectively. The highest proportions of positive partial correlations are between GPP and rainfall in all periods in Fig. 7, which clearly shows the importance of rainfall. The results in this section further enhance the credibility of the response characteristics of GPP to El Niño in eastern China being mostly due to its close relationship with rainfall.

Figure 6.  Spatial patterns of partial correlations between the growing season (May–September) BMA GPP and (a) rainfall (RAIN), (b) downward shortwave radiation (DWSR), and (c) temperature (TEMP) in China in 1957–78. (d), (e), and (f) are the same as (a), (b), and (c), respectively, but in 1979–92. (g), (h), and (i) are also the same as (a), (b), and (c), respectively, but in 1993–2010. Gray regions indicate insignificant partial correlations (P > 0.05). The target zones (28°–32°N, 111°–122°E; 32°–38°N, 111°–122°E) are shown by rectangles. The GPP and meteorological factors were detrended in advance.

Figure 7.  Proportions of positive (+) and negative (–) partial correlations between the growing season (May–September) BMA GPP and rainfall (P), temperature (T), and downward shortwave radiation (R) in (a) northern China (32°–38°N, 111°–122°E) and (b) the Yangtze River valley (28°–32°N, 111°–122°E) in strong and weak East Asian Summer Monsoon periods (1957–78, 1979–92, and 1993–2010). Gray bars indicate the proportions of insignificant partial correlations. All data were detrended in advance.

The mismatch in spatial scale between the observed GPP at the flux sites and the modeled GPP from MsTMIP in the corresponding 0.5° × 0.5° grids may lead to statistical errors. The monthly GPP from the MODerate-resolution Imaging Spectroradiometer (MODIS) in 0.05° × 0.05° grids and machine learning Model Tree Ensemble (MTE) in 0.5° × 0.5°grids (see Appendix A.5) were chosen to validate the quality of the BMA GPP with such errors. From the scatterplots averaged from northern China and the Yangtze River valley (see Fig. 8), we found that although the amplitude of the BMA GPP anomaly was less than the MODIS GPP anomaly, their R² values all passed 99% confidence tests. Compared with the single model in MsTMIP, BMA GPP is more similar to the multi-year averaged distribution of MTE GPP, especially in eastern China (see Fig. 9). This spatiotemporal consistency between BMA GPP and validation data in Fig. 8 and Fig. 9 not only proves the effectiveness of validating simulations in coarse grids with site observations, but also shows the rationality of simply classifying PFTs as trees, grasses, and crops when using BMA. There are no rainforests or alpine meadows in eastern China, so the uncertainties from these ecosystems can be neglected. We conclude that the quality of the BMA GPP in eastern China is reliable.

Figure 8.  Anomaly scatterplots of the monthly MODIS GPP against the monthly BMA GPP averaged in (a) northern China (32°–38°N, 111°–122°E) and (b) the Yangtze River valley (28°–32°N, 111°–122°E) in the time period 2000–10. All of the GPP data were detrended in advance.

Figure 9.  Spatial patterns of the annual GPP in China averaged from MTE, BMA, and the various MsTMIP model contributors during 1982−2010.

The general consensus is that GPP in the tropics is usually reduced due to the warm and dry conditions caused by El Niño (Cavaleri et al., 2017; Yang et al., 2018; Qian et al., 2019; Yan et al., 2019). By contrast, the current study was carried out in the temperate zone (28°–38°N, 111°–122°E) of China. We found that GPP was not simply increased or decreased, but showed opposite north–south anomalies (divided by 32°N) when El Niño occurred under different phases of the PDO. However, the spatial anomalies of the BMA GPP in eastern China for 1957, 1968, 1976, and 1994 did not follow the general pattern seen in most selected El Niño years. This inconsistency suggests that using a combination of El Niño and the phase of the PDO to predict GPP is not completely reliable. The misplacement of the BMA GPP anomalies in eastern China for 1968 (minor El Niño) compared with 2006 (minor El Niño), 1994 (moderate El Niño) compared with 2004 (moderate El Niño), and 1957 (major El Niño) compared with 1991 (major El Niño) also indicates that a classification of the El Niño event as “strong” or “weak” is not sufficient to judge the inter-annual variations in the GPP. The diversities in how GPP responds to different El Niño episodes further complement and deepen our understanding of the nonlinear climatic effects of sea surface temperature on local plants.

The inter-annual variabilities of the GPP derived from eddy covariance network observations are controlled by rainfall in temperate regions (Jung et al., 2011), which is consistent with the results in this study. The different responses of GPP to El Niño in tropical and temperate regions are seen not only for the enzyme-driven land surface models in this paper, but are also reported for the model driven by remote sensing data (Zhang et al., 2019b). Specific terrestrial biosphere models (e.g., DLEM, LPJ, VEGAS, VISIT, CLM, ISAM, and ORCHIDEE) with a proven ability to reproduce the deviations in carbon flux caused by El Niño (Wang et al., 2016; Chang et al., 2017; Bastos et al., 2018) also took part in MsTMIP and contributed to the BMA GPP. Therefore, our conclusions based on the BMA GPP are credible.

The GPP in eastern China is difficult to describe quantitatively as a result of the uncertainties in the modeled GPP and the variety of climatic effects from different El Niño events (Schwalm et al., 2011), therefore the responses of GPP to El Niño are only qualitatively shown here. This study only focused on the inter-annual relationships between different climate variables (rainfall, temperature, and downward shortwave radiation) and GPP, so the roles of specific climate changes (e.g., an increase/decrease in the intensity or frequency of rainfall) were neglected. Given the complexities of the Earth’s climate system, the physical mechanisms in dry episodes associated with El Niño and the resulting state of the vegetation in eastern China still require further research. TBMs simulated the physiological processes of crops more poorly than other PFTs that were less affected by human activities (see Appendix A.2). A benefit from the advantages of the BMA method is that this defect in eastern China where crops occupy a large proportion of the area can be minimized by weighting different models according to their performance. Additionally, the intensity of agricultural management (especially irrigation) is difficult to evaluate at a regional level, so its impacts on the conclusions found in this study also need further investigation.

We selected BMA, a method that can integrate the advantages of a number of different models into the final result, to produce a set of long time-series GPP data in China. When training BMA at ChinaFLUX sites and applying weight coefficients in spatial scale, PFTs were simply classified as trees, grasses, and crops. The validation of BMA GPP proved this kind of up-scaling method is reliable, especially in eastern China.

In contrast to the tropics, where it is generally believed that local plant physiology can be significantly affected by El Niño, this paper focused on the neglected eastern monsoon region of China. The spatial anomalies of BMA GPP averaged from the multiple El Niño years in different phases of PDO produced the following hypothesis: El Niño during the cool phase of the PDO led to greater GPP in northern China (32°–38°N, 111°–122°E) and less GPP in the Yangtze River valley (28°–32°N, 111°–122°E); but El Niño during the warm phase of PDO reversed the GPP anomalies in these two regions. Even though there were some discrepancies, most individual El Niño years supported the hypothesis mentioned above, and therefore the conclusions have high applicability. It should also be noted that the effects of El Niño are more significant in northern China than in the Yangtze River valley.

The synchronized spatiotemporal distribution in selected El Niño years and high partial correlation both in strong and weak periods of the EASM all demonstrated that rainfall linked El Niño and GPP in eastern China. Although EASM systems are complex (Ding and Chan, 2005) and the movements of the EASM rain belt are difficult to forecast (Gao et al., 2011, 2014; Fan et al., 2012), this study showed that we can still predict how GPP will change in eastern China (28°–38°N, 111°–122°E) when El Niño occurs if we can identify the phase of the PDO. The results of this study also provide a reference for future research on the impact of El Niño on the carbon cycle of terrestrial ecosystems in China.

Acknowledgments. This study was jointly supported by the National Key Research and Development Program of China (Grant Nos. 2016YFA0602501 and 2018YFA0606004) and the Strategic Priority Research Program of the Chinese Academy of Sciences (Grant Nos. XDA20040301 and XDA20020201).

The MsTMIP outputs (https://nacp.ornl.gov/mstmipdata/) and drivers (https://daac.ornl.gov/NACP/guides/NACP_MsTMIP_Model_Driver.html) can be downloaded freely. The observed GPP from ChinaFLUX can be obtained from www.cnern.org.cn/data/initDRsearch. The MOD17A2 products are available from http://files.ntsg.umt.edu/data/NTSG_Products/MOD17/GeoTIFF/Monthly_MOD17A2/GeoTIFF_0.05degree/. The BMA GPP for China are available in Zenodo with the identifier http://doi.org/10.5281/zenodo.3899605.

The monthly observed GPP data from eight typical flux measurement sites in the ChinaFLUX network (Yu et al., 2006) were used to validate the modeled GPP. The selected eddy covariance tower sites cover all the main ecosystem types in China and have been successfully applied to validate the performance of different models (Li et al., 2013; Wang et al., 2015; Zhang et al., 2016). Detailed information about these sites is given in Table A1.

The monthly MsTMIP GPP values in the grid cells corresponding to the location of each flux site were selected for comparison with the observed GPP. The statistical properties obtained in this way only had a credible significance when the dominant plant function types (PFTs) in the selected grid cells were the same as the observed PFTs at the flux sites and their proportions were >50% (Peng et al., 2015). The land use and land cover change maps used by MsTMIP were prescribed by the SYNergetic land cover MAP (SYNMAP) (Jung et al., 2006). Fig. A1 displays the SYNMAP biome type fractions in the grid cells corresponding to the ChinaFLUX sites (see Table A1). Even if we accounted for the number of mixed trees, the proportion of the target PFT in the grid cell corresponding to Dinghushan was still <50% and did not meet the requirements mentioned above. The grid cell containing the Haibei site had similar problems. Evergreen broad leaf trees (EBT) and shrubs (SHB) are widely distributed in China. Unfortunately, Dinghushan and Haibei are the only two sites in China that share observations, respectively, for subtropical EBT and alpine SHB in the time period 2003–2010, and therefore cannot be replaced.

Figure A1.  Biome type compositions in the SYNMAP grids corresponding to the ChinaFLUX sites.

Taylor diagrams were used to indicate how well the modeled results matched the observed values in terms of three statistics: the Pearson correlation coefficient (R), the root-mean-square difference (RMSD), and the standard deviation (STD) (Taylor, 2001). The STD and RMSD are proportional to the radial distance from the origin to the observed point. R is defined as the cosine of the azimuth angle. In general, the closer a modeled point is to the observed point, the better the model performs. The Taylor diagrams in Fig. A2 show how closely the GPP from the 14 MsTMIP models matched the observations from 2003 to 2010. The GPP simulations had large uncertainties, especially at YuCheng, where the ground cover was cropland. This means that simple equal-weighted averaging cannot correct the bias when most models have a poor performance for specific PFTs. The diversity in the multi-model results in Fig. A2 shows the need to use a multi-model averaging method that can weight each model according to its simulation skill.

Figure A2.  Taylor diagrams for GPP simulations from 14 MsTMIP models at eight ChinaFLUX sites. The red dashed lines are the observed standard deviations at the sites.

When use is made of the observed GPP from ChinaFLUX to train BMA, the dominant PFTs in the qualified grid cells cannot cover all of the PFTs in the SYNMAP legend. To reduce the errors caused by mixed PFTs in the grid cells and the limited number of ChinaFLUX sites, the PFTs at ChinaFLUX sites and used by MsTMIP were categorized as trees, grasses, and crops when BMA was applied. The detailed classifications of the PFTs at ChinaFLUX sites and SYNMAP legend are listed in Table A1 and Table A2, respectively. After being trained by the observed GPP at the grouped ChinaFLUX sites and the MsTMIP GPP in the corresponding grid cells during 2003–08, the weight coefficients of trees, grasses, and crops for different MsTMIP models were calculated by BMA (see Table A3).

The Nash–Sutcliffe efficiency index (E_f; Nash and Sutcliffe, 1970) was used to evaluate the performance of the GPP from BMA and individual model in MsTMIP. Its mathematical expression is:

where ${o_i}$ is the ith observed data point, $\overline o $ is the mean of the observed data, and ${m_i}$ is the ith modeled data point. The E_f values range between −∞ and 1. E_f values <0 indicate that the mean value of the observed time series is a better predictor than the model. E_f values between 0 and 1 generally indicate that the performance of the model is at an acceptable level. E_f values equal to 1 indicate a perfect fit.

The monthly MODerate-resolution Imaging Spectroradiometer (MODIS) Version MOD17A2 GPP product (Zhao et al., 2005) was used to evaluate the temporal evolution of the BMA GPP in the target zones. To reduce the disadvantages caused by the difference in spatial resolution (MODIS: 0.05° × 0.05°; MsTMIP: 0.5° × 0.5°) and land cover classification system (MODIS: MOD12Q1; MsTMIP: SYNMAP), we only compared their spatially averaged anomalies from 2000 to 2010.

The monthly Model Tree Ensemble (MTE) GPP upscaled from FLUXNET observations to global 0.5° × 0.5° grids during 1982−2010 (Jung et al., 2011) was selected to estimate the magnitude of BMA GPP. Considering that the empirical relationships in machine learning algorithm are time-sensitive, we weakened the time dimension and only compared the spatial patterns of multi-year averaged GPP from MTE, BMA, and MsTMIP to show the rationality of BMA GPP distribution.