Interannual Variation in Carbon Sequestration Depends Mainly on the Carbon Uptake Period in Two Croplands on the North China Plain

Interannual variation in plant phenology can lead to major modifications in the interannual variation of net ecosystem production (NEP) and net biome production (NBP) as a result of recent climate change in croplands. Continuous measurements of carbon flux using the eddy covariance technique were conducted in two winter wheat and summer maize double-cropped croplands during 2003–2012 in Yucheng and during 2007–2012 in Luancheng on the North China Plain. Our results showed that the difference between the NEP and the NBP, i.e., the crop economic yield, was conservative even though the NEP and the NBP for both sites exhibited marked fluctuations during the years of observation. A significant and positive relationship was found between the annual carbon uptake period (CUP) and the NEP as well as the NBP. The NEP and the NBP would increase by 14.8±5.2 and 14.7±6.6 g C m−2 yr−1, respectively, if one CUP-day was extended. A positive relationship also existed between the CUP and the NEP as well as the NBP for winter wheat and summer maize, respectively. The annual air temperature, through its negative effect on the start date of the CUP, determined the length of the CUP. The spring temperature was the main indirect factor controlling the annual carbon sequestration when a one-season crop (winter wheat) was considered. Thus, global warming can be expected to extend the length of the CUP and thus increase carbon sequestration in croplands.


