Compensatory Effects Between CO2, Nitrogen Deposition, and Nitrogen Fertilization in Terrestrial Biosphere Models Without Nitrogen Compromise Projections of the Future Terrestrial Carbon Sink

Although terrestrial biosphere models (TBMs) with and without nitrogen cycling successfully reproduce the historical terrestrial carbon sink, the influence of nitrogen cycling under interacting and intensifying global change drivers in the future is unclear. Here, we compare TBM projections with and without nitrogen cycling over alternative future scenarios (the Shared Socioeconomic Pathways) to examine how representing nitrogen cycling influences CO2 fertilization as well as the effects of a comprehensive group of physical and socioeconomic global change drivers. Because elevated nitrogen deposition and nitrogen fertilization have stimulated terrestrial carbon sequestration over the historical period, a model without nitrogen cycling must exaggerate the strength of CO2 fertilization to compensate for these unrepresented nitrogen processes and to reproduce the historical terrestrial carbon sink. As a result, it cannot realistically project the future terrestrial carbon sink, overestimating CO2 fertilization as the trajectories of CO2, nitrogen deposition and nitrogen fertilization diverge in future scenarios.


Introduction
The terrestrial C sink has increased over recent decades driven primarily by CO 2 fertilization and it currently sequesters approximately 30% of anthropogenic CO 2 emissions Walker et al., 2020). The persistence of the terrestrial C sink over the 21st century is uncertain due to the combined influences of multiple global change drivers-rising CO 2 alongside rising temperature, varying precipitation, and land use change. In particular, nitrogen (N) is an essential limiting nutrient (Elser et al., 2007;Fernández-Martínez et al., 2014;Greaver et al., 2016;LeBauer & Treseder, 2008;Wright et al., 2018) and constrains CO 2 fertilization (Terrer et al., 2019;S. Wang et al., 2020). However, intensifying N fertilizer use and fossil fuel burning cause elevated N deposition which could alleviate N limitation (O'Sullivan et al., 2019;R. Wang et al., 2017). Elevated temperature drives soil organic matter decomposition which releases N, that is, N mineralization, and this could further alleviate N limitation (Liu et al., 2017). Consequently, the extent to which N limitation will constrain the future terrestrial C sink under this cast of interacting and intensifying global change drivers is unresolved.
Abstract Although terrestrial biosphere models (TBMs) with and without nitrogen cycling successfully reproduce the historical terrestrial carbon sink, the influence of nitrogen cycling under interacting and intensifying global change drivers in the future is unclear. Here, we compare TBM projections with and without nitrogen cycling over alternative future scenarios (the Shared Socioeconomic Pathways) to examine how representing nitrogen cycling influences CO 2 fertilization as well as the effects of a comprehensive group of physical and socioeconomic global change drivers. Because elevated nitrogen deposition and nitrogen fertilization have stimulated terrestrial carbon sequestration over the historical period, a model without nitrogen cycling must exaggerate the strength of CO 2 fertilization to compensate for these unrepresented nitrogen processes and to reproduce the historical terrestrial carbon sink. As a result, it cannot realistically project the future terrestrial carbon sink, overestimating CO 2 fertilization as the trajectories of CO 2 , nitrogen deposition and nitrogen fertilization diverge in future scenarios.

Plain Language Summary
Terrestrial biosphere models simulate the terrestrial carbon sink (in plant and soil biomass), which takes up a third of anthropogenic CO 2 emissions. However, the influence of nitrogen cycling and nitrogen limitation of plant growth in the Shared Socioeconomic Pathways (SSPs) (recently adopted alternative scenarios for the 21st century) is unclear. Here, we compare a model with and without nitrogen cycling in comprehensive simulations of the SSPs. We find that a model without nitrogen must exaggerate the influence of elevated CO 2 on plant growth to compensate for unrepresented nitrogen cycling processes in order to correctly simulate the historical terrestrial carbon sink. Specifically, it does not represent how elevated atmospheric nitrogen input (due to intensive agriculture and fossil fuel burning) and elevated nitrogen fertilization have increased plant growth over the historical period. As a result, models without nitrogen cannot realistically project the future terrestrial carbon sink because they are calibrated to reproduce the historical terrestrial carbon sink but the trajectories of CO 2 , atmospheric nitrogen input, and nitrogen fertilization diverge from their historical trajectories in future scenarios. This leads to an overestimation of the future terrestrial carbon sink by models without nitrogen with implications for climate change projections and policy.
Terrestrial biosphere models (TBMs) are the principal tool for simulating the terrestrial C sink, serving as the land components of Earth System Models (ESMs). In TBMs, CO 2 fertilization is suggested to be overestimated without a representation of N limitation: When TBM projections of the future terrestrial C sink under a high atmospheric CO 2 scenario were constrained with observations of N inputs and N stoichiometry, terrestrial C sequestration was reduced by 20% (Wieder et al., 2015). However, this study only examined estimated constraints applied post hoc to simulations of TBMs without N cycling rather than including an explicit representation of N cycling in TBMs. N limitation of CO 2 fertilization has motivated an increasing number of TBMs to include a representation of coupled C and N cycling. However, only around half of TBMs include N cycling: 10 out of 17 TBMs represented coupled C and N cycling in the most recent Global Carbon Project Kou-Giesbrecht et al., 2023).
Despite their major structural differences, intercomparisons between TBMs with and without N cycling have found that TBMs with N cycling and TBMs without N cycling perform similarly in reproducing the historical terrestrial C sink Kou-Giesbrecht et al., 2023;Seiler et al., 2022). That both TBMs with and without N cycling reproduce the historical terrestrial C sink suggests that TBMs are implicitly (or explicitly) tuned to reproduce the historical terrestrial C sink and as such even TBMs without N cycling intrinsically account for N limitation given historical N inputs. Importantly however, intercomparisons between TBMs with and without N cycling have found that their responses to individual global change drivers acting over the historical period diverge (Canadell et al., 2021;Huntzinger et al., 2017). In particular, the strength of the CO 2 fertilization effect over the historical period simulated by TBMs without N cycling was found to be over twice that simulated by TBMs with N cycling. Elevated N deposition increased terrestrial C sequestration by approximately 20% over the historical period as simulated by TBMs with N cycling (but was not represented in TBMs without N cycling). TBMs with N cycling simulated lower CO 2 emissions from decomposition than TBMs without N cycling in response to historical climate change. These divergent responses to individual global change drivers acting over the historical period between TBMs with and without N cycling have unclear consequences for projecting the future terrestrial C sink. Few studies have comprehensively examined N cycling in plausible future scenarios that encompass all global change drivers, either examining historical simulations or focusing solely on CO 2 and/or N deposition for a single (often "worst case") future scenario (Goll et al., 2012;Smith et al., 2014;Sokolov et al., 2008;Thornton et al., 2009;Y. P. Wang et al., 2015;Wårlind et al., 2014;Zaehle et al., 2010).
Here we use the Canadian Land Surface Scheme including Biogeochemical Cycles (CLASSIC), the land component of the Canadian Earth System Model (CanESM5), to examine how representing coupled C and N cycling influences the response of terrestrial C sequestration to global change by comparing simulations of CLASSIC with and without coupled C and N cycling over the 21st century. By examining a single TBM with and without coupled C and N cycling, we can isolate the impact of explicitly representing coupled C and N cycling whereas intercomparisons across TBMs with and without N cycling are confounded by other structural and parametric differences between TBMs that may obscure the effects of N cycling. CLASSIC represents both flexible vegetation C:N stoichiometry and the upregulation of symbiotic biological N fixation under N limitation (described and evaluated in Asaadi and Arora (2021) and Kou-Giesbrecht and Arora (2022)) thereby presenting an advanced representation of coupled C and N cycling in a TBM. We evaluate the role of N cycling under individual and combined contributions of a comprehensive group of physical and socioeconomic global change drivers: CO 2 , climate, N deposition, and land use (crop area and N fertilization). We simulate the historical period and three alternative future scenarios that are based on the Shared Socioeconomic Pathways (SSPs) (Riahi et al., 2017), which are the most recent scenarios adopted by the Intergovernmental Panel on Climate Change (IPCC): SSP126 ("sustainability") has low greenhouse gas emissions, SSP370 ("regional rivalry") has high greenhouse gas emissions, and SSP585 ("fossil-fueled development") has very high greenhouse gas emissions.

