Digital Mapping of Soil Organic Carbon Contents and Stocks in Denmark

Estimation of carbon contents and stocks are important for carbon sequestration, greenhouse gas emissions and national carbon balance inventories. For Denmark, we modeled the vertical distribution of soil organic carbon (SOC) and bulk density, and mapped its spatial distribution at five standard soil depth intervals (0−5, 5−15, 15−30, 30−60 and 60−100 cm) using 18 environmental variables as predictors. SOC distribution was influenced by precipitation, land use, soil type, wetland, elevation, wetness index, and multi-resolution index of valley bottom flatness. The highest average SOC content of 20 g kg−1 was reported for 0−5 cm soil, whereas there was on average 2.2 g SOC kg−1 at 60−100 cm depth. For SOC and bulk density prediction precision decreased with soil depth, and a standard error of 2.8 g kg−1 was found at 60−100 cm soil depth. Average SOC stock for 0−30 cm was 72 t ha−1 and in the top 1 m there was 120 t SOC ha−1. In total, the soils stored approximately 570 Tg C within the top 1 m. The soils under agriculture had the highest amount of carbon (444 Tg) followed by forest and semi-natural vegetation that contributed 11% of the total SOC stock. More than 60% of the total SOC stock was present in Podzols and Luvisols. Compared to previous estimates, our approach is more reliable as we adopted a robust quantification technique and mapped the spatial distribution of SOC stock and prediction uncertainty. The estimation was validated using common statistical indices and the data and high-resolution maps could be used for future soil carbon assessment and inventories.


Introduction
Digital Soil Mapping uses statistical tools to quantify the spatial relationship between soil property values to its environmental covariates [1]. The digital mapping of soil organic carbon (SOC) at fine resolution is a challenging task [2] and the mapping is also a high priority for SOC assessment and monitoring [3]. Spatial models on SOC prediction has a long history (e.g., [4]). A range of techniques have been used to predict and map SOC from landscape to national or continental levels and Minasny et al. [5] provided a comprehensive review. Several researchers applied splines to model the vertical distribution of SOC in the soil profiles and predicted SOC at a landscape scale using data-mining tools and environmental variables as predictors [6][7][8]. Mishra et al. [9] calculated SOC pool in each soil horizon and applied geographic weighted regression to map SOC pool at a regional scale in midwest USA. Odgers et al. [10] used splines to derive SOC content from soil map units and predicted SOC at six standard soil depths for the entire USA. Arrouays et al. [11], Bou Kheir et al. [12], Chaplot et al. [13], Doblas-Miranda et al. [14], mapped SOC at national level using different statistical tools ranging from statistical aggregation to advanced regression tools such as regression trees.
The distribution of SOC changes across the landscape and it also varies by depth. In most soils, SOC is higher in the surface horizons and it decreases with depth. Such depth-wise variability is mostly continuous [15][16][17] except in soils with a strong human impact like some soils in the Netherlands [18]. Although most SOC studies and inventories are confined to 30 cm soil depth [19][20][21], the amount of SOC stored below 30 cm is of relevant in many ecosystems [22,23]. For accurate quantification of SOC stocks, a depth function needs to be modeled. Several tools exist: spline function [15,16], exponential decay function [24] or soil-type specific or profile depth functions [18,25]. For modeling SOC with depth, the equal-area spline function has been proven to be useful in several studies [7,8,10]. Spline predicted SOC values with depth act as a geo-referenced point data to which environmental variables are joined and prediction models are generated using digital soil mapping techniques.
In Denmark, studies on SOC dynamics and its quantification has started after a national wide soil database was established between the years 1975-1985. Based on simple statistical scalingup techniques, Krogh et al. [26] calculated a total stock of about 579 Tg and reported that 69% of it was stored in the soils under agriculture. Greve et al. [27] estimated topsoil SOC contents (g 100 g 21 ) for the whole country but have not assessed the SOC stocks. Bou Kheir et al. [12] predicted the spatial extent of organic soils across Denmark. Similarly, Olesen [28] and Taghizadeh-Toosi et al. [29] estimated the stock but did not map its distribution. Most of these studies have estimated the SOC content and stock but they have not explored the spatial distribution of SOC stock nor quantified the uncertainty of the SOC predictions.
We applied digital soil mapping techniques to quantify the SOC content and stocks for Denmark. The main objectives of this study were: to model the vertical distribution of SOC content and bulk density in soil profiles, to predict and map their spatial distribution using environmental variables, to identify major environmental variables responsible for SOC distribution, to estimate the SOC stock for the soils of Denmark, and to assess the uncertainties in the SOC predictions.

