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

Yueyue LI ,Li DAN ,Jing PENG ,Junbang WANG ,Fuqiang YANG ,Dongdong GAO ,Xiujing YANG,and Qiang YU

1Key Laboratory of Water Cycle and Related Land Surface Processes,Institute of Geographic Sciences and Natural Resources Research,Chinese Academy of Sciences,Beijing 100101,China

2College of Resources and Environment,University of Chinese Academy of Science,Beijing 100049,China

3Key Laboratory of Regional Climate-Environment for Temperate East Asia,Institute of Atmospheric Physics,Chinese Academy of Sciences,Beijing 100029,China

4National Ecosystem Science Data Center,Key Laboratory of Ecosystem Network Observation and Modeling,Institute of Geographic Sciences and Natural Resources Research,Chinese Academy of Sciences,Beijing 100101,China

5State Key Laboratory of Soil Erosion and Dryland Farming on the Loess Plateau,Northwest A & F University,Shaanxi 712100,China

6School of Life Sciences,University of Technology Sydney,Sydney 2007,NSW,Australia

ABSTRACT Gross primary production (GPP) plays a crucial part in the carbon cycle of terrestrial ecosystems.A set of validated monthly GPP data from 1957 to 2010 in 0.5° × 0.5° grids of China was weighted from the Multi-scale Terrestrial Model Intercomparison Project using Bayesian model averaging (BMA).The spatial anomalies of detrended BMA GPP during the growing seasons of typical El Ni?o years indicated that GPP response to El Ni?o varies with Pacific Decadal Oscillation(PDO) phases:when the PDO was in the cool phase,it was likely that GPP was greater in northern China (32°–38°N,111°–122°E) and less in the Yangtze River valley (28°–32°N,111°–122°E);in contrast,when PDO was in the warm phase,the GPP anomalies were usually reversed in these two regions.The consistent spatiotemporal pattern and high partial correlation revealed that rainfall dominated this phenomenon.The previously published findings on how El Ni?o during different phases of PDO affecting rainfall in eastern China make the statistical relationship between GPP and El Ni?o in this study theoretically credible.This paper not only introduces an effective way to use BMA in grids that have mixed plant function types,but also makes it possible to evaluate the carbon cycle in eastern China based on the prediction of El Ni?o and PDO.

Key words:East China,Bayesian model averaging,Gross primary production,El Ni?o,Pacific Decadal Oscillation,Monsoon rainfall

1.Introduction

Gross primary production (GPP) is the total amount of carbon dioxide "fixed" by land plants per unit time through the photosynthetic reduction of COinto organic compounds (

Gough,2011

).The capacity of terrestrial GPP to partly offset anthropogenic COemissions is especially important for China,where the emissions of fossil fuel COinto 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 COinduced 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.

2.Materials and Methods

2.1.MsTMIP GPP and CRU–NCEP Datasets

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 COconcentration (CO) and nitrogen deposition (N).The timevarying 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-CTEMN model was excluded because its GPP estimates for China are much lower than is reasonable (

Shao et al.,2016

).

Table 1.Characteristics of the 16 El Ni?o years in this study.

The available climatology datasets from the Climate Research Unit (CRU) (

Harris et al.,2014

) and the NationalCenters 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).

2.2.Bayesian Model Averaging

The optimum of the likelihood function is:

After implementation of training with data

y

,βfrom the solutions of Eq.(2) should satisfy the conditions:

where

w

=βfrac+βfrac+βfrac.The BMA GPP in the target grid cell was given as:

2.3.Data Preprocessing

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.

3.Results

3.1.Validation of the BMA GPP at the ChinaFLUX Sites

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 Eindices (see Appendix A.4).The Eindices 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 Eindices 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 Eindexes suggested that the manner in which BMA was used in this study was unable to reduce the uncertainties from these two specific PFTs.

Table 2.Nash–Sutcliffe efficiency index values for GPP simulated by 14 terrestrial biosphere models and Bayesian model averaging at eight ChinaFLUX sites.

3.2.The Response of GPP to El Ni?o in Different Phases of the PDO

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).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 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.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.

3.3.The Solid Relationship between GPP and Rainfall in Eastern China

The inter-annual relationships between meteorological factors and vegetation activity can be changed in different dec-ades (

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).

Fig.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.

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.

4.Discussion

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 Rvalues 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.

Fig.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.

Fig.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.

Fig.4.Latitude–time cross-sections of the (a) rainfall (mm yr?1) and (b) BMA GPP (gC m?2 yr?1) 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.

Fig.5.Latitude–time cross-sections of the (a) rainfall (mm yr?1) and (b) BMA GPP (gC m?2 yr?1) 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.

Fig.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.

Fig.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.

Fig.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.

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 interannual 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.

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

5.Conclusions

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.

.

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).

Data availability statement

The MsTMIP outputs (

https://nacp.ornl.gov/mstmipdata/

) and drivers (

https://daac.ornl.gov/NACP/guides/NACP_MsTMIP_Mod el_Driver.html

) can be downloaded freely.The observed GPP from ChinaFLUX can be obtained from www.cnern.org.cn/data/init-DRsearch.The MOD17A2 products are available from

http://files.ntsg.umt.edu/data/NTSG_Products/MOD17/GeoTIFF/Monthly_M OD17A2/GeoTIFF_0.05degree/

.The BMA GPP for China are available in Zenodo with the identifier

http://doi.org/10.5281/zenodo.3899605

.

APPENDIX A Supplementary Information for Materials and Methods

APPENDIX A.1 Observed GPP from ChinaFLUX

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

.

APPENDIX A.2 Uncertainties in the MsTMIP GPP

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 weaccounted 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.

Table A1. Information for the eight ChinaFLUX sites used in this study.

Fig.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.

APPENDIX A.3 Merging of PFTs when using BMA

When use is made of the observed GPP from China-FLUX 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

).

Table A2. Summary of SYNMAP land cover types used for Bayesian model averaging.

Table A3. Weight coefficients for 14 different models calculated by Bayesian model averaging.

Fig.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.

APPENDIX A.4 Nash–Sutcliffe efficiency index

The Nash–Sutcliffe efficiency index (E;

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

is the ith observed data point,

o

is the mean of the observed data,and

m

is the ith modeled data point.The Evalues range between ?∞ and 1.Evalues <0 indicate that the mean value of the observed time series is a better predictor than the model.Evalues between 0 and 1 generally indicate that the performance of the model is at an acceptable level.Evalues equal to 1 indicate a perfect fit.

APPENDIX A.5 MODIS GPP and MTE GPP

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) GPPupscaled 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.