Agricultural fertilization significantly enhances amplitude of land-atmosphere CO₂ exchange

Increased seasonal amplitude of atmospheric CO₂ measured in remote locations signal that anthropogenic activities have changed terrestrial carbon exchange with the atmosphere¹. Observations illustrate that the seasonal amplitude of CO₂ increased in the Northern Hemisphere, with an increase of 15% over the past 50 years at Mauna Loa, Hawaii, and more than 30% at Point Barrow, Alaska, where the seasonal cycles are much stronger^2,3. These enhancements have largely been attributed to changes in high latitude ecosystems associated with climate change^2,4,5,6 and agriculture^7,8.

Several mechanisms may be responsible for observed changes in the amplitude of seasonal CO₂ fluxes. For example, warmer summers and longer growing seasons across ecosystems in the high latitudes may accelerate rates of photosynthesis and accentuate draw-down of atmospheric CO₂ in the boreal summer^{2,4,5,6,9,10,11}. Forest regrowth in mid-latitudes and forest expansion across high latitudes may similarly increase photosynthetic drawdown of growing season CO₂², as does the CO₂ fertilization of photosynthesis^4,5. Although higher respiratory CO₂ release during the fall may also increase CO₂ amplitude, increased photosynthesis is the dominant driver of the increased seasonal CO₂ amplitude^2,6. All these mechanisms – longer growing seasons, forest expansion, and CO₂ fertilization of photosynthesis – emphasize increases in plant productivity in high latitude ecosystems, a greening trend that has been observed from space¹². However, warmer temperatures and higher CO₂ concentrations cannot fully explain the increased CO₂ seasonal amplitude^2,7,8.

The large increase in agricultural productivity associated with intensified management has also contributed to the increased Northern Hemisphere seasonal CO₂ amplitude^7,8 and may help to explain the missing contribution to increased CO₂ seasonal amplitude. Agricultural practices have intensified over the past half century, initiated by the Green Revolution, increasing overall yield per area primarily through hybridization, irrigation, and fertilization¹³. Atmospheric circulation transports CO₂ from temperate regions into high latitudes¹⁴, so agricultural intensification in the mid-latitudes may contribute to the large increases in seasonal CO₂ amplitude observed at high latitude sites. While the combination of agricultural management and climate may explain the enhanced seasonal amplitude of the atmospheric CO₂, research has yet to explicitly quantify the relative contributions of each. Existing efforts use a synthesis of literature and expert opinion¹⁵, which is limited by the lack of explicit agricultural management representation in models and therefore quantitative attribution. Additionally, no previous work has disentangled specific contributions of individual agricultural management practices like irrigation and fertilization. Here we investigate changes in the seasonal amplitude of atmospheric CO₂ exchange that are simulated by ESMs and investigate the potential role of explicit agricultural management and climate change in driving changes in land-atmosphere CO₂ exchange over the observational record. We illustrate that nitrogen fertilization contributes more to the increased CO₂ amplitude than irrigation or climate change.

Changes in CO₂ amplitude simulated by Earth system models

The Coupled Model Intercomparison Project version 6 (CMIP6) ensemble of models brackets the amplitude of the seasonal cycle of atmospheric CO₂ observations at Point Barrow Alaska, whereas nearly all the CMIP5 ensemble of models underestimate the amplitude (Fig. 1). Earth system models are often used to attribute the underlying mechanisms contributing to the observed carbon cycle changes because large-scale fluxes of carbon between the atmosphere and biosphere are not directly observable and inversion models cannot discern process contributions. However, the capability to attribute changes in the magnitude and trend of CO₂ seasonality has been limited because older versions of ESMs used in CMIP5 underestimate the observed amplitude and the change in amplitude over the observational record at Point Barrow (Fig. 1a) and throughout the Northern Hemisphere high latitudes³. A variety of reasons for model underestimation of CO₂ amplitude and trends have been proposed, including poor representation of high latitude vegetation and soil processes, inaccurate responses to CO₂ fertilization, and an oversimplification in representing key terrestrial processes like nutrient limitation¹¹; these biases have been used as the basis of an emergent constraint on the response of northern ecosystems to elevated CO₂^11,16. However, the robustness of emergent constraints may be affected by the absence of atmospheric transport or representation of other important terrestrial processes in most models, like changes in crop productivity which exert a strong influence on the CO₂ seasonality^7,8. Since the seasonal CO₂ amplitude is an integrated signal of change in terrestrial ecosystem processes, improving process representation to more accurately capture ecosystem dynamics may increase confidence in future carbon cycle projections, including terrestrial carbon sequestration and carbon-climate feedbacks. Newer iterations of models used for CMIP6 updated carbon cycle processes and better capture the amplitude in the seasonal cycle of CO₂ at Point Barrow (Fig. 1b). Given these broad improvements in the CMIP6 models, we take a deeper look at atmospheric CO₂ dynamics that are simulated by successive generations of the Community Earth System Model (CESM).