Introduction
The potential of cropland to mitigate the greenhouse effect has received much attention in recent years [1]. Cropland can serve as an accumulator of carbon and can sequester approximately 50-100 Tg C each year [2]. The North China Plain, which is the largest agricultural production center in China, has an area of 3610 5 km 2 and supports the production of food for over 300 million people [3]. This plain is a large area planted in winter wheat and summer maize double-cropping systems and supplies more than 50% of China's wheat and 33% of its maize [4]. Owing to its large area, this agricultural center has great potential to store carbon and reduce the emissions of greenhouse gases. Describing the strength of the carbon sink or source in this area will aid in the assessment of the regional greenhouse gas balance, in the development and validation of models that simulate a crop's carbon cycle and in the evaluation of remotely sensed data [5]. Long-term continuous measurements of carbon flux using eddy covariance technology is a practical method to assess the carbon balance and can sufficiently reflect the interannual variations in carbon exchange, which would provide scientists with valuable information for understanding the response of ecosystem processes to climate change [6].
There is an on-going debate regarding whether cropland acts as a carbon sink or as a carbon source [7]. According to published studies, the long-term variation in the ecosystem carbon exchange is usually described by the standard deviation of annual flux sum or the coefficient of variations (i.e. the ratio of standard deviation and mean value of annual flux sum) [8,9]. In agricultural systems, annual carbon sequestration is represented by the net biome carbon production (NBP) but not the net ecosystem carbon production (NEP) [10]. The NBP in cropland is mainly consisted of three parts: the carbon absorbed by plants as NEP, carbon inputs such as the application of fertilizer or pesticide and carbon outputs resulting from harvesting and machinery [11]. Compared with the amount of carbon lost through harvesting, other forms of carbon input or emission constitute a small fraction of the NBP [12]. Owing to the large amount of carbon stored in grain, the NBP has been reported as being close to zero (i.e., carbon neutral) [13][14][15][16][17], greater than 100 g C m 22 y 21 (i.e., a carbon sink) [4,12,18,19] or significantly less than zero (i.e., a carbon source) [12].
Although most previous studies that examined the seasonal variations mechanisms in agroecosystem carbon exchange were based on less than 3 years of datasets [1,[20][21][22][23][24], more recent researches have been conducted on crops to interpret the interannual variations in carbon sequestration using long-term observations of eddy covariance data [17,25,26]. These studies indicated that the interannual variations in agroecosystem carbon fluxes can easily be explained based on the annual variations in environmental factors such as temperature [17] and precipitation [14]. Moreover, climate factors are strong indicators of spatial patterns in terrestrial ecosystem carbon fluxes [27,28].
Vegetation phenology, which aims to study the timing and length of the growing season in terrestrial ecosystems and their relationship with the climate, has become a meaningful earth systems science [29]. Plant phenology has been shown to be affected by climate changes and, in turn, influences ecosystem processes such as the carbon cycle [30]. As an intermediary of ecosystems and climate change, variations in plant phenology can be used to interpret the seasonal and interannual variations in carbon fluxes and their relationship with climate change [31]. The effects of variation in plant phenology on the annual net ecosystem carbon exchange have been reported in different ecosystems at local [32] and continental scales [33]. These studies indicated that the annual net ecosystem carbon storage was significantly related to the growing season length (GSL) [34][35][36][37][38][39][40], the CUP [37,41] and the leaf area index (LAI) [31]. According to published studies, the impact of the interaction between climate change and the abovementioned typical land surface phenological variables on annual ecosystem carbon sequestration is variable [9,31,32,42,43]. Most studies have focused on one season per year for forest and grassland ecosystems; it has been suggested that climate factors, especially air temperature, would affect the starting [44] or the ending date of the growing season [43] or the time interval in days between the gross ecosystem production (GEP) and the CUP [42] and would ultimately influence the annual ecosystem carbon budget.
Many studies have been conducted in the critical area of the North China Plain to explore temporal carbon budget variations based on continuous eddy covariance datasets. Li et al. (2006) [19] interpreted the mechanisms underlying the seasonal variations in crop carbon exchange using 2 years of flux data. Shen et al. (2013) [31] quantified the carbon fluxes of crops under different water conditions based on 4 years of data. Tong et al. (2014) [45] studied the light response characteristics of crop carbon exchange using 4 years of flux data. In the present study, continuous measurements of carbon flux using eddy covariance were conducted in two winter wheat/summer maize double-cropped croplands from 2003 to 2012 in Yucheng and from 2007 to 2012 in Luancheng. Our objectives were to (1) quantify the interannual variability in the NEP and the NBP and (2) investigate the environmental and biotic controls over the interannual variations in the NEP and the NBP.

Ethics Statement
The two sites in this study are maintained by the Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences and the Institute of Genetics and Developmental Biology, Chinese Academy of Sciences. Both areas are practice bases for researchers at the Chinese Academy of Sciences. All necessary permits were obtained for the fields study. The field study did not involve endangered or protected species. Nitrogen fertilizer was applied to the farmlands after the planting date. The crops were irrigated using well water during the reviving and jointing stages of winter wheat and during planting or the jointing stage of summer maize. The crops were irrigated with 100-150 mm of water on each irrigation date [31]. In Yucheng, the ground water table varies from 1.5 to 3.5 m with an average of 2.5 m, the soil in the depth of 1,20 cm is consist of clay loam, silt loam and sand loam with the portion of 22.1%, 65.1%, 12.8%, respectively. In Luancheng, the soil in the depth of 1,25 cm is predominated by sand loam. The bulk density of the soil is 1.3 gcm 23 [56] and 1.4 gcm 23 in Yucheng and Luancheng, respectively [57].

Flux and meteorological measurements
Similar monitoring instrumentation and data collection methods were used in Yucheng and Luancheng. The eddy covariance system consisted of a three-dimensional sonic anemometer (Model CSAT 3, Campbell Scientific Inc., USA) to monitor fluctuations in vertical wind velocity and an open-path and fast-response infrared gas analyzer (Model LI-7500, Li-Cor Inc., USA) to monitor the fluctuations in the concentrations of CO 2 and water vapor. The flux towers were located in the center of the crop fields at both sites, and the fetch was greater than 200 m. The eddy covariance (EC) systems were set on the tower at a height of 2.8 m in Yucheng and a height of 3.5 m in Luancheng. All of the raw data were collected continuously at 10 Hz using a data logger (Model CR5000, Campell Scientific Inc., USA), and the 30 min mean data were outputted.

Flux calculation, gap filling and partitioning
The half-hourly turbulent CO 2 flux (F c , mmol CO 2 m 22 s 21 ) above the crop canopy, i.e., the net ecosystem carbon exchange (NEE, mmol C m 22 s 21 ), can be calculated from the covariance between the vertical wind velocity (w, m s 21 ) fluctuation and the CO 2 density fluctuation ( r c , mmol CO 2 m 23 ) according to Reynolds' decomposition rule: where the primes denote the turbulent fluctuations (departures from the mean) and the overbar indicates a time-averaged mean (30 min). Several procedures have been performed to correct the 30 min mean output data. First, to satisfy the hypothesis of the eddy covariance technology, a tilt correction was necessary for the systems. The tilt error caused by the nonparallel fixation of the anemometers with the field surface was corrected by double rotation to align u with the mean wind direction, forcing v to 0 and the mean vertical velocity (w) to 0 [47]. Second, as with all measuring instruments, the eddy covariance system is not able to completely capture the true turbulence when sufficiently high and low frequencies occur, which results in the loss of information compared with ideal conditions. These missing flux data can arise for several reasons, such as an inadequate sensor frequency response, separation of the instruments (particularly the sonic anemometer and the infrared gas analyzer), line averaging and distributed sampling [48]. Therefore, the spectral correction method of Eugster and Senn (1995) [48] was indispensable in this study to compensate for the missing raw covariances. Finally, the Webb-Pearman-Leuning (WPL) correction was applied to correct the error caused by the transfer of heat and water vapor because fluctuations in temperature and humidity can cause variation in the density of trace gases [49]. During monitoring, there were occurrences of sensor malfunction, rain events, instrument maintenance and power failures, which resulted in missing or anomalous values of the flux and the meteorological data. For the meteorological data, we deleted the data over the normal range, and the data which has deviated more than 1.96 standard deviation of a time series from the mean value. Then, we applied a linear interpolation method to fill in the missing data that had a time interval less than 2 h and the mean diurnal variation (MDV) method to fill in the missing data that had a time interval longer than 2 h [50]. In addition, missing flux data resulted from the exclusion of these data under lowturbulence conditions at night when the wind friction velocity (u * ) was less than the threshold value. The threshold value was determined from the relationship between u * and the nighttime NEE. The total missing nighttime flux data is about 36%63% (the mean and standard deviation of ten years) and 43%62% (the mean and standard deviation of five years) of annual total flux data in Yucheng and Lucheng, respectively. The total missing daytime flux data was about 14%62% and 24%67% of annual total flux data in Yucheng and Luancheng, respectively. Linear interpolation was also used to fill the flux data gaps shorter than 2 h. The marginal distribution sampling (MDS) method was used to fill the flux data gaps longer than 2 h, by which the missing flux data could be ''looked up'' on the basis of environmental factors such as the photosynthetic photon flux density and the temperature correlations with the missing data [51].
The NEE is the sum of the gross ecosystem carbon exchange (GEE) and the ecosystem respiration (RE). The daytime GEE could not be obtained directly from the EC systems, but it could be obtained as an estimate of RE. In general, the nighttime RE is equal to the nighttime NEE because the nighttime GEE is zero. By using the nighttime flux data under high turbulence, a regression model (the LT model) described by Lloyd and Taylor (1994) [52] was used to fit the model parameters R ref and E 0 , and then the new model was used to estimate the daytime RE. The Lloyd and Taylor model is as follows: where R ref denotes the ecosystem respiration rate at a reference temperature (T ref = 10uC), E 0 is a parameter associated with the activation energy and determines the temperature sensitivity of RE and T 0 is a constant temperature parameter (246.02uC) [52]. Because the crop LAI changed dramatically, the RE was affected by the temperature and the crop phenology. To consider the effects of phenology on the calculation of RE, the long-term dataset was divided into a series of short subperiods (the window was 15 days), which were shifted by 5 days, i.e., the overlap between adjacent windows was 10 days. The regression of the LT model was performed separately for each subperiod [51].
In our analysis, we used the net ecosystem production (NEP) rather than the NEE (NEP = -NEE) and the gross ecosystem production (GEP) rather than the GEE (GEP = -GEE). The sign of the NEP is positive when CO 2 is transported from the atmosphere to the ecosystem and is negative when CO 2 is transported in the opposite direction.

Calculations of the plant phenological factors
In farmland, the growing season length (GSL) is usually defined as the total number of days from the sowing date to the harvest date [53]. The carbon uptake period (CUP) could be calculated only by using the eddy covariance data because the CUP was defined as the number of days when the ecosystem acts a carbon sink [32]. We determined the CUP as the number of days when the agroecosystem is a net carbon sink during the growing season [13], and we calculated this value from the smoothed moving average of the NEP [32]. The leaf area index (LAI, m 2 m 22 ) of the crop was measured using an area meter (LI-3100, Li-Cor, NE, USA) 1-4 times per month during the growing season.

Calculation of the net biome production (NBP)
In agroecosystems, the net biome production (NBP) can be obtained from the following equation: where C gr refers to the carbon content of the grains, which is a product of the harvest, and C ag refers to the other forms of carbon transfer, including diesel combustion, the application of fertilizer and pesticide and the transport of inputs and grains. C gr was estimated from the following: where W gr is the grain water content, i.e., 0.140 for wheat and 0.155 for maize; f c is the fraction of the carbon in the grain, i.e., 0.45 for wheat and 0.447 for maize; and Y is the crop yield [20]. It is worth noting that Bernacchi et al. (2005) [7] indicated that the C ag accounted for approximately 7% and 10% of the NBP in a soybean crop and a corn crop, respectively, when all of the carbon transformation events were accounted for. In fact, according to the experimental records of the Chinese Ecosystem Research Network (CERN), the workers in Yucheng and Luancheng manually weed the fields rather than apply pesticide. In addition, they also use limited nitrogen fertilizer (only 1-2 times per growing season), thereby releasing very little CO 2 to the atmosphere. The carbon released by fuel combustion and machinery transportation only occurs during the harvest. As a result, we inferred that the C ag would account for a small fraction of the NBP and we ignored it when calculating the NBP in this study.

Data statistics and analysis
To determine the dominant factor controlling the interannual variations in carbon exchange, several procedures were adopted to analyze the dataset in this study. The initial variables used to interpret the variation in the mechanisms of the carbon flux were environmental factors that included the air temperature, the soil temperature at a depth of 5 cm, the precipitation, the photosynthetic photon flux density (PPFD) and phenological factors such as the GSL, the LAI and the CUP. First, we applied the stepwise method in the multiple regression by PASW statistic 18 software (2010, ver18.0; SPSS Inc) to select the significant independent variables (p,0.05). The independent variable(s), which entered multiple regression models, would have significant effect on the dependent variables. And then, the path analysis was performed to interpret the relationship between the independent variables and their direct and indirect effects on the interannual variations in the NEP and NBP. The path analysis were conducted by Amos (2003, ver5.0; SPSS Inc).  Table 1 also lists the annual values of these variables. The seasonal variation in the daily mean air temperature exhibited a single-peak curve, and the maximum and minimum values occurred from July to August and from November to the next January, respectively. The annual air temperature ranged from 11.5 to 13.9uC (12.960.7uC, mean6standard deviation) and 8.0 to 12.7uC (10.462.2uC) in Yucheng and Luancheng, respectively. The seasonal trend in the soil temperature was similar to the seasonal trend in the air temperature, but the daily soil and air temperature were not equal to each other. The annual soil temperature was close to the annual air temperature in Yucheng, the annual soil temperature was higher than the air temperature (13.160.4uC) in Luancheng. Because the straw residue of the summer maize remained in the field after harvesting, the rate of the decrease in the air temperature was more dramatic than the soil temperature when extreme cold weather occurred during the years of 2008-2009 and 2011-2012 in Luancheng. The asynchronous changes in the air and the soil temperature resulted in different decline amplitudes in Luancheng. Generally, the precipitation frequency in the second half of the year was higher than that in the first half of the year. The annual precipitation was 528.36197.2 mm and 399.56170.0 mm in Yucheng and Luancheng, respectively. The seasonal variation in the PPFD also exhibited a single peak, with the strongest radiation in April or May each year. The total annual PPFD was 8250.66535.6 mol m 22 in Luancheng, which was higher than that in Yucheng (7670.16369.0 mol m 22 ). The LAI is considered to be an indicator of crop development (Fig. 1c). The LAI was much lower during the germination stage and then gradually increased to reach a maximum during the bloom stage. The mean annual LAI ranged from 0.

Interannual variations in carbon fluxes and net biome production (NBP)
The seasonal and interannual variations in the daily NEP, GEP and RE in Yucheng and Luancheng are presented in Fig. 2. The 10-day running mean of the NEP, GEP and RE are also shown. In winter, the carbon fluxes were close to zero from the sowing stage in October to the reviving stage the next March due to the slow growth of winter wheat. The leaves photosynthesized at a low rate because of the low LAI, and there were low rates of respiration because of the low temperatures. Moreover, there were not sufficient photosynthetic substrates for autotrophic respiration. However, it is worth noting that during the winter from 2008 to 2009 in Luancheng, a portion of the wheat seed have died as the result of the low air temperature, and the reduction in photosynthesis converted the crop ecosystem to a carbon source. With an increase in the temperature in spring, the winter wheat began to grow rapidly, and the dynamics of the GEP and the RE corresponded to the changes in the LAI. The NEP became positive, and the crops started to absorb carbon. The NEP, GEP and RE had almost reached their maximum value simultaneously in late April when the winter wheat LAI was at its maximum. However, the maximum value of the RE of the winter wheat was delayed compared with the GEP in 2010 in Luancheng. The lag may be attributed to both environmental and phenological factors. The carbon fluxes decreased because of plant senescence. The seasonal dynamics of the summer maize carbon fluxes also corresponded with the seasonal trends in the LAI. Because the germination stage of summer maize is shorter than that of winter wheat, the summer maize began to absorb carbon in mid-July (approximately one month after sowing). The NEP of the summer maize increased dramatically after the germination stage and usually reached its maximum value in mid-August.
The interannual variations in annual and seasonal NEP and NBP are showed in Fig.3 Fig. 4 shows the dependence of the annual NEP and NBP on the mean annual air temperature (MAT) and the annual carbon uptake period (CUP). Table 1 also lists the annual values of the growing season length (GSL), the CUP and the LAI. The interannual variations in the NEP and NBP were significantly correlated with both the MAT (R 2 = 0.78 for the NEP and 0.72 for the NBP, p,0.001) (Fig. 4a) and the annual total CUP (R 2 = 0.74 for the NEP and 0.64 for the NBP, p,0.001) (Fig. 4b) when the two sites were combined. The starting date of CUP (CUP start-date ) had a strongly effect on annual CUP (R 2 = 0.78, p,0.0001 for two sites combined; R 2 = 0.78, p,0.01 in Yucheng; R 2 = 0.67, p, 0.0001 in Lucheng) (Fig. 5), and the CUP start-date was significantly affected by the MAT (R 2 = 0.49, p = 0.04 for sites combined; R 2 = 0.49, p = 0.025 in Yucheng; R 2 = 0.047, p = 0.21 in Lucheng) (Fig. 6). No relationship was observed between the annual NEP or the NBP and the GSL or between the NEP or the NBP and the mean annual LAI. The annual total GSL was very close to the length of an entire year because of the short duration (only several days) of the fallow stage (this refers to the interval stage between the two crops). Compared with the GSL, the crop CUP experienced a large fluctuation between years. The CUP of the winter wheat ranged from 78 to 108 days (88.9610.6 days) in Yucheng and from 64 to 97 days (82.0611.8 days) in Luancheng. The CUP of the summer maize was shorter than that of the winter wheat, which was 74.066.0 days and 60.0611.6 days in Yucheng and Luancheng, respectively. As a result, the annual total CUP, which was the sum of the CUP of the winter wheat and the CUP of the summer maize, increased from 147 to 186 days with a longterm mean of 162.3612.5 days in Yucheng and from 116 to 161 days with a mean of 142.2617.9 days in Luancheng. The linear regression models suggested that the NEP would increase by 144.1647.2 g C m 22 yr 21 per additional degree Celsius and by 14.865.2 g C m 22 yr 21 per additional CUP day. The NBP also increased by a similar amount when one additional unit temperature was added and the CUP was extended (Fig. 4). Despite both the NEP and NBP being affected by the CUP, the difference between them, i.e., the crop yield was not obviously influenced by the CUP and was conservative over a long time scale (Fig. 4c).

Factors affecting the interannual variation in carbon sequestration
Given that the annual carbon sequestration of the crops was affected by both the MAT and the annual CUP, it was necessary to distinguish which was the main driving factor. Fig. 7 demonstrates the direct and indirect effect of the environmental and phenological factors on the control of the interannual variations in the NEP and the NBP in croplands. The annual NEP and NBP were controlled mainly by the annual CUP and had path coefficients of 0.87 and 0.80, respectively. Because of the significant relationship between the MAT and the annual CUP (the path coefficient was 0.88), the MAT had an indirect effect on the NEP and the NBP through its effect on the interannual variation in the CUP (Fig. 7a, d). This finding indicated that the annual CUP was the primary driving factor. In Yucheng, the CUP explained 40% (p = 0.049) of the interannual variation in the NEP and 19% (p = 0.21) of the variation in the NBP. In Luancheng, the values were 96% (p = 0.009) and 85% (p = 0.03), respectively. No relationship was observed between the MAT and the NEP or the NBP in Yucheng, but the MAT in Luancheng explained 95% (p = 0.004) and 99% (p,0.0001) of the interannual variation in the NEP and the NBP, respectively. The better linear relationship between MAT or CUP and NEP or NBP in Luancheng than Yucheng may result from the different soil texture and water related soil properties in both sites.
The factors controlling the interannual variation in the NEP of the winter wheat were the 5-cm depth soil temperature (Ts) of winter wheat and the winter wheat CUP, both of which were statistically significant (p,0.05). No relationship was found between the winter wheat soil temperature and the winter wheat CUP (R 2 = 0.03, p.0.05); hence, the path coefficient was small (2 0.17; the negative sign represents a negative relationship). Therefore, no indirect effect of the winter wheat soil temperature on the winter wheat NEP through the winter wheat CUP was possible. Because the path coefficient of the relationship between the winter wheat CUP and the NEP (0.61) was larger than the path coefficient of the relationship between the winter wheat Ts and the NEP (20.48) (Fig. 7b), the CUP was the predominant factor controlling the interannual variation in the winter wheat NEP. Only the winter wheat CUP was adopted in the path analysis and this was considered the single factor that controlled the interannual variation in the winter wheat NBP (Fig. 7e). Similarly, the summer maize CUP was the only factor that controlled the summer maize NEP and the NBP (Fig. 7c and f).

Annual carbon budget and comparison with other terrestrial ecosystem
Our study showed that the annual NEP was 475.66159.  There are wide fluctuations in the annual carbon budget between different agroecosystems in other regions ( Table 2). A series of studies have suggested that the spatial pattern of carbon flux in terrestrial ecosystems (mainly refer to forest and grassland) is driven by climate factors such as temperature or precipitation [12,28]. Similarly, the spatial pattern of the NEP for different crop sites was also clearly explained by the mean air temperature (R 2 = 0.40 p = 0.027) (Fig. 8). However, in the present study, the NEP in Luancheng deviated from the regression line because the crops suffered from abnormal climate conditions in late 2008 and at the beginning of 2009, as noted above, and much carbon was lost in both years. A similar result was not observed for the NBP. The combined effect of the differences in grain removal, field management and environmental conditions influenced the variations in the NBP between sites [54].

Arguments about uncertainty
There was a systemtematic uncertainty associated with choosing any algorithm for gap-filling. A series of studies have discussed the effects of gap-filling method on the sum flux [58,59]. Moffat et al. (2007) [59] found that in most cases, the algorithms being used were approaching the noise limit (uncertainty) of the measurement. However, highly empirical approaches, such as marginal distribution sampling (MDS), would performed the best across a range of forested European sites. At annual scale, difference among methods were generally modest, because most calculated annual NEE integrals were within 625 gCm 22 y 21 of the mean. According to those analyses, it is possible to suggest that the uncertainty of preference of gap-filling method is relatively small when one of the highly empirical approaches is used if the gap is  not too long and if the data set is large enough and have good quality to get model parameters.
Leaving crop residues (e.g., stems, leaves and roots) in the field is an important farmland management practice. The impact of crop residues on the change in seasonal carbon fluxes can be observed in Fig. 2. After the winter wheat harvest in mid-June, the RE did not decrease to zero but exhibited peaks. This pattern occurred because the suppression of the autotrophic respiration of the winter wheat resulting from the harvest was compensated for the increase in the heterotrophic respiration caused by leaving the residues in the field. Moreover, the high temperatures and the rain events in summer might also contribute to this increase in the RE. In this study, we did not estimate the carbon released from the decomposing crop residue. A series of studies have evaluated the contribution of crop residue carbon emissions to total ecosystem respiration after complete decomposition. These studies suggested that the additional emissions from the crop residues may account for a considerable part of the total carbon release. Moureaux et al. (2008) [23] indicated that the carbon content of winter wheat residue was 290 g C m 22 yr 21 , which corresponds to approximately 30% of the RE during the cultivation period. In a threecrop rotation, the carbon content of sugar beet residue contributed approximately 40% to the next crop ecosystem RE [10]. Therefore, it is necessary to consider the effect of this important field management technique when evaluating the annual carbon sequestration of a crop and its response to climate change in future research. By not calculating the carbon emissions from the crop residue in the present study, we may have introduced an important source of uncertainty in the assessment of carbon sequestration and its response to environmental or biotic factors.
Another source of uncertainty in the analysis of carbon budget in this study may come from the neglect of C ag . As aforementioned above, C ag usually includes the carbon emission from applying agrochemicals (fertilizer, herbicides, insecticides and fungicides), from the use of machinery and from the combustion of fuel. Many published researches have suggested the indirect carbon emission from C ag [60][61][62]. Maraseni et al. [60] have indicated that the amount of greenhouse gas emission because of the three main farm input was large in cropping land over 30 years, and the emission amount was 17094 kg carbon dioxide equivalents (CO 2 e)  for use of agrochemicals, 13272 kg CO 2 e for the use of fuel and 1910 kg CO 2 e for the farm machinery, respectively. When the site-scale research was conducted, this part of carbon emission to atmosphere would have some potential influence on assessment of annual carbon budget. Furthermore, the North China Plain is a so large area that is predominated by wheat-maize rotation system that this part of carbon emission would not be low.

The increase of carbon sequestration per extended CUP day and comparison with other ecosystems
The role of the CUP in affecting carbon sequestration has been well studied in forest and grassland systems, and these studies indicated that the annual NEP was strongly correlated with the critical phenological factor [3,41]. Our study also suggested that the interannual variations in both the NEP and the NBP were controlled mainly by the annual CUP ( Fig. 4 and Fig. 7a and d) when the two sites were combined. The annual NEP and NBP would increase by 14.865.2 and 14.766.6 g C m 22 yr 21 , respectively, if one CUP-day was extended. The increment amplitudes differed from those of other ecosystems, according to published studies. In a temperate deciduous forest, the NEP increased by 5.9 g C m 22 yr 21 per extended growing season day [37]. In a boreal coniferous forest, the increase was 6.9 g C m 22 yr 21 [31]. The NEP value may be lower for grassland systems. For example, Ma et al. (2007) [40] reported that the NEP would increase by 2.0 g C m 22 yr 21 and 4.0 g C m 22 yr 21 per extended growing season day for a savanna and an open grassland, respectively. The increase of carbon sequestration per extended growing season day in cropland in present study is much higher than forest and grassland ecosystem. The cropland is an artificial ecosystem and usually have simple ecosystem constituent to guarantee crop economic yield, and the cropland CUP is not as flexible as forest and grassland, which have multiple plant species

Interpreting the effect of MAT on the CUP
After path analysis, we have found a clearly relationship between MAT and CUP (Fig.7) and it can be interpreted by Fig. 5 and Fig. 6. Our result shows that the length of CUP was strongly controlled by the beginning date of CUP (CUP start-date) (Fig. 5), an early starting of CUP can lead to a longer CUP. This important date would be affected by MAT significantly (Fig. 6), a higher MAT can result an early beginning of CUP. As a result, the higher MAT, through its positive effect on the CUP start-date , prolonged the length of the CUP and thus increased the carbon sequestration each year in the croplands.
The conclusion that a higher mean annual air temperature would result in an earlier beginning of the CUP and prolong the duration of the CUP in this study contradicts other literature. For example, Black et al. (2000) [46] and Welp et al. (2007) [55] indicated that it was a higher spring air temperature rather than a higher mean annual air temperature that led to an earlier leaf-out and an earlier start of the growing season, which ultimately increased the annual NEP. To analyze the effect of a warm spring on the interannual variation in carbon sequestration, we did not define the 'spring period' based on the meteorological standard division (i.e., from March to May) but rather as February to March according to Welp et al. (2007) [55], who noted that 'spring' was the period during which plants would capture significant ecosystem carbon fluxes during the early growing season. In this study, the winter wheat was sown in October but began to grow rapidly during the next February to March after a period of dormancy to acclimate to the cold environment. Therefore, we defined the spring temperature (Tsp) as the mean temperature of February and March. Our results indicated that through its significant effect on the CUP start-date (R 2 = 0.45 for Tsp vs CUP start-date ; R 2 = 0.49 for MAT vs CUP start-date ), the Tsp and the MAT could explain 83% and 51%, respectively, of the interannual variation in the winter wheat NEP. Moreover, the Tsp and the MAT could explain 42% and 78%, respectively, of the interannual variation in the annual NEP. Therefore, the Tsp was the main factor controlling the annual NEP in a one-season crop (winter wheat), but it was the MAT in a two-season rotation crop that determined the annual NEP. Consequently, our results are in agreement with previous studies when we only consider winter wheat. However, when the two rotation crops are taken into account, our findings contradict the results of previous studies.
There were also studies focused on the effect of ending date of carbon uptake period on annual carbon sequestration. For example, Piao et al. (2008) [43] suggested that the ending date of CUP (CUP end-date ) would be advanced because of warming autumn and could therefore lead to loss carbon. Wu et al. (2012) stated that the interval between the ending date of the GEP and the ending date of the CUP (i.e., the autumn interval) in forest ecosystem would affect annual carbon sequestration. However, the current study has neither found the relationship between CUP enddate and annual NEP or NBP nor the relationship between the autumn interval and annual carbon sequestration, partly because the differences scales in those researches and other differences in types of ecosystems and thus in the photosynthetic pathways, in climate conditions, in ecosystem management and uncertainties in the flux data.

Outlook
Irrigation is an important management in cropland and can impact vegetation growth, eco-physiological characteristic of plant, ecosystem ecological processes and crop production obviously [63,64]. However, the effects of irrigation on dynamic of carbon sequestration and its relationship between environmental or biological factor can not be detected in this study because of the well and similar irrigation conditions in the two sites. There were many studies have reported the influence of water conditions on the carbon exchange in cropland. Suyker et al. (2012) found that the annual NEP was 552.0673.0 gCm 22 y 21 and 471.0652.0 gCm 22 y 21 in irrigate and rain-fed maize crop, respectively, who has also indicated that the peak value of weekly ensemble average of NEP were 2.5 mgCm 22 s 21 and 2.1 mgCm 22 s 21 in irrigate and rain-fed maize. Thus it can be seen that irrigation may be a source of variations of carbon flux, and, in future, how the temperature or carbon uptake period influence carbon sequestration in rain-fed field would become an interesting topic in North China Plain.
The current study further confirms the fact that the growing season length is a strong indicator of annual carbon sequestration not only in forests and grasslands but also in crop systems. However, it seems that different climate factors can affect the growing season length. Because of different physiological and ecological processes in plant leaves and canopies and their required climate conditions, it is not reasonable to state that one dominant factor that affects the annual carbon absorption in one ecosystem is the same factor in another ecosystem. Therefore, further research is needed on the interannual variation in ecosystem carbon flux in different sites and regions.

Conclusions
In this study, the interannual variations in the carbon fluxes were presented based on eddy covariance data collected in two winter wheat/summer maize double-cropping systems from 2003 to 2012 in Yucheng and from 2007 to 2012 in Luancheng. The annual NEP was 475.66159.4 g C m 22 yr 21 and 13.06272.5 g C m 22 yr 21 for Yucheng and Luancheng, respectively. After accounting for the grain harvest, the annual carbon sequestration (NBP) was 276.16174.5 g C m 22 yr 21 for Yucheng and 2 564.06272.2 g C m 22 yr 21 for Luancheng. The carbon uptake period (CUP) had a strong effect on the annual NEP and NBP. The crop could gain 14.865.2 g C m 22 per CUP day on an annual basis. Our study also showed that farmland will sequester more carbon as a result of a higher MAT. The path analysis indicated that there was a significant positive correlation between the MAT and the annual NEP or NBP, which was attributed to the strong effect of the MAT on the CUP. A higher MAT could result in an earlier start date of the carbon uptake season, which would significantly affect the length of the CUP and thereby increase carbon sequestration in croplands. Moreover, the spring temperature had a strong effect on annual carbon sequestration of winter wheat, whereas it was not the determined factor control on annual carbon balance of a double-cropping agroecosystem.

Author Contributions
Conceived and designed the experiments: XB XW XS. Performed the experiments: XB XW XS FZ YW. Analyzed the data: XB XW XS. Contributed reagents/materials/analysis tools: XB XW XS. Wrote the paper: XB XW XS.