Reducing uncertainty in projections of terrestrial carbon uptake

Carbon uptake by the oceans and terrestrial biosphere regulates atmospheric carbon dioxide concentration and affects Earth’s climate, yet global carbon cycle projections over the next century are highly uncertain. Here, we quantify and isolate the sources of projection uncertainty in cumulative ocean and terrestrial carbon uptake over 2006–2100 by performing an analysis of variance on output from an ensemble of 12 Earth System Models. Whereas uncertainty in projections of global ocean carbon accumulation by 2100 is <100 Pg C and driven primarily by emission scenario, uncertainty in projections of global terrestrial carbon accumulation by 2100 is >160 Pg C and driven primarily by model structure. To statistically reduce uncertainty in terrestrial carbon projections, we devise schemes to weight the models based on their ability to represent the observed change in carbon accumulation over 1959–2005. The weighting schemes incrementally reduce uncertainty to a minimum value of 125 Pg C in 2100, but this reduction requires an impractical observational constraint. We suggest that a focus on reducing multi-model spread may not make terrestrial carbon cycle projections more reliable, and instead advocate for accurate observations, improved process understanding, and a multitude of modeling approaches.


Introduction
The accumulation of carbon in the ocean and terrestrial biosphere reduces the atmospheric carbon dioxide (CO 2 ) burden and thus the influence of anthropogenic carbon emissions on global climate (Ciais and Sabine 2013). Observations and models suggest that the ocean has accumulated 170 Pg C since 1750 (Le Quéré et al 2015). This uptake of CO 2 by the ocean is primarily driven by the fast response of air-sea gas exchange to increasing CO 2 in the atmosphere and the slow transport of this CO 2 from the surface into the ocean interior (Graven et al 2012). In contrast, the terrestrial biosphere is estimated to have lost 25 Pg C since 1750, the sum of 190 Pg C emitted from land use change and 165 Pg C accumulated by terrestrial ecosystems (Le Quéré et al 2015). Models and observations suggest that the accumulation or loss of carbon by these reservoirs is changing with time. For example, from the decade of the 1960s to the most recent decade (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014), the rate of ocean carbon uptake increased from 1.1 to 2.6 Pg C yr À1 , land use change emissions decreased from 1.5 to 0.9 Pg C yr À1 , and the residual land sink grew from 1.7 to 3.0 Pg C yr À1 (Le Quéré et al 2015).
Future changes in the accumulation of carbon in the ocean and terrestrial biosphere will affect atmospheric CO 2 concentration and thus climate, yet many studies report high uncertainty in projections of air-sea CO 2 flux, land use change emissions, and net terrestrial biospheric production over the coming century (Friedlingstein et al 2006, Arora et al 2013, Hoffman et al 2014, Hewitt et al 2016. Quantifying the relative importance of the sources of uncertainty in these projections, namely internal variability, emission scenario, and model structure, is a necessary step to realizing reductions in projection uncertainty (Hawkins and Sutton 2009 climate processes (e.g. El Niño-Southern Oscillation). Emission scenario uncertainty is caused by unknowns in the future behavior of society, while model structural uncertainty is a product of different representations of the physical climate system and the biology of the biosphere across a range of Earth System Models.
Here, we use output from models participating in the 5th Coupled Model Intercomparison Project (CMIP5) to quantify and assess the relative importance of the sources of projection uncertainty in globally-integrated, cumulative uptake of carbon by the ocean and terrestrial biosphere over the next century. Results from this analysis are used to inform a strategy for reductions in terrestrial carbon uptake projection uncertainty.

CMIP5 models
We analyze output from 12 CMIP5 earth system models that simulated the historical period  and the future period (2006-2100) under a common set of anthropogenic forcings, including the historical atmospheric CO 2 concentration and the projected atmospheric CO 2 concentration from a collection of 4 emission scenarios (Representative Concentration Pathways, or RCPs). Table S1 (available at stacks.iop. org/ERL/12/044020/mmedia) shows the CMIP5 models analyzed in this study, and lists the RCP simulations conducted for each model.

Analysis of variance
We analyze future projections of ocean and terrestrial carbon uptake from the CMIP5 models to quantify their uncertainty and to partition this uncertainty into three sources: internal variability, emission scenario, and model structure.
The globally-integrated, cumulative carbon uptake by land or ocean since 2006, T(m, s, t), is a function of model m, emission scenario s, and time t. Each individual prediction, T(m, s, t) was fit with a 4th order polynomial over the years 2006-2100, F(m, s, t), resulting in a timeseries of residuals R(m, s, t), Rðm; s; tÞ ¼ Tðm; s; tÞ À Fðm; s; tÞ: The projection uncertainty is quantified on the basis of ensemble spread, or the standard deviation of all the projections in a given year, where var m,s is the variance across all models and scenarios.
We use the method outlined in Hawkins and Sutton (2009) to quantify the fractional contribution of the three sources of projection uncertainty. The total variance, U(t) 2 is equal to the sum of the variance due to internal variability, U 2 V , the variance due to emission scenario, U S (t) 2 , and the variance due to model structure, U M (t) 2 , The variance due to internal climate variability is calculated as where var t is the temporal variance, N em is the total number of ensemble members in the CMIP5 suite (N em ¼ 40), and W m is a normalized weight for each model (see next section). This method for calculating internal variance assumes that variance is appropriately captured by the residuals from a polynomial fit to the projections, as described in Hawkins and Sutton (2009). Other studies have estimated internal variability in CO 2 uptake based on model ensembles (Lombardozzi et al 2014. Our method further assumes that the internal variance does not change with time and is unaffected by emission scenario. Finally, by averaging the internal variance across all CMIP5 models, our method masks the subtle differences in internal variance across a range of model structures that may be prevalent in multi-century preindustrial control simulations of the same models (Resplandy et al 2015).
The variance due to emission scenario is the crossscenario variance (var s ) of the weighted multi-model mean forced signal: where N m is the number of ensemble members that simulated a given emission scenario (N RCP2:6 The variance due to model structure is the multiscenario mean of the weighted intermodel variance in the forced signal: where N s is the number of scenarios in the CMIP5 suite (N s ¼ 4), and var W m is the weighted variance. The fractional variance is then 2 , and U M ðtÞ 2 U ðtÞ 2 for internal variability, emission scenario, and model structure, respectively. In this simple statistical framework, the sum of three uncertainty terms equals the total uncertainty, since we assume that the modelscenario interaction term is negligible, as in Hawkins and Sutton (2009).

Weighting schemes
We introduce weighting schemes to our analysis of variance to evaluate the sensitivity of terrestrial carbon Environ. Res. Lett. 12 (2017) 044020 projection uncertainty to model structure. Both unweighted and weighted uncertainty estimates of terrestrial carbon uptake are presented. We develop 7 weighting schemes that weight the models by their ability to simulate observed changes in terrestrial carbon uptake over 1959-2005. Models that successfully reproduce past changes are given higher weights in all schemes. Schemes are based on the probability density of normally-distributed, linear trends in cumulative carbon uptake over 1959-2005, where m and s represent the mean and standard deviation, respectively, of the modeled trends in cumulative terrestrial carbon uptake over 1959-2005. Each weighting scheme (WS) uses different values of m and s, derived from trends from a subset of models in the corresponding model tier(s) (table S1). For example, WS 1 uses m and s values derived from Tier 1 models, WS 2 uses m and s values derived from Tier 1 and 2 models, and so on. WS 7 uses m and s values derived from all of the models (Tiers 1 through 7). Models were categorized into tiers according to the absolute value of the difference between the modeled trend and the observationally-based trend over this period (table S1). Each model is given a weight, w m , according to the resulting probability density ( f (x)) of the modeled trends in cumulative land uptake for a given weighting scheme. These weights are expressed as normalized weights, W m , in our analysis of variance,

Results and discussion
The time series of globally-integrated, cumulative ocean carbon uptake from 2006 to 2100 reveals a growing ocean carbon sink for all CMIP5 ensemble members, with higher uptake corresponding to higher CO 2 concentration pathways ( figure 1(a) Environ. Res. Lett. 12 (2017) 044020 to our focus on global and cumulative, rather than regional mean or annual-mean CO 2 fluxes.
The evolution of globally-integrated, cumulative terrestrial carbon uptake over the next century ( figure 1(b)) depends on the evolution of land use change emissions and terrestrial ecosystem production. While some CMIP5 ensemble members show an accumulation of carbon in the terrestrial biosphere, others show a loss of carbon over 2006-2100 ( figure 1(b)). The uncertainty in these projections increases linearly from 0 Pg C in 2006 to 163 Pg C in 2100 ( figure 1(d )), reflecting greater divergence than the ocean projections for all prediction lead times. Model structure accounts for ∼80% of the projection uncertainty in terrestrial carbon uptake for all prediction lead times, with emission scenario and internal variability playing much smaller roles (figure 1(f )), consistent with a previous study of decadal-mean terrestrial carbon fluxes (Hewitt et al 2016). Thus, in the terrestrial biosphere, we have observed higher overall projection uncertainty and an important role for model structural uncertainty that were not evident in our analysis of global ocean carbon accumulation.
What steps can we take to reduce uncertainty in projections of ocean and land carbon accumulation? In the global ocean, our analysis points to a clear role for emission scenario uncertainty in the latter half of the century, when overall uncertainty is high. Thus, a significant reduction in uncertainty here is principally attainable by narrowing the uncertainty in future emission trajectories. While Dunne (2016) argues that recently proposed climate stabilization targets will lead to natural reductions in scenario uncertainty, this source of uncertainty is nevertheless dependent on future societal behavior and technological advancements and largely outside the realm of physical science. In the terrestrial biosphere, however, uncertainty is dominated by model structure, a source of uncertainty that is potentially reducible through advancements in ecological theory and modeling. Hoffman et al (2014) suggest that significant decreases in model structural uncertainty may be achievable through (1) closer coordination among modeling centers, and (2) systematic evaluation of models through comparison with observations. In their study, they use historical observations of atmospheric CO 2 from Mauna Loa to constrain CMIP5 model predictions of atmospheric CO 2 , yielding a considerably narrowed distribution of potential atmospheric CO 2 concentrations by the end of the century. Here, we attempt a similar exercise for terrestrial biosphere carbon accumulation.
Each year, the Global Carbon Project (GCP) publishes a plausible history of land use change and terrestrial carbon fluxes that is based on observations of emissions, atmospheric growth rate, and ocean uptake and can be used to address model fidelity  2(a)), with a linear trend of 0.49 Pg C yr À1 (r 2 ¼ 0.91, figure 2(b)). This observational estimate is bounded by the cumulative uptake estimated from the CMIP5 models over the same period ( figure 2(a)), such that the GCP trend falls close to the middle of a nearnormal distribution of linear trends in the CMIP5 models ( figure 2(b)). These historical trends form the basis of schemes that we devise to weight the CMIP5 models based on their ability to represent past observed changes (see Methods, figure 2(c)). The most extreme weighting scheme (WS 1) gives all the weight to a few of the best-performing models, while the least extreme weighting scheme (WS 7) gives near-equal weight to all the models ( figure 2(c)). The most extreme weighting schemes are thus likely to yield the largest reductions in model structural uncertainty and overall projection uncertainty in terrestrial carbon uptake.
Terrestrial carbon uptake projection uncertainty decreases incrementally from 163 Pg C in 2100 with no weighting scheme to 125 Pg C in 2100 when we apply the most extreme weighting scheme (WS 1) to our analysis of variance ( figure 3(a)), suggesting that a single observational constraint of past model performance can be used to inform future projections and lower overall projection uncertainty. We note, however, that terrestrial carbon uptake projection uncertainty under the most extreme model weighting scheme (WS 1) is >25% larger than that of the global ocean at the end of the century. The primary source of projection uncertainty on land is model structure for all but the most extreme weighting scheme (WS 1), where emission scenario becomes the primary source of uncertainty after ∼2040 ( figure 3(b)). Further statistical reductions in model structural uncertainty beyond the weighting schemes presented here require unrealistic measures, such as excluding all but one model from the analysis of variance (not shown).
Our analysis indicates that a meaningful reduction in model structural and overall uncertainty in terrestrial carbon uptake projections is obtainable only with a fairly impractical observational constraint. For example, in order to allow emission scenario to be the dominant end-of-century source of uncertainty, we had to design a model weighting scheme (WS 1) that excluded all but 2 models (HadGEM2-CC and HadGEM2-ES) from the analysis of variance. In fact, however, these two models differ only in the inclusion of tropospheric chemistry in HadGEM2-ES, which has little impact on the simulated carbon cycle (Martin et al 2011). Inclusion of just two additional models (WS 2) causes model structural uncertainty to exceed scenario uncertainty for all prediction lead times ( figure 3(c)).
Observational constraints have previously been advocated as a means to reduce model uncertainty in the terrestrial biosphere (Randerson et al 2009, Luo et al 2012, Cox et al 2013, Schimel et al 2015, however, one needs to consider the uncertainty in the observational constraint itself. Here, we use the GCP estimate of Environ. Res. Lett. 12 (2017) 044020 net land carbon flux with a reported uncertainty of 0.9 Pg C yr À1 in recent decades (Ciais andSabine 2013, Schimel et al 2015). Only one CMIP5 model (IPSL-CM5B-LR; table S1) falls outside the range of uncertainty in the observed trend in cumulative land uptake over the 1959-2005 period (0.49 ± 0.9 Pg C yr À1 ). Constraining models to an uncertain observation may not increase their predictive skill.
Improvements to model structure require deep knowledge of the real-world processes being represented by the models. Net land carbon flux is a complex quantity in an earth system model, as it depends on many different parameters that together determine land use change, drought, fire, and CO 2 fertilization fluxes and their sensitivities to changes in climate. Thus, even models that strongly agree with historical observations can exhibit divergent projections in the future (Knutti and Hegerl 2008).
Current research suggests that the large structural uncertainty in projections of terrestrial carbon fluxes from Earth System Models is driven both by differences in simulated climate, and by the challenge of representing life in these models, with the rich diversity of lifeforms and complexity of ecological systems. Even when forced with the same climate, the current generation of models give divergent depictions of the carbon cycle (Sitch et al 2015). Previous studies have demonstrated large spread within a single model structure arising from parameter uncertainty (Booth et al 2012). Yet there is the added difficulty of Environ. Res. Lett. 12 (2017) 044020 mathematically representing ecological processes. Even when models utilize the same theoretical underpinnings of photosynthesis, their varying numerical implementation of that theory leads to divergent simulations (Rogers et al 2017). Other studies have focused on structural uncertainty within a single model and have identified photosynthetic and respiratory temperature acclimation, biological nitrogen fixation, and microbially-based soil organic matter dynamics as key uncertainties in carbon cycle projections over the twenty-first century (Lombardozzi et al 2015, Wieder et al 2015a, 2015b. Nitrogen limitation of terrestrial productivity is a critical determinant of carbon cycle projections, but our understanding of how to model biogeochemical processes is poor (Zaehle et al 2014, Medlyn et al 2015, as is our ability to model land use emissions (Lawrence et al 2016) and wildfire (Hantson et al 2016). Ecological complexity necessitates a multitude of modeling approaches to capture the range of possible outcomes. The focus on reducing multimodel spread does not necessarily make carbon cycle projections more reliable and may, in fact, limit scientific progress. Environ. Res. Lett. 12 (2017) 044020