Continued increase in atmospheric CO 2 seasonal amplitude in the 21 st century projected by the CMIP 5 Earth system models

In the Northern Hemisphere, atmospheric CO2 concentration declines in spring and summer, and rises in fall and winter. Ground-based and aircraft-based observation records indicate that the amplitude of this seasonal cycle has increased in the past. Will this trend continue in the future? In this paper, we analyzed simulations for historical (1850–2005) and future (RCP8.5, 2006–2100) periods produced by 10 Earth system models participating in the fifth phase of the Coupled Model Intercomparison Project (CMIP5). Our results present a model consensus that the increase of CO2 seasonal amplitude continues throughout the 21st century. Multi-model ensemble relative amplitude of detrended global mean CO2 seasonal cycle increases by 62±19 % in 2081–2090, compared to 1961–1970. This amplitude increase corresponds to a 68± 25 % increase in net biosphere production (NBP). The results show that the increase of NBP amplitude mainly comes from enhanced ecosystem uptake during Northern Hemisphere growing season under future CO2 and temperature conditions. Separate analyses on net primary production (NPP) and respiration reveal that enhanced ecosystem carbon uptake contributes about 75 % of the amplitude increase. Stimulated by higher CO2 concentration and high-latitude warming, enhanced NPP likely outcompetes increased respiration at higher temperature, resulting in a higher net uptake during the northern growing season. The zonal distribution and spatial pattern of NBP change suggest that regions north of 45 N dominate the amplitude increase. Models that simulate a stronger carbon uptake also tend to show a larger increase of NBP seasonal amplitude, and the cross-model correlation is significant (R = 0.73, p < 0.05).


Introduction
Modern measurements at Mauna Loa, Hawaii (19.5 • N, 155.6 • W, 3400 m altitude) have shown an increase in atmospheric CO 2 concentration from < 320 ppm in 1958 to 400 ppm in 2013.There is also a mean seasonal cycle that is characterized with a 5-month decrease (minimum in October) and a 7-month increase (maximum in May).The peakto-trough amplitude of this seasonal cycle is approximately 6.5 ppm, which represents a close average of a large portion of the Northern Hemisphere (NH) biosphere (Kaminski et al., 1996) where the amplitude ranges from about 3 ppm near the Equator to 17 ppm at Point Barrow, Alaska (71 • N).The seasonal variation of Mauna Loa (MLO) CO 2 reflects the imbalance of growth and decay of the NH biosphere.Early studies have speculated that global primary production would decrease because of global changes such as acid rain and deforestation (Reiners, 1973;Whittaker and Likens, 1973).If this is the case, assuming changes in respiration are similar at peak and trough of the CO 2 seasonal cycle, we might observe a reduction of CO 2 seasonal amplitude.However, Hall et al. (1975) found no evidence of long-term amplitude change from 15 years of MLO CO 2 record (1958)(1959)(1960)(1961)(1962)(1963)(1964)(1965)(1966)(1967)(1968)(1969)(1970)(1971)(1972).They concluded that either the biosphere is too big to be affected yet or the degradation of biosphere is balanced by enhanced CO 2 fertilization and increased use of fertilizers in agriculture.
In the 1970s through the 1980s, the metabolic activity of the biosphere seems to be getting stronger, as indicated by Published by Copernicus Publications on behalf of the European Geosciences Union.
F. Zhao and N. Zeng: Continued increase in atmospheric CO 2 seasonal amplitude in the 21st century rapid increase in MLO CO 2 amplitude (Pearman and Hyson, 1981;Cleveland et al., 1983;Bacastow et al., 1985).Enhanced CO 2 fertilization was considered as a major factor, and climate change a possible cause (Bacastow et al., 1985).Keeling et al. (1996) linked the amplitude increase with climate change by showing the 2-year phase lag relationship between trends of CO 2 amplitude and 30-80 • N mean land temperature.Unlike CO 2 fertilization, the combined effect of climate (temperature, precipitation, etc.) introduces strong interannual variability to the CO 2 amplitude change.In the early 1990s, despite of the continuing rise of 30-80 • N mean land temperature, CO 2 seasonal amplitude at MLO has declined.Buermann et al. (2007) attributed this decline to the severe drought in North America during 1998-2003.In late 1990s, the increasing trend resumed at MLO.The latest analysis shows a 0.32 % yr −1 increase in MLO amplitude and a 0.60 % yr −1 increase in Point Barrow (BRW) amplitude (Fig. 1a, Graven et al., 2013).This trend (over 50 years) corresponds to an increase of 16 % in MLO, and 30 % in BRW CO 2 seasonal amplitude, respectively.Graven et al. (2013) also compared aircraft measurements taken at 500 and 700 hPa heights in 1958-1961 and 2009-2011, suggesting an even larger (∼ 50 %) increase of CO 2 seasonal amplitude north of 45 • N. Furthermore, to infer the modelsimulated CO 2 amplitude increase at 500 hPa, they applied two atmospheric transport models to monthly net ecosystem production (NEP) from the historical simulation (Exp3.2) results of eight CMIP5 models.Compared with aircraft data, they found the CMIP5 models simulated a much lower amplitude increase.
Surface CO 2 monitoring stations have two major limitations.First, they are sparse.For several decades, the Global Monitoring Division of NOAA/Earth System Research Laboratory (ESRL) has measured CO 2 from more than 100 surface monitoring sites (Conway et al., 1994).Only some have over 30 years of record.Similarly, Randerson et al. (1997) determined the CO 2 amplitude trend north of 55 • N by averaging flask data from five stations.Second, the surface CO 2 stations do not measure carbon exchange between the land/ocean and atmosphere directly.Atmospheric inversion models are capable of providing surface carbon fluxes with global coverage.However, the resolution and accuracy of these models are inherently limited due to a small number of stations used, and errors in atmospheric transport (Peylin et al., 2013).
Process-based terrestrial biosphere models (TBMs) can generate surface fluxes over the past for longer period, usually with a spatial resolution of half to three degrees.Thus, they offer opportunities to understand the mechanisms of CO 2 amplitude increase better.McGuire et al. (2001) calculated amplitude trends of total land-atmosphere carbon flux (north of 30 • N) from four TBMs.Compared to Mauna Loa CO 2 , they found the trend was overestimated by one of the four models and underestimated by the other three.They suggest the observed trend may be a consequence of the com-bined effects of rising CO 2 , climate variability and land use changes, which have also been recognized in previous studies (Kohlmaier et al., 1989;Keeling et al., 1995Keeling et al., , 1996;;Randerson et al., 1997Randerson et al., , 1999;;Zimov et al., 1999).Models show varied extent of amplitude increase, possibly due to their different sensitivities to CO 2 concentration and climate.Interestingly, Graven et al. (2013) found that CMIP5 models underestimate the CO 2 amplitude increase in the mid-troposphere at latitudes north of 45 • N.However, previous observations indicated that the models might overestimate CO 2 fertilization effect (Piao et al., 2013), suggesting that our understanding of the amplitude trend is still limited.
In the future, we do not know if the CO 2 amplitude will increase or decrease.With temperature rise and CO 2 increase, we may see a further increase of CO 2 amplitude.On the other hand, the frequency and/or duration of heat waves are very likely to increase over most land areas, and the increases in intensity and/or duration of drought and flood are likely (IPCC, 2013).As a result, the ecosystem productivity may decrease, which may reduce the CO 2 amplitude.In this study, we analyzed the fully coupled CMIP5 Earth system model runs as part of the Fifth Assessment Report (AR5) of the United Nations' Intergovernmental Panel on Climate Change (IPCC).Specifically, we looked into the emissiondriven simulations, which include many of the aforementioned processes and feedbacks.Our specific questions are the following.(1) How do the CMIP5 models predict the amplitude and phase changes of CO 2 seasonal cycle in the future?(2) Are the changes mostly driven by changes in ecosystem production or respiration?(3) Where do the models predict the largest CO 2 amplitude changes will occur?
Section 2 describes the CMIP5 experiments, models used and our analyzing method; Sect. 3 presents the major results of our multi-model analyses; Sect. 4 discusses amplitude increases at individual stations, physical mechanisms and uncertainties; and Sect. 5 concludes our main findings.