The CMIP5 (a) and CMIP6 (b) values plotted for each model are the change in the mean CO₂ amplitude over the last ten years of the CMIP5 historical simulations (1995–2005) versus the change in amplitude at the end of the simulation relative to the 1970–1980 mean. Observations are averaged for the same timeframe as model simulations. X-axis error bars represent the standard deviation in the CO₂ amplitude across years.

Full size image

The magnitude of the seasonal CO₂ cycle has improved in successive generations of the CESM (Fig. 2). The CESM1¹⁷ was one of two CMIP5 models that included nitrogen limitation on plant growth. However, like most models participating in the CMIP5 experiment, CESM1 underestimated the mean CO₂ amplitude for 1995–2005 compared to observations³ at Kumukahi, Mauna Loa, Niwot Ridge, and Point Barrow (Fig. 2a–d), and underestimated trends in CO₂ amplitude compared to Point Barrow (Fig. 2h). The magnitude of CO₂ amplitude in CESM version 2¹⁸ is more comparable to observations than CESM1 (Figs. 2; also see Fig. S1 and Table S1), although CESM2 does not clearly improve CO₂ amplitude trends over time. Both model versions underestimate the CO₂ amplitude and the change in amplitude through time at high elevation sites Mauna Loa and Niwot Ridge (Fig. 2b, c, f, g), suggesting some inaccuracies in vertical atmospheric transport^14,19. For example, CESM2 underestimates observed amplitude at higher altitude sites compared to nearby sites at sea level (e.g., comparing Kumukahi, Fig. 2a, and Mauna Loa, Fig. 2b). Analysis verifies that the seasonality of atmospheric CO₂ increases throughout the free troposphere in CESM2, although the seasonality is weaker compared to observations as altitude increases (Fig. S2) highlighting the need for improvements to the representation of vertical CO₂ transport. There have been numerous developments between CESM1 and CESM2 that have improved carbon cycle fluxes (e.g., photosynthetic capacity, soil decomposition, plant nitrogen use) and may contribute to differences in CO₂ amplitude^18,20, with the representation of agriculture²¹ one of the key contributions to changes in CO₂ amplitude (Fig. S3).

Plots illustrate the annual cycle (a–d), averaged for 1995-2005 and trends (e–h) at Cape Kumukahi, HI (19.7^oN,155.0^oW; altitude = 0.3 m), Mauna Loa, HI (19.5^oN, 155.6^oW, altitude = 3397.0 m), Niwot Ridge, CO (40.1^oN, 105.6^oW, altitude = 3523.0 m), and Point Barrow, AK (71.3^oN, 156.61^oW, altitude = 11.0 m) from NOAA’s CO₂ flask network (black), CESM1 (orange), and CESM2 (purple) simulations. Shaded area in a-d illustrates the standard deviation across years. Simulation data are from the model grid cell and atmospheric level nearest the ___location of observation.

Full size image

Agricultural N fertilization and the amplitude of atmospheric CO₂

The role of human activities in emissions-driven ESM simulations is typically represented through external model inputs of fossil fuel emissions and land use change, and few ESMs represent land management processes associated with agriculture²². Agriculture is an important contributor to seasonal carbon fluxes in the extratropics^7,8,23, which underscores the potential importance of representing agriculture in ESMs²⁴. However, in the CMIP5 ensemble few ESMs included a representation of agricultural management, and none included agricultural N fertilization (only two represent N limitation). In contrast, five of the fourteen CMIP6 models contributing ‘esm-hist’ simulations reported cropland N fertilizer fluxes (UKESM-1.0, CESM2, MIROC-ES2L, EC-Earth3-CC, NorESM).

