Decadal variation in CO2 fluxes and its budget in a wheat and maize rotation cropland over the North China Plain

Carbon sequestration in agroecosystems has great potential to mitigate global greenhouse gas emissions. To assess the decadal trend of CO2 fluxes of an irrigated wheat– maize rotation cropland over the North China Plain, the net ecosystem exchange (NEE) with the atmosphere was measured by using an eddy covariance system from 2005 to 2016. To evaluate the detailed CO2 budget components of this representative cropland, a comprehensive experiment was conducted in the full 2010–2011 wheat–maize rotation cycle by combining the eddy covariance NEE measurements, plant carbon storage samples, and a soil respiration experiment that differentiated between heterotrophic and below-ground autotrophic respirations. Over the past decade (from 2005 to 2016), the cropland exhibited a statistically nonsignificant decreasing carbon sequestration capacity; the average of total NEE, gross primary productivity (GPP), and ecosystem respiration (ER), respectively, were −364, 1174, and 810 gC m−2 for wheat and −136, 1008, and 872 gC m−2 for maize. The multiple regression revealed that air temperature and groundwater depth showed pronounced correlations with the CO2 fluxes for wheat. However, in the maize season, incoming shortwave radiation and groundwater depth showed pronounced correlations with CO2 fluxes. For the full 2010–2011 agricultural cycle, the CO2 fluxes for wheat and maize were as follows: for NEE they were −438 and −239 gC m−2, for GPP 1078 and 780 gC m−2, for ER 640 and 541 gC m−2, for soil heterotrophic respiration 377 and 292 gC m−2, for below-ground autotrophic respiration 136 and 115 gC m−2, and for above-ground autotrophic respiration 128 and 133 gC m−2. The net biome productivity was 59 gC m−2 for wheat and 5 gC m−2 for maize, indicating that wheat was a weak CO2 sink and maize was close to CO2 neutral to the atmosphere for this agricultural cycle. However, when considering the total CO2 loss in the fallow period, the net biome productivity was −40 gC m−2 yr−1 for the full 2010–2011 cycle, implying that the cropland was a weak CO2 source. The investigations of this study showed that taking cropland as a climate change mitigation tool is challenging and that further studies are required for the CO2 sequestration potential of croplands.

Abstract. Carbon sequestration in agroecosystems has great potential to mitigate global greenhouse gas emissions. To assess the decadal trend of CO 2 fluxes of an irrigated wheatmaize rotation cropland over the North China Plain, the net ecosystem exchange (NEE) with the atmosphere was measured by using an eddy covariance system from 2005 to 2016. To evaluate the detailed CO 2 budget components of this representative cropland, a comprehensive experiment was conducted in the full 2010-2011 wheat-maize rotation cycle by combining the eddy covariance NEE measurements, plant carbon storage samples, and a soil respiration experiment that differentiated between heterotrophic and below-ground autotrophic respirations. Over the past decade (from 2005 to 2016), the cropland exhibited a statistically nonsignificant decreasing carbon sequestration capacity; the average of total NEE, gross primary productivity (GPP), and ecosystem respiration (ER), respectively, were −364, 1174, and 810 gC m −2 for wheat and −136, 1008, and 872 gC m −2 for maize. The multiple regression revealed that air temperature and groundwater depth showed pronounced correlations with the CO 2 fluxes for wheat. However, in the maize season, incoming shortwave radiation and groundwater depth showed pronounced correlations with CO 2 fluxes. For the full 2010-2011 agricultural cycle, the CO 2 fluxes for wheat and maize were as follows: for NEE they were −438 and −239 gC m −2 , for GPP 1078 and 780 gC m −2 , for ER 640 and 541 gC m −2 , for soil heterotrophic respiration 377 and 292 gC m −2 , for below-ground autotrophic respiration 136 and 115 gC m −2 , and for above-ground autotrophic respira-tion 128 and 133 gC m −2 . The net biome productivity was 59 gC m −2 for wheat and 5 gC m −2 for maize, indicating that wheat was a weak CO 2 sink and maize was close to CO 2 neutral to the atmosphere for this agricultural cycle. However, when considering the total CO 2 loss in the fallow period, the net biome productivity was −40 gC m −2 yr −1 for the full 2010-2011 cycle, implying that the cropland was a weak CO 2 source. The investigations of this study showed that taking cropland as a climate change mitigation tool is challenging and that further studies are required for the CO 2 sequestration potential of croplands.
The widely used eddy covariance technique has fostered our understanding of the integrated fluxes of NEE, GPP and ER but cannot provide detailed CO 2 budget components, which consist of carbon assimilation (i.e., GPP), soil heterotrophic respiration (R H ), above-ground autotrophic respiration (R AA ), below-ground autotrophic respiration (R AB ), lateral carbon export at harvest, and import at sowing or through organic fertilization . These different CO 2 components result from different biological and biophysical processes (Moureaux et al., 2008) that may respond differently to climatic conditions, environmental factors and management strategies (Ekblad et al., 2005;. Differentiating among these components is a prerequisite for understanding the response of terrestrial ecosystems to changing environment (Heimann and Reichstein, 2008), thus the carbon budget evaluations have been reported for a few croplands (e.g., Moureaux et al., 2008;Ceschia et al., 2010;Wang et al., 2015;Demyan et al., 2016;Gao et al., 2017). In particular, to account for the literal carbon export, the net biome productivity (NBP) is often estimated by combining the eddy covariance technique and field carbon measurements associated with harvests and residue treatments Kutsch et al., 2010). As a detailed CO 2 budget might facilitate better predictions of agroecosystems' responses to climate change, CO 2 budget evaluations in different croplands remain necessary.
The North China Plain (NCP) is one of the most important food production regions in China, and it guarantees national food security by providing more than 50 % and 33 % of the nation's wheat and maize, respectively (Kendy et al., 2003). Irrigation by diverting water from the Yellow River is common to alleviate water stress during spring in the NCP, resulting in a very shallow groundwater depth (usually range from 2 to 4 m) along the Yellow River (Cao et al., 2016) (Fig. 1). Wang et al. (2015) suggested that groundwater-fed cropland in the NCP had been losing carbon, and other studies also reported croplands in this region as carbon sources (e.g., Li et al., 2006;Luo et al., 2008). However, the long-term variations (e.g., > 10 years) of the CO 2 fluxes over the NCP remain lacking, leaving the trend of carbon sequestration capacity of this region unknown.
To this end, this study is designed to assess the long-term variation in CO 2 fluxes and its budget of the representative wheat-maize rotation cropland in the NCP. The eddy covariance system was used to measure the CO 2 exchange from 2005 to 2016. For the full 2010-2011 agricultural cycle, we measured soil respiration and sampled crops to quantify the detailed CO 2 budget components. These measurements allow us to (1) investigate the decadal CO 2 flux (NEE, GPP and ER) trend over this cropland; (2) provide detailed CO 2 budget components; and (3) estimate the net primary productivity (NPP), net ecosystem productivity (NEP), and NBP.