Study area
Denmark is situated in Northern Europe covering an area of approximately 43,000 km 2 . The temperate climate is characterized by a mild winters with annual mean winter and summer temperatures of about 0uC and 16uC [30]. Precipitation is well distributed throughout the year with an average annual rainfall of 800 mm in the west to 500 mm to the east. The country is relatively flat with a mean elevation above mean sea level is about 31 m. Denmark is divided into 10 physiographic regions2referred in this paper as geo-regions2based on geographical, climatic and soil-formation criteria.
The soils are coarse sandy to clayey to heavy clays as defined in the Danish Soil Classification System [31]. Soil in the western part of the country are developed in non-glaciated sandy parent material along the glacial flood-plains and Saalian moraine, whereas the soils from the central and eastern region have been developed on relatively young basal till high in finer materials [31].
Most of the soils in the north have been formed in sand mixed with uplifted marine sediments. More than 66% of the soils are classified as Podzols (Spodosols) and Luvisols (Alfisols). Podzols occupy a major portion in the west and Luvisols and Cambisols (Inceptisols, Entisols) in the central and eastern part of the country [32]. Peat soils (Histosols) occur in poorly drained basins throughout the country. About 66% of the total land area is used for agriculture with grain and potato as the main crops. Forest areas include spruce, pine and beach and these cover more than 12%.

Point data
Point (pedon) data on soil organic carbon (SOC) (g kg 21 ) and bulk density (D b ) (Mg m 23 ) were derived from two databases: Danish Soil Classification database (DSC) and Danish Soil Profile database (DSP). DSC consists of about 36,000 point observations from the topsoil (0220 cm) sampled randomly from agricultural fields in the period 197521980. From about 6,000 same pedons, soils from the subsoil (35255 cm) was also sampled. SOC content was determined by dry combustion using a LECO IR-12 furnace. DSP consists of a grid based data (767 km spacing) established during the 1990s for improve fertilizer recommendation in Denmark [33]. At each 850 grid intersection, soil samples were collected based on genetic horizons and were analyzed for SOC by dry combustion. For some selected horizons, samples were taken to determine D b . In addition, soil data from about 1100 profiles were used, and these were collected during the establishment of main gas pipeline system and other research activities across Denmark [34].
In total 40,250 topsoil point samples and 1,994 soil profiles with horizon based SOC data were used in this study. About 1,133 soil profiles included D b measurements.

Environmental covariates
The environmental covariates data used in this study are terrain parameters from the Digital Elevation Model (DEM) of Denmark derived using Light Detection and Ranging (LiDAR) technology. The LiDAR points were interpolated using Delaunay Triangulation and the output Triangular Irregular Network (TIN) surface was converted to a grid DEM with a 1.661.6 m spacing. This study resampled the original 1.6 m DEM to 30.4 m using simple aggregation considering the mean value. This 30.4 m grid size was also used in the previous studies in Denmark (e.g., [32,35]). Before aggregating to a coarser grid, the DEM was corrected by removing the pits and peaks of about 50 cm dimensions in order to ensure a regular flow on the surface. Once the DEM was processed, a number of terrain parameters (e.g., slope aspect, slope gradient, elevation, SAGA wetness index, multi-resolution index of valley bottom flatness (MrVBF), altitude above channel network, slopelength factor, over-land flow distance, and valley depth) were derived. Tools and algorithms incorporated in SAGA GIS [36] and Arc GIS V10.2 [37] were used to process the DEM and to derive these parameters.
Other covariates used were six choropleth maps which were compiled at different cartographic scales including: soil map, landscape types, geo-regions, geology, land use, and wetlands -see Table 1 and for more detail Adhikari et al. [35].
Key environmental variables affecting the spatial distribution of SOC and D b in Denmark were identified. The relative usage of the environmental variables used during SOC and D b prediction was calculated for each depth and their importance was expressed in percentage. This was done with Cubist software which determines the relative importance of variables based on their usage in the prediction rules. The prediction method adopted in this study is summarized in Figure 1.  Table 2. Average soil bulk density (Mg m 23 ) for different soil organic carbon levels (g 100 g 21 ) within the central wetlands [Source: [44]].