For agricultural systems, which are designed to maximize plant productivity, the addition of industrial N fertilizer is a fundamental aspect of modern crop management and has sustained human population growth^{24,25,26,27,28}. It is arguably the most important contribution to increased crop yields over the past century^21,29. Nitrogen is a key nutrient required for plant growth, and its availability limits terrestrial carbon gain through primary production^30,31. Using the new representation of agricultural management^20,21 in CESM2, we disentangle the contributions of increasing industrial N fertilizer use by comparing simulations with and without industrial N fertilizer in agricultural regions and illustrate that fertilizer contributes ~3 ppm of the Northern Hemisphere seasonal amplitude of CO₂ (Figs. 3a, Fig. S3) and up to 9 ppm in key agricultural regions for 1996–2015 (Fig. 3a). Considering the large impact of agricultural N fertilization on photosynthesis and growth^21,23, its large, local effect on seasonal carbon fluxes is perhaps not surprising. However, the large contribution of industrial N fertilization to the seasonal CO₂ amplitude has not been previously identified due to difficulty in disentangling the signal in observational data and inversion models, as well as the limited representation in process models.

Panel a illustrates the change in mean amplitude of surface CO₂ in CESM2 relative to a simulation without agricultural N fertilization for 1996–2015 (a). Panels b-e compare control CLM5 simulations (fully transient forcings and management) to simulations with individual forcings turned off to illustrate the change in land-only CLM5 simulation amplitude of net biome production (NBP), the net flux of carbon from terrestrial ecosystems to the atmosphere for 1996–2015 due to agricultural N fertilization (b), CO₂ fertilization (c), irrigation (d), and climate change (e). The CO₂ fertilization experiment uses 1850 CO₂ concentrations; and the climate change experiment cycles over 1900–1920 climate. All differences illustrated here are the mean of 1996-2015.

Full size image

Other contributions to seasonal carbon flux amplitude

Agricultural N fertilization causes similar spatial patterns of change in net biome production (NBP) from land-only simulations (Fig. 3b) and CO₂ amplitude from coupled simulations using CESM2 (Fig. 3a). Net biome production represents the net land-air CO₂ flux that determines terrestrial contributions to atmospheric CO₂ concentrations in the CESM2-esm-hist simulations. The seasonal cycle of NBP, therefore, can be used as a proxy for CO₂ amplitude³². Net biome production, however, does not capture the effects of atmospheric transport, as evident in the large changes in CO₂ amplitude north of agricultural regions in East Asia from the CESM2-esm-hist simulations (Fig. 3a) and smaller changes in non-agricultural regions throughout the Northern Hemisphere, compared to the local distribution of NBP fluxes in Fig. 3b. Additionally, while ocean CO₂ fluxes are an important component of the global carbon cycle, they are only a small contribution to the seasonal cycle of CO₂³². Given the dominant role of NBP in driving patterns of change in the CO₂ annual cycle, we leverage data available from the Land Use Model Intercomparison Project (LUMIP) to identify the relative contributions of forcings (CO₂ and climate) and agricultural management (fertilization and irrigation) on terrestrial carbon fluxes to the atmosphere. Analysis focuses on the Community Land Model version 5 (CLM5, the terrestrial component of CESM2) since it is the only model to contribute simulations with and without fertilizer and irrigation.

Agriculture is a key component of the seasonal CO₂ amplitude^7,8, and in our simulations, industrial N fertilizer has the largest impact on the seasonal amplitude of terrestrial carbon fluxes (0.54 g C m⁻² d⁻¹ in 2010–2015, contributing to a 45 ± 2.5% enhancement from the century prior (1910–1915) compared to the enhancement over the same time in the control CLM5 simulation; Figs. 3b, 4). Until now, increasing atmospheric CO₂ concentrations was considered the dominant contribution to the increased amplitude in high latitude carbon fluxes^4,32. Our results highlight the importance of elevated CO₂ in increasing the seasonal amplitude of carbon fluxes throughout the Northern Hemisphere (0.39 g C m⁻² d⁻¹ in 2010–2015, a 40 ± 2.0% increase; Figs. 3c, 4). However, the impact emerges later (Fig. 4b) and is smaller in magnitude (peaks ≤  2 g C m⁻² d⁻¹) than agricultural N fertilization (peaks >  2 g C m⁻² d⁻¹) in most regions (Fig. 3b, c). While the carbon cycle response to CO₂ fertilization is a large uncertainty in previous generations of ESMs³³, the response of CLM5 to elevated CO₂ is within the observed range of responses in Free-Air CO₂ Enrichment (FACE) experiments³⁴.