CLASSIC Overview
The CLASSIC (Melton et al., 2020;Seiler et al., 2021) is the land component in the family of the Canadian Earth System Models (CanESM) (Swart et al., 2019). CLASSIC simulates land-atmosphere fluxes of energy, momentum, water, carbon (C), and nitrogen (N). The physical component of CLASSIC simulates fluxes of energy, momentum, and water (Verseghy, 1991;Verseghy et al., 1993). The biogeochemical component of CLAS-SIC simulates the land-atmosphere exchange of C via photosynthesis, autotrophic respiration, heterotrophic respiration, land use change, and fire (Arora & Boer, 2005a). For biogeochemical processes, vegetation is partitioned into nine plant functional types (PFTs): needleleaf evergreen trees, needleleaf deciduous trees, broadleaf evergreen trees, broadleaf cold deciduous trees, broadleaf drought deciduous trees, C 3 crops, C 4 crops, C 3 grasses, and C 4 grasses. CLASSIC prognostically simulates the amount of C in vegetation, litter, and soil organic matter pools for each PFT and over the bare soil fraction in each grid cell. Vegetation C is represented by leaf, stem, and root, each of which consists of structural and non-structural carbohydrate pools. CLASSIC simulates the land-atmosphere exchange of N via biological N fixation (free-living and symbiotic), specified N deposition and N fertilizer application, nitric oxide (NO) emissions, nitrous oxide (N 2 O) emissions, N 2 emissions, ammonia (NH 3 ) volatilization, N leaching, and land use change (Asaadi & Arora, 2021;Kou-Giesbrecht & Arora, 2022). CLASSIC prognostically simulates the amount of N in vegetation, litter, soil organic matter, and inorganic soil N (ammonium (NH 4 + ) and nitrate (NO 3 − )) pools for each PFT and over the bare soil fraction in each grid cell. Vegetation N is represented by leaf, stem, and root, each of which consists of structural and non-structural N pools. See Text S1 in Supporting Information S1 for a more detailed description of the physical and biogeochemical components of CLASSIC.
In CLASSIC-CN, photosynthesis is dependent on leaf N such that, when leaf N is low, photosynthesis is downregulated and, when leaf N is high, photosynthesis is upregulated (Asaadi & Arora, 2021). Additionally, vegetation exhibits a dynamic response to N limitation of vegetation growth. First, vegetation upregulates and downregulates symbiotic biological N fixation in response to weak N limitation and strong N limitation, respectively (Kou-Giesbrecht & Arora, 2022). Second, vegetation has flexible stoichiometry and thus the vegetation C:N ratio responds to changing N limitation (Asaadi & Arora, 2021). In CLASSIC-C, N cycling is turned off and the downregulation of photosynthesis under increasing CO 2 is controlled by a parameter as explained in Arora et al. (2009). Briefly, this parameter, which ranges between 0 and 0.9, determines the rate of increase of photosynthesis with increasing CO 2 . When it is set to 0, photosynthesis does not increase with increasing CO 2 . When it is set to 0.9, photosynthesis increases with increasing CO 2 at an unconstrained rate. When it is set to 0.35, CLASSIC-C simulations estimate a global net atmosphere-land CO 2 flux that lies within uncertainty range of estimates from the Global Carbon Project .

