Global distribution of soil organic carbon – Part 2 : Certainty of changes related to land use and climate

Global biosphere models vary greatly in their projections of future changes of global soil organic carbon (SOC) stocks and aggregated global SOC masses in response to climate change. We estimated the certainty (likelihood) and quantity of increases and decreases on a half-degree grid. We assessed the effect of changes in controlling factors, including net primary productivity (NPP), litter quality, soil acidity, water saturation, depth of permafrost, land use, temperature, and aridity associated with probabilities (Bayesian network) on an embedded, temporally discrete, three-pool decomposition model. In principle, controlling factors were discretized into classes, where each class was associated with a probability and linked to an output variable. This creates a network of links that are ultimately linked to a set of equations for carbon (C) input and output to and from soil C pools. The probability-weighted results show that, globally, climate effects on NPP had the strongest impact on SOC stocks and the certainty of change after 75 years. Actual land use had the greatest effect locally because the assumed certainty of land use change per unit area was small. The probability-weighted contribution of climate to decomposition was greatest in the humid tropics because of greater absolute effects on decomposition fractions at higher temperatures. In contrast, climate effects on decomposition fractions were small in cold regions. Differences in decomposition rates between contemporary and future climate were greatest in arid subtropical regions because of projected strong increases in precipitation. Warming in boreal and arctic regions increased NPP, balancing or outweighing potential losses from thawing of permafrost. Across contrasting NPP scenarios, tropical mountain forests were identified as hotspots of future highly certain C losses. Global soil C mass will increase by 1 % with a certainty of 75 % if NPP increases due to carbon dioxide fertilization. At a certainty level of 75 %, soil C mass will not change if CO2-induced increase of NPP is limited by nutrients.


Introduction
Soil organic carbon (SOC) represents about three-quarters to four-fifths of the terrestrial organic carbon (C) mass (Prentice et al., 2001).The mean turnover rate of SOC is slower than that of any other terrestrial organic pool (Reeburgh, 1997).Due to its size, small relative changes in the SOC mass can have large effects on atmospheric CO 2 concentration and hence on climate change.
Global SOC mass in five general circulation models (GCM) was projected to change by between −46 and +51 Pg (Schaphoff et al., 2006) by the end of the century.Eleven earth system models even showed a range of projected changes between −72 and +253 Pg for a high-CO 2 scenario (Todd-Brown et al., 2014).Projections also differ in where changes occur (Sitch et al., 2008).The large variation in expected future changes is due to the balance of, on one hand, different expected increases of C input from net primary productivity (NPP) by CO 2 fertilization and higher temperatures and, on the other hand, faster decomposition accelerated by higher temperatures (Davidson and Janssens, 2006;Smith et Published by Copernicus Publications on behalf of the European Geosciences Union. al., 2008).The point of balance may vary over the course of time (Jones et al., 2005).Furthermore, although NPP might increase in the future because of increasing concentrations of CO 2 in the atmosphere (CO 2 fertilization), productivity may still be limited by the availability of nitrogen or other resources (Gedalof and Berg, 2010;Norby et al., 2010;Todd-Brown et al., 2014).
With this wide range of projected changes in the global mass of SOC, one may wonder how likely increases or decreases of SOC stocks (mass of organic C per volume of soil) are in response to potential changes in C input and climate across the world.One way to address the certainty of projections is to obtain a frequency distribution of the ensemble output of several models as has been done for changes in climate (e.g.Power et al., 2011).This approach, however, cannot formally address the uncertainty in the parameters of the models.Rather, this requires consideration of the frequency distribution of the values of many potentially controlling factors of the organic C cycle.
The effect of the frequency distribution of parameters on SOC stocks, the global SOC mass, and SOC changes was assessed for one model (Hararuk et al., 2014) using a Markov Chain-Monte Carlo approach for model calibration.The distribution parameters of the model caused SOC losses ranging between 15 and 100 Pg for a scenario implying high greenhouse-gas emissions and mean global temperature increases ("RCP8.5").This study focused on the approach and a global perspective and left out other important impacts that affect SOC stocks and masses regionally and globally, e.g.fire, insect outbreaks, erosion, landslides, windthrow, and flooding.
For a richer picture of the change of global SOC mass, wetland soils, including peatlands, containing at least 6-12 % of the global SOC mass in the upper 1 m (depending on the definition and estimated area of wetland, Köchy et al., 2015), and permafrost regions containing about 40 % of the global SOC mass in the upper 1 m (Köchy et al., 2015) must be considered as well.Furthermore, SOC stocks are not only affected by climate change but probably even more so by change in land use (Brovkin et al., 2013).This is especially true for organic soils because they contain > 15 % of the global SOC mass in the top 1 m (Köchy et al., 2015).SOC losses from organic soils in Scotland, for example, are expected to be ca.3.5 times greater than losses from C-poor soils (Smith et al., 2010).In C-rich wet or waterlogged soils, decomposition of organic matter is slow because of a lack of oxygen (Armentano and Mengeo, 1985;Mitra et al., 2005;Smith et al., 2008).Draining of wetlands exposes C to oxygen.Outside wetlands, conversion of natural forests to grassland or cropland causes drastic losses in temperate (Poeplau et al., 2011) and tropical regions (Holmes et al., 2006;Don et al., 2011).Organic matter in the soil can also be physically protected from microbial decomposition by adsorption to soil particles (Six et al., 2002;Davidson and Janssens, 2006) or permafrost.
The relative impacts of climate change, land use change, and thawing of permafrost on SOC stocks of mineral and organic soils have been studied with great detail only at small or regional scales (e.g.Grosse et al., 2011).There is a lack of a comprehensive global assessment of the certainty of soil C changes (Vesterdal and Leifeld, 2007).In the present study we quantify the uncertainty of changes of present-day SOC stocks (ca.2010) due to projected changes in climate and land use by aggregating the uncertainty or variability in controlling variables and their effects on SOC stocks.In addition, we identify where soils are most likely to be vulnerable at the global scale and with relevance to the global C cycle.