Site description and field management
The experiment was conducted in a rectangular-shaped (460 m × 280 m) field of the representative cropland over the NCP (36 • 39 N, 116 • 03 E, Weishan site of Tsinghua University, Fig. 1). The soil is silt loam with a field capacity of 0.33 m 3 m −3 and saturation point of 0.45 m 3 m −3 for the top 5 cm of the soil. The mean annual precipitation is 532 mm and the mean air temperature is +13.3 • C. The winter wheat-summer maize rotation system is the representative cropping style in this region. On average, the winter wheat is sown around 17 October and harvested around 16 June of the following year with crop residues left on the field; summer maize is sown following the wheat harvest around 17 June and harvested around 16 October. Prior to sowing wheat of the next season, the field is thoroughly plowed to fully incorporate maize residues into the top 20 cm of the soil. The canopies of both wheat and maize are very uniform across the whole season. Nitrogen fertilizer is commonly applied at this site with the amount being 35 gN m −2 for wheat and 20 gN m −2 for maize. The crop density is 775 plants per square meter for wheat with a ridge spacing of 0.26 m and 4.9 plants per square meter for maize with a ridge spac- ing of 0.63 m on average. Wheat is commonly irrigated with water diverted from the Yellow River and the irrigation is about 150 mm every year; maize is rarely irrigated because of the high precipitation in the summer. During the 2010-2011 agricultural cycle, when CO 2 budget components were evaluated, winter wheat was sown on 23 October 2010 and subsequently harvested on 10 June 2011 and summer maize was sown on 23 June 2011 and harvested on 30 September 2011. The entire year from 23 October 2010 to 22 October 2011 was studied for the annual CO 2 budget evaluation.