Simulations
We use CLASSIC-C and CLASSIC-CN to simulate energy, momentum, water, C, and N fluxes at the global scale over the historical period (1851-2014) and over the future period (2015-2100) for three SSPs; SSP126, SSP370, and SSP585. For the historical period, we conducted simulations following the TRENDY protocol (for contributions to the Global Carbon Project ). We also conducted historical simulations following the Inter-Sectoral Impact Model Intercomparison Project (ISIMIP) protocol (Warszawski et al., 2014) in order to launch future simulations following the ISIMIP protocol. Forcings are described in Table S1 in Supporting Information S1. For both the historical period and future period, we conducted simulations with all global change drivers acting concurrently as well as four separate simulation experiments to disentangle the contributions of CO 2 , climate, N deposition, and land use (which includes changes to both crop area and to N fertilization) to the global net atmosphere-land CO 2 flux. We did not isolate the influence of population density and CH 4 because these forcings regulate fire C emissions and soil CH 4 fluxes respectively, and have minimal influence on the global net atmosphere-land CO 2 flux in comparison to CO 2 , climate, land use, and N deposition. Simulations are described in detail in Text S2 and in Tables S2 and S3 in Supporting Information S1.
Over the historical period, we compared the global net atmosphere-land CO 2 flux simulated by CLASSIC-C and CLASSIC-CN to estimates from the Global Carbon Project . CLASSIC-C and CLASSIC-CN simulations over the historical period have been validated previously in other studies (Asaadi & Arora, 2021;Kou-Giesbrecht & Arora, 2022;Melton et al., 2020;Seiler et al., 2021).