Plots illustrate the 1996-2015 average annual cycle (a) and historical trends in amplitude (b) of NBP from land-only CESM2 (e.g., CLM5) simulations, averaged north of 30^oN. The CLM5 simulation (purple) includes transient atmospheric forcing and active crop management (irrigation and fertilization). Additional perturbation simulations show the impacts of removing agricultural N fertilization (blue), not irrigating (yellow), preindustrial climate (orange), and preindustrial CO₂ (green). Shaded area in a illustrates the standard deviation across years.

Full size image

Although climate change has a smaller influence on Northern Hemisphere seasonal amplitude of NBP fluxes (18 ± 1.4% by 2010-2015; Fig. 3c, d), it has a strong latitudinal gradient and may act to further increase the CO₂ amplitude in high latitudes (Fig. 3e)^7,32. The effect of climate change across much of the boreal region is similar in magnitude to the effect of CO₂ fertilization, but with little or no changes in mid and low latitudes (Fig. 3e). Warmer temperatures are often associated with a longer growing season in high latitudes and can increase the CO₂ amplitude through increased ecosystem productivity that draws down CO₂, and CO₂ releases due to increased respiration through higher decomposition rates⁶. Previous work indeed highlights that the increasing Northern Hemisphere seasonal carbon amplitude is primarily due to increased gross primary productivity (GPP) that sequesters carbon during the growing season, with smaller contributions from ecosystem respiration^2,6. Additionally, disentangling climate drivers highlights that warmer temperatures contribute to the increased amplitude more than precipitation in most high latitude regions³². In low- and mid-latitudes, however, other studies^4,35,36 found that warmer temperatures decreased the seasonal amplitude of carbon fluxes due to increased water stress.

While other agricultural management practices are critical for increased crop production and can influence the carbon cycle³⁷, they contribute less to the seasonal amplitude of carbon fluxes than N fertilization. For example, although irrigation is an important aspect of agricultural management, we find that its effect on Northern Hemisphere seasonal carbon fluxes is small (3.7 ± 0.46%) due to the relatively small proportion of crop area that is irrigated (Figs. 3d, 4). Irrigation contributions to local carbon fluxes, for example in the Indo-Gangetic Plain, are similar in magnitude to fertilization (Fig. 3d), suggesting that future expansion of irrigation³⁸ may have increasingly important impacts on carbon fluxes. Land conversion from forest to cropland is another important component of agricultural management. Indeed, land cover change significantly impacts the carbon cycle^39,40 and accounts for 18% of the total anthropogenic carbon emissions over the last 50 years⁴¹. However, it does not strongly contribute to changes in the seasonal cycle of Northern Hemisphere carbon fluxes in CLM5 or other models^4,32, but it may contribute to the trends in amplitude change. The design of the LUMIP simulations does not allow for disentangling interactions among the drivers of change, although these are likely to be small relative to the direct effects isolated in this analysis. For example, previous work illustrates that interaction between irrigation and fertilization contributes ~10% of their combined effect on crop yields in CLM²¹.

The increased use of irrigation, fertilizers, and high-yield cultivars over the past fifty years has tripled crop production with important carbon cycle consequences^8,42. Our work illustrates that industrial N fertilization has a large impact on seasonal CO₂ amplitude that has not been previously isolated, accounting for a large fraction of the agricultural contribution to the amplitude^7,8,24. While CO₂ fertilization and climate are also important drivers of changes in seasonal carbon flux amplitude, especially in high latitudes^2,4,32, industrial N fertilization has a larger impact throughout the Northern Hemisphere, contributing to contemporary mean amplitude increases as large as 9 ppm CO₂ in some regions, and increasing land-atmosphere carbon fluxes by 45% across the entire Northern Hemisphere. The magnitude of change due to industrial N fertilization on land-atmosphere carbon fluxes in 2010–2015 (0.54 g C m⁻² d⁻¹) is similar to that of CO₂ and climate combined (0.55 g C m⁻² d⁻¹).