Modeling SOC and D b distribution
The vertical distribution of SOC and D b in the soil profiles was modeled with mass preserving equal-area quadratic splines [15] in R [38]. The mathematical derivation of the spline has been described by Malone et al. [8]. As the fitting quality of splines to profile attribute data depends on a smoothing parameter2lambda (l), we tested seven l values (0.00001, 0.0001, 0.001, 0.01, 0.1, 1, 10) for SOC and D b data from all the profiles and selected a l value that showed the best fit for all the profiles using the root mean square. With increasing l value, the fit becomes more rough. During the fit, we pegged the spline by introducing a 1 cm thick slice with a same SOC value at the topmost layer to prevent unnecessary extrapolation on the surface horizon.
Once the depth function of SOC and D b were modeled, a weighted-average value of these properties were derived for five soil depths (025, 5215, 15230, 30260, and 602100 cm) based on the GlobalSoilMap specifications [39]. To these new values for the standard depths from all soil profiles, environmental variables were intersected and used for statistical analysis and spatial prediction.

Mapping to the spatial domain
Spatial prediction of SOC content and D b at five depths was based on Regression kriging (RK) [40]. RK assumed that the spatial prediction function consists of a deterministic model formed by a regression, and the residuals of the regression (unexplained variation) are spatially correlated. The general principle of RK includes (1) regression, and (2) simple kriging of the residuals from the regression, where outputs from these two steps are added to obtain the final prediction. For the regression step, we adopted Regression-rules (RR) derived using Cubist software [41]. This tool generates a set of condition2based rules where each rule comes with a multiple regression prediction function that only operates once the conditions specified by the rule are met [42].
To build the SOC and D b prediction model in Cubist, the data set was split randomly into two sets: 75% for calibration and 25% for model validation. Before the data split, SOC content from the topsoil observations (0220 cm) were joined to the spline predicted SOC content from 025 and 5215 cm soil depths. Similarly, SOC content from subsoil observations (35255 cm) were attached to the spline predicted SOC at 30260 cm soil depth. This way, a larger number of point SOC observations were incorporated to the splined SOC data from the 767 km grid profiles. The Cubist tool was run for log transformed SOC [log SOC g kg 21 ] and D b data from each depth interval and the output was converted to a regular grid map using a program written in FORTRAN. For each calibration location, the difference between measured and RR predicted value was calculated and its spatial distribution over the study area was generated using local variogram and point kriging in VESPER program [43]. This continuous residual surface was added to the corresponding RR output to get a final prediction of SOC and D b for all five depths. Together with the prediction, a map showing the uncertainty of the prediction was generated. Both SOC prediction and uncertainty maps at each depth were then back-transformed to SOC in g kg 21 .
Bulk density in the peat areas Peat lands are mostly present along the central part of the wetlands across Denmark. Bulk density in those areas were adjusted according to Greve et al. [44]. For the three surface layers (i.e., 025, 5215, and 15230 cm soil depths), D b from 0230 cm was used, and for the 30260 and 602100 cm, D b from 302 60 cm, and 60290 cm were used ( Table 2). Table 3. Descriptive statistics of soil organic carbon content (g kg 21

SOC stocks
SOC stock for each soil depth was calculated according to Eq. (1) using SOC content and D b data. SOC stocks from the five layers were summed to obtain a SOC stock to 1 m soil depth. D b in Eq. (1) was corrected for gravel content.
Where SOC stock is the SOC stock (t ha 21 ), SOC content the SOC content (g kg 21 ), D b the soil bulk density (Mg m 23 ) and D the given soil layer thickness (cm).

Model validation
Model performance in predicting SOC content and D b was evaluated on 25% of the point data. The following three indices that were calculated: RMSE~ffi Where obs and pred are observed and predicted SOC and D b values from n number of observations at i th locations, ME is mean error, and RMSE is the root mean square error. Predicted SOC stock was determined for all five depth intervals and then aggregated to 0230 cm and 02100 cm soil depth. These two C stock maps were stratified based on soil and land use types. The soil class map of Denmark based on the FAO legend [32] and a land use map (CORINE data) were used. The CORINE database for Denmark has 31 classes [45] but in this study the legend was reduced to five major land uses types: artificial surface (urban areas, industrial areas, road network and ports etc.), agricultural areas, forest and semi-natural areas, wetlands, and others (e.g., coastal lagoon, estuaries etc.) for ease of comparison to other studies (e.g., [11,22,26]). We also stratified our stock maps based on Danish geo-regions. We first derived the area of each soil, land use and geo-regions class based on the number of predicted pixels within those classes and then calculated an average and total SOC stock for each class at 0230 and 02 100 cm soil depth. The 0230 cm depth represents the plough depth from agriculture areas and estimation of carbon for this depth is of interest to farm management. The top 1 m soil depth mostly represents a rooting depth of many field crops and may act as an important soil depth section for carbon balance and accounting studies.

Summary of SOC and bulk density
SOC content was highly variable and ranged from 0 to 562 g kg 21 for the topsoil (0220 cm) and from 0 to 569 g kg 21 in the subsoil (Table 3). Mean SOC decreased with soil depth and SOC at 602100 cm was about four times lower than the SOC in the 025 cm layer. With depth, the coefficient of variation (CV) of the SOC content increased. The CV at 025 cm was 184% and that for the 602100 cm was about 466%. The SOC data was positively skewed at all soil depths with a maximum skewness coefficient at 602100 cm. The equal-area splines modeled the depth-wise distribution and generated a continuous SOC profile to 1 m depth. The best l value to fit all soil profiles for both SOC and D b data was 0.1. Also average D b was found to be increased with soil depth. Up to 30 cm depth the D b was on average 1.44 Mg m 23 , whereas it increased to 1.52 Mg m 23 below 60 cm depth. Bulk density appeared to be less variable with depth (Table 3). Figure 2 shows a measured and spline predicted SOC and D b for a coarse sandy soil under agriculture area from the Saalian moraine soilscape in western Denmark (West Jutland). Measured SOC from different horizons decreased with depth except at 532 73 cm where it increased due to podsolization. Although a spline should pass through the mid-point of each measured horizon for this soil profile, in this case the spline slightly extrapolated the SOC value at 35255 cm due to the selected l at 0.1.

Prediction rules for SOC and D b
Depending on the soil depth, 17 to 54 condition-based regression rules were generated while predicting SOC and D b .  Table 4. Relative usage (%) of the environmental variables to predict soil organic carbon at different soil depths in Denmark. This paper only included one of the 54 rules produced during SOC prediction at 025 cm soil depth, as an example. This rule used elevation, SAGA wetness index, MrVBF, precipitation and slope gradient to predict SOC in the areas where the MrVBF index is higher than 6.6 and consisted of fine sandy soils from a Moraine landscapes from the geo-regions (e.g. Himmerland). This rule was only valid for 400 locations where the mean SOC was 3.1 [log SOC g kg 21 ] and the prediction error was about 0.2 [log SOC g kg 21 ].

Identifying predictors
Variables were identified for SOC and D b prediction based on their relative usages in the model (Table 4). In all models, precipitation appeared to be the most dominant variable followed by altitude above channel network and SAGA wetness index to predict SOC. As an example, to predict SOC at 5215 cm, precipitation had a usage of 98% for both rule setting and developing a linear prediction model followed by SAGA wetness index which had a relative usage of 67% and 94%. For this model, insolation and slope aspect had the lowest contribution. Geology became robust with increasing soil depth, whereas land use was important for SOC prediction of the surface layers. Geology, soil map, wetlands, land use, precipitation, MrVBF, SAGA wetness index, elevation, slope gradient, slope-length factor, and altitude above channel network were among the predictors that had a relative importance of .60%. Similarly, land use, soil map, geology, slope gradient, SAGA wetness index, MrVBF, and elevation appeared to be the most important variables for predicting D b at all soil depths.

Predicted maps
Predicted maps of SOC content ( Figure 3) and D b (map not shown) at five soil depths were produced at a resolution of 30.4630.4 m. The highest mean SOC content was in the 025 cm layer (mean 20 g kg 21 ; sd 11 g kg 21 ). Predicted SOC content decreased with soil depth and at 602100 cm, it was on average 2.2 g kg 21 . The soils of western Denmark have relatively higher SOC content than the rest of the country. The northern part has a moderate amount of SOC with two large raised bogs with high SOC contents. Along the coastline, especially in the west, soils with lower SOC concentration were mapped.
The prediction errors were higher towards the west and along the coastline. The prediction error increased with soil depth. For the 025 cm layer, the mean error was 1.1 g kg 21 and it increased to 1.8 g kg 21 at 602100 cm soil depth.
Predicted SOC content for the FAO2UNESCO soil groups is shown in Table 5. Average SOC content ranged between 11.8 to 52.6 g kg 21 at 0230 cm and between 1.9 to 37.5 g kg 21 at 302 100 cm. The highest SOC was observed in Histosols and the lowest in Arenosols at most soil depths. Model validation SOC prediction models performance is summarized in Table 6. The best prediction was found at 5215 cm soil depth for both training and test data sets. The model performance at 602100 cm was relatively poor compared to the rest of the soil depths.
Negative ME values suggested that almost all the prediction models were negatively biased suggesting some under prediction of the mapping of SOC distribution.

SOC stocks
SOC stock maps were made for two soil depths (Figure 4). For 0230 cm, average SOC stocks were about 72 t ha 21 and for the top 1 m depth, it was about 120 t SOC ha 21 . Most of the western and northern parts of the country have more than 80 t SOC ha 21 in the top 30 cm whereas the average stock in the eastern part of the country was less than 80 t SOC ha 21 . Total SOC stock was calculated for each geo-region ( Figure 5) and Himmerland and West Jutland had an average stock of about 135 t SOC ha 21 followed by North Jutland and Thy both having a mean stock of . 120 t SOC ha 21 . The soils of West Jutland and East Denmark contain almost 50% of the total SOC stock in Denmark.
Luvisols and Podzols contain about 60% of the total SOC stock ( Table 7). Other soil groups that contain significant amounts of SOC were Cambisols (6%), Gleysols (9%), and Arenosols (9%). Although Histosols had SOC stock of 176 t ha 21 , its total content was 20 Tg SOC. For all soil groups, more than 58% of the total stock was in the top 30 cm. Unmapped areas in the soil map representing major Danish cities that may also contain a substantial amount of carbon [46].
Of the total stock of 570 Tg SOC about 59% is in the upper 30 cm. Soils under agriculture have an average stock of 121 t ha 21 and contain about 444 Tg which is almost 78% of the total estimated stock (Table 8). Another large fraction of SOC stock is found in the soils of the forest and semi-natural areas, and these had a stock of 39 Tg in the top 30 cm and about 67 Tg up to 1 m soil depth. Wetland areas contain large amounts of SOC, and average SOC stock within 1 m soil depth was about 152 t ha 21 which is nearly 2% of the total stock. Almost 90% of the total SOC stock within the top 1 m soil depth is found in the soils under agriculture, forest and natural areas.

Discussion
Here we have predicted the distribution of SOC contents and stocks across Denmark including an assessment of the prediction error. We have also estimated the SOC contents and stocks for different land uses and soil types. This discussion focuses on the  prediction model, the importance of the variables, the uncertainty, and the SOC contents and stocks.

Prediction model
The equal-area spline fit the discrete horizon data but it also harmonized the profile by disaggregating the horizon bulk data and generated a continuous function of SOC distribution. Several other researchers have advocated the usefulness of such splines in depth-wise mapping of SOC in different parts of the world applied from local to regional extents (e.g., [7,8,10]). Pegging of spline by introducing an artificial horizon on the surface benefitted our splines that restricted biased extrapolation of SOC on the surface horizon.
The spatial prediction method (i.e. rule-based regression using the Cubist software) was capable of exploring the relationship of SOC to its environment predictors. The prediction rules were conditioned to a given environmental settings such that each rule is valid only to that specific boundary within which SOC distribution presumably less heterogeneous compared to the areas outside where other conditions and rules prevail. For example, SOC content in forests or clayey soils might be different than the SOC from agriculture or sandy soils hence different prediction models operated in these two specific areas. Such a beneficial and advanced function has been used by Lacoste et al. [7] who found regression-rules combined with the spline depth function useful for producing a detailed pseudo-3D map of SOC content in heterogeneous agricultural landscape in France. Minasny and McBratney [42] found that regression rules provided a greater accuracy compared to partial least squares while predicting total carbon content. Several other studies have applied this tool in digital soil mapping (e.g., [35,47,48]).
Variable importance. The environmental variables used to predict SOC content showed a varying level of importance in the model. There was a large influence of precipitation, land use, soil type and some terrain parameters such as elevation, slope gradient, SAGA wetness index, and MrVBF on the spatial distribution of SOC content. The influence of topographic parameters on SOC distribution has been documented in other studies [6,8,9,[49][50][51]. Similarly, land use, precipitation, soil types, wetlands were found important while mapping SOC distribution [11,14,18,[52][53][54][55]. The influence of terrain parameters on SOC composition and distribution can be linked to its behavior on soil re-distribution through erosion and deposition, in the maintenance of vegetation cover and rooting depths, and in soil drainage that affects SOC decomposition as well as vegetation. In Denmark, the influence of elevation, soil types, geology, and slope gradient was also documented by Bou Kheir et al. [12] when predicting SOC in wet cultivated lands.
The categorical variables such as soil types, geo-regions, and land use were used in defining the rule conditions and continuous variables (terrain parameters) in regression functions. Some terrain parameters like elevation, SAGA wetness index, MrVBF, altitude above channel network etc. were utilized in setting rule conditions. This combined approach in defining the conditions partitioned the study area into several possible strata where SOC distribution was supposed to be more homogeneous and in each stratum several terrain parameters were again used to make sure that the withinstratum SOC variability was well captured by the model. This has made the model robust in predicting SOC content.
Uncertainty. Based on validation indices, the prediction models showed a higher performance (i.e., higher R 2 and a lower RMSE) in calibration data (75%) than in validation data (25%). The uncertainty of the SOC prediction increased with depth. The R 2 value ranged between 0.23-0.63 -the higher values for the surface layers. It suggests that our prediction was able to capture up to 63% of total SOC variability. The range of R 2 values was comparable to similar SOC mapping studies where internal validation was applied [24,50,56,57]. Values higher than 0.7 are unusual and values ,0.5 are quite common in soil attribute predictions [58]. The difference in R 2 value between the two datasets could probably be attributed to the different data density used in prediction. Minasny et al. [5] reviewed several previous SOC mapping studies, and reported that with increasing data density the R 2 of prediction was larger. Similarly, SOC maps from the upper soil depths were associated with a low prediction error compared to the maps with depth. This could be linked to the terrain parameters used as predictors because most of these parameters explain soil surface phenomena and the uncertainty increases with depth [24]. Moreover, higher data density from the surface layers (e.g., 025 and 525 cm) might have a positive influence on prediction performance.
SOC content and stock. SOC data from two main sources were used in this study: high-density (about 1 observation per km 2 ) topsoil and subsoil samples, and from soil profiles at each 7-km grid-intersections together with the data from profiles along the pipe lines across Denmark. The grid and point data covered the entire country.
The SOC data were all from dry combustion analysis and from the period 1975285 so we have modeled and predicted here the SOC contents and stocks for that period. Assuming a steady-state condition at that period, i.e. the carbon levels represent the equilibrium with its physical environment and landuse. Likely, SOC contents have changed since that time, thus the map can be used as a baseline to indicate spatio-temporal changes.
A higher SOC content was found in surface soils, and the lowest average SOC content of 3 g kg 21 was found at 602100 cm soil depth. The western part of the country mostly consists of glacial floodplains and old Saalian moraine landscapes where soils are predominantly sandy. To increase the soil fertility in those areas, farmers has been adding for decades large amounts of manure in the form of pig slurry. This may have led to increased SOC content and gave the highest SOC stock in those areas compared to other geo-regions in Denmark as also reported in Taghizadeh-Toosi et al. [29].
Our results are in line with the previous SOC estimations in Denmark. According to Krogh et al. [26], total stock in Denmark to 1 m soil depth ranged between 5632598 Tg with 579 Tg as the average which was comparable to our prediction of 570 Tg. Similarly, SOC stock from the topsoil (0228 cm) was about 230 Tg while our prediction showed about 266 Tg. Small differences  in stocks estimates between the two studies were noticed. For example, SOC stock in agricultural areas as predicted by Krogh et al. [26] was 404 Tg, whereas the estimate in the current study is 444 Tg. This could be due to a difference in areas estimated for agricultural lands in two studies. Our prediction used CORINE data that suggested an agricultural area of nearly 32,942 km 2 , whereas this area was 28,883 km 2 in the previous study based on AIS (Area Information System). In a separate study, Olesen [28] reported a total stock of 604 Tg from 0260 cm soil depth calculated for an area of 34,000 km 2 based on AIS. Unlike in our study where bulk density from the peat lands were adjusted, Olesen applied a standard bulk density of 1 Mg m 23 for all peats or organic soils and that might have increased the stock leading to over estimation. Previous estimates of SOC stocks had not quantified the spatial distribution of SOC stock nor validated their prediction. Our approach seemed to be more reliable and the data generated could be useful for future SOC content and stock assessments.
A large amount of SOC was present in the soils under agriculture. Our estimation of 121 t ha 21 is in a comparable to the SOC stock of 111 t ha 21 estimated for arable lands in Scotland [59] but it was slightly higher than the arable stock in Southeast Germany [60]. Similarly, our predicted average SOC stock of 121 t ha 21 from the agricultural areas was lower than the findings of Taghizadeh-Toosi et al. [29]. The later study used a soil-type dependent standard D b established for 767 km grid and perhaps overestimated the SOC stock.
The estimation of average SOC stock in different soil groups within 0230 cm depth was comparable to Arrouays et al. [11], but the total stock was approximately nine times less than reported for France (which is almost 13 times larger than Denmark). Contrary to our study, the soils under forest soils contained more SOC than the arable soils in France. Similar results were also reported by Chaplot et al. [13] from Laos. A slightly lower average stock in the soils under forest in our study might appear due to the inclusion of non-forested areas (semi-natural areas) in the same class and exclusion of litter layer that possibly lowered its value. Decade long intensive soil management practices such as addition of large amount of manure to the agricultural soils might have increased SOC contents and consequently the SOC stock from those areas [29]. Likewise, our estimated average stock for 1 m soil depth from different soil groups, for example, Cambisols (110 t ha 21 ), and Gleysols (140 t ha 21 ) were in a agreement with the global estimate of Batjes [22] using the same soil depth. We also observed .50% stock stored within the top 30 cm soil depth for almost all soil groups which also corresponds to the finding of Batjes [22].
Although the SOC stock in Denmark might have changed substantially over the past few decades [29], our estimation based on the available data has provided a baseline SOC. Together with the estimation of uncertainty, the maps are more reliable and could be useful in environmental research in Denmark. It could support national carbon accounting and carbon balance studies and also act as a back ground information for monitoring temporal SOC changes in Denmark.

Conclusions
This study predicted the spatial distribution of SOC content (g kg 21 ) at five soil depths intervals (025, 5215, 15230, 30260, and 602100 cm) and quantified its stock (t ha 21 ) to 1 m soil depth for Denmark. DEM derivatives and soil maps were used as predictors where condition-based regression rules were applied to quantify the spatial relationship between measured SOC and D b with the predictors. The following can be concluded from this study: N Equal area spline modeled the continuous depth function of SOC and D b data from discrete soil horizons in soil profiles from Denmark.
N The most important variables that influenced SOC distribution across Denmark were precipitation, wetlands, land use, soil types, elevation, and saga wetness index.  N This article is an example for a national level SOC mapping based on GlobalSoilMap procedures and the methods applied can be tested and used in other part of the world.