Historical Net Biome Productivity
The net biome productivity (NBP) is the global net atmosphere-land CO 2 flux and is the sum of CO 2 sequestration in growing vegetation (net primary productivity), CO 2 emissions from decomposing soil organic matter (heterotrophic respiration), and CO 2 emissions from disturbance. NBP quantifies the terrestrial C sink (or source). NBP simulated by CLASSIC-C (which does not represent N cycling) and CLASSIC-CN (which represents coupled C and N cycling) are similar over the historical period (Figure 1). Positive NBP values since the 1960s indicate a terrestrial C sink over the historical period. Between 2000 and 2010, CLASSIC-C and CLASSIC-CN simulate a terrestrial C sink of 1.3 Pg C yr −1 and 1.2 Pg C yr −1 , respectively, which lie within the uncertainty range of estimates from the Global Carbon Project (1.2 ± 0.3 Pg C yr −1 ). NBP is driven by contributions from a comprehensive group of physical and socioeconomic global change drivers: CO 2 , climate, N deposition, and land use (crop area and N fertilization) (Figure 2). Despite both CLASSIC-C and CLASSIC-CN exhibiting a similar NBP when all global change drivers act concurrently, the cumulative NBP contribution from each global change driver over the historical period differ between CLASSIC-C and CLASSIC-CN (Figure 3a and Figure S1 in Supporting Information S1). CLASSIC-C exhibits a significantly stronger NBP increase driven by elevated CO 2 over the historical period than CLASSIC-CN (due to higher vegetation growth; Figure S2 in Supporting Information S1). Climate change over the historical period decreases NBP due to increased decomposition (associated with elevated temperature). In CLASSIC-C there is more standing litter and soil C which results in higher decomposition and thus a slightly stronger NBP decrease in response to climate change than in CLASSIC-CN ( Figure S3 in Supporting Information S1). While climate change enhances decomposition at lower latitudes, it enhances vegetation growth at higher latitudes (Figures 3c and 3e). CLASSIC-C does not represent the effects of N deposition or N fertilization. In CLASSIC-CN, both elevated N deposition and N fertilization stimulate NBP over the historical period (due to increased vegetation growth; Figure S2 in Supporting Information S1). Land use change over the historical period decreases NBP due to CO 2 emissions associated with the conversion of natural vegetation to crops. This NBP decrease is weaker for CLASSIC-CN than for CLASSIC-C because it is offset by increased vegetation growth associated with N fertilization. In CLASSIC-CN, the contribution of N deposition is stronger at higher latitudes (Figures 3c and 3e), which are more N-limited (Kou-Giesbrecht & Arora, 2022).
Overall, CLASSIC-C and other TBMs without N cycling do not represent the stimulation of terrestrial C sequestration by N deposition or N fertilization over the historical period (Huntzinger et al. (2017) and Figure S1 in Supporting Information S1). Therefore, the stimulation of terrestrial C sequestration by elevated CO 2 over the historical period must be exaggerated to compensate for these unrepresented N processes in order to reproduce the historical terrestrial C sink. However, this introduces compensatory effects. We now explore the consequences of these compensatory effects in projections of the future terrestrial C sink.