Model descriptions
We analyzed historical and future emission-driven simulation results from 10 CMIP5 Earth System Models (ESMs).The historical simulations, referred to as experiment 5.2 or ESM historical 1850-2005 run (Taylor et al., 2012), were forced with gridded CO 2 emissions reconstructed from fossil fuel consumption estimates (Andres et al., 2011).The future simulations, referred to as experiment 5.3 or ESM RCP8.5 2006-2100 run, were forced with projected CO 2 emissions, following only one scenario-RCP8.5 (Moss et al., 2010).We chose the emission-driven runs because the fully coupled ESMs in these runs have interactive carbon cycle component.Global atmospheric CO 2 concentrations are simulated prognostically, therefore they reflect the total effect of all the physical, chemical, and biological processes on Earth, and their interactions and feedbacks with the climate system.We obtained model output primarily from the Earth System Grid Federation (ESGF), an international network of distributed climate data servers (Williams et al., 2011).
For the GFDL model, we retrieved results from its data portal (http://nomads.gfdl.noaa.gov:8080/DataPortal/cmip5.jsp).The main characteristics of the 10 models are listed in Table 1.

Analysis procedure
We first analyzed the monthly output of prognostic atmospheric CO 2 concentrations to evaluate the change of CO 2 seasonal amplitude (defined as maximum minus minimum of detrended seasonal cycle) from 1961 to 2099.Atmospheric CO 2 was obtained primarily as the area-and pressureweighted mean of CO 2 across all vertical levels -a better representation of atmospheric carbon content than surface CO 2 .The INM-CM4 model does not provide CO 2 concentration, so we converted its total atmospheric mass of CO 2 to mole fraction.We excluded the IPSL model from analyses in Sects.3.1 and 3.2 because its CO 2 output is not available.Only CanESM2 provides three different realizations for both historical and future runs, and we simply use its first realization in our comparison.We believe this choice would lead to a more representative result than including all realizations of CanESM2 in multi-model averaging.
To extract the CO 2 seasonal cycle from the monthly records, we applied the curve-fitting procedures using the CCGCRV software developed at the National Oceanic and Atmospheric Administration Climate Monitoring and Diagnostics Laboratory (Thoning et al., 1989; http://www.esrl.noaa.gov/gmd/ccgg/mbl/crvfit/crvfit.html).This algorithm first fits the long-term variations and the seasonal component in the monthly CO 2 record with a combination of a trend function and a series of annual harmonics.Then the residuals are filtered with fast Fourier transform and transformed back to the real domain.Specifically, we followed the default setup of a quadratic polynomial for the trend function, 4-yearly harmonics for the seasonal component, and long/short-term cutoff values of 667 days/80 days for the filtering in our analyses.We examined the phase change of CO 2 detrended seasonal cycle by counting how frequent the maxima and minima occur in different months.We used two definitions of seasonal amplitude in our analyses that yield similar results: one directly comes from the CCGCRV package, and another definition is simply maximum minus minimum of detrended seasonal cycle in each year.For each model's monthly global mean CO 2 , we first computed the detrended CO 2 seasonal cycle as the annual harmonic part plus the filtered residue using the short-term cutoff value.Then we started to investigate the global carbon budget in each model: The left term is the change of CO 2 concentration (or CO 2 growth rate), which we simply computed as the difference between the current month and previous month's concentration -this  emission (FFE), net biosphere production (NBP, or net terrestrial-atmosphere carbon exchange, positive if land is a carbon sink), and net ocean-atmosphere flux (FOA, positive if ocean releases carbon).For each model, we checked and ensured that the sum of individual flux terms on the RHS of Eq. ( 1) equals to the CO 2 growth rate.Previous studies have limited the impact of FFE and FOA on trends in CO 2 amplitude to less than a few percent change (Graven et al., 2013).Therefore we focused on examining the seasonal cycle of NBP in this study.To investigate whether the NBP amplitude change is mostly due to enhanced production or respiration, we inspected the seasonal cycle of NPP and respiration separately.The INM model does not provide NPP output, so it is excluded in this part of analyses.For respiration, one complication is that even though NBP represents the net terrestrial-atmosphere carbon exchange in all models (thus allowing model comparison), its further breakdown varies.For example, the GFDL-ESM2M model's NBP has component fluxes including NPP, heterotrophic respiration (R h ), land use change (fLuc), fire (fFire), harvest (fHarvest) and grazing (fGrazing).In contrast, NBP approximately equals to NPP minus R h in CanESM2.Instead of directly adding all flux components such as R h and fLuc for each model (which would be unnecessary and difficult since not all fluxes are provided), we defined Additionally, we analyzed the spatial patterns of NBP change between future (2081-2090) and historical (1961)(1962)(1963)(1964)(1965)(1966)(1967)(1968)(1969)(1970) period.We approximated NBP amplitude change as the difference between the peak seasons of carbon uptake and release by the biosphere, namely May-July and October-December averages, respectively.We chose 3-month averages for multimodel ensemble, because not all models simulate peak uptake in June and peak release in October.Monthly output of NBP, NPP and R * h (derived from NBP and NPP) from all models were first resampled to 2 × 2 • grids.Then the spatial and zonal means for both May-July and October-December were computed.

Changes of CO 2 and NBP seasonal amplitude
The CMIP5 models project that the increase of CO 2 seasonal amplitude continues in the future.Figure 1a shows detrended and globally averaged monthly column atmospheric CO 2 from 1961 to 2099, averaged over nine models (no IPSL).The models project an increase of CO 2 seasonal amplitude (defined as maximum minus minimum in each year) by about 70 % over 120 years, from 1.6 ppm during 1961-1970 to 2.7 ppm in 2081-2090.The increase is faster in the future than in the historical period.Another feature is that the trend of minima (−0.63 ppm century −1 ) has a larger magnitude than the trend of maxima (0.41 ppm century −1 ), suggesting that enhanced vegetation growth contributes more to the amplitude increase than respiration increase.Gurney and Eckels (2011) found the trend of net flux in dormant season is larger than that of growing season.However, they applied a very different definition for amplitude, considering all months instead of maxima and minima, to analyze the atmospheric CO 2 inversion results from 1980-2008.Specifically, they defined growing season net flux (dormant season net flux) as the total of any month for which the net carbon flux is negative (positive), and amplitude as the difference of the two net fluxes.It is no surprise they reached a conclusion that seems to contradict ours, since growing season is much shorter than dormant season at global scale.averaged over 1961-1970 and 2081-2090 for the nine models, and their multi-model ensemble (MME) and standard deviation (SD).
In addition, we notice a phase advance of maxima and minima by counting their time of occurrence (data not shown).
Excluding models above one standard deviation from the ensemble mean yields similar results.
To illustrate how well the models reproduce the seasonal variations of CO 2 , we compared the multi-model ensemble global CO 2 at the lowest model level -not equivalent to the height of surface CO 2 measurement, but relatively close -with ESRL's global mean CO 2 over 1981-2005 (Fig. 1d).The surface CO 2 seasonal amplitude estimated by the model ensemble is lower than that of ESRL's global CO 2 estimate (Ed Dlugokencky and Pieter Tans, NOAA/ESRL, www.esrl.noaa.gov/gmd/ccgg/trends/),however the amplitude increases are similar (Table 3).This surface station-based global CO 2 estimate also indicates that the amplitude increase is dominated by the trend of minima.
We further calculated the change of relative amplitude (relative to [1961][1962][1963][1964][1965][1966][1967][1968][1969][1970] for each model.The amplitude here is computed by the CCGCRV package.As illustrated in Fig. 2, all nine models show an increase in both global mean CO 2 and total NBP seasonal amplitude.CO 2 seasonal amplitude has increased by 62 ± 19 % in 2081-2090, compared to 1961-1970; whereas NBP seasonal amplitude has increased by 68 ± 25 % over the same period (see Table 4 for details of individual models).The trend of increase is much higher in the future (CO 2 /NBP: 0.70 %/0.73 % yr −1 during 2006-2099) than in the historical period (0.25 % and 0.28 % yr −1 during 1961-2005 for CO 2 and NBP), albeit the model spread also becomes larger in the future.When we applied the same procedure to the Northern Hemisphere (25-90  (excluding INM-CM4 which only has global CO 2 mass), we saw a higher amplitude increase and larger model spread: 81 ± 46 and 77 ± 43 % for CO 2 and NBP, respectively.

Production vs. respiration
Our next major question is whether the amplitude increase of NBP is largely driven by NPP or respiration.We computed the mean seasonal cycle of detrended CO 2 growth rate, −NBP, −NPP (reverse signs so that negative values always indicate carbon uptake) and R * h in two periods: 1961-1970 (black) and 2081-2090 (red), for the nine models (for this and following analyses, we excluded INM which does not provide NPP, and included the IPSL model except for CO 2 growth rate).The seasonal cycle of −NBP resembles that of detrended CO 2 growth rate (Fig. 3a-d), confirming that the activities of land ecosystem dominate the CO 2 seasonal cycle and its amplitude increase in the model simulations.Except for CanESM2 (also noted in Anav et al., 2013), and BNU-ESM (which simulates a second peak carbon uptake around November) to some extent, most models can reproduce the net uptake of carbon during spring and summer (when increasing NPP overcomes respiration) and the net carbon release during fall and winter at global scale: net carbon uptake peaks in June (five models) or July (three models) for the historical period, and exclusively in June for the future period.However, the model spread on amplitude is large: CESM1-BGC and NorESM1-ME, which has the same land model (CLM4) that features an interactive nitrogen cycle, are characterized by a small seasonal amplitude of −NBP -merely 30 % of those on the high end of the models (IPSL-CM5A-LR and MPI-ESM-LR).The seasonal amplitude of multi-model ensemble NBP, computed as maximum minus minimum (June-October), has increased from 2.7 to 4.7 PgC month −1 (Fig. 3d).
Earth Syst.Dynam., 5, 423-439, 2014 This 2 PgC month −1 amplitude increase is the sum of enhanced net carbon uptake in June and higher net release in October, and the enhancement in uptake (1.4 PgC month −1 ) is nearly 3 times as large as the release increase (0.5 PgC month −1 ).
We then investigate the June and October changes of −NPP and R * h , respectively.By definition, their sum should equal to the amplitude change of −NBP.NPP has increased in all months (Fig. 3e and f), with much larger changes during the NH growing season.The amplitude of multi-model ensemble NPP has increased from 4.8 to 7.1 PgC month −1 , and an increase from 2.7 to 4.3 PgC month −1 is found for R * h .In June, NPP increase (4.5 PgC month −1 ) is larger than that of R * h (3.1 PgC month −1 ), resulting in enhanced net uptake.In October, NPP increase (1.9 PgC month −1 ) is smaller www.earth-syst-dynam.net/5/423/2014/Earth Syst.Dynam., 5, 423-439, 2014 than that of R * h (2.4 PgC month −1 ), leading to enhanced net release.These results are consistent with trends of maxima and minima in Fig. 1.The models also indicate a shift in peak NPP from July to June, consistent with the shift of CO 2 minima.

Spatial and latitudinal contributions
To further investigate the regional contribution to NBP amplitude increase, we plotted the 10-model mean −NBP changes (Fig. 4) over peak NH growing season (May-July, Fig. 4a) and dormant season (October-December, Fig. 4b).Because the models disagree on the time of maximum and minimum NBP (Fig. 3), our choice of doing seasonal averages would be more representative of the models than averaging over 1 month.Note that the difference between the two seasonal averages is smaller than the peak-to-trough amplitude, but here we are only concerned with the spatial pattern.We saw a stronger net carbon uptake in May-July almost everywhere north of 45 • N, and also over the Tibetan Plateau and some places near the equator.Net carbon uptake weakens over western United States and Central America, South and Southeast Asia and central South America.The change of net carbon release in October-December generally shows an opposite spatial pattern, with a noticeably smaller magnitude north of 45 • N.
In addition, we calculated the corresponding zonal averages (Fig. 3c).The area-weighted totals of the zonal mean curves correspond to the future minus historical averages of global total −NBP (Fig. 3d), averaged over May-July and October-December, respectively.These two curves do not account for phase difference; instead, they approximate latitudinal contribution to the amplitude increase of global total −NBP.It is apparent that this increase is dominated by regions north of 45 • N with a weak contribution from the Southern Hemisphere tropics (0-25 • S).The northern subtropical region and Southern Hemisphere (10-30 • N, 35-55 • S) partly offset the amplitude increase.It is also clear that the amplitude increase is dominated by changes in peak growing season (the green shade is larger on the left of the zero line than on its right), consistent with our findings in the previous section.
Analogous to the cold-warm seasonality in the temperate/boreal region, the tropics has distinctive dry and wet seasons, and recently Wang et al. (2014) suggested the tropical ecosystem is becoming more sensitive to climate change.
In our analyses on the multi-model ensemble patterns, the tropical region exhibits a small negative contribution to the seasonal amplitude increase of global total −NBP.This does not mean the net carbon flux in the tropics, which has a different seasonal cycle phase, would experience an amplitude decrease in the future.To illustrate the seasonal amplitude change at different latitudes, we show the zonal amplitude of NBP in the historical (black) and future (red) periods for all models (Fig. 5).At every 2 • band, we first calculated a 10-year mean seasonal cycle, then compute its amplitude (maximum minus minimum).Most models predict an increase in NBP seasonal amplitude at almost every latitude under the RCP85 emission scenario.Only two of the models, CanESM2 and MIROC-ESM, predict decreased seasonality for parts of the tropics and subtropics.Unlike in Fig. 4c, an area-weighted integral cannot be performed due to different phases zonally.The Southern Hemisphere has an opposite phase from its northern counterpart, but its magnitude is small due to its small land area.The two subtropical maxima around 10 • N and 10-15 • S reflect the wet-dry seasonal shift in the Intertropical Convergence Zone (ITCZ) and monsoon movement.They are comparable to the NH maxima in terms of both amplitude and amplitude increase for about a third of the models, however they are out of phase and largely cancel each other out.
To further illustrate this cancelation effect, we aggregated monthly −NBP over six large regions: the globe (90 • S-90 • N), northern boreal (50-90 • N), northern temperate (25-50 • N), northern tropics (0-25 • N), southern tropics (25 • S-0 • ) and Southern Hemisphere (90-25 • S; Fig. S1).It is clear that the changes of global −NBP seasonal cycle mostly come from the northern boreal region; it partly comes from the northern temperature region in a few models.The seasonal cycle of the northern tropics is characterized by spring maxima and fall minima, and prominent increases of its seasonal amplitude are found for BNU-ESM, GFDL-ESM2M and IPSL-CM5A-LR.However, they are largely counterbalanced by the southern tropics.For GFDL-ESM2M, changes in the southern tropics are larger than its northern counterpart, but even so, the net contribution of tropical regions to its global −NBP seasonal amplitude (September maxima minus June minima) increase is limited to about 25 %, the largest of all models.

Mechanisms for amplitude increase
As discussed in Sect. 1, two major mechanisms for amplitude increase identified in previous literature are CO 2 fertilization effect and high latitudes "greening" in a warmer climate.Both mechanisms lead to enhanced ecosystem productivity during peak growing season, and consequently more biomass to decompose in dormant season, therefore increasing the amplitude of NBP seasonal cycle.Because models have different climate and CO 2 sensitivity (Arora et al., 2013), their relative importance may vary.In the case of  1961-1970 (black) and 2081-2090 (red).For each model, NBP is first regridded to a 2 × 2 • common grid.Monthly zonal totals are then computed for every 2-degree band, which determine the amplitude (maximum minus minimum) at every band.The Southern Hemisphere has an opposite phase from its northern counterpart, but its magnitude is small due to its small land area.The two subtropical maxima around 10 • N and 10-15 • S reflect the wet-dry seasonal shift in the Intertropical Convergence Zone (ITCZ) and monsoon movement.They have similar magnitude as the Northern Hemisphere maxima in about a third of the models, however their net contribution to global total NBP seasonal amplitude is small, because they are out of phase and largely cancel each other out.CMIP5 ESMs, two additional sensitivity experiments are recommended: fixed feedback 2 (esmFdbk2) and fixed climate 2 (esmFixClim2).The former keeps CO 2 concentration fixed but allows physical climate change responding to increasing historical and future (RCP4.5)concentrations; the latter keeps climate fixed under preindustrial CO 2 condition but allows the carbon cycle to respond to historical and future (RCP4.5)CO 2 increase.This setup does not permit quantifying the contribution of CO 2 increase and climate change to NBP amplitude increase: one major difference is the use of RCP4.5 concentrations instead of RCP8.5 emissions.How-ever, we can still make qualitative assessments by examining the spatial patterns.We will focus on the high latitude regions, which contribute most to amplitude increase of global total NBP.
Of the 10 models we studied, only CanESM2, GFDL-ESM2M and IPSL-CM5A-LR have submitted NBP output for these two experiments (MIROC submitted output for esmFixClim2 only).Here we display the spatial patterns of −NBP changes for GFDL-ESM2M (Fig. 6) and IPSL-CM5A-LR (Fig. 7).CanESM2 results are not shown because it does not correctly reproduce the phase of global total NBP Earth Syst.Dynam., 5, 423-439, 2014 www.earth-syst-dynam.net/5/423/2014/seasonal cycle.The changes of −NBP for both models during peak growing season are clearly dominated by CO 2 fertilization effect (right panels).In contrast, climate change under fixed CO 2 fertilization conditions has mixed effects on high latitude regions.Northern high latitude net carbon release in October-December is increased both under climate change (Fig. 6c) and elevated CO 2 conditions (Fig. 6d) for GFDL-ESM2M, but over different regions.For IPSL-CM5A-LR however, net carbon release increase in regions north of 45 • N is only obvious under elevated CO 2 condition.
Our results only indicate CO 2 fertilization effect is the dominant factor for NBP seasonal amplitude increase in some models.For models with strong carbon-climate feedbacks and weak/moderate water constraints in northern high latitude regions, climate change may be more important.However, we cannot find a clear example due to data availability.MIROC-ESM is known to have strong carbonclimate feedback (Arora et al., 2013).From its simulation under fixed climate (figure not shown), we found no obvious patterns of widespread net carbon release increase in dormant season, suggesting climate change may play a bigger role for this model.The HadGEM model is another possible candidate; it is also a particularly interesting model to analyze since one of its historical simulations represented the largest increase in CO 2 amplitude in Graven et al. (2013).Unfortunately, for the ESM simulations, both CO 2 and NBP from HadGEM are not available on the ESGF servers.

Relationship with mean carbon sink
Our analyses above suggest CO 2 fertilization effect is a major mechanism causing the amplitude increase in some models.If it is important in most models, we expect to see models with a larger change in mean carbon sink simulate a higher increase in seasonal amplitude.By plotting the −NBP change against NBP seasonal amplitude increase for all 10 models (Fig. 8), we found there is indeed a negative crossmodel correlation (R = −0.73,p < 0.05), indicating models with a stronger net carbon uptake are likely to simulate a larger increase in NBP seasonal amplitude.Note that this result is based on the 10 models we analyzed; it is subject to large uncertainty and may change substantially with inclusion or exclusion of certain model(s).Again all models show an increase in NBP seasonal amplitude, even though they disagree on the direction of future NBP change.While our study hint at a possible relationship between mean carbon sink and NBP seasonal amplitude, it is beyond our scope to discuss further, or comment on why models show such different mean sink estimate.Interested readers may refer to the insightful discussion on this issue in Friedlingstein et al. (2013).-2090 and 1961-1970 for 10 CMIP5 ESMs.The negative cross-model correlation (R = −0.73,p < 0.05) suggests that a model with a larger net carbon sink increase is likely to simulate a higher increase in NBP seasonal amplitude.

Discussions
We have primarily focused on model ensembles of aggregated quantities.Ensemble patterns are sometimes dominated by only a few models due to large seasonality variations among the models.However, the close examination of each individual model show that the spatial patterns of −NBP change during peak growing season (May-July) are all dominated by high latitude regions (approximately north of 45 • N).In CESM1-BGC and NorESM1-ME models, enhanced net carbon uptake are confined to some of the high latitude regions (Fig. S2 in the Supplement).Models differ on finer details.For example, about half of the models predict an obvious increase of net carbon uptake for the Tibetan Plateau.It is worth mentioning that the esmFixClim2 experiment of MIROC-ESM shows little change in NBP for this region under elevated CO 2 alone.High latitude regions also dominate the increase of net carbon release in October-December for most models (Fig. S3).One exception is INM-CM4, which displays very small change in the dormant season, and most of its NBP amplitude increase comes from enhanced carbon uptake during peak growing season.BNU-ESM and CanESM2 have some limitations in reproducing the correct phase of global −NBP seasonal cycle.Exclusion of these two models from ensemble mean calculation exhibits very similar spatial and zonal patterns as shown in Fig. 4. Another caveat is the assumption of 1961-1970 as the historical condition and 2081-2090 as future condition.This choice is valid if the selected variables have roughly monotonic trends, and 10 years is long enough to smooth out most of the interannual variability.Figure 2  sumption is quite reasonable for model ensembles, and acceptable for individual models.
We presented aggregated quantities due to large model uncertainty in space.We have largely omitted model evaluation against observations (due to limited observation during [1961][1962][1963][1964][1965][1966][1967][1968][1969][1970].However, this step can be helpful in model evaluation studies (Anav et al., 2013;Peng et al., 2014).One concern is to examine whether the models can reproduce observed CO 2 seasonal amplitude increase at the two stations with longest observation records -Mauna Loa, Hawaii and Point Barrow, Alaska.To address this issue, we extracted simulated CO 2 concentration from eight models at their model grid that is closest to Mauna Loa in the threedimensional space (similar procedure for Point Barrow).The results of this comparison at one model grid can reflect multiple sources of model uncertainties (such as uncertainties in the atmospheric tracer transport and mixing simulations).For example, GFDL-ESM2M is known to simulate a damped CO 2 gradient (Dunne et al., 2013) which has long been identified as a deficit in models of the atmospheric CO 2 cycle (Fung et al., 1987).
Figure 9 (and Fig. S4 for more details) presents the changes of CO 2 seasonal amplitude at Mauna Loa for the models and observation.CO 2 seasonal amplitude is underestimated by a factor of 2 in three-quarters of the models.However, the amplitude increase from ensemble model estimate (0.36 ± 0.24 % yr −1 , error range represents one standard deviation model spread) is much closer to observation (0.34±0.07 % yr −1 , error range represents one standard error of the least-squared trend calculation).MPI-ESM-LR reproduces both the magnitude and trend of Mauna Loa CO 2 seasonal amplitude reasonably well.For Point Barrow (Figs. 10  and S5), MPI-ESM-LR also simulates a similar amplitude increase to observation, but the magnitude of amplitude is much larger (almost twice).All other models underestimate the amplitude, but for the amplitude increase, the model ensemble (0.46 ± 0.21 % yr −1 ) again is similar to observation (0.43 ± 0.10 % yr −1 ).MRI-ESM1 is found to reproduce both the magnitude and increase of Point Barrow CO 2 amplitude quite well.Graven et al. (2013) found the CMIP5 models substantially underestimate the amplitude increase of CO 2 north of 45 • N at altitude of 3 to 6 km.However, we did not find that the models underestimate Point Barrow CO 2 amplitude increase at surface level.One big difference is the observational data used for comparison.During the 1974-2005 period, CO 2 seasonal amplitude increases by 0.43 % yr −1 , or 21.5 % over 50 years at the Point Barrow station.This is much lower than the ∼ 50 % amplitude increase found between the two aircraft campaigns during 1958-1961and 2009-2011(Graven et al., 2013)).This difference might be attributed to mechanisms controlling the vertical profile of CO 2 concentration.It is also not clear to what extent the large interannual variability of CO 2 seasonal amplitude affects the trend estimation of observed CO 2 amplitude increase.
Under the RCP8.5 emission scenario, CMIP5 showed a 62 ± 19 % increase of CO 2 seasonal cycle globally from 1961-1970 to 2081-2090.The increase is 85 ± 48 % at Mauna Loa (range indicates one standard deviation model spread), and 110 ± 42 % at Point Barrow.Even though the CMIP5 models are able to reproduce the increase of CO 2 seasonal amplitude at the two locations, some of the models rely heavily on the CO 2 fertilization mechanism, which may be too strong compared to observational evidence.Previous research suggest it should explain no more than 25 % of the observation at a high fertilization effect permitted by lab experiments (Kohlmaier et al., 1989).Similarly, Randerson et al. (1997) found the linear factor of CO 2 fertilization has to be 4 to 6 times greater than the mean of the experimental values, in order to explain the 0.66 % yr −1 amplitude increase (north of 55 • N) during 1981-1995.Recent studies have indicated that some important mechanisms, such as changes in ecosystem structure and distribution (Graven et al., 2013) and land use intensification (Zeng et al., 2014), are missing in the current CMIP5 models.Yet another main source of uncertainty is future CO 2 emissions.The RCP8.5 scenario used to drive the ESMs is on the high side of future scenarios.Also, the emission-driven runs simulate higher CO 2 than observed over the historical period, and such biases are likely to accumulate over time as the increase of atmospheric CO 2 growth rate accelerates (Hoffman et al., 2014).
The models do not have the same strength of carbonclimate feedback, but even if they do, their response to climate change may vary significantly simply because they simulate very different climate change.To briefly address this issue, we present soil moisture (Figs.S6 and S7) and nearsurface temperature (Figs.S8 and S9) changes for all models.All the models show temperature increase, but in different ranges.The more prominent difference was observed in the spatial pattern of soil moisture changes predicted by models.The combined effect of soil moisture regimes, temperature change and plant functional type (PFT) specifications could cause diverse behaviors of models over same regions.Such www.earth-syst-dynam.net/5/423/2014/Earth Syst.Dynam., 5, 423-439, 2014 are important caveats that highlight the importance of sensitivity experiments and warrant more in-depth future studies.The combined effect of climate and CO 2 changes not only alters the balance between production and respiration for existing ecosystems, but also leads to changes of ecosystem types.For example, Fig. 11 shows that the tree fraction has increased over wide areas of the northern high latitude regions for MPI-ESM-LR and INM-CM4.Figure 12 reveals notable natural grass increase over the northern high latitude regions for BNU-ESM.Such widespread vegetation change has not been observed during the satellite era, and it is possibly yet another highly uncertain mechanism contributing to amplitude increase in some CMIP5 models.
The major crops are characterized by high productivity in a short growing season, and they tend to have larger NBP seasonal amplitude compared to the natural vegetation they replace (usually natural grass).An increase in cropland fraction over high latitude regions could contribute to the seasonal ampltiude increase of NBP.As far as we know, no CMIP5 model has accounted for agricultural intensification, and only some models have implemented a conversion matrix (Brovkin et al., 2013).Therefore, the most important change implemented in the CMIP5 models is fractional land cover change based on Hurtt et al. (2011).In Fig. S10 we present the change of crop fraction, available from five mod-els.It is apparent that crop area has increased mostly in the tropics, while regions north of 30 • N have actually seen a decrease (due to a variety of factors: cropland abandonment, reforestation, urbanization, etc.).Therefore, crop fractional cover change alone may decrease the NBP seasonal amplitude in CMIP5 simulations.A better representation of land use change, especially the agricultural intensification, is needed in CMIP5 models to represent the CO 2 and NBP seasonal cycle better.On a side note, the other major part of land cover change -pasture (often treated as natural grass in ESMs, Brovkin et al., 2013) fraction change is unlikely to have a significant effect on NBP seasonal amplitude in the CMIP5 simulations.

Conclusions
Under the RCP8.5 emission scenario, all models examined in this study project an increase in seasonal amplitude of both CO 2 and NBP.The models' results indicate an earlier onset and peak of Northern Hemisphere biosphere growth and decay under future climate and CO 2 conditions.The amplitude increase is dominated by changes in net primary productivity, and changes in regions north of 45 • N. Our results suggest the models simulating a larger mean carbon sink increase are likely to project a larger increase in NBP seasonal ampli-tude.Considerable model spread is found, likely due to different model setup and complexity, different climate conditions simulated by the models, sensitivity to CO 2 and climate and their combined effects, and strength of feedbacks.Our findings indicate factors including enhanced CO 2 fertilization and lengthening of growing season in high-latitude regions outcompetes possible severe drought and forest degradation (leading to loss of biosphere productivity) in the future.
Despite of the model consensus in global CO 2 and NBP seasonal amplitude increase, and a reasonable representation of CO 2 seasonal amplitude increase at Mauna Loa and Point Barrow compared to surface in situ observations, the mechanisms contributing to these changes are debatable.CO 2 fertilization may be too strong, and factors like ecosystem change and agricultural intensification are underrepresented or missing in the CMIP5 ESMs.Future modelintercomparison projects should encourage models to participate in consistent and comprehensive sensitivity experiments.
The Supplement related to this article is available online at doi:10.5194/esd-5-423-2014-supplement.
Figure  1b and cpresent detrended global mean CO 2 growth rate (1 ppm = 2.12 PgC month −1 for unit conversion) and global total −NBP, two quantities showing very similar character-

Figure 2 .
Figure 2. Time series of the relative seasonal amplitude (relative to 1961-1970 mean) of (a) global mean atmospheric CO 2 ; and (b) global total NBP from 1961 to 2099.Thick black line represents multi-model ensemble, and one standard deviation model spread is indicated by light grey shade.

Figure 3 .
Figure 3. Seasonal cycle of detrended global mean CO 2 growth rate (a and b), global total −NBP (c and d), global total −NPP (e and f), and global total R * h (g and h, computed as NPP minus NBP), averaged over1961-1970 and 2081-2090  for the CMIP5 models (excluding INM, also excluding IPSL for CO 2 growth rate).Seasonal cycles of individual models are presented in the left panel(dashed for 1961-1970,  and solid for 2081-2090).Ensemble mean and one standard deviation model spread (black/grey for 1961-1970, red/pink for 2081-2090) are displayed in the right panels.Blue arrows mark the changes in June and October (NBP maxima and minima), except for CO 2 growth rate and −NPP, where arrows also indicate phase shifts of minima between the two periods.We show −NBP and −NPP so that the negative values represent carbon uptake by the biosphere, and positive values indicate carbon release from the biosphere.Note that −NBP and its two components −NPP and R * h are not detrended, so that the sum of (f) and (h) equals to (d).Detrended −NBP seasonal cycle (not shown) looks very similar to (d), as its trend is small compared to the magnitude of seasonal cycle.

Figure 4 .
Figure 4. Spatial patterns and latitudinal distributions of 10-model mean −NBP (gC m −2 day −1 ) changes between 2081-2090 and 1961-1970, during mean (a) peak growing season (May-July) and (b) dormant season (October-December).(c) Aggregates the spatial patterns in (a) and (b) zonally, where the black curve corresponds to the −NBP changes in May-July (a), and the red curve corresponds to the −NBP changes in October-December (b).Further reduction of −NBP in peak growing season -where the black curve falls on the left of the zero line, and increase of −NBP in dormant season -where the red curve is on the right of the zero line, both contribute to amplitude increase.We shade those instances in green, and shade the reversed case (contribute negatively to global total −NBP amplitude increase) in yellow.It is clear that the amplitude increase is dominated by the boreal regions, and by changes in peak growing season.

Figure 5 .
Figure5.Zonal amplitude of NBP from the 10 CMIP5 models (PgC month −1 per 2-degree band), averaged over1961-1970  (black) and 2081-2090 (red).For each model, NBP is first regridded to a 2 × 2 • common grid.Monthly zonal totals are then computed for every 2-degree band, which determine the amplitude (maximum minus minimum) at every band.The Southern Hemisphere has an opposite phase from its northern counterpart, but its magnitude is small due to its small land area.The two subtropical maxima around 10 • N and 10-15 • S reflect the wet-dry seasonal shift in the Intertropical Convergence Zone (ITCZ) and monsoon movement.They have similar magnitude as the Northern Hemisphere maxima in about a third of the models, however their net contribution to global total NBP seasonal amplitude is small, because they are out of phase and largely cancel each other out.

Figure 6 .
Figure 6.Spatial patterns of GFDL-ESM2M −NBP (gC m −2 day −1 ) changes from 2081 to 2090 and from 1961 to 1970, during mean peak growing season (May-July, first row) and dormant season (October-December, second row) for the esmFdbk2 (first column, constant CO 2 fertilization and changing climate) and esmFixClim2 (second column, constant climate and rising CO 2 ) experiments.The northern high latitude regions show mixed response to climate change during peak growing season (a), and most of the northern temperate and boreal regions see enhanced carbon uptake under elevated CO 2 (b).Net carbon release is increased both under climate change (c) and elevated CO 2 conditions (d), however they have different spatial patterns.

Figure 7 .
Figure 7. Same as Fig. 6, but for IPSL-CM5A-LR.Both the carbon uptake in peak growing season and net carbon release in dormant season are clearly dominated by changes in atmospheric CO 2 rather than climate for this model.

Figure 8 .
Figure 8. Relationship between −NBP change and increase of NBP seasonal amplitude, calculated as the differences between 2081-2090 and 1961-1970 for 10 CMIP5 ESMs.The negative cross-model correlation (R = −0.73,p < 0.05) suggests that a model with a larger net carbon sink increase is likely to simulate a higher increase in NBP seasonal amplitude.

Figure 10 .
Figure 10.CO 2 mean seasonal amplitude (ppm) during 2001-2005 and increase in CO 2 seasonal amplitude at Pt. Barrow during 1974-2005 (% yr −1 , linear trend) from eight CMIP5 ESMs and observation.The big black circle represent surface CO 2 observation at Point Barrow, Alaska (71.3 • N, 156.5 • W; 11 m above sea level).The colored squares represent the CO 2 output at lowest model level (four models at 1000 hPa, and four at 925 hPa) at the original grid that covers Point Barrow from each of the eight models.Error bars indicate ±1 standard error in the trend calculation.Compared to the surface observation, only MPI-ESM-LR overestimate the CO 2 mean seasonal amplitude at Point Barrow, while the other models underestimate this amplitude.Models split between overestimating and underestimating the CO 2 seasonal amplitude increase at Point Barrow.

Figure 11 .
Figure11.Changes of tree cover fractions between future (2081-2090) and historical(1961)(1962)(1963)(1964)(1965)(1966)(1967)(1968)(1969)(1970) periods from six CMIP5 ESMs.The values represent fractional cover changes relative to the whole grid cell, instead of relative change of tree cover.For MPI-ESM-LR and INM-CM4, tree fraction has increased over wide areas of the northern high latitude regions.For MIROC-ESM, tree fraction has generally decreased over the same regions, possibly in response to a hotter and drier climate condition.

Table 1 .
List of Models used and their characteristics.

Table 2 .
Amplitude (maximum minus minimum) of global mean column atmospheric CO 2 , CO 2 growth rate (CO 2 g) and global total NBP, The multi-model ensemble (MME) here is a simple average over the nine models in the table.The values are slightly larger than given in text because of averaging method (in the main text, multi-model averaging of detrended variables are done first, then their amplitude are computed and mean amplitude changes are derived).

Table 3 .
Amplitude increase (ppm) and trends of maxima/minima of surface CO 2 from eight models, their multi-model ensemble (MME), and ESRL's Global mean CO 2 (CO 2 GL).
istics as expected.All models simulate an increase in amplitude, although considerable model spread is found (Table

Table 4 .
Column atmospheric CO 2 and NBP amplitude (computed by CCGCRV, slightly different from max minus min) Increases of nine models by 2081-2090 relative to their 1961-1970 values and their multi-model ensemble (MME).
suggests that this aswww.earth-syst-dynam.net/5/423/2014/Earth Syst.Dynam., 5, 423-439, 2014 Figure 9. CO 2 mean seasonal amplitude (ppm) during 2001-2005 and increase in CO 2 seasonal amplitude at Mauna Loa during 1959-2005 (% yr −1 , linear trend) from eight CMIP5 models and observation.The big black circle represent surface CO 2 observation at Mauna Loa, Hawaii (19.5 • N, 155.6 • W; 3400 m above sea level).The colored squares represent the 700 hPa (close to the altitude of Mauna Loa station surface) CO 2 output at the original grid that covers Mauna Loa from each of the eight models.Error bars indicate ±1 standard error in the trend calculation.Compared to the surface observation, only MPI-ESM-LR and GFDL-ESM2M overestimate CO 2 mean seasonal amplitude at Mauna Loa, while the other models underestimate this amplitude.Models split between overestimating and underestimating the CO 2 seasonal amplitude increase at Mauna Loa.