Methods
We assess the effects of climate (temperature, aridity), soil (acidity, permafrost, aerobicity, C adsorption), vegetation (vegetation type, litter quality), and land use (via NPP, harvest factor, and litter quality) on SOC stocks at a spatial resolution of pixels with 0.5 • latitude by 0.5 • longitude.We acknowledge that variation in soil, vegetation, environmental, land use, and other factors controlling SOC decomposition exists within a pixel.This variation is partly included in our analysis but we cannot quantify its contribution to overall uncertainty.
To quantify the certainty of SOC changes over a period of 75 years (2010-2085), we apply a quasi-steady-state threepool model of SOC decomposition to historic SOC stocks within an environmental framework (Fig. 1) and compare differences between future reference and target conditions.Historic conditions correspond to SOC stocks of ca.1950-2000 of the top 1 m under current (ca. 1980-2010) climate and land use.Future conditions are characterized by projected climate and land use of 2075-2100.Reference conditions imply no change in environment, whereas target conditions imply changes in climate, vegetation, and land use.For comparisons of reference and target conditions, we keep local (within pixel) settings of soil pH, CEC, and constraints of O 2 availability (unless caused by a change of land use to or from wetland) constant.
SOC stocks in reference and target conditions are prescribed by decomposition under constant ranges of monthly temperature, monthly precipitation, litter input, and land use for 75 years starting with the historic SOC stock (we account for the expected gradual change in climate and land use in an additional step).We compare reference with target SOC stocks after 75 years instead of steady states with t → ∞ in order to compare the same points of time across all soil conditions of the world in the spatial analysis and for compatibility with the point of time of projected NPP and climate values.Examination of our decomposition formulas in a spreadsheet indicated that a steady state is reached within 75 years in most mineral soils and non-extreme environments.In other words, by considering a defined, limited time period we com- pare two possible outcomes (reference vs. target conditions) for the same year.This facilitates comparisons with projections of C stocks for the end of the century by other authors.
In the following sections we describe first the core decomposition model and then how environmental factors affect the values of the decomposition parameters.