Future Net Biome Productivity
Despite the agreement between NBP simulated by CLASSIC-C and CLASSIC-CN over the historical period and their agreement with NBP estimates from the Global Carbon Project, there are major differences between NBP projected by CLASSIC-C and CLASSIC-CN over the 21st century, especially in future scenarios characterized Figure 1. Net biome productivity (NBP) over the historical period and for SSP126 ("sustainability"), SSP370 ("regional rivalry"), and SSP585 ("fossil-fueled development") simulated by CLASSIC-C (which does not represent N cycling) and CLASSIC-CN (which represents coupled C and N cycling). The yellow boxes indicate NBP estimates from the Global Carbon Project. by high atmospheric CO 2 . Between 2090 and 2100, projections of NBP by CLASSIC-C and CLASSIC-CN differ by 0.5 Pg C yr −1 for SSP126 ("sustainability"), by 1.4 Pg C yr −1 for SSP370 ("regional rivalry"), and by 3.4 Pg C yr −1 for SSP585 ("fossil-fueled development") ( Figure 1).
Over the historical period, CO 2 , N deposition, and N fertilization all increase simultaneously. In SSP126, CO 2 stabilizes then decreases after 2050 while N deposition decreases and N fertilization increases ( Figure 2). Thus, CLASSIC-C and CLASSIC-CN project decreasing NBP over the 21st century and the terrestrial C sink transitions to a terrestrial C source because vegetation growth plateaus and is surpassed by decomposition which persists due to its longer timescale ( Figure S4 in Supporting Information S1). In SSP370, although CO 2 , N deposition, and N fertilization increase simultaneously as in the historical period, CO 2 increases at a faster rate than N deposition and N fertilization. Under SSP585, CO 2 increases whereas N deposition and N fertilization peak then decrease (Figures 2f  and 2g). Over the historical period, the contribution of N deposition change to cumulative NBP was 20% (19.6 Pg C of 98.4 Pg C) (Figure 3a and Table S4 in Supporting Information S1). The contribution of N deposition change to cumulative NBP is much weaker in future scenarios than over the historical period. For SSP370 and SSP585, the contribution of N deposition change to cumulative NBP is only 2% (6.6 Pg C of 324.5 Pg C) and 1% (4.4 Pg C of 425.1 Pg C), respectively (Figure 3b and Table S4 in Supporting Information S1). While it is not possible to disentangle the effects of crop area from those of N fertilization, the contribution of land use change to cumulative NBP is much weaker in future scenarios than over the historical period. Whereas over the historical period the contribution of land use change to cumulative NBP is 25% (24.3 Pg C of 98.4 Pg C), the contribution of land use change to cumulative NBP are only 12% (38.4 Pg C of 324.5 Pg C) and 8% (33.0 Pg C of 425.1 Pg C) for SSP370 and SSP585, respectively, for CLASSIC-CN (Figure 3b and Table S4 in Supporting Information S1). Additionally, as temperature increases in future scenarios, N mineralization is stimulated in CLASSIC-CN. N mineralization alleviates N limitation as indicated by higher vegetation growth in CLASSIC-CN than in CLASSIC-C ( Figure S2 in Supporting Information S1), leading to a weaker NBP decrease from climate change in CLASSIC-CN than in CLASSIC-C (Figure 3b). Overall, N deposition and N fertilization contribute less to terrestrial C sequestration in future scenarios than over the historical period, whereas CO 2 contributes substantially more to terrestrial C sequestration in future scenarios than over the historical period. Because the CO 2 contribution relative to the N deposition and N fertilization contributions differs in future scenarios from the historical period, the difference between CLASSIC-C and CLASSIC-CN is considerable in future scenarios. For SSP370, CLASSIC-CN projects a terrestrial C sink that is 1.4 Pg C yr −1 lower (58% lower) than that projected by CLASSIC-C at the end of the 21st century (Figure 1). For SSP585, the discrepancy between CLASSIC-C and CLASSIC-CN is substantial: CLASSIC-CN projects a terrestrial C sink that is 3.4 Pg C yr −1 lower (64% lower) than that projected by CLASSIC-C at the end of the 21st century (Figure 1).