Although the role of agriculture has historically been debated, ESMs that are typically used to disentangle the drivers of changing CO₂ amplitude do not include detailed representations of agriculture or management and therefore suggest that agriculture has a marginal contribution. Other studies^2,7,8 that include agriculture illustrate that intensified agricultural production contributes to 17–45% of the increase in CO₂ amplitude, with lower estimates from studies that do not represent agricultural N fertilizer² or identify impacts from four key crop types⁷ rather than all crops. Our results align with the high-end estimate⁸ that accounts for global agricultural intensification and illustrate that agricultural N fertilization is an important component of seasonal carbon cycle fluxes that is currently overlooked by most ESMs. Given the key role of agricultural N fertilization in regulating carbon fluxes, we strongly encourage that models used for future carbon cycle simulations include a representation of agricultural N fertilization.

Additionally, emergent constraints on photosynthesis^11,16 typically assume that climate change is the primary driver of increased seasonal CO₂ amplitude. Our findings underscore the need to reconsider this assumption, given the important role of agricultural fertilization on seasonal CO₂ amplitude changes. There are, of course, other contributions to increased CO₂ amplitude, like shifts towards younger forests in high latitudes⁴³ or other agricultural management practices⁴⁴ that are not incorporated in these scenarios and might further improve model fidelity with observations. Ultimately, this work highlights that agricultural N fertilization is a profound impact of human activities on seasonal CO₂ amplitude.

Coupled model intercomparison project