Eddy covariance measurements
A flux tower was set up at the center of the experiment field in 2005 (Lei and Yang, 2010;. The NEE was measured at 3.7 m above ground with an eddy covariance system consisting of an infrared gas analyzer (LI-7500, LI-COR Inc., Lincoln, NE, USA) and a three-dimensional sonic anemometer (CSAT3, Campbell Scientific Inc., Logan, UT, USA). The 30 min averaged NEE was calculated from the 10 Hz raw measurements with TK2 (Mauder and Foken, 2004) from 2005 to 2012 and TK3 software pack-age (Mauder and Foken, 2011) from 2013 to 2016. The storage flux was calculated by assuming a constant CO 2 concentration profile. Nighttime measurements under stable atmospheric conditions with a friction velocity lower than 0.1 m s −1 were removed from the analysis (Lei and Yang, 2010). In the gap-filling procedure, gaps less than 2 h were filled using linear regression, while other short gaps were filled using the Mean Diurnal Variation (MDV) method ; gaps longer than 4 weeks were not filled. NEE was further partitioned to derive GPP and ER using the nighttime method (Reichstein et al., 2005;Lei and Yang, 2010), which assumes that daytime and nighttime ER follow the same temperature response, which thereby estimates the daytime ER using the regression model derived from the nighttime measurements. In particular, this study adopted the method proposed by Reichstein et al. (2005) to quantify the short-term temperature sensitivity of ER from nighttime measurements as described by the van 't Hoff equation, where T S is soil temperature, ER ref is the reference respiration at 0 • C and b is a parameter associated with the com-Q. Zhang et al.: Carbon budget over wheat-maize cropland monly used temperature sensitivity coefficient Q 10 , The long-term temperature sensitivity b of the season (either wheat or maize) was determined by averaging all the estimated short-term b in each of the 4 d windows with the inverse of the standard error as a weighing factor. The longterm temperature sensitivity b was then used to estimate the ER ref parameter in each of the 4 d windows by fitting Eq. (1). Following this, ER ref of each day was estimated by using the least-squares spline approximation (Lei and Yang, 2010).
To quantify the contribution of source areas to the CO 2 flux measurement of the eddy covariance, we used an analytical footprint model (Hsieh et al., 2000), where D = 0.28 and P = 0.59 are similarity constants for unstable condition (Hsieh et al., 2000), κ = 0.4 is von Karman constant, χ represents the horizontal coordinate, L represents the Obukhov length, z m represents the measurement height, and z u represents the length scale expressed as follows: where z 0 represents the roughness height set to be 0.1 H c (canopy height). Note that the eddy covariance system failed from 23 October 2010 to 1 April 2011 during the wheat dormant season. To evaluate the seasonal CO 2 budget of this rotation cycle, the flux gap of this period was filled by using the machine learning Support Vector Regression (SVR) algorithm (Cristianini and Shave- Taylor, 2000), which has been proved to be an appropriate tool for flux gap filling (e.g., Kang et al., 2019;Kim et al., 2019) (see Appendix A).

Meteorological and environmental condition measurements
The meteorological variables were measured at 30 min intervals by a standard meteorological station on the tower. Among these variables were the air temperature (T a ) and relative humidity (RH) (HMP45C, Vaisala Inc, Helsinki, Finland) at a height of 1.6 m and precipitation (P ) (TE525MM, Campbell Scientific Inc), incoming shortwave radiation (R si ) (CRN1, Kipp & Zonen, Delft, Netherlands), and photosynthetic photon flux density (PPFD) (LI-190SA, LI-COR Inc) at a height of 3.7 m. The 30 min interval edaphic measurements included soil temperature (T S ) (109-L, Campbell Scientific Inc.) and volumetric soil moisture (θ ) (CS616-L, Campbell Scientific Inc.) for the top 5 cm of the soil; soil matric potential (ψ) (257-L, Campbell Scientific Inc.) has been measured since 2010 at the same depth. The groundwater depth (WD) (CS420-L, Campbell Scientific Inc.) was measured at a location close to flux tower in 30 min intervals.

Biometric measurements and crop samples
To trace crop development and carbon storage, we measured canopy height (H c ), leaf area index (LAI), crop dry matter (DM) and carbon content of crop organs at an interval of 7-10 d in the footprint of eddy covariance. Due to inclement weather, measurement intervals were occasionally extended to 2 weeks or longer. The H c was measured with a ruler, and LAI was measured with LAI-2000 (LI-COR Inc.) at 10 locations randomly distributed in the field. For crop samples, four locations were randomly selected at the start of the growing season, and crop samples were then collected close to these four locations throughout the experimental period. At each location, 10 crop samples were collected for wheat and 3 crop samples were collected from maize. To reduce the sample uncertainty at harvest, 200 crops and 5 crops were collected in each location for wheat and maize, respectively. The crop organs were separated and oven-dried at 105 • C for kill-enzyme torrefaction for 30 min and then oven-dried at 75 • C until a constant weight. The crop samples were used to estimate the average field biomass (Dry Matter). The carbon content was analyzed using the combustion-oxidationtitration method (National Standards of Environmental Protection of the People's Republic of China, 2013) to estimate carbon storage. The crop samples provided a direct estimate of the NPP.

Soil respiration measurements
Soil respiration was measured every day in the footprint of the eddy covariance between 13:00 and 15:00 UTC+8 from March to September 2011 using a portable soil respiration system LI-8100 (LI-COR Inc.). Below-ground autotrophic respiration and heterotrophic respiration were differentiated using the root exclusion method . The total soil respiration (R S ) and R H were measured at treatments with and without roots, respectively, and the corresponding difference is R AB . To reduce the uncertainty associated with spatial variability, we set three replicate pairs of comparative treatments (i.e., with root and without root) randomly in the field. The uniform field condition contributes to reducing the measurement uncertainty associated with the spatial variability (see . To assess the seasonal variations and total amount of soil respirations, the seasonal continuous R H was constructed using the Q 10 model by incorporating soil moisture as follows : where θ f is the field capacity. The parameters were inferred by fitting the R H and T S measurements by using the least-squares method (see , where A = 1.16 µmol m −2 s −1 , B = 0.0503 • C −1 and a = −44.9 (unitless) (see . Note that the plant biomass was negligible before 14 March, during which R H was set to equal to the ecosystem respiration and the R AB was assumed to be 0. R AB of other periods was estimated based on the R H measurement and the ratio of R AB to R S estimated previously , and the continuous R AB /R S ratio was interpolated from the daily records (Fig. 2). This estimation method is robust because the R AB /R S ratio is nearly constant around its diurnal average (Zhang et al., 2015b).

Synthesis of the CO 2 budget components
The CO 2 budget components were derived by combining the eddy covariance measurements, soil respiration experiments and crop samples. Eddy-covariance-measured NEE is the difference between carbon assimilation (i.e., GPP) and carbon release (i.e., ER). The ER consists of R H , R AB (i.e., root respiration) and above-ground autotrophic respiration (R AA ). The total soil respiration is the sum of R H and R AB , The total autotrophic respiration (R A ) is the difference between the eddy-covariance-derived ER and R H , The above-ground autotrophic respiration (R AA ) is the difference between the eddy-covariance-derived ER and R S in Eq. (6), NPP is plant biomass carbon storage and can be quantified as the difference between GPP and R A , where the subscript "EC" represents that the NPP is estimated from the eddy-covariance-derived GPP. In parallel, NPP can also be directly inferred from biomass samples as follows: where the subscript "CS" indicates that NPP is based on crop samples and C cro is the plant biomass carbon storage at harvest. We used the average of the two independent NPPs as the measurement for this site. NEP is commonly estimated by the NEE measurement (NEP EC = −NEE). In this study, the crop samples and soil respiration measurements also provided an independent estimate as follows: We used the average of the two NEPs as the measurement for this site. At this site, there were no fire and insect disturbances and no manure fertilizer application. The carbon input from seeds was negligible, and all crop residues were returned to the field. Thus, NBP can be quantified as the difference between NEP and grain export carbon loss (C gra ), 3 Results

Meteorological conditions and crop development
The interannual variations in major meteorological variables are shown in Fig. 3, and they showed no clear trend for both wheat and maize seasons. For the full 2010-2011 cycle with comprehensive experiments, the average R si and T a were very close to other years; however, the P during maize season was a little higher than other years (Fig. 3c), leading to a shallow WD in the maize season (Fig. 3d). The intra-annual variations in field microclimates for the full 2010-2011 cycle are shown in Fig. 4. The seasonal maximum and minimum T a occurred in July and January, respectively, and the variations in vapor pressure deficit (VPD) followed the T a well. The WD mainly followed the irrigation events in winter and spring but followed P in summer and autumn. In particular, the WD varied from 0 to 3 m throughout the year. The wet soil conditions prohibited the field from experiencing water stress (Fig. 4d) because even the lowest soil matric potential (−187.6 kPa) remained a lot higher than the permanent wilting point of crops (around −1500.0 kPa). Figure 5 shows  the wheat season, the stages of regreening, jointing, booting, heading and maturity started approximately on 1 March, 20 April, 1 May, 7 May and 5 June, respectively. The seasonal variations in DM agreed well with the crop stages (Fig. 6), and the wheat biomass mainly accumulated in April and May, while maize biomass mainly accumulated in July and August. The total DM was 1718 g m −2 for wheat and 1262 g m −2 for maize at harvest. Upon harvest, the wheat DM was distributed as 3 % root, 43 % stem, 9 % leaf and 45 % grain, while the maize DM was distributed as 2 % root, 29 % stem, 7 % green leaf, 5 % dead leaf, 4 % bracket, 7 % cob and 46 % grain. The seasonal average carbon contents of the root, stem, green leaf, dead leaf, and grain were 410, 439, 486, 452, and 457 gC kg −1 DM for wheat and 408, 438, 477, 457, and 456 gC kg −1 DM for maize (see Table 1 for the seasonal variation).

The interannual variations in the NEE, GPP and ER
For the period from 2005 to 2016, if grain export was not considered, wheat was a consistent CO 2 sink, as the seasonal total NEEs were consistently negative, and maize was a CO 2 sink in most years, except for 2012 and 2013 when NEE was positive (Fig. 7a). NEEs of both wheat and maize fields became less negative during the past decade (though not in a statistically significant way), implying a progressive decline of the carbon sequestration potential of this cropland. The GPPs of both wheat and maize showed an increasing trend, though they were not statistically significant (Fig. 7b). The ERs of both wheat and maize also showed an increasing trend in these years, but only the trend of maize  was significant (Fig. 7c). The decadal average of NEE, GPP, and ER were −364 (SD ± 98), 1174 (SD ± 189), and 810 (SD ± 161) gC m −2 for wheat and −136 (SD ± 168), 1008 (SD ± 297), and 872 (SD ± 284) gC m −2 for maize. The NEE, GPP and ER for both wheat and maize were correlated with the three main environmental variables of R si , T a and WD using the multiple regression (see Appendix B for Figure 6. Seasonal variations in the total dry biomass (DM) and its major components of root, stem, green leaf and grain. The error bars of total biomass denote 1 standard deviation of the four sample points. details). In the wheat season, T a showed its relatively great importance (compared to R si and WD) to all three of the CO 2 fluxes with a higher T a increasing both GPP and ER and also enhancing NEE (more negative) (Fig. 8a). WD correlated negatively with GPP, thereby reducing net carbon uptake (less negative NEE). WD exhibited almost no effect on ER. R si exhibited almost no effect on all three CO 2 fluxes. Therefore, T a explained most of the interannual variations in NEE, GPP and ER, followed by WD. In the maize season, WD had good correlations with all three fluxes of GPP, ER and NEE, where a deeper WD contributed to lower both GPP and ER and also drove higher net carbon uptake (more negative NEE). T a showed almost no effect on all three CO 2 fluxes. R si had a positive correlation with ER but almost no correlation with GPP (Fig. 8b). Ultimately, higher R si in the maize season lowered the net carbon uptake (more positive NEE). Overall, R si and WD showed their great importance in influencing the interannual variation in maize NEE, with R si having a positive correlation and WD having a comparable negative correlation (Fig. 8b).

Intra-annual variations in the NEE, GPP and ER
The intra-annual variations in NEE, GPP and ER exhibited a bimodal curve corresponding with the two crop seasons (Fig. 9). All three CO 2 fluxes were almost in phase, with peaks appearing at the start of May during the wheat season and in the middle of August during the maize season. During some of the winter season, the field still sequestered a small amount of CO 2 because of the weak photosynthesis, which was confirmed by leaf level gas exchange measurement (data not shown). Net carbon emission happened during the fallow periods, in addition to the start of the maize season when the plant was small and high temperatures enhanced heterotrophic respiration. During the wheat season, two evident spikes appeared on 21 April and 8 May with positive NEE values (i.e., net carbon release). These spikes resulted from the radiation decline during the inclement weather (Fig. 4b), which suppressed the photosynthesis rate; similar phenomena also appeared during the maize season. Figure 10 shows the variations in ER and its components. During the wheat season, the variation in ER closely followed crop development and temperature, but there were two evident declines at the end of April and the start of May due to low temperatures associated with the inclement weather. During the early growing stage of maize, R H was the main component of ER. When waterlogging conditions occurred in late August and early September, both R H and R AB were suppressed to zero.

CO 2 budget synthesis in the 2010-2011 agricultural cycle
CO 2 budget analysis showed that this wheat-maize rotation cropland has the potential to uptake carbon from the atmosphere (Fig. 11). In the full 2010-2011 cycle, the total NEE, GPP, and ER values were −438, 1078, and 640 gC m −2 for wheat and −239, 780, and 541 gC m −2 for maize. The NPP values were 750 and 815 gC m −2 for wheat based on crop samples and the eddy covariance and that complemented with soil respiration measurements, respectively, and were 592 and 532 gC m −2 for maize based on the two methods. We used the average of these two methods for NPP measurements, which were 783 (SD ± 46) gC m −2 for wheat and 562 (SD ± 43) gC m −2 for maize. We also used the average of NEP from the two independent methods for the measurement, and the NEP was 406 gC m −2 for wheat and 269 gC m −2 for maize. Furthermore, when considering the carbon loss associated with the grain export, the NBP values were 59 gC m −2 for wheat and 5 gC m −2 for maize, respectively. Considering the net CO 2 loss of −104 gC m −2 during the two fallow periods, NBP of the whole wheat-maize crop cycle was −40 gC m −2 yr −1 , suggesting that the cropland was a weak carbon source to the atmosphere under these specific climatic conditions and field management practices.

Discussion
This study investigated the decadal variations in the NEE, GPP and ER for an irrigated wheat-maize rotation cropland over the North China Plain, and the results exhibited a decreasing trend of the CO 2 sink capacity during the past decade. The interannual variations in the carbon fluxes of wheat showed close dependence on temperature and groundwater depth, while those of maize were mostly regulated by solar radiation and groundwater depth. Furthermore, the detailed CO 2 budget components were quantified for a full wheat-maize agricultural cycle. Investigating the decadal trend of the CO 2 fluxes and quantifying the detailed CO 2 budget components for this representative cropland will provide useful knowledge for regional greenhouse gas emission evaluation over the North China Plain. Figure 10. Seasonal variations in the components of ecosystem respiration (ER), total soil respiration (R S ) and soil heterotrophic respiration (R H ). The difference between ER and R S denotes aboveground autotrophic respiration (R AA ), and the difference between R S and R H denotes below-ground autotrophic respiration (R AB ).

Comparison with other croplands
The cropland has been reported as carbon neutral to the atmosphere (e.g., Ciais et al., 2010), as a carbon source (e.g., Anthoni et al., 2004a;Verma et al., 2005;Kutsch et al., 2010;Wang et al., 2015;Eichelmann et al., 2016) and also as a carbon sink (e.g., Kutsch et al., 2010). Such inconsistency probably results from the different crop types and management practices (residue removal, the use of organic manure, etc.), in addition to variations in the climatic conditions (Béziat et al., 2009;Smith et al., 2014) and fallow period length (Dold et al., 2017). Our results show that the fully irrigated wheatmaize rotation cropland with a shallow groundwater depth was a weak CO 2 sink during both the wheat and maize seasons in the full 2010-2011 cycle, but the CO 2 loss during Figure 11. Carbon budget of wheat (a), maize (b) and the full wheat-maize rotation cycle, with fallow periods included (c). Note that the absolute value of NEE is shown here; NBPs of wheat and maize are the average of two independent methods (i.e, an eddy covariance-based method and a crop-sample-based method).
the fallow period reversed the cropland from a sink into a weak source with an NBP of −40 gC m −2 yr −1 . These results are consistent with previous studies that reported the wheatmaize rotation cropland as a carbon source (Li et al., 2006;Wang et al., 2015). However, the net CO 2 loss was much lower at our site, most likely due to the shallow groundwater depth. Field measurements of the long-term CO 2 fluxes over croplands remain lacking, and we found the carbon sequestration capacity of this cropland has been progressively decreasing, though it was not statistically significant. The cropland has been widely suggested as a climate change mitigation tool (e.g., Lal, 2001), but the potential in the future is challenging. However, since cropland management greatly impacts the carbon balance of cropland (Béziat et al., 2009;Ceschia et al., 2010), it remains required investigating if the management adjustment can foster the cropland carbon sink capacity over the long term.
The annual total NPP of 1345 gC m −2 yr −1 at our site is approximately twice the average of the model-estimated NPP for Chinese croplands (714 gC m −2 yr −1 ) with a rotation index of 2 (i.e., two crop cycles within 1 year) (Huang et al., 2007), more than 3 times the value estimated by MODIS (400 gC m −2 yr −1 ) (Zhao et al., 2005) and slightly higher than the value of the same crop rotation at the Luancheng site (1144 gC m −2 yr −1 ) (Wang et al., 2015). The higher NPP at our site may partially result from the sufficient irrigation and fertilization (Huang et al., 2007;Smith et al., 2014).
The contrasting respiration partitioning of the same crop in different regions (Table 2) indicate that the respiration processes may also be subject to climatic conditions and management practices. Though the ratio of ecosystem respiration to GPP at our site is comparable to other studies, the ratio of autotrophic respiration to GPP is much lower at our site, and while the ratio of heterotrophic respiration to ecosystem respiration is greater at our site, these findings are different from those at the other sites with similar crop variety (Moureaux et al., 2008;Aubinet et al., 2009;Suleau et al., 2011;Wang et al., 2015;Demyan et al., 2016), as they showed that ecosystem respiration is usually dominated by below-ground and above-ground autotrophic respirations. The higher soil heterotrophic respiration at our site probably results from the full irrigation and shallow groundwater, which both alleviate soil water stress.

The effects of groundwater on carbon fluxes
The groundwater table at our site is much closer to the surface because of the irrigation by water diverted from the Yellow River. In contrast, the nearby Luancheng site (Wang et al., 2015) is groundwater-fed with a very deep groundwater depth (approximately 42 m) (Shen et al., 2013), and their CO 2 budget components had some differences with our study. Comparing the net CO 2 exchange of wheat, the GPP at our site is a little higher than the Luancheng site, implying the irrigation at our site may better sustain the photosynthesis rate for wheat; ER at our site is also a little higher than the Luancheng site. For maize, both sites are not irrigated due to the high summer precipitation. GPP and ER at our site were comparable to Luancheng site, implying that the irrigation method prior to the maize season had no discernible effect on the integrated CO 2 fluxes for maize. However, the three components of ER in our study showed pronounced differences from the Luancheng site, where they reported the R AA was 411 gC m −2 for wheat and 428 gC m −2 for maize, 3 times the results of our study (128 gC m −2 for wheat and 133 gC m −2 for maize). However, their R AB for wheat (36 gC m −2 ) and maize (16 gC m −2 ) were less than a quarter of our results (136 gC m −2 for wheat and 115 gC m −2 for maize). Their R H of wheat (245 gC m −2 ) was less than our estimate (377 gC m −2 ), but R H of maize (397 gC m −2 ) was greater than our result (292 gC m −2 ). In general, the above-ground crop parts in our site respired more carbon than the Luancheng site, possibly because the shallow groundwater depth at our site increased the above-ground biomass allocation but lowered the root biomass allocation (Poorter et al., 2012). These independent cross-site comparisons demonstrate that carbon budget components may be subject to the specific groundwater depth influenced by the irrigation type, and even the same crop under similar climatic conditions can behave differently in carbon consumption. Our site experienced a short period of waterlogging during the 2010-2011 cycle due to the combined effects of full irrigation and the high precipitation during the summer. This distinct field condition reduced soil carbon losses in the maize season, potentially maintaining the CO 2 captured by the cropland. Waterlogging events were occasionally reported in upland croplands. For example, Terazawa et al. (1992) and Iwasaki et al. (2010) suggested that waterlogging causes damage to plants, resulting in a decline in GPP as reported by Dold et al. (2017) and our study. Our study further shows that waterlogging reduces ER to a greater degree than GPP, possibly because of the low soil oxygen conditions, and thereby reduces the overall cropland CO 2 loss. However, the CH 4 released over the short term may be pronounced in waterlogged soils. As CH 4 emission in this kind of cropping system over the North China Plain cropland remains lacking, additional field experiments are required to understand how irrigation and water saturation field condition impact the overall carbon budget.

Uncertainty in the estimation and limitation of this study
In the comprehensive experiment period for the full 2010-2011 agricultural cycle, the NEE of the wheat season from 23 October 2010 to 1 April 2011 was calculated using a calibrated SVR model. The SVR model performs well for predicting GPP and ER with very high R 2 of 0.95 and 0.97 and an acceptable uncertainty level of 22.9 % and 15.2 % for GPP and ER, respectively. Hence, these estimates should have a negligible effect on the seasonal total carbon evaluation. The footprint analysis showed that 90 % of the measured eddy flux comes from the nearest 420 and 166 m in wheat and maize crops under unstable conditions, respectively, confirming that both soil respiration experiments and crop samples paired well with the EC measurements . Root biomass was difficult to measure, but the uncertainty should be low because the root ratio (the ratio of the root weight to the total biomass weight) accounts for 15 %-16 % of the crop for wheat and maize (Wolf et al., 2015), and our measurements are very close to these values; i.e., the averaged seasonal root ratio was 15 % for wheat and 10 % for maize at our site. However, the relatively low root ratios (3 % for wheat and 2 % for maize) at harvest probably result from the root decay associated with plant senescence. The estimates of annual soil respiration are based on the Q 10 model validated by the field measurements that may generate some uncertainty in the soil respiration budget due to the hysteresis response of soil respiration to temperature (Phillips et al., 2011;Zhang et al., 2015aZhang et al., , 2018. However, the Q 10 model remains robust in soil respiration estimations if it is well validated (Tian et al., 1999;Latimer and Risk, 2016), allowing for confidence in the estimates. During the wheat season, the cumulative curves of NPP EC and NPP CS were not perfectly consistent in the main growing season, as clear differences emerged during the dormant season of wheat from 15 December 2010 to 8 March 2011 (Fig. 12). These differences may result from the small wheat sample number. However, the sample number at harvest was sufficiently big, and no discernible difference was found between the two NPPs at harvest. These two independent estimates of NPP were similar throughout the maize season (Fig. 12).
This study provides a comprehensive quantification of the CO 2 budget components of the cropland, but it remains limited to a relatively wet year (see Fig. 3c and d). The integrated carbon fluxes (NEE, GPP and ER) have pronounced interannual variations, also suggesting further investigations are required on the interannual variations in the carbon budget components.

Conclusion
Based on the decadal measurements of CO 2 fluxes over an irrigated wheat-maize rotation cropland over the North China Plain, we found the cropland was a strong CO 2 sink if grain export was not considered. When considering the grain export, the cropland was a weak CO 2 source, with an NBP of −40 gC m −2 yr −1 in the full 2010-2011 agricultural cycle. The net CO 2 exchange during the past decade from 2005 to 2016 showed a statistically nonsignificant decreasing trend, implying a decreasing carbon sequestration capacity of this cropland, discouraging the potential of taking agroecosystems as the mitigation tool of climate change. In the wheat season, air temperature showed the best correlation with the CO 2 fluxes followed by the groundwater depth, whereas in the maize season both shortwave radiation and groundwater depth showed good correlation with the CO 2 fluxes. The comprehensive investigation showed most of the carbon sequestration occurred during the wheat season, while maize was close to being CO 2 neutral. Soil heterotrophic respiration in this cropland contributes substantially to CO 2 loss in both the wheat and maize seasons. This study provides detailed knowledge for estimating regional carbon emissions over the North China Plain. The Support Vector Regression (SVR) method is a machine learning technique-based regression, which transforms regression from nonlinear into linear by mapping the original low-dimensional input space to higher-dimensional space (Cristianini and Shave-Taylor, 2000). The SVR method has two advantages: (1) the model training always converges to global optimal solution, with only a few free parameters to adjust, and no experimentation is needed to determine the architecture of SVR. (2) The SVR method is robust to small errors in the training data (Ueyama et al., 2013). The Support Vector Machine (SVM) software package obtained from LIBSVM (Chang and Lin, 2005) is used in this study.

A2 Data processing and selection of explanatory variables
Gross primary productivity (GPP) is influenced by several edaphic, atmospheric and physiological variables, among which air temperature (T a ), relative humidity (RH), leaf area index (LAI), net photosynthetically active radiation (PAR) and soil moisture (θ ) are the dominant factors. Hence, we select T a , RH, LAI, PAR and θ as explanatory variables of GPP. Ecosystem respiration (ER) consists of total soil respiration and above-ground autotrophic respiration. The total soil respiration is largely influenced by soil temperature and soil moisture, while above-ground autotrophic respiration is largely influenced by air temperature and aboveground biomass. Therefore, we select T a , soil temperature at 5 cm (T S ), θ and LAI as explanatory variables of ER. LAI is estimated from the Wide Dynamic Range Vegetation Index derived from the MOD09Q1 reflectance data (250 m, 8 d average, https://ladsweb.modaps.eosdis.nasa. gov/missions-and-measurements/products/MOD09Q1/, last access: 15 May 2016; see Lei et al., 2013). The three wheat seasons of 2005-2006, 2009-2010 and 2010-2011 are selected for model training, and the original half-hourly measurements of GPP and ER, together with the explanatory variables, are averaged to the daily scale, but we remove days missing more than 25 % of half-hourly data. We have GPP available from 466 d and ER from 483 d for model training. The explanatory variables for the equipment failure are also averaged into daily scale, which will be used to calculate GPP and ER with the trained model described in the following section.

A3 SVR model training and flux calculation
In order to eliminate the impact of variables with different absolute magnitudes, we rescale all the variables in the training data set to the [0, 1] range prior to SVR model training. In the training process, the radial basis function (RBF, a kernel function of SVR) is used and the width of insensitive error band is set as 0.01. The SVR model training follows these steps: 1. All training data samples are randomly divided into five nonoverlapping subsets, and four of them are selected as the training sets (also calibration set); the remaining subset is treated as the test set (or validation set). This process is repeated five times to ensure that every subset has a chance to be the test set.
2. For the selected training set, the SVR parameters (cost of errors c and kernel parameter σ ) are determined using a grid search with a five-fold cross-validation training process. In this approach, the training set is further randomly divided into five nonoverlapping subsets. Training is performed on each of the four subsets within this training set, with the remaining subset reserved for calculating the root-mean-square error (RMSE), and model parameters (c and σ ) yielding the minimum RMSE value are selected.
3. The SVR model is trained based on the training set from step (1) and initialized by the parameters (c and σ ) derived from step (2).
4. The test set from the step (1) is used to evaluate the model obtained from the step (3) by using the coefficient of determination (R 2 ) and RMSE. 5. The model is trained with all of the available samples that achieved good performance, as R 2 are 0.95 and 0.97 for GPP and ER, respectively, and the mean RMSE is 1.28 and 0.44 gC m −2 d −1 . The RMSE can be further used as a metric quantifying uncertainty, which accounts for 22.9 % and 15.2 % for the averaged GPP and ER, respectively. GPP and ER during the equipment failure period are then calculated with the trained model complemented with the observed explanatory variables, and NEE is derived as the difference of GPP and ER.

Appendix B: Multiple regression for NEE, GPP and ER with microclimate variables
The flux of NEE, GPP or ER is correlated with incoming shortwave radiation (R si ), air temperature (T a ) and groundwater depth (WD) as flux = aR si + bT a + cWD + d, where flux is NEE, GPP or ER; a, b, c and d are regression parameters. All the variables are normalized to derive their z score before the regression, where the z score is calculated by subtracting the mean from the data and dividing the result by the standard deviation. The coefficient of each variable represents the relative importance of the corresponding variable in contributing to the dependent variable.
Q. Zhang et al.: Carbon budget over wheat-maize cropland Data availability. The data of this study are available to the public by request to the corresponding author (Huimin Lei).
Author contributions. QZ and HL designed the study and methodology with substantial input from all co-authors. QZ conducted the field experiment. BF conducted the SVR calculation for gap filling. All authors contributed to interpretation of the results. QZ drafted the manuscript, and all authors edited and approved the final manuscript.