Discussion
We show here that, while a TBM that does not represent coupled C and N cycling can reproduce the historical terrestrial C sink, it cannot realistically project the future terrestrial C sink. This is because, in a TBM that does not represent coupled C and N cycling, calibrating terrestrial C sequestration to reproduce the historical  , and (f) in extratropical regions for future scenarios (2015-2100), simulated by CLASSIC-C (which does not represent N cycling) and CLASSIC-CN (which represents coupled C and N cycling). The black dots indicate the cumulative NBP (i.e., the sum of the contributions of all drivers to cumulative NBP). The red line indicates the cumulative NBP estimate over the historical period  from the Global Carbon Project. 7 of 10 terrestrial C sink in the absence of N cycling introduces compensatory effects: the stimulation of terrestrial C sequestration by elevated CO 2 over the historical period (i.e., the CO 2 fertilization effect) must be exaggerated to compensate for the absence of N cycling, that is, the stimulation of terrestrial C sequestration by elevated N deposition and N fertilization over the historical period. The result of these compensatory effects is that a TBM that does not represent coupled C and N cycling but reproduces the historical terrestrial C sink correctly cannot realistically project the future terrestrial C sink as global change drivers follow divergent trajectories and occur in unprecedented combinations. Specifically, it will overestimate the CO 2 fertilization effect as CO 2 increases faster than N deposition and N fertilization (which also decrease in some future scenarios) as well as temperature which drives increasing N mineralization. This is supported by our simulations of the terrestrial C sink in future scenarios characterized by high atmospheric CO 2 : the TBM used here with coupled C and N cycling projects a global net atmosphere-land CO 2 flux that is between 1.4 and 3.4 Pg C yr −1 lower (58%-64% lower) than that projected by the TBM used here without coupled C and N cycling at the end of the 21st century.
Numerous lines of evidence have suggested the importance of N limitation in constraining CO 2 fertilization, including meta-analyses of elevated CO 2 experiments (Terrer et al., 2019) and temporal analyses of satellite-based estimates of terrestrial photosynthesis paired with foliar N observations (S. Wang et al., 2020). Consistent with these studies, we show that explicitly representing N limitation in a TBM reduces the CO 2 fertilization effect. This has been previously recognized in TBM development (Thomas et al., 2015) and has motivated the incorporation of N cycling in TBMs (Davies- Barnard et al., 2020). Of novel importance, we show that N limitation of the CO 2 fertilization effect operates within the framework of N cycling as a whole and that the influences of human-induced N deposition, N fertilization, and N mineralization must be considered as they operate in tandem with N limitation of CO 2 fertilization.
These results extend to climate change projections by ESMs, which also include the feedbacks between climate and the terrestrial biosphere. CLASSIC serves as the land component of CanESM5, which contributed to the sixth phase of the Coupled Model Intercomparison Project (CMIP6) and the IPCC Assessment Reports. Within the C 4 MIP-CMIP6 ensemble, ESMs with terrestrial N cycling exhibit a 25%-30% lower CO 2 fertilization effect than ESMs without terrestrial N cycling (Canadell et al., 2021). For CanESM5 with the older version of CLASSIC without N cycling, the projected global net atmosphere-land CO 2 flux reaches a staggering 12.0 Pg C yr −1 at the end of the 21st century for the future scenario with the highest CO 2 (SSP585). This estimate was the highest among participating ESMs, despite CanESM5's ability to reproduce several aspects of the historical global C budget (Arora & Scinocca, 2016), and was closely followed by estimates from two other ESMs without N cycling Koven et al., 2022). Despite the critical importance of N cycling and the complexity of its interactions with C processes, only 6 out of 11 ESMs include coupled C and N cycling within C 4 MIP-CMIP6 and the sixth IPCC Assessment Report (Canadell et al., 2021).

Conclusions
Our analyses show that reproduction of the historical terrestrial C sink, which is achieved successfully by most TBMs regardless of whether they represent N cycling , cannot be considered an indicator for the reliability of their projections of the future terrestrial C sink. Our findings show that a TBM that does not represent coupled C and N cycling cannot realistically represent the combined influences of multiple global change drivers, overestimating CO 2 fertilization as the trajectories of CO 2 , N deposition, and N fertilization diverge over the 21st century. Scaling fundamental ecological understanding of C and N interactions to the global scale through the explicit representation of physical and biological processes rather than calibration to reproduce the historical terrestrial C sink is key for realistically projecting the future terrestrial C sink under global change with TBMs and ultimately climate change with ESMs.

Data Availability Statement
The source code for CLASSIC is available on the CLASSIC community Zenodo page (https://doi.org/10.5281/ zenodo.6499554).