Characterization of the core decomposition model
We consider three pools of SOC -fast (C fast ), slow (C slow ), and inaccessible -that differ in their maximum annual decomposition fractions under optimal conditions.At the beginning of each year, above-and below-ground coarse and fine litter (equal to NPP in the long term) is added to C fast .Removal of C from the ecosystem by disturbances (e.g.fire) or harvest is taken into account in relation to land use and is described below.The C fast pool is reduced by the annual decomposition fraction (F f ) and a fraction moving from the fast to the slow pool (to slow ).Litter quality controls to slow .The maximum decomposition fraction is constrained by temperature, soil humidity, soil acidity, and oxygen availability.The constraints are jointly expressed as a fraction modifying factor (fmf) ranging between 0 and 1 (Fig. 2).The maximum decomposition fraction of the slow pool (F s ) is constrained in the same way as F f .The inaccessible pool is the fraction of SOC that is frozen, submerged in water, or adsorbed to soil matter and whose decomposition fraction we assume to be negligible under extant conditions within the time perspective of our study.The fraction of C in the inaccessible pool, however, may differ between reference and target conditions when the comparison reflects changes in water logging or frozen soil.We use a minimum function that considers the constraint of decomposition by a high water table after the thawing of permafrost in wetlands.
The total amount of C in the fast and slow pool after 1 year are C fast,t+1 = (C fast,t , with C fast,t=0 = 0 and C slow,t=0 = C 0 • af or the accessible fraction (af) of the initial total C stock (C 0 ).NPP is supplied by external models and described below in Sect.2.3.After several decades, the sizes of the fast and slow pools depend mostly on the amount of annually added C, the decomposition fractions, and the distribution of matter between the fast and slow pools but little on the initial amount of accessible C for not-too-small values of fmf (fmf < 0.1, see Supplement 3 for a summary equation).Limited substrate availability could reduce decomposition rates (Davidson and Janssens, 2006;Kirschbaum, 2006) but contributed very little to the prediction of existing SOC stocks in an earlier version of the model.Therefore, substrate availability was not included in the final version.
The values of the parameters and variables of the decomposition model are controlled by variables in the environmental framework characterizing the physical and biotic environment in a particular location: soil, vegetation, climate, land use.

The environmental framework
The environmental variables and their causal relationships are described by a probability network (Spiegelhalter et al., 1993) using Netica (version 5,Norsys,Vancouver,Canada).
www.soil-journal.net/1/367/2015/SOIL, 1, 367-380, 2015 Probability networks, also known as Bayesian networks, associate classes (e.g.levels) of each variable, with probability distributions of their occurrence contingent on the occurrence of classes (levels) within other variables.These probabilities can be interpreted as certainties of potential outcomes.Joint probabilities of classes are calculated according to the laws of probabilities or, in the case of continuous variables, sampled by Monte Carlo techniques.Probability networks allow the inclusion of the uncertainty or variability of variables and expert knowledge.The networks (one for each of two contrasting NPP submodels in two NPP scenarios) with all probability tables and class borders are supplied in Supplements 1 and 2. The widths of classes can be thought of as encompassing sub-pixel variation of soil, vegetation, and environmental variables, but we cannot quantify the contribution of sub-pixel variation to overall certainty of SOC changes.

NPP scenarios in the environmental framework
The strong effect of NPP via litter input on SOC stocks in models is well established and founded on theory (Todd-Brown et al., 2014).However, the long-term net effect of climate change on NPP and, consequently, litter input is still unclear.NPP might increase because of increasing concentrations of CO 2 in the atmosphere (CO 2 fertilization) but productivity may be limited by the availability of nitrogen or other resources (Gedalof and Berg, 2010;Norby et al., 2010;Todd-Brown et al., 2014).To present the range of the effects of nutrient limitation and CO 2 fertilization we used two contrasting scenarios of future NPP together with climate conforming to the A1B emissions scenario (IPCC, 2000).The first scenario, "limited NPP", represents a change in productivity caused by temperature and precipitation alone, i.e. without CO 2 fertilization.This limitation of potential increases could be similar to the net effect of CO 2 fertilization and nutrient-constrained growth.This NPP is based on the empirical NCEAS model (Del Grosso et al., 2008, an extension of the Miami or Lieth model), a function of mean annual temperature, mean annual precipitation, and vegetation type (Supplement 3: Data processing and sources).The second scenario, "enhanced NPP", represents an increase of productivity due to CO 2 fertilization in addition to changes in temperature and precipitation without limitation of the additional growth by nutrients.This NPP is derived from LPJ, a process-based dynamic global vegetation model (DGVM) (Sitch et al., 2003;Gerten et al., 2004).In comparisons between the NPP scenarios we keep functional relations (e.g.decomposition fractions, fraction modifying factors, decomposability) and plant type composition within each vegetation type constant.

Environmental framework structure and parameterization
The amount of C added to the fast pool is set equal to the NPP in natural ecosystems and consists of leaf litter, fine and coarse woody debris, and fine and coarse dead roots.C removal with harvested products and higher NPP input by agricultural fertilization must be accounted for in land use effects.We calculated the harvest factor, the mean ratio of NPP after harvest including agricultural fertilization (NPP t ) to NPP of the potential zonal vegetation (NPP 0 ), using a global database (Haberl et al., 2007) for each combination of 13 climatic vegetation zones and six land use classes (zonal, built up, herbaceous crops, pasture, woody crops/plantations, wetlands; Supplement 3, Table S4.2 in Supplement 4) for present environmental conditions."Wetlands", a special vegetation type and land cover, is included in this list for convenience.The harvest factor allowed us to use existing models of NPP 0 and apply them to future conditions.We note that using this procedure glosses over regional differences within vegetation zones.
The fmf aggregates the effects of temperature, soil moisture, oxygen availability, and soil reaction by multiplication similar to the rate-modifying factors in other decomposition models (e.g.Roth-C, Coleman and Jenkins, 1999).The probability distribution of fmf temperature (Fig. 2a) is a discretization of the equation exp(−2.5 + 0.07• temperature) on laboratory incubation data from several sources (Fig. 2 in Paul et al., 2002).We used laboratory data because we were interested in the effect of temperature, isolated from other variables, on the maximum decomposition fraction under optimal conditions.Our discretization of the temperature effect (Fig. 2a) encompasses many other empirical temperature functions that are used in established soil C models like Roth-C, APSIM, or Century (Paul, 2001).In situ microbial communities might respond to increased constant tempera-tures with acclimation (Allison et al., 2010).Our approach considers monthly variation of temperatures and aridity, and it is unclear from the current literature how strong acclimation is relative to this variability and how much current latitudinal patterns are caused by climate.We assume here that the class widths used in the parameterization of fmf (Fig. 2a) and distinction of monthly temperatures encompass sub-monthly effects of acclimation.
The classification of the moisture effect (Fig. 2b) is associated with Walter and Lieth's (1967) climatic aridity classes.These are arid (MMP/2 < MMT), dry (MMP/3 < MMT), moist (MMP > 100), and mesic (the remainder), where MMT is mean monthly air temperature ( • C) and MMP is mean monthly precipitation (mm).In addition, "wet" is used for wetland soils.The shape of the probability distribution associated with the aridity classes corresponds to the probability distribution of the 0.75-1.0quantile range of laboratory decomposition studies (Paul et al., 2002) and of subsamples of the moisture functions in the decomposition models AP-SIM (Probert et al., 1998), ED-RAMS (Ise et al., 2008), andECOSSE (SEERAD, 2007).We selected these models because they included relationships for water-saturated soils.To include seasonal effects, fmfs of temperature and aridity were calculated by month, multiplied with each other for each month and averaged (geometric mean) per quarter and then per year.Averaging was necessary to reduce the complexity of the probability network to a level that could be calculated with Netica.Oxygen availability was taken from the Harmonized World Soil Database (HWSD) supplementary data and is based on soil drainage and takes into account soil type, soil texture, soil phases, and topographic position (Fischer et al., 2008).Changes in soil moisture are reflected in the decomposition fraction and the accessible pool fraction.The probability distribution of oxygen availability was set so that completely anaerobic conditions reduce the maximum decomposition fraction to one-seventh (Freeman et al., 2001).
Following the discussion of ECOSSE (SEERAD, 2007), we specified the probability distribution of the soil acidity effect so that decomposition fraction increases from a medium level within the aluminum-buffer pH range (acidic) to optimal within the carbonate-buffer pH range (neutral) and decreases to low as soils become alkaline (pH > 8.5; Fig. 2d).
Organic matter in the inaccessible C pool decomposes extremely slowly or not at all due to lack of oxygen in waterlogged soils, permafrost, or adsorption to soil particles (Six et al., 2002).It can have a history of millennia or, in the case of permafrost, can date back to the last interglacial period.We defined this quasi-constant pool as the inaccessible fraction of SOC in the top 1 m.It was calculated as the maximum of the fraction of the volume of soil that is permanently frozen, plus the fraction of the remaining C protected by soil particle adsorption.We set the fraction of waterlogged soil to 0.9-1 with 80 % probability for wetlands and to 0-0.1 with > 94 % probability for other land cover types (Ta-ble S4.6 in Supplement 4).The fraction of permafrost is calculated as a function of the number of days in a year with a mean daily temperature > 0 • C (degree days, DD0), 1-0.02 • √ DD0 (Anisimov et al., 2002).The fraction of SOC protected by particle adsorption is calculated for mineral soils, roughly following the SOCRATES model (Grace et al., 2006), as (0.14• CEC + 40)/100 if CEC < 100 mmol kg −1 or else (0.04• CEC + 50)/100 with a maximum of 1.The CEC of soils with > 20 % organic C content was set to 0.
The decomposability of litter that controls the fractions of the slow and fast C pools depends on plant type and indirectly on vegetation type and land use (Tables S4.4 and S4.5 in Supplement 4).The classification of plant types and their association with decomposability is our interpretation of metaanalyses of leaf litter (Cornwell et al., 2008) and wood decomposition (Weedon et al., 2009).We estimated, based on studies in many types of forests (Rodin and Bazilevich, 1967;Laiho and Prescott, 2004;Rice et al., 2004;Steinaker and Wilson, 2005), that above-and below-ground fine (leaf and fine roots respectively) and above-and below-ground coarse (coarse woody debris and coarse roots respectively) litter contribute, on average, equal proportions to litter input entering forest soils.About 10 % of the C allocated to roots may be lost as easily decomposable exudates (van Hees et al., 2005).This is implicitly reflected in the fraction of total litter attributed to the high and very high decomposability class and its fraction (not) going to the slow C pool (Table S4.5 in Supplement 4).The proportions of plant types within vegetation types (Table S4.2 in Supplement 4) are based on those reported in Sterling and Ducharne (2008) proportions of vegetation classes in the Global Land Cover Characterization (Loveland et al., 2000) global ecosystem legend (plantations), USGS legend (tundra, wetlands), and personal experience for non-vascular plants.
Zonal vegetation of target conditions are linked to the reference zonal vegetation by transition probabilities for each vegetation type (Gonzalez et al., 2010).Probability distributions of target land use are contingent on reference land use, target temperature, and target aridity.We compared crop and pasture land use maps for 2000 and 2075 (de Noblet-Ducoudré and Peterschmitt, 2008) to calculate probabilities of changes among land use classes for each crossclassification of reference temperature and aridity.We modified these probabilities and assigned probabilities for land use changes among other classes according to a set of rules based on our experience (Table S4.3 in Supplement 4).These land use changes are thus hypothetical and constrained mostly by climate.Wetlands, however, are set to remain wetlands with a very high probability (≥ 97 %).

Calibration
The decomposition model was calibrated at the global scale (using half of the data points for training and the other half for validation) by systematically varying the decomposition www.soil-journal.net/1/367/2015/SOIL, 1, 367-380, 2015 fractions F f and F s in steps of 0.05.Smaller steps produced hardly perceivable changes to the output.We selected the combination of decomposition fractions so that the HWSD-SOC stock classes across all pixels were predicted most often by the most probable reference SOC stock class.This resulted in maximum decomposition fractions F f = 0.75 -corresponding to a decomposition rate k = −log e (1 − F f ) = 1.4 -and F s = 0.35 (k = 0.4).Using these decomposition fractions and NPP 0 (Haberl et al., 2007), the model correctly predicted 78 % of all HWSD-C stock classes (Fig. S3.6 in Supplement 3).Further details are provided in Supplement 3.

Presentation of results
All variables in the model were associated with a frequency distribution that affected the certainty of changes in SOC stocks.This paper focuses on the certainty of changes in SOC stocks.Therefore, we report and discuss the probabilityweighted means of their distributions.SOC stocks presented in this paper are averages across SOC stock classes (class borders: 0, 2, 5, 10, 15, 20, 30, 60, 80, 100 kg m −2 ) weighted by class probability.Changes in SOC stocks are expressed as half the mean of the frequency distribution of the differences between the reference and target SOC stock distributions throughout this paper.Using the half puts the numerical value of changes more in line with the gradual change between reference and target conditions as projected by GCM-DGVM combinations.We express the certainty of gains or losses as the certainty of the difference between reference and target SOC stock being > 0. We call changes associated with a certainty P > 0.67 "fairly certain" and changes with P > 0.75 "highly certain".In this context, 50 % certainty would imply that losses and gains are equally certain.These certainties are not expressions of statistical significance but represent the likelihood that a certain outcome (or distribution of outcomes) might be or become true given the assumptions and width of classes of the contributing variables.All reported SOC stocks from our simulations are standardized to the SOC stocks of the HWSD as processed by Köchy et al. (2015) (Supplement 3, Fig. S3.3).
For examination of specific land use changes we examined individual pixels (Table S5.1, Fig. S5.1 in Supplement 5) including those used by Jones et al. (2005) and Schaphoff et al. (2006).

Agreement of steady state with HWSD stocks
We compared the SOC stocks for reference conditions to the SOC stocks calculated from HWSD.We did not use observed time series of SOC for validation because our study is aimed at assessing (un)certainties of effects and assumes constant conditions until a defined endpoint.Our SOC stocks are the average of the class mean weighted by class proba- in 89 % of all cases.The greatest sensitivity of the total reference SOC stock was, in decreasing order, to NPP t and NPP 0 , followed by vegetation zone, HWSD-SOC stock, and fmf (Supplement 3).

Effect of climate change on environmental factors
Climate change until the end of the century (i.e.changes in monthly temperatures and aridity) is reflected in changes of the fmf and the probabilities of shifts in vegetation zone and land use with secondary effects via plant type composition, decomposability, and proportion of C input going directly to the slow pool.In general, the differences in fmf between reference and target conditions were small (Fig. 3b).Strong increases were confined mostly to relatively arid tropical regions where increases in winter precipitation were strong.
Increases in the depth of the active layer of permafrost were < 20 cm in most locations.Increased thawing depth of 20-30 cm was projected only for the central Asian mountain ranges.Although permafrost soils were projected to thaw deeper, not all thawing regions were exposed to decomposition of more organic matter because C, at least in mineral soils, was stabilized by CEC.
Decomposability of litter varied between around 0.3 and 0.6 relative units, but changes in decomposability between reference and target conditions were comparatively small (Supplement 4, Fig. S4.1).In the limited NPP scenario, the low increase of NPP had little effect on the plant type distribution so that changes of decomposability were within ±0.05 units.In the enhanced NPP scenario, plant type distributions shifted more strongly so that changes of decomposability ranged between −0.2 and +0.1 units.The greatest decreases (causing higher SOC stocks) occurred in tundra where the proportion of woody plants was projected to increase by external models.The greatest increases of decomposability occurred (1) in the southern boreal forest where deciduous trees replaced evergreen trees and (2) in those parts of the tropical forests with a high probability that pristine forest would be converted to cropland or perennial plantations.

Effect of climate and land use change on SOC stocks under contrasting NPP scenarios
In the limited NPP scenario, NPP t changed little by the end of the century (Fig. 4a).Most notable were extensive increases in the boreal forest and scattered losses in tropical regions.Under limited NPP conditions, global SOC mass might decrease (P > 0.5) by a net mean of 21 Pg due to the effect of projected climate and land use change (Fig. 5a).At midlatitudes (35-65 • N), gains were greater than losses, whereas losses were greater than gains between 30 • N and 30 • S. In most locations the certainty of changes was low (P < 0.67) (Fig. 6a).If we consider only SOC changes at a certainty level ≥ 0.75, the global change of SOC mass is 0 Pg (Fig. 5a), i.e. gains and losses were balanced.Patches with highly certain losses occurred in high elevations of the northern Andes, New Guinea, and eastern central Africa.Locations with highly certain gains were high elevations in the southwestern US, the highlands of southern Africa, and the central Asian mountain ranges.
In the enhanced NPP scenario, NPP t increased in almost all regions.The increases were greatest in the humid tropics and higher altitudes (Fig. 4b).The global SOC mass (0-1 m) might increase (P > 0.5) by a net mean of 55 Pg due to climate change, land use change, and CO 2 fertilization (Fig. 5b).If we consider only changes at a certainty level ≥ 0.75 (Fig. 6b), the global net mean gain of SOC is 11 Pg (Fig. 5b).Regions with highly certain losses and means of 2-10 kg m −2 were few and small, representing high elevations in the northern Andes, eastern Central Africa, and New Guinea.Highly certain gains with means of 2-10 kg m −2 occurred in arid higher elevations of southeastern North America, southern Africa, central Asia, and scattered in other highlands.
We located "hotspots" of SOC vulnerability where losses occurred across both NPP scenarios.There was, however, little overlap across both NPP scenarios.In both NPP scenarios, SOC stocks are fairly certain to increase in northern Labrador (Canada) and the Chukotsky peninsula (most western tip of Siberia), whereas SOC stocks are fairly certain to decrease in parts of the mountain ranges of the northern Andes, in the Ethiopian and eastern-central African highlands, and in the mountain range of New Guinea (Fig. 5c).There was no discernable pattern of coincidence with current land use.In the tropics, 23 % of the 169 pixels with losses were located below 500 m altitude, 19 % within 500-1000 m alti- tude, 58 % within 1000-2500 m altitude, and 84 % on ridged terrain (Fig. 5c).Ten of the 23 pixels with gains in the tropics (43 %) were located at altitudes > 2500 m.Pixels outside the tropics showed no discernable patterns with topography.

Effect of climate change with enhanced NPP by vegetation zone and by land use
In each vegetation zone the median net change of SOC stocks across all pixels was positive (Fig. S4.3 in Supplement 4).
For each vegetation zone as a whole, the direction of net change is uncertain due to the great environmental heterogeneity within each zone.Most changes between land use types across vegetation zones were similarly uncertain because of the great heterogeneity involved.Conversion of cropland to wetland, however, has a high certainty of SOC stock gains with a mean of 3 kg C m −2 .Similarly, conversion of wetlands to crops is highly certain to incur losses of SOC  with a mean of 5 kg C m −2 across vegetation zones.Conversion of wetlands to pastures or plantations are also fairly certain to incur losses.
In order to reduce the heterogeneity in the assessment, we also considered land use changes within vegetation zones.The same patterns as reported above occurred.In addition, conversions from zonal vegetation to croplands were found to fairly certainly incur losses in temperate vegetation zones, tropical forests, and woodlands.In contrast, conversion from zonal vegetation to cropland in hot deserts is fairly certain to incur gains with a mean of 4 kg C m −2 .This effect arises because the imposed land use change assumes strong increases in SOC input in a region with natural low productivity.Similarly, reconversion of cropland to zonal vegetation is fairly certain to incur gains in C stocks in temperate forests and tropical deciduous and evergreen forests with means ranging between 2 and 3 kg C m −2 .Obviously, the opposite was fairly certain for hot deserts.
By inspection of individual pixels, we reduced the spatial heterogeneity as far as possible in our approach (Table S5.2 in Supplement 5).A general picture emerged.Changes with fairly or high certainty were associated with strong increases in NPP t .
Almost all non-wetland tundra is currently not subject to land use.Under target conditions, 15 % of the tundra is projected to be used as boreal forest plantations.Changes in fmf were small and increase in NPP 0 balanced or overcompensated losses of thawed fossil SOC.Increases in NPP 0 and consequently SOC stocks were projected to be higher and gains fairly certain in alpine tundras of central Asia.In many arctic locations, SOC stocks would increase fairly certainly if the permafrost soil turned into wetlands.
Boreal regions would be increasingly used as timber plantations and probably as arable land.Timber harvesting of formerly pristine boreal forests is not likely to decrease SOC stocks.Tree removal for cropping, however, reduces SOC stocks fairly certainly in some pixels with mean losses of up to 5 kg m −2 compared to reference conditions.Creation of wetlands in boreal forests, e.g. by thermokarst processes, would fairly certainly increase SOC stocks with a mean gain of 3 kg m −2 .
In some pixels of temperate grasslands (steppes) in eastern Asia, the projected increase in NPP was high and increase in SOC stock in grazed steppes (up to 6 kg m −2 ) was fairly certain.Where temperate grasslands are currently used as arable land and were used for pasture or left to return to zonal vegetation, SOC stocks would increase with fair or higher certainty.Conversion of cultivated steppe (perhaps former wet depressions, e.g.prairie potholes) to wetlands would fairly certainly increase SOC stocks by up to 5-6 kg m −2 .
SOC stocks in pixels with temperate forests also echoed the general picture that higher inputs were associated with fairly or highly certain gains, so that conversions to cropland or pastures often incurred losses (e.g. 5 kg C m −2 ) and land use changes to woody vegetation incurred fairly or highly certain gains of similar size.Pixels where changes in fmf were rather high for the vegetation zone showed that this was not sufficient for making decreases fairly certain if the increase of NPP t was not also low.
In dry tropical and subtropical regions encompassing desert, grassland, shrubland, and savanna, projected increases in productivity of zonal vegetation were high, generally entailing highly certain increases in SOC stocks with means in the range of 2-6 kg m −2 .In several shrubland pixels, zonal vegetation was projected to suffer from strong decreases in NPP, resulting in fairly certain decreases, especially when the land use simultaneously changed to crops or pasture.
Pixels in tropical forest zones were characterized by projected changes in zonal NPP between −0.07 and +0.72 and increases in fmf between 0 and 0.22.In locations where fmf increased strongly, reductions in NPP input due to land use change acerbated the loss of SOC from soil.Losses due to cropping or conversion to pasture were highly certain in most cases, with mean losses between 2 and 10 kg C m −2 .

Discussion
We assessed the certainty of changes in SOC stocks due to climate and land use change using a framework that explicitly considers the frequency distribution of values of controlling variables of decomposition processes.The greatest changes globally in mean SOC stocks after 75 years were due to absolute changes in NPP and thus C input into the soil.The effects of two contrasting global models of NPP on changes in soil C mass differed in sign (Fig. 5), showing the great importance of improving projections of NPP, especially with respect to N limitation and CO 2 fertilization.Direct climate effects on mean decomposition fractions via temperature and moisture were greatest at the drier edge of the tropics (Fig. 3).The reason is that increased precipitation throughout the tropics (as projected by ensemble GCM results) had a greater effect on decomposition fractions in the drier regions than in the humid tropics.The global patterns of SOC stock changes in our study and the enhanced NPP scenario agree with results obtained by mechanistic earth system models (Todd-Brown et al., 2014) and underline the validity of our modelling approach.At high latitudes the increase in temperatures were considerable, but the absolute effect of the increase on decomposition fractions remained small.The slight negative effect on SOC stocks was matched or surpassed by greater input from NPP.In warm regions, increases in temperature entailed relatively great decomposition fractions whose effects on SOC stocks were not necessarily matched by greater input from NPP.
Our results based on the limited-NPP scenario (Fig. 5a) may be more likely if future increases of NPP are limited by nutrient availability despite CO 2 fertilization (Norby et al., 2010), a concern also raised in a recent review of earthsystem models (Todd-Brown et al., 2014).The results based on enhanced NPP may become more likely if anthropogenic N deposition increases and alleviates nutrient limitations of growth (Hyvönen et al., 2007).The probability-weighted mean change of global SOC of 11 Pg (at a level of 75 % certainty) in the elevated-NPP scenario in our study, with consideration of wetlands, permafrost, and land use, corresponds to the lower end of the range of 15-100 Pg indicated by another study using frequency distributions for parameters (Hararuk et al., 2014) when both results are expressed as a fraction of the total SOC mass derived from the reference database used in each study (1061 and 1567 Pg).
Tropical mountain forests showed losses across both NPP scenarios and emerged as hotspots of SOC vulnerability.Most of these hotspots were characterized as zonal vegetation.Total modelled SOC loss in tropical hotspots (1 Pg, both NPP scenarios) comprised a small portion of the whole tropwww.soil-journal.net/1/367/2015/SOIL, 1, 367-380, 2015 ical loss (11 Pg, limited-NPP scenario; 6 Pg, enhanced-NPP scenario) on 2 % of the tropical land area.SOC losses in identified hotspots were highly certain.Pixels with gains were less consistent among NPP scenarios.Only the easternmost tip of Siberia emerged as a hotspot of certain gains.
The low certainty of SOC changes shown in this study reflects the limited certainty of the key variables in the terrestrial C cycle and their future changes, especially current SOC stocks, NPP t , and decomposition fractions.The uncertainty is unfortunately greatest for the location, spatial extent, and actual stock of soils rich in organic C (Köchy et al., 2015).Our study includes only part of the uncertainty in NPP 0 expressed by the variability among different global climate models linked to different vegetation models (Schaphoff et al., 2006;Sitch et al., 2008).This uncertainty is exacerbated by the need to estimate the harvest factor for future conditions.Currently, human activity reduces global NPP 0 by 24 % (Haberl et al., 2007).This fraction can be expected to increase in the future with increasing global populations.The uncertainty regarding the effect of temperature and moisture on decomposition fractions at the global scale and at long time scales is still great (Kirschbaum, 2006).Reported large-scale correlations between temperature or moisture and decomposition fractions (or rates) have been suggested to be spurious due to differences in litter quality and moisture along latitudinal gradients (Giardina and Ryan, 2000).This view has been contested (Davidson et al., 2000).On the one hand, our model framework supports a strong positive association between litter quality and aridity across vegetation zones.On the other hand, it shows only a weak association between litter quality and mean annual temperature and thus latitude.At the global scale, temperature sensitivity may be lower (Q 10 = 1.4-1.9,Hararuk et al., 2014;Ise and Moorcroft, 2006) than generally assumed from short-term incubation studies (Q 10 ≥ 2) and maximum decomposition occurs at greater soil moisture than assumed in most models (Ise and Moorcroft, 2006).If temperature sensitivity at the global scale is indeed lower than generally assumed, decomposition fractions in warm regions would not increase as much as projected, and tropical SOC stocks would decrease less strongly.
Differences in global climate models and vegetation models cause spatial variation in projected distributions of plant types (Alo and Wang, 2008;Gonzalez et al., 2010).In addition, global vegetation models do not account for adaptation of species in response to climate change.Our sensitivity assessment showed that plant types had an overall small effect on decomposition fractions via litter quality.This was mostly due to the fact that the woody fraction of NPP is the strongest predictor of litter quality at the global scale but did not change drastically in most places.The moderate sensitivity of SOC stocks to plant types suggests that our conclusions about changes of SOC stocks and their certainty are robust to uncertainties in the global distribution of plant types.
Changes in SOC stocks after 75 years due to changes in land use were on average small because the average probabilities of land use change per vegetation zone that we applied to each pixel were low.This effect is also exemplified by a retrospective study where drastic deforestation of Amazonian forest only caused a net decrease of 0.5 % of SOC stocks across the whole study area (Holmes et al., 2006).Quantitative, spatially explicit projections of the change of land use towards the end of this century require socio-economic models linked to vegetation models, which adds another layer of complexity and uncertainty.Instead, we examined different prescribed land use changes in individual locations.For tree-dominated ecosystems, SOC stocks decreased in the order zonal > plantation > pasture > annual crops, emphasizing that C input is one of the most important driving variables for predicting the net rate of SOC change; the greater the remainder after harvest, the greater the increase (or the smaller the losses).Draining wetlands is exposing SOC and was linked to considerable and highly certain SOC losses in all vegetation zones.Our model suggests that agroforestry producing sufficient C input or wetland cropping (rice paddies) may be ways to conserve SOC and extract food at the same time.Our modelling results agree in trend with reviews on land use change effects in temperate and tropical regions (Don et al., 2011;Poeplau et al., 2011).In contrast, our results do not show the great increases (> 100 %) in SOC stocks after conversion from cropland to grassland or forest in temperate zones and our projected increases for land use changes in tropical regions are greater (by a factor of 1.5-2) than those reported by Don et al. (2011).
The goal of our assessment of SOC stocks was to identify patterns of fairly and highly certain changes of SOC at the global scale.The SOC stocks and their changes simulated by our model were within the wide range of outcomes produced by different combinations of global circulation models (Schaphoff et al., 2006), global vegetation models (Sitch et al., 2008), and different soil modules (Yurova et al., 2010) and several earth system models (Todd-Brown et al., 2014).The delicate balance between higher NPP and higher decomposition rates as the main control of the gain or loss of SOC stocks that emerged from the aforementioned global studies was confirmed by our assessment.In addition, our study showed that the certainty for strong changes in SOC stocks is rather low.
Increasing aridity with climate change could substantially reduce SOC stocks in tropical peatlands (Li et al., 2007).Our framework allows attaching a probability to this.Assuming an SOC stock of 30 kg m −2 in evergreen tropical forests with 0.8-1.0kg yr −1 NPP, a water table < 20 cm below ground, and average annual temperature and aridity, our model projects a distribution of SOC stock with a mean of 24 kg m −2 .If precipitation decreased so that the water table decreases to > 20 cm, the potential SOC loss due to higher decomposition fractions would not be balanced by the higher NPP of 1.0-1.2kg yr −1 projected for future conditions.The distribution of potential net losses has a mean of 3 kg m −2 , with a 68 % certainty that the losses are > 2 kg m −2 .
In boreal and arctic regions we found few hotspots of vulnerability because increases in decomposition fractions and exposed frozen soil were smaller than increases in productivity.Our assessment agrees with ensemble results of 10 global climate models (Qian et al., 2010) incorporating CO 2 fertilization of NPP.The ensemble results, however, indicated that increases in NPP would level off towards the end of the 21st century while C emissions would continue increasing, causing a net source in northern high latitudes in the following century.Loss of SOC stocks could be greater than simulated by our model because of additional heat produced by microbial activity (Khvorostyanov et al., 2008) and if permafrost thawing conforms more to the algorithm of Khvorostyanov et al. (2008) from a global circulation model (Poutou et al., 2004) than to the algorithm used by us (Anisimov et al., 2002).Another study found higher C emissions in some locations in the permafrost region that were not matched by increased NPP.It is, however, unclear if this local observation is representative of the heterogeneous thermokarst landscape (Kuhry et al., 2010), where ponds from meltwater might limit C emissions (Moore and Knowles, 1989;Rouse et al., 1997).Spatially integrating measurements using eddy flux methodology over a ca. 10 ha area in a Siberian tundra site (70 • 50 N, 147 • 30 E; Parmentier et al., 2011) suggested that ecosystem C uptake would be low under climate warming, which is also projected by our assessment in both scenarios.
The accumulation and preservation of high SOC stocks is to a large extent due to conditions of low oxygen availability in wetlands.Existing global wetland maps show great heterogeneity in wetland extent and classification (Lehner and Döll, 2004;Köchy and Freibauer, 2010).For the purpose of soil C modelling, a (global) map indicating water table depth (and its variability) within the top 1 m may be more suitable than the combination of three indices -the index of O 2 constraint, the "wet" moisture class, and water-table estimates per land use class -that we used in the absence of such a map.Such a map might be achievable by assimilation of remote sensing (Finn et al., 2011) in wetland models (Fan and Miguez-Macho, 2011) and could provide the base to incorporate effects of aridity in wetlands (Ise et al., 2008).Changes in aridity and sea level rise may create new and vast wetland areas.In arctic regions, wetland areas might increase due to thawing of permafrost soil, snowmelt, and flooding (Rouse et al., 1997).These wetland dynamics are commonly not taken into account in global C models.The probability tables used in our model can be a starting point for such global modelling activities.

Conclusions
Our assessment showed, in a spatially explicit way, the great uncertainties associated with potential global SOC stock changes.Changes of SOC had a high certainty only in a few locations.In general, the strength of changes did not correlate with the certainty of changes.Therefore, conclusions about local and global changes would differ depending on what level of certainty one accepts for accounting changes.This aspect has been considered only via ensemble runs of models but not within models so far.Assessments of uncertainties within models could direct future research.The fertilization effect of CO 2 on NPP in the long term and its variation at the global scale is one of the major uncertainties.Global maps of current SOC stocks are derived from soil surveys dating as far back as the 1920s and complete information for soils rich in SOC comes from few soil profiles.Up-to-date records of up-to-date location, extent, and water table variation of wetlands at the global scale are incomplete.The same holds for permafrost and the active layer.Monitoring of these variables is crucial for decomposition models, assessments of global SOC stocks, and certainty of changes.

Figure 2 .
Figure 2. Probability distribution of four fraction modifying factors (fmf).Box width within each class of decomposition fraction controlling variable (x axis) is proportional to probability of the class of the fraction modifying factor.O 2 constraint labels: no, moderate, severe, very severe, rock or bare, permafrost, wetland.Note that the y axes are not linear.

Figure 3 .
Figure 3. Decomposition fraction modifying factor (fmf) under reference climatic conditions (a) and the difference in fmf between reference and target conditions (b).fmf summarizes the effects of temperature, humidity, soil reaction, and oxygen availability on the decomposition fraction.

Figure 4 .
Figure 4. Difference in NPP t between target (future climate and land use) and reference (current) conditions using (a) the NCEAS model (limited NPP) or (b) the LPJ model (enhanced NPP) for NPP 0 .

Figure 5 .
Figure 5. Changes of soil organic carbon stocks and masses using (a) limited NPP and (b) enhanced NPP.(c) Consensus: locations where changes ≥ |1 kg m −2| in both scenarios, overlaid on topography.

Figure 6 .
Figure 6.Certainty of changes of SOC stocks using (a) limited NPP and (b) enhanced NPP.(c) Consensus: locations where certainties of both scenarios are > 0.66 (P limited > 0.66 and P enhanced > 0.66).
Network design: fmf is the decomposition fraction modifying factor and F is the decomposition fraction.Boxes with thick edges indicate directly supplied information.Boxes with dashed edges (on the target side) indicate values copied from the reference side.Shaded boxes indicate dependent variables.