Analysis of CMIP simulations uses esm-hist simulations, which are emissions-driven transient historical simulations, where CO₂ is calculated and anthropogenic sources are prescribed⁴⁵. Data from esm-hist simulations used in our analysis includes all available CO₂ data (variable ‘co₂’) from the monthly atmospheric data (Amon) contributed to the Coupled Model Intercomparison Project (CMIP), versions 5 and 6. Data are available through the Earth System Grid Federation (ESGF; https://esgf-node.llnl.gov/projects/esgf-llnl/).

Observations for model evaluation

Model simulations were evaluated against observational data of CO₂ concentrations, which were obtained from the National Oceanic and Atmospheric Administration’s CO₂ flask network⁴⁶ (https://www.esrl.noaa.gov/gmd/ccgg/). Data from Point Barrow, Alaska (71.3^oN, 156.61^oW, altitude = 11.0 m); Niwot Ridge, Colorado (40.1^oN, 105.6^oW, altitude = 3523.0 m); Cape Kumukahi, Hawaii (19.7^oN,155.0^oW; altitude = 0.3 m); and Mauna Loa, Hawaii (19.5^oN, 155.6^oW, altitude = 3397.0 m) were extracted for use in this analysis and compared to the nearest model grid points and atmospheric layers. Data for 24 additional sites were used for extended model evaluation (see Fig. S1).

Simulations of NEE from CLM5 were evaluated at six agricultural flux tower sites (Fig. S4). Sites were selected by combing three flux networks (Ameriflux, FluxNet, and National Ecological Observatory Network (NEON) flux tower networks) to determine locations where observations were collected in agricultural ecosystems. Sites were further filtered to include only those that grew crops simulated by CLM (corn, wheat, soy, cotton, rice, sugarcane) and included more than three years of data with no major gaps, which limited analysis to sites growing corn. Meteorological forcing and NEE evaluation data from the towers were used from the NCAR-NEON system⁴⁷, the Plumbing of Land Surface Models (PLUMBER) project⁴⁸, or downloaded from Ameriflux (https://ameriflux.lbl.gov/sites/site-search/) and FluxNet (https://fluxnet.org/data/fluxnet2015-dataset/). The data from these networks have limited information on management, such as crop rotations, planting dates, and rates and timing of fertilizer applications. Thus, simulations at each site use a single crop type, prognostic planting dates based on growing degree days, and fertilizer rates from the dataset used in global simulations²¹. Since management has strong influence on crop growth and productivity, mismatches between simulations and observations are expected. Nevertheless, our simulations show overall much better agreement with NEE observations when N fertilizer is accounted for than otherwise (Fig. S4), supporting our hypothesis at the site scale.

Community Earth System Model

The Community Earth System Model is an open-source community model that couples atmosphere, ocean, sea ice, land, land ice, and river transport components of the Earth system to simulate biogeophysical and biogeochemical processes for past, present, and future climates¹⁸. The CESM2 includes numerous improvements to previous CESM versions, including low-latitude precipitation, shortwave cloud forcing, and global land carbon accumulation¹⁸ which is affected by improvements to the representation of agriculture, photosynthetic capacity, plant nitrogen use, and vertically-resolved soil decomposition. Of particular relevance for this work is the new capability of the Community Land Model version 5 (CLM5)²⁰, the land component of CESM2, to explicitly simulate agriculture.

The agriculture module simulates global transient agricultural expansion as six major crop types (corn, soy, wheat, cotton, rice, sugarcane) and represents key land management practices (irrigation, fertilization), all of which can impact carbon, nitrogen, water, and energy fluxes between the land and atmosphere. Specifically, simulated crops use the same physiology as natural vegetation (e.g, photosynthesis, stomatal conductance, etc.) but with crop-specific parameters and crop-specific temperature-triggered phenological phases (planting, leaf emergence, grain fill, and harvest), as described in Lombardozzi et al²¹. Industrial fertilizer is prescribed by crop type, year, and country based on fertilization rates from the observationally derived Land Use Harmonization Version 2 (LUH2)⁴⁹ and is applied to all existing and new crop area each year, while irrigation is applied to all irrigated crop fractions if soil moisture is less than a specified threshold. Compared to observations, CLM5 underestimates global crop yields over the past 20 years due to not representing processes related to agricultural intensification (e.g, increased planting density, improved cultivars, etc.), which may contribute to an underestimation of simulated NBP and CO₂ amplitude. Additional details about the crop model can be found in Lombardozzi et al²¹.

Both CESM2 and CLM5 have been extensively evaluated against available observations, which notably illustrate that the cumulative land carbon sink in CLM5 and CESM2 simulations are similar to observationally based estimates^18,50. Specifically, Lawrence et al²⁰. use the International Land Benchmarking system (ILAMB)⁵¹ to evaluate CLM5 net ecosystem exchange (NEE) and GPP against FLUXNET2015⁵² data and also evaluate GPP using FLUXCOM⁵³ and WECANN⁵⁴ datasets. Additionally, Lombardozzi et al.⁴⁷ evaluated CLM5 NEE and GPP at 48 National Ecological Observatory Network sites. Evaluation of CLM5 at six agricultural flux tower sites is included in Fig. S4 and illustrates that CLM5 simulations of NEE with fertilizer are similar in magnitude to observations and reduces model error compared to not including fertilizer, although the model does not capture the interannual variability in NEE amplitude observed at the sites, likely due to limited site-level information on agricultural management (crop rotations, timing and rates of fertilizer applications, planting dates, etc.) and limitations in model representation of agricultural management. Global net ecosystem carbon balance is evaluated against the Global Carbon Project⁴¹ and Hoffman et al.⁵⁵ datasets, and response to elevated CO₂ and nitrogen additions is evaluated against meta-analyzes of FACE and nitrogen addition experiments³⁴. Specific to the Arctic-boreal zone, Birch et al.⁵⁶ found that the timing of leaf onset and offset in CLM5 is typically later than observed. Additionally, CLM5 overestimates GPP and ecosystem respiration in this region due to inaccuracies in the representation of phenology, allocation, and photosynthetic capacity parameterizations, potentially contributing to inaccuracies in the simulated magnitude of NBP and CO₂ amplitude.

Simulations from CESM2¹⁸ are a large focus of the analysis because this is the only model contributing simulations to CMIP6 and LUMIP projects that includes a comprehensive representation of agriculture, including spatiotemporally changing distributions (based on the Land Use Harmonization version 2, LUH2)^20,49, as well as explicit crop types and phenological phases, and fertilization and irrigation management (described above). Within CESM, the transient land use datasets are produced using the CLM Land Use Data Tool (http://www.cgd.ucar.edu/iam/projects/thesis/thesis-clm-landuse-tool.html), which merges present-day land cover distribution with historical and future transitions described in LUH2⁴⁹, translating these into changing distributions of the plant categories and management information used in CESM²⁰. All simulations used in this analysis were run from 1850 through 2014 at nominal 1^o horizontal resolution.

The impact of individual forcings on terrestrial carbon fluxes were identified using additional CESM2 simulations. To understand the impact of agricultural management, two simulations were branched from the CESM2-esm-hist CMIP6 simulation in 1970 (data available at the UCAR/NCAR Climate Data Gateway doi:10.5065/q6br-jr78), when the impact of agriculture is apparent on carbon cycle processes²¹. The first of these turned off the explicit representation of agriculture so that all crop areas are represented as generic C₃ crops, where crop phenology is simulated as C₃ grasses and do not include irrigation or fertilization (referred to as ‘generic crop’; this is the CESM1 representation of agriculture). The second uses the explicit representation of agriculture but removes industrial N fertilization from all existing and new crop area (referred to as ‘no fertilization’). To ensure that changes in CO₂ fluxes were minimally impacted by model drift, each simulation equilibrated carbon fluxes in 1970 by cycling over a single year of forcing for ten years. The CESM2 simulated these alternative representations of agriculture in a CO₂ emissions-forced historical scenario following the esm-hist experimental protocol.

Land Use Model Intercomparison Project

We additionally leveraged land-only CLM5 simulations available through the CMIP6 Land Use Model Intercomparison Project⁵⁷ (LUMIP; data available through ESGF: https://esgf-node.llnl.gov/projects/cmip6/) to understand the impacts of anthropogenic forcings (CO₂ and climate), agricultural management (irrigation and fertilization), and land cover change on carbon fluxes. We used data from land-hist (fully transient historical, including transient CO₂, N-deposition, aerosol deposition, etc.), and other simulations that were the same as land-hist except with a single factor held constant: land-cCO2 (constant CO₂), land-cClim (constant climate), land-crop-noIrrig (irrigated area held at 1850 levels), and land-crop-noFert (fertilization rates and area held at 1850 levels/distribution; see Table 2 in Lawrence et al⁵⁷. for detailed description of LUMIP and these simulations). The LUMIP simulations use the land component of ESMs run with GSWP3v1 climate forcing data (http://hydro.iis.u-tokyo.ac.pf/GSWP3), and we focus analysis on CLM5, the land model used within CESM2, because it is the only LUMIP model that submitted simulations with agricultural fertilization and irrigation. Because CO₂ concentrations are not available from land-only simulations, we use NBP, the net carbon flux from the land surface (the carbon flux variable that is passed to the atmosphere in coupled simulations), as a proxy for the atmospheric CO₂ concentration dynamics. We recognize that land-only simulations cannot account for the impacts of circulation and oceanic CO₂ uptake on the distribution and seasonal cycling of CO₂.

Amplitude and trend calculations

Seasonality in terrestrial carbon fluxes is due to many factors, including GPP, ecosystem respiration, and carbon losses from fires and land use change. Our analysis focuses on the seasonal amplitude of CO₂ from coupled simulations (CMIP and CESM) and NBP fluxes from land-only simulations (CLM5). Net biome production is the terrestrial carbon flux from CLM5 that is passed to the atmosphere in coupled simulations. We quantify the seasonal amplitude of CO₂ and NBP as the absolute difference between the annual maximum and minimum monthly averaged fluxes, representing the net source or sink of carbon fluxes to the atmosphere. Note that maximum NBP values indicate the peak net carbon sink. The linear trend in amplitude through time was removed from the calculations of the mean annual cycle. All amplitude analysis focuses on the Northern Hemisphere, > 30 ^oN, unless otherwise specified. The Northern Hemisphere has larger land area and more dramatic changes in the seasonal cycle of carbon fluxes.

The trends in amplitude for CESM1 and CESM2 were evaluated against trends in observational data (see description above) from the National Oceanic and Atmospheric Administration’s CO₂ flask network⁴⁶ (https://www.esrl.noaa.gov/gmd/ccgg/). The model data were available from 1970 through 2020 by combining historical (through 2005 for CESM1 and 2014 for CESM2) and future (RCP8.5 for CESM1 and SSP5-85 for CESM2) scenarios. Each station has a different record length within this time, so we truncate the model timeseries to have the same start and end year as the observational record for each station. Some stations have missing data for individual years within the middle of their time series, and our approach effectively assumes zero anomalies from climatology for these years. Regression slopes (b) and intercepts (a) were calculated using a least squares linear regression approach. Uncertainties on the slope (b) were calculated by resampling the residuals from the original regression slope. The original data were given by CO₂(y)  =  a  +  by + e(y), where y is year and e represents the residuals from the original regression slope for each year. For a timeseries of length N, we randomly sample, with replacement, N values from the distribution of e(y), generate a new time series CO₂’(y) = a + by + e’(y) where e’(y) are the randomly resampled residuals, and compute the least squared regression slope of this new time series. We repeat this 1000 times to generate a distribution of 1000 estimates of the regression slope, centered on the original, and estimate the uncertainty on the regression slope as the 2.5^th and 97.5^th percentile range of these bootstrapped slopes, encompassing 95% of the resampled slopes.