Wetlands (2023) 43:65 https://doi.org/10.1007/s13157-023-01705-3 WETLAND SOILS Modeling the Spatial Distribution of Soil Organic Carbon and Carbon Stocks in the Casanare Flooded Savannas of the Colombian Llanos Javier M. Martín‑López1  · Louis V. Verchot1 · Christopher Martius2 · Mayesse da Silva1 Received: 8 March 2022 / Accepted: 4 June 2023 / Published online: 12 July 2023 © The Author(s) 2023 Abstract Flooded savannas are valuable and extensive ecosystems in South America, but not widely studied. In this study, we quantify the spatial distribution of soil organic carbon (SOC) content and stocks in the Casanare flooded savannas. We sampled 80 sites at two soil-depth intervals (0-10 and 10-30 cm), where SOC values ranged from 0.41% in the surface and 0.23% in the sub-surface of drier soils to over 14.50% and 7.51%, in soils that experienced seasonal flooding. Spatial predictions of SOC were done through two digital soil mapping (DSM) approaches: Expert-Knowledge (EK) and Random-Forest (RF). Although both approaches performed well, EK was slightly superior at predicting SOC. Covariates derived from vegetation cover, topography, and soil properties were identified as key drivers in controlling its distribution. Total SOC stocks were 55.07 Mt with a mean density of 83.1±24.3 t·ha-1 in the first 30 cm of soil, with 12.3% of this located in areas that experience long periods of flooding (semi-seasonal savannas) , which represented only 7.9% of the study area (664,752 ha). Although the study area represents only 15% of the total area of the Casanare department, the intensive pressure of human development could result in the reduction of its SOC stocks and the release of important amounts of greenhouse gases into the atmosphere. At regional level, the impact of a large-scale land use conversions of the flooded Llanos del Orinoco ecosystem area (15 Mha) could transform this area in a future source of important global emissions if correct decisions are not taken regarding the land management of the region. Keywords Digital soil mapping · Random Forest · Expert knowledge · Carbon sequestration · Wetlands · Driving factors Resumen Las sábanas inundadas son ecosistemas extensos y valiosos en América del Sur, pero no han sido ampliamente estudiados. En este estudio, cuantificamos la distribución espacial del contenido y las reservas de carbono orgánico del suelo (COS) en las sábanas inundadas de Casanare. Se muestrearon 80 sitios a dos diferentes profundidades del suelo (0-10 y 10-30 cm), donde los valores de COS variaron de 0,41% en la superficie y 0,23% en la subsuperficie de los suelos más secos, a 14,50% y 7,51% en suelos que experimentan inundaciones estacionales. Las predicciones espaciales de COS se estimaron a través de dos enfoques de mapeo digital de suelos (MDS): Expert-Knowledge (EK) y Random-Forest (RF). Aunque ambos enfoques presentan un buen desempeño, EK fue ligeramente superior en la predicción del COS. Las covariables derivadas a partir de la cobertura vegetal, la topografía y las propiedades del suelo se identificaron como variables clave en explicar su distribución. Las reservas totales de COS fueron de 55,07 Mt, con una densidad media de 83,1±24,3 t·ha-1 en los primeros 30 cm del suelo, de los cuales el 12,3% se localizan en zonas que experimentan largos periodos de inundación (sabanas * Javier M. Martín-López j.m.martin@cgiar.org * Mayesse da Silva m.a.dasilva@cgiar.org 1 International Center for Tropical Agriculture (CIAT), Multifunctional Landscapes, Km 17 Recta Cali-Palmira, Valle del Cauca Z.C. 763537 - A.A. 6713 Palmira, Colombia 2 Center for International Forestry Research (CIFOR) Germany gGmbH, Bonn, Germany Vol.:(012 3456789) 65 Page 2 of 20 Wetlands (2023) 43:65 semiestacionales), lo que representa solo 7,9% del área de estudio (664.752 ha). Aunque el área de estudio representa solo el 15% del área total del departamento de Casanare, la intensa presión producto del desarrollo humano podría resultar en la reducción de sus reservas de COS y la liberación de importantes cantidades de gases de efecto invernadero a la atmósfera. A nivel regional, el impacto de una conversión del uso del suelo a gran escala en los ecosistemas inundables de los Llanos del Orinoco (15 Mha) podría transformar esta área en una fuente futura de importantes emisiones globales si no se toman decisiones correctas en cuanto al manejo del suelo de la región. Palabras clave Mapeo digital de suelos · Random forest · Expert Knowledge · Secuestro de carbono · Humedales · Factores de formación Introduction Köchy et al. (2015) estimated that they cover around 9% of the tropical land area and contain around 40 Pg of soil organic Tropical wetlands provide important ecosystem services, carbon (SOC), mostly in floodplains and marshes. In addition including global regulating roles in biogeochemical cycles to uncertainty about the extent of tropical wetlands, the carbon and biosphere-atmosphere interactions. They store large density or carbon stock per unit area is also uncertain. Assess- amounts of carbon, are a sink for atmospheric C O2, and ments of SOC storage vary widely due to both definition and make up around 60-80% of the natural atmospheric meth- methodological differences, from about 120 Gt to 535 Pg of ane (CH4) source (Chmura et al. 2003; Köchy et al. 2015; SOC (Mitra et al. 2005). This is an important constraint in Rice et al. 2016; Saunois et al. 2020; Whiting and Chanton earth systems modelling, where lack of data on the carbon 2001; Zhu et al. 2017). In addition to being part of natural balance and gas fluxes in tropical wetlands limits accurate pre- hydrology systems and to ensuring water quality and regu- dictions about future climate and the importance of wetlands lating its flow (IPBES 2019), they are also key reservoirs in climate change (Melton et al. 2013; Sjögersten et al. 2014). of biodiversity, they support fisheries, they are important South America has the largest extent of tropical wetlands, sources of food, filter and retain sediments, and are cultur- most of them found in floodplains where they are subject to ally and aesthetically important to many communities and seasonal fluctuations in water levels (Gumbricht et al. 2017; cultures (Mitsch et al. 2015; Shiel et al. 2016). One study Junk et al. 2013). The continent has large extents of season- estimated that 25% of the global value of ecosystem services ally flooded savannas, including the Pantanal, Bananal Island, is provided by wetlands (Costanza et al. 2014). and the Lavrado in Brazil; the Iberá in Argentina, the western Wetlands soils, under regular anaerobic conditions, Guianan flooded savannas in, Guiana, Venezuela, and Suri- reduce the decomposition of organic plant materials, result- name; the Llanos de Moxos in Bolivia; and the Llanos del ing in the accumulation of organic matter. This enables Orinoco (Barbosa et al. 2012; Hamilton et al. 2002), which is wetlands to effectively store substantial amounts of soil the subject of this study. The Llanos del Orinoco cover about organic carbon (SOC) (Nahlik and Fennessy 2016). Despite 38.31 Mha (Barreto and Armenteras 2020) in the Orinoco wetlands only occupying about 5-8% of the earth's surface, River Basin of Colombia and Venezuela. The Llanos del Ori- (Mitsch and Gosselink 2007), point out that they store more noco landscape is a mosaic of ecosystems that include forests, SOC than all types of vegetation on earth. Wetlands are extensive areas of Moriche (Mauritia flexuosa) palms, upland believed to store between 20-30% of the world's soil organic savannas, and flooded savannas. The flooded area in this region carbon (SOC), and in some cases, the SOC concentration shows significant interannual variation. Hamilton et al. (2004) can exceed 40% (Vepraskas and Craft 2016). This character- found that the long-term mean inundation area was about 3.47 istic makes wetlands a key ecosystem for regulating water, Mha, with a median area of 2.53 Mha, which makes it the sec- climate, and biodiversity. ond largest flooded savanna wetland in South America after the In this order, the development of accurate maps of wet- Brazilian Pantanal. For this study we focused on the flooded lands and soil organic carbon (SOC) is crucial, as it enables savannas of the Casanare in Colombia, which drains into the an adequate inventory and monitoring of these carbon-rich Meta River, one of the major tributaries of the Orinoco River ecosystems, leading to better quantifications of SOC stocks that originates in the Andes. Few data exist for this region, and greenhouse gas (GHG) emissions (Page et al. 2011). As but the area is undergoing rapid change with the introduction stated by Minasny et al. (2020), a better and rapid estima- of exotic pasture grasses and extensive areas of monoculture tion of peatlands distribution allows stakeholders to make cropping. Additionally, fire has been used in this landscape informed decisions regarding the management, conserva- for centuries, although its duration and intensity in the flooded tion, and restoration of these ecosystems. Nonetheless, the savannas typically is low due to the presence of rocky out- extent and state of tropical wetlands and the quantification of crops, sand dunes, and low biomass on sandy patches, which soil organic carbon stocks are still uncertain (Xu et al. 2018). are spread throughout the landscape (Lasso et al. 2010). 1 3 Wetlands (2023) 43:65 Page 3 of 20 65 SOC accumulation in flooded savannas appears to be wet season. Rainfall distribution follows a unimodal regime related to length of flooding; short-term measurements in with total annual average precipitation of 2,684 mm, char- the flooded Colombian Llanos del Orinoco suggest stable acterized by a rainy season from April to November, during to slightly accumulating SOC stocks (Vega et al. 2014). The which most of the area is flooded, and a short dry season objectives of the present study were to map the spatial dis- from December to March when the water drains and many of tribution of SOC content and stocks throughout two digital the flooded wetlands dry out (Buriticá Mejia 2016; Castillo- soil mapping (DSM) approaches and determine the factors Figueroa et al. 2019; Lasso et al. 2010). that influence their variability. To do this, we focused on the The land cover is mainly composed of native eastern Casanare Department in Colombia, where the largest savanna, native forest, and permanent water bodies. extents of flooded savannas in the Orinoco Basin are located Native savanna is sub-classified in this study as sea- (Gumbricht et al. 2017; Hamilton et al. 2004). sonal savanna, hyperseasonal savanna, and semisea- sonal savanna, based on its vegetation and soil-water conditions (Cabrera-Amaya et al. 2020; Mora-Fernández Material and Methods and Peñuela-Recio 2013; Peñuela et al. 2014; Romero Duque et al. 2018; Sarmiento 1983). Seasonal savanna is Study Area a tree savanna, characterized by discontinuous grasses, dicotyledons and herbaceous vegetation, with the pres- The study area (Fig. 1), known as the “Sabanas Inunda- ence of scattered trees and shrubs. It develops under bles de Casanare” (Casanare flooded savannas) is part of a bi-seasonal regime marked by a dry season (4 to 6 an extensive floodplain ecosystem located in the “Llanos” months) and a wet season, in which it remains f lood ecoregion in the Orinoco River Basin (Olson et al. 2001). free. In the dry season, the vegetation acquires a wilted With a total area of 664,752 ha, the study area is located appearance, and large numbers of plants die off, increas- between the latitudes 6.30°N and 5.52°N and longitudes ing fire frequency. This ecosystem is located mainly on 71.07°W and 69.84°W. It has an elevation ranging from 78 plain-convex banks, continental dunes and sandy depos- to 146 m.a.s.l. and an average slope of 1.55 degrees. The its, where coarse-textured soils, compaction, and low annual average temperature in the region is approximately water availability create difficult conditions for vegeta- 26°C, with an average maximum temperature of 33°C during tion growth. Hyperseasonal savanna is the most widely the dry season and an average minimum of 21°C during the spread ecosystem in the region. It is mostly a treeless Fig. 1 Study location in the Casanare Flooded Savannas, Colombia. a) location within the tropical savannas ecoregion “Llanos” (Olson et al. 2001) b) position within the Casanare department (dark red line) and c) local land cover map of the study area 1 3 65 Page 4 of 20 Wetlands (2023) 43:65 savanna, where perennial grasses and annual dicotyle- wetlands, resulting in the maintenance of water through- dons predominate and sporadic woody plants, shrubs, out the region for most of the year (Sarmiento 1984). and palms are found. This ecosystem is temporarily flooded during the rainy season, where water table level Data Collection and Soil Characterization can reach up to 20 cm above the surface, but it dries out quickly after the rain ends. It experiences four seasons: A total of 80 sites were sampled in February 2019 at two drought, onset of soil saturation, flooding, and decline soil depth intervals (0-10 cm and 10-30 cm), resulting in 160 in soil saturation. Soil texture is mostly represented by soil samples (disturbed and non-disturbed). The sites were sandy-loam and clay-silty loam. Semiseasonal savanna sampled following two strategies: 43 sites located in four is a treeless savanna, with very rare presence of woody transects along the peat depth distribution predicted by the plants or palms. Vegetation (grasses, sedges, march and Global Wetlands Map (Gumbricht et al. 2017) and 37 sites aquatic plants) varies in function of flooding duration distributed across the study area based on the conditioned and water table level, which could reach up to 1 m above Latin hypercube sampling (cLHS) strategy (Fig. 2). The the surface. This ecosystem is found in topographic cLHS method is a stratified random procedure for sampling depressions with poorly drained conditions, and clay variables from a multidimensional distribution of covari- soils. Semiseasonal savannas experience long periods ates, which defines sampling locations that capture the great- of flooding (between 8 and 11 months) and never pre- est variability among covariates (Minasny and McBratney sent water deficit, even during the dry season, serving 2006). Maps of covariates used in the cLHS to define the as an important food and water source in the region. sampling site locations included: the Orinoquía land cover Similarly, the native forest includes Moriche palm and map at 30 m resolution ; geology and lithology at a scale of gallery forests. Moriche palm corresponds to plant com- 1:100.000 (IGAC 2014); peat depth and wetland predictions munities that are dominated by Mauritia flexuosa palms with a spatial resolution of 123 m (Gumbricht et al. 2017); and numerous herbaceous and shrubby species. They landforms (Jasiewicz and Stepinski 2013) with a spatial are present in flooded areas or in soils with medium- resolution of 12 m (derived from the DEM); and the cost fine textures that do not present water deficits during distance, which represents the cost of reaching each point the dry season. Gallery forest, in turn, is dominated by in the landscape from the road network based on topography a great variety of woody trees, palms, and shrubs and (slope at 12 m resolution) and distance from roads (Aitken are developed close to permanent water bodies, such as 1977; Roudier et al. 2012). The Orinoquía land cover map of riverbank levees, alluvial, and overflow plains of rivers, Colombia was developed from Landsat 5-8 images at a reso- channels, and ravines. lution of 30 m, the landforms were mapped using the DEM For mapping purposes, a local land cover map derived ALOS/PALSAR scaled up at a resolution of 12 m based on from remote sensing images (Fig. 1) was created through the geomorphons approach (Jasiewicz and Stepinski 2013), a supervised image classification algorithm (further infor- and cost distance was mapped according to Douglas, (1994). mation in “Environmental Covariates” section). After Due to the adverse field conditions, some sampling sites classification, five categories were defined: seasonal needed to be readjusted to locations with similar conditions, savanna, hyperseasonal savanna, semiseasonal savanna, keeping original sampling characteristics (Fig. 2). forest, and water. Gallery forest and Moriche palm were The soil samples were air-dried, sieved, and analyzed at combined within a unique forest category since their the International Center for Tropical Agriculture (CIAT) lab- spectral signatures were not clearly differentiated. oratory for determining bulk density (BD) and soil texture Soils of Casanare flooded savannas are predominantly (clay, sand, silt). BD and soil texture were performed by the Inceptisols and Entisols formed from alluvial, fluvio- volumetric ring (Blake 1965) and hydrometer (Bouyoucos torrential, and fluvio-glacial sediments deposited from 1936) methods, respectively. SOC, stable 13C isotope (δ13C), the upper Pleistocene to the lower Holocene, upon which total nitrogen (Nt), and stable 15N isotope (δ15N) analyses Aeolian materials subsequently accumulated (IGAC were performed using the dry combustion method and mass 2014). Relief is mainly composed of large extensions of spectrometry at the UC Davis Laboratory. Carbon to nitro- aeolian and alluvial terraces, dissected by depressions and gen ratios (C:N ratios) were calculated from SOC and Nt. channels covered by fine materials (floodplains), which Other characteristics such as land cover (seasonal remain flooded for a great part of the year, and some ele- savanna, hyperseasonal savanna, semiseasonal savanna, vated sandy areas called “Medanos” or continental dunes, gallery forest and Moriche palm); geomorphic classifica- which have been formed longitudinally by the northeast- tion (rise: slightly elevated areas (1-3% slope), talf: very erly trade winds (IGAC 2014). The combination of the low slope gradients (0-1%), dip: depressions, and continen- regional geomorphology and the high annual precipitation tal dunes: aeolian mounds of sand (> 3% slope)) defined creates the hydrological connectivity between rivers and according to Haskins et al. (1999) and Soil Science Division 1 3 Wetlands (2023) 43:65 Page 5 of 20 65 Fig. 2 Spatial distribution of the original cLHS sampling sites (purple points), field sampled sites (blue points) and peat depth (m) distribution predicted by Gumbricht, et al. (2017) Staff (2017); and generalized soil texture (sandy – coarser comparisons. In both cases, significant tests were carried out textures: sand, loamy sand, sandy loam; loamy – medium using the “PMCMRplus” R-package including a Bonferroni textures: loam, silt loam, silt, sandy clay loam, clay loam, p-adjustment method. silty clay loam; and clayey – finer textures: sandy clay, silty clay, clay) simplified from the soil texture taxonomy classi- Environmental Covariates fication (Soil Science Division Staff 2017) at each sampling site were reported. To map the distribution of SOC content at 0-10 cm and 10-30 cm using digital soil mapping (DSM) (McBratney et al. 2003) a set of 32 environmental covariates was considered (see Statistical Analysis Table S1 in the supplementary material). Terrain attributes were derived from the DEM – ALOS/PALSAR scaled up at Descriptive statistics of chemical and physical soil proper- 12 m spatial resolution (Shimada et al. 2014). Geomorphology ties (0-10 cm and 10-30 cm) were done using SAS Studio. maps at a scale of 1:100.000 were obtained from the geoportal An ANOVA on log-transformed data, followed by the Stu- of the Instituto Geográfico Agustin Codazzi (IGAC). Because dent-Newman-Keuls test was carried out for normally dis- of lack of soil texture differences at the two depths, soil texture tributed data, while the Kruskal-Wallis non-parametric test was generated by using the DSM approach random forest at followed by a pair-wise separation using the Dwass, Steel, depth of 0-30 cm. Vegetation indices (VI) were calculated in Critchlow-Flinger method was applied for data that could Google Earth Engine from SENTINEL-2 spectral images for not be normalized. 3 distinct periods (dry: December (2018) – February (2019), SOC content at each depth was analyzed according to wet: May (2018) – October (2018), transition dry-wet: Novem- different environmental conditions: land cover, geomorphic ber (2018) – December (2018)). Climate data (temperature and feature, and generalized soil texture. Statistical analyses precipitation from 2010 to 2018) was obtained from CHELSA considered were data normality, homogeneity of variances, and MODIS spectral images. The local land cover (Fig. 1) map and design balance analysis tested using Shapiro-Wilk was created through a supervised classification method com- (“stats” R package), Levene (“car” R package), and the bal- bining field data and SWIR-1 (B11), near-infrared (B8), and ance test (“nlme” R-package) functions. A Kruskal-Wallis blue (B2) bands from Sentinel-2 images. non-parametric test followed by the Dunn post-hoc test was performed to identify significant differences among SOC 1 3 65 Page 6 of 20 Wetlands (2023) 43:65 Digital Soil Mapping Where: Vij is the estimated potential value at location i, j ; S k ij are the fuzzy membership values at each location (similarity DSM combines different environmental covariates accord- values); vk is a typical value of soil type k. ing to the soil-forming factors defined by Jenny (1941) Random Forest (RF) is a machine learning method used with field soil point data in a mathematical model to to predict continuous and categorical variables from a large produce high-resolution soil maps (raster/pixel-based) collection of non-correlated decision trees (Breiman 2001). that are continuum and variable in the space (McBratney The construction of each tree depends on the classification et al. 2003; Moore et al. 1993; Zhu et al. 2001). Two DSM of the environmental covariates from the training dataset. To approaches, Expert Knowledge (EK) and Random Forest select the most appropriate covariates to enter the model, a (RF), were used to predict and map SOC content distribu- Variance Inflation Factor (VIF) diagnosis was carried out tion at 0-10 cm and 10-30 cm of soil depth. Both models and covariates with high multi-collinearity were identified were calibrated using a random set of 70% of the field (Alin 2010). VIF values close to zero represent null or scarce data (training dataset; same points for each model) and the collinearity, while values above 5 indicate high collinearity other 30% of the field data (validation dataset) were used between variables (Allam et al. 2020; Vu et al. 2015). In this for model assessment. study, the VIF assessment was done independently for each The EK model was based on the Soil Inference Engine depth, removing variables with VIF greater than 5 in order (SIE) system described by Shi et al. (2004). This model to reduce overfitting and improve models’ performances combines expert knowledge of soil scientists with fuzzy (Valbuena et al. 2017). logic by creating instance rules based on the spatial correla- The RF was performed with the randomForest R-package tion between soil properties and environmental covariates (Liaw and Wiener 2003), using the training dataset and ntree: (Ashtekar and Owens 2013; Heitkamp et al. 2020; Zhu 1997; 1,000. In order to identify which covariates mainly influ- Zhu et al. 2001). To perform the model, four units were ence the prediction of the SOC, we computed the variable defined according to the most dominant local land covers importance assessment mean decrease in accuracy (%Inc- (Fig. 1) and their respective instance rules were set to create MSE), which is the percentage increase in mean squared a Membership Function (MF) map (Zhu 1999, 1997) using errors (Zhou et al. 2019). The most important variables were the training dataset. The MFs were created using the FRBS defined as those that resulted in a significant reduction of and sets R-packages in which two-piece normal and sig- predictive power when removed from the model (Behrens moid functions were identified as more frequent descriptors et al. 2014). The most relevant covariates for predicting SOC of the soil-landscape relationship within every unit. After content based on the RF approach in the Casanare flooded individual environmental covariate’s fuzzification was done, savannas were established from the top-five most important a weighted average was computed between the MF within covariates used at 0-10 cm and 10-30 cm soil depth models. every instance in order to get the overall optimality value for the unit. The covariates used for mapping SOC with EK were DSM Assessment chosen according to the SOC-landscape relationship observed during the sampling campaign, prior knowledge Model performance for EK and RF approaches was assessed from the literature, information collected during the field based on three  statistical  indices  applied to the valida- trip, and the distribution of soil samples data across the tion dataset: Root Mean Squared Error (RMSE), Mean environmental covariates. Scatter plots were created for Absolute Error (MAE), and the coefficient of determina- each covariate using SOC content values as dependent tion (R2). Lower values of RMSE and MAE indicated better variable to understand and create SOC-landscape relation- model performance, while higher values of R 2 indicated that ships. Categorical covariates which did not provide clear the model explained a higher proportion of the variation. differentiation by grouping many sites into each category were discarded. SOC Stocks Spatial Distribution The MF maps were used to map the continuous SOC content from the weighted average of the typical property SOC stocks (t ha-1) at 0-30 cm soil depth were calculated values sampled in the field for each soil unit according to the for each sampling point by multiplying the SOC content and following equation (Zhu 1997): bulk density using equation (2) and subsequently compared ∑n S kvk under different geomorphic features and land cover condi- k=1 ij Vij = (1) tions through a mean comparison test (“Statistical Analysis” ∑n S k k=1 ij section). 1 3 Wetlands (2023) 43:65 Page 7 of 20 65 SOC Stock (ton∕ha) = SOC ⋅ BD ⋅ Dep ( e t th ⋅ 1 − u G∕100) 10 (2) R s l s ⋅ Where: SOC is the soil organic carbon content (%) ; Soil Properties Characterization ( ) BD is the soil b ulk density Mg∕m3 ; Depth is the depth of the soil layer (m) ; G is the volume (%) occupied with coarse Descriptive statistics of the measured physicochemical soil fragments (gravel or stones) and 10 is the conversion factor properties are presented in Table 1. The highest mean values to t/ha. Note that sampled soils in the area did not contain of SOC and Nt contents were found in the semiseasonal coarse fragments. savannas and the lowest values were found in the seasonal The spatial distribution of SOC stocks was modeled using savanna. Isotopic signatures of δ13C indicate the relative the best-performing DSM approach for SOC content at 0-10 proportion of carbon originating from trees and shrubs with cm and 10-30 cm and the respective BD maps developed a C3 photosynthetic pathway or grasses with a C4 photosyn- previously for the study area with a random forest model. thetic pathway. As expected, depleted values in the forest Basic statistics of predicted SOC stocks per land cover were stands showed carbon from predominantly C 3 sources, while calculated . both seasonal savanna and hyperseasonal savanna showed predominantly C4 sources. Depletion was intermediate in the moriche palm and semiseasonal savanna showing mixed sources of SOC in these systems. The measurement of δ15N provides an integrated index of the N cycle and all sites Table 1 Descriptive statistics of chemical and physical soil proper- cm depths. Values are mean ± SD and values for each parameter ties for the different land cover in the flooded savanna landscape of followed by the same letter are not significantly different from each Casanare. Values are presented for two soil layers 0-10 cm and 10-30 other (P < 0.05) SOC (%) δ13C (‰) Nt (%) δ15N (‰) C:N 0 – 10 cm   Semiseasonal savanna (n = 13) 8.71 ± 0.94a -18.36 ± 1.05b 0.86 ± 0.11a 2.42 ± 0.27bc 10.79 ± 0.55ab   Forest (n = 10) 2.09 ± 0.49c -26.18 ± 0.97c 0.17 ± 0.04c 3.68 ± 0.41ab 11.93 ± 0.77ab   Moriche palm (n = 5) 4.51 ± 0.93b -21.01 ± 0.60bc 0.45 ± 0.10b 2.05 ± 0.17c 10.27 ± 0.24b   Hyperseasonal savanna (n = 27) 3.73 ± 0.38b -15.62 ± 0.31a 0.31 ± 0.04b 3.37 ± 0.27b 12.70 ± 0.29a   Seasonal savanna (n = 25) 1.51 ± 0.20c -15.00 ± 0.33a 0.12 ± 0.01c 4.74 ± 0.26a 12.39 ± 0.36a 10 – 30 cm   Semiseasonal savanna 3.03 ± 0.52a -17.47 ± 1.04b 0.29 ± 0.06a 4.42 ± 0.39b 11.08 ± 0.63   Forest 1.50 ± 0.56bc -23.35 ± 1.37c 0.14 ± 0.06bc 5.16 ± 0.48ab 11.21 ± 0.83   Moriche palm 1.92 ± 0.29b -20.32 ± 0.61bc 0.19 ± 0.03ab 3.33 ± 0.43b 10.69 ± 0.51   Hyperseasonal savanna 1.84 ± 0.20ab -14.67 ± 0.36ab 0.15 ± 0.02b 5.10 ± 0.30b 12.66 ± 0.31   Seasonal savanna 0.81 ± 0.10c -13.66 ± 0.30a 0.06 ± 0.01c 6.98 ± 0.38a 12.50 ± 0.45 BD (g cm-3) Sand (%) Silt (%) Clay (%) SOC Stock (t ha-1) 0 – 10 cm   Semiseasonal savanna 0.74 ± 0.06c 30.8 ± 2.9 38.5 ± 3.1 30.7 ± 2.3 59.06 ± 3.62a   Forest 1.10 ± 0.09b 34.4 ± 5.2 42.1 ± 4.6 23.5 ± 2.3 19.42 ± 2.74c   Moriche palm 1.06 ± 0.12b 24.4 ± 1.8 43.9 ± 4.0 31.6 ± 4.7 44.08 ± 7.21b   Hyperseasonal savanna 1.16 ± 0.05b 28.2 ± 2.2 43.2 ± 1.9 28.6 ± 2.0 39.08 ± 2.78b   Seasonal savanna 1.42 ± 0.02a 42.8 ± 5.1 34.9 ± 4.0 22.3 ± 1.9 20.59 ± 2.47c 10 – 30 cm   Semiseasonal savanna 1.02 ± 0.04c 23.3 ± 3.0 40.1 ± 2.8 36.5 ± 3.9 58.54 ± 8.02a   Forest 1.25 ± 0.08b 33.2 ± 6.1 41.8 ± 4.2 24.9 ± 3.7 29.86 ± 6.70b   Moriche palm 1.19 ± 0.11b 28.6 ± 4.7 42.5 ± 3.7 28.9 ± 3.9 43.81 ± 5.38a   Hyperseasonal savanna 1.27 ± 0.03b 27.1 ± 2.4 45.0 ± 2.2 27.8 ± 1.7 44.13 ± 3.73a   Seasonal savanna 1.46 ± 0.02a 40.6 ± 5.6 35.3 ± 4.2 24.1 ± 2.7 23.23 ± 2.83b SOC Soil organic carbon; Nt Total nitrogen; δ13C Stable 13C isotope; δ15N Stable 15N isotope; C:N ratio Carbon to nitrogen ratio; BD Bulk density 1 3 65 Page 8 of 20 Wetlands (2023) 43:65 showed δ15N enrichment indicating an open or “leaky” N with talfs and continental dunes, while these last two dif- cycle (Craine et al. 2015). Enrichment was lower in the fered significantly at both soil depths. The highest SOC Moriche palm and semiseasonal savannas compared to the content of the study area (14.5%) was found at soil depth other land covers, suggesting that these ecosystems were 0-10 cm in the dip feature, while the highest value in talf, less leaky. Carbon to nitrogen ratios were similar across land rise, and dune for the same depth were 6.3%, 5.5%, 1.3%, covers. respectively. The largest mean reduction of SOC content Soil bulk density ranged from 0.74 to 1.42 g cm-3 (0-10 between soil depths occurred in the dip feature, where cm) and from 1.02 to 1.46 g c m-3 (10-30 cm). Higher values mean SOC content decreased from 7.2% at 0-10 cm to of BD were measured in the seasonal savanna, which are 2.82% at 10-30 cm, followed by talf (2.48% » 1.28%), rise characterized by a high sand content, while lower BD values (2.01% » 1.02%) and continental dune (0.61% » 0.38%). were found in the semiseasonal savannas as a result of higher Differences in SOC content were also observed among soil organic matter and clay contents. SOC stocks were sig- the different land covers. For soil depth 0-10 cm, SOC con- nificantly higher in clay soils and semiseasonal savannas tent was higher in semiseasonal savannas, with values rang- compared to seasonal savanna and forest. ing from 4.2% to 14.5%, and lower in seasonal savanna, with values ranging from 0.4% to 3.6%. For soil depth 10-30 cm SOC Distribution Across the Different Environmental the same pattern was observed, but the differences in SOC Conditions content among land covers were lower than those observed in the surface layer. Measured SOC content at soil depths of 0-10 cm and 10-30 Soil texture also played an important role in SOC accu- cm and its distribution across the environmental conditions mulation for this region. Although loamy soils presented of geomorphic features, soil textures, and land cover are pre- SOC content at the surface up to 14.5%, the average was sented in Fig. 3. In general, SOC content was around 1.64% 3.8% and this was much lower than the minimum SOC higher on average at soil depth 0-10 cm than at 10-30 cm in content (4.2%) found in clayey soils. Coarse textured soils all environmental conditions. had the lowest SOC content, with values inferior to 3.6% Among the different geomorphic features, dips had (0-10 cm). At 10-30 cm soil depth, higher values of SOC significantly higher SOC content than talfs, rises, and were found on clayey and loamy soils and lower content continental dunes. Rise features did not show differences in coarse soils. Fig. 3 Boxplot diagrams showing soil organic carbon (SOC) content blue crosses represent outliers, vertical black lines (whiskers) denote (%) distribution at soil depth 0-10 cm and 10-30 cm under different the maximum and minimum values (excluding outliers). Uppercase environmental conditions in the Casanare flooded savannas, Colom- and lowercase letters denote statistically significant differences at P < bia. Grey dots represent single values of sampled soils, the black 0.05 using the Dunn post-hoc test at 0-10 cm and 10-30 cm, respec- horizontal line in the middle of the box represents the median value, tively lower and upper boundaries indicate the 25th and 75th percentiles, 1 3 Wetlands (2023) 43:65 Page 9 of 20 65 Table 2 Model performance of soil organic carbon (SOC) content (%) calculation. For the RF model, an independent set of covariates was selected for each depth after the elimina- Soil depth Method RMSE MAE R2 tion of highly collinear covariates through the VIF assess- ment (see Table S1 in supplementary material). Figure 5 0–10 cm RF 1.76 1.48 0.49 shows the top-five most relevant covariates for predicting EK 1.81 1.45 0.42 SOC ranked by the %IncMSE for the RF model at each 10–30 cm RF 0.74 0.59 0.32 depth. The variable importance ranking identified local EK 0.70 0.53 0.36 land cover, normalized difference vegetation index for the EK Expert knowledge; RF Random forest; RMSE Root mean squared wet period – NDVI (wet), sand, and SAGA wetness index error; MAE Mean absolute error; R2 Correlation coefficient as common covariates for predicting SOC at both depths. Local land cover was found to be the most important covariate at 0-10 cm, contrary to the 10-30 cm depth in Assessment of DSM Approaches which it was placed at the bottom of the ranking. Normal- ized difference moisture index for the transition: wet-dry The performance assessment using the validation dataset period – NDMI (trans), in turn, was the most important showed small differences of RMSE, MAE, and R 2 between variable at 10-30 cm. NDVI (wet) and sand remained at the EK and RF approaches at both soil depths, with smaller the same position of the ranking at both depths, while the prediction errors (RMSE and MAE) at the 10-30 cm depth SAGA wetness index was more important for the 10-30 cm (Table 2). The RF presented lower values of RMSE and depth. Several covariates, such as local land cover, SAGA higher values of MAE and R 2 compared to the EK at 0-10 wetness index and sand, were utilized in both EK and RF cm, while at 10-30 cm EK presented lower values of RMSE approaches. Additionally, covariates regarding the vegeta- and MAE and higher R2 values. Figure 4 shows the scatter tion properties were the most frequently used in both DSM plots of the relationship between the dataset used for SOC approaches. (%) model validation using RF and EK approaches and the Maps of SOC content with a spatial resolution of 12 m measured data. developed using both the EK and RF methods at soil depths of 0-10 cm and 10-30 cm are shown in Fig. 6. In Digital Soil Mapping of SOC Content general, RF overestimated the minimum and underesti- mated the maximum SOC contents while EK underesti- Prior to implementing the EK and RF models, environ- mated both minimum and maximum (Table 3). Comparing mental covariates for each approach were selected. Covari- between minimum values of the predicted SOC content ates selected for implementing the EK at 0-10 cm and maps (EK and RF) and the measured data, we observed 10-30 cm were land cover, moisture stress index for the smaller differences for EK than RF at both depths, while dry period – MSI (dry), green-red vegetation index for RF had smaller differences for predicting the maximum the dry and transition: wet-dry period – GRVI (dry and values. At 0-10 cm, both approaches presented similar trans), relative slope position – RSP, clay, silt, sand, and meanw SOC content of 4.20% (RF) and 4.22% (EK), rep- the SAGA wetness Index, which is a 'Topographic Wet- resenting a higher m eanw prediction in comparison to the ness Index' (TWI), based on a modified catchment area measured meanw value (3.63%). At soil depth 10-30 cm, Fig. 4 Scatter plots showing the relationship between the dataset used for SOC (%) model valida- tion using RF (blue dots) and EK (red dots) approaches and data measured in the Casanare flooded savannas, Colombia for 0-10 cm and 10-30. Dots repre- sent the values for each variable while the lines represent the best fit for the correlation between them 1 3 65 Page 10 of 20 Wetlands (2023) 43:65 Fig. 5 Importance of environ- mental covariates used in the RF models for the study area according to the mean %Inc- MSE (percentage increase in mean square error). High values represent a high importance of the covariates in predicting soil organic carbon (SOC) content (%). BSI: Bare Soil Index, NDVI: Normalized Difference Vegetation Index, NDMI: Nor- malized Difference Moisture Index, wet: wet period -> May (2018) - October (2018), trans: transition dry-wet: November (2018) - December (2018) RF overestimated the m eanw value compared to the meas- that experienced longer periods of seasonal flooding (darker ured value (difference of 0.28%), while EK slightly under- areas) and the lowest values were predicted in the continental estimated the meanw values (difference of 0.06%). dunes with sandy soils and no seasonal flooding. However, there In general, the two approaches presented a similar pattern are differences in the predicted SOC values in some areas such regarding the spatial distribution of SOC content (Fig. 6). The as the large gallery forests located in the floodplains of the Ari- highest values of predicted SOC content were found at locations poro and Meta Rivers where RF estimated higher values of SOC Fig. 6 Spatial distribution of SOC (%) using two different DSM approaches (EK and RF) in the Casanare flooded savannas, Colombia at depth intervals: 0-10 cm and 10-30 cm 1 3 Wetlands (2023) 43:65 Page 11 of 20 65 Table 3 Descriptive statistics of SOC (%) predicted using EK and to semiseasonal savanna land cover and dip geomorphic fea- RF approaches and measured data collected in the Casanare flooded tures, where lower values of BD and higher values of SOC con- savannas, Colombia for 0-10 cm and 10-30 cm tents were identified (Table 1). Lower SOC stock values were Depth Method Min Max Meanw sdw attributed to seasonal savanna land cover and continental dunes, SOC (%) (0-10 cm) Field Data 0.41 14.52 3.63 1.90 which were characterized by lower SOC contents and high BD. Map (RF) 0.85 11.92 4.20 0.85 Due to the slightly better capacity of EK to describe and Map (EK) 0.30 11.30 4.22 0.51 map the spatial distribution of SOC content for the study area, SOC (%) (10-30 cm) Field Data 0.23 7.51 1.65 1.08 the map of SOC stocks at 0-30 cm soil depth was developed Map (RF) 0.48 5.09 1.93 0.45 using this approach and selected covariates described in “Dig- Map (EK) 0.10 4.81 1.59 0.19 ital Soil Mapping” section. The corresponding model perfor- mance parameters were RMSE: 26 -, MAE: 23, and R2: 0.69. Meanw: weighted arithmetic mean and sdw: weighted arithmetic (Circles indicate the larger areas of higher SOC stock standard deviation based on the proportion of data per land cover accumulation) Predicted SOC stocks ranged from 6 to 210 t h a-1, with for both soil depths (0-10 cm and 10-30 cm). Predicted values of average stocks of 83.13 ± 24.32 t C ha-1 (Fig. 8). Simi- SOC in the hyperseasonal savannas tend to be lower for the RF lar to SOC content, high (150-180 t h a-1) and very high than for the EK model at 0-10 cm; while at 10-30 cm, predicted (> 180 t ha-1) SOC stocks were found in areas charac- SOC values across the hyperseasonal savannas are similar. terized by their proximity to water bodies, their location in floodplains with clay-loam soils, and the presence of SOC Stock Distribution long-term flooded conditions (circles indicate the larger areas of higher SOC stock accumulation). Medium-high SOC stocks over the entire 30 cm profile determined using field values (120-150 t h a-1) were found adjacent to areas with data and scaling by geomorphic feature and land cover type are higher SOC stocks. The hyperseasonal savanna land cover, summarized in Fig. 7. Comparisons between measured SOC which occupies the larger proportion of the study area, stocks across the different geomorphic features and land covers is characterized as a mix mainly between medium-high showed significant differences (P < 0.05) between dip (111 ± and medium SOC stocks (90-120 t h a-1) in dip features 30 t ha-1) and the other features, while no significant differences and the medium and medium-low (60-90 t ha-1) stocks were identified between talf (63 ± 24 t ha-1), rise (49 ± 28 t in well-drained talf features. Low values (30-60 t ha-1) of ha-1), and continental dune (20 ± 8 t h a-1). On an ecosystem SOC stocks were largely predicted in riversides along the basis, semiseasonal savannas (118 ± 34 t ha-1) differed sig- Ariporo River and the Meta River as well as along other nificantly from forests (48 ± 26 t ha-1) and seasonal savannas rivers, where coarser soil textures predominate. Very low (44 ± 25 t h a-1) and did not show significant differences with values (< 30 t ha-1), in turn, were characteristic of areas Moriche palm (88 ± 28 t h a-1) or hyperseasonal savanna (83 ± where continental dunes (northeasterly elongated geomor- 29 t ha-1). Overall, higher mean measured values were attributed phic features) and sandy patches predominate. Fig. 7 Boxplot diagrams showing SOC stocks (t h a-1) dis- tribution at soil depth 0-30 cm under different environmental conditions. Grey dots represent single values of sampled soils, the black horizontal line in the middle of the box represents the median value, lower and upper boundaries indicate the 25th and 75th percentiles, blue crosses represent outliers, vertical black lines (whiskers) denote the maximum and minimum values (excluding outliers). Lowercase letters denote statistically sig- nificant differences at P < 0.05 using the Dunn post-hoc test 1 3 65 Page 12 of 20 Wetlands (2023) 43:65 Fig. 8 Map of SOC stocks (t ha-1) predicted from EK model at 0-30 cm soil depth in the Casanare flooded savannas, Colombia. Using the land cover map (Fig. 1) as a mask, descriptive sta- is often associated with topography that concentrates water tistics of the predicted SOC stock map per land cover were cal- and creates anaerobic conditions that slow decomposition culated (Table 4). We found that the highest mean predicted SOC (Gumbricht et al. 2017; Paul 2016; Sjögersten et al. 2014). stocks values were in the semiseasonal savannas, followed by the Carbon accumulation is also regulated by the types and hyperseasonal savanna, seasonal savanna, and forest, which com- quality of vegetation, particularly the proportion of plant- bines riparian forest and Moriche palm. Total soil organic carbon derived C and N that is incorporated into soil organic matter, stock predicted for the first 30 cm of soil depth in the entire study and by soil matrix interactions that control its stabilization area (Fig. 8) was equivalent to 55.07 Mt of C. Hyperseasonal (Cotrufo et al. 2013; Paul 2016). Clay and silt fractions play savannas contributed 65% (42.73 Mt) of the total stocks, while important roles in stabilizing soil organic matter through forest contributed with 22% (3.55 Mt), semiseasonal savanna both physical and chemical protection. We found strong with 7% (6.79 Mt), and seasonal savanna with 6% (2.00 Mt). associations between geomorphology, plant communities, and soil texture that affected the SOC concentrations in the flooded Llanos del Orinoco and in the surrounding ecosys- Discussion tems. There are few data from South American floodplain wetlands for comparison and available data do not represent SOC Characterization systematic sampling, but we found that mean SOC contents are higher in the Casanare flooded savannas compared to the Carbon accumulation in tropical wetlands is typically the levels reported for flooded savannas in the Brazilian Panta- result of small differences between input and outputs, and nal (Vega et al. 2014), the Brazilian Cerrado (Wantzen et al. Table 4 Basic statistics of Land cover Range SOC stock Mean ± SD SOC stock Total SOC Landscape predicted SOC stocks (t h a-1) (t ha-1) (t ha-1) stock (Mt) area occupied per land cover within the study (%) area Water - - - 0.5 Seasonal savanna 6.27 – 119.83 47.13 ± 14.16 2.00 6.39 Semiseasonal savanna 25.68 – 210.52 129.59 ± 11.43 6.79 7.89 Forest 22.88 – 108.38 42.64 ± 8.19 3.55 12.51 Hyperseasonal savanna 36.61 – 138.02 88.39 ± 10.33 42.73 72.71 Total 55.07 100 1 3 Wetlands (2023) 43:65 Page 13 of 20 65 2012) or the Bolivian Llanos do Moxos (Boixadera 2003), independent trees created in the model (Callens et al. 2020), and they were similar to flooded savannas in Roraima, Brazil and generates an approximation to the field values based (Barbosa et al. 2012). on averages of most similar areas, which limits its ability We found high SOC contents in locations characterized to describe local extreme values (Szabó et al. 2019). This by clayey soils with high grass productivity (semiseasonal situation probably caused the misprediction along the rivers, savannas), and depressions (dip features). SOC accumula- which may have been influenced by other areas with simi- tion in these areas have been favored by frequent and exces- lar environmental conditions. RF has been used effectively sive precipitation, high concentrations of sediments, and in several studies for estimating SOC at landscape scale long periods of water accumulation (Barthelmes and Joos- (Grimm et al. 2008; Wang et al. 2020; Zeraatpisheh et al. ten 2018; Sarmiento 1984; Mora-Fernández et al. 2015). 2019; Zhou et al. 2019). Its prediction power has been shown Sediment transportation processes occurred during the Ter- to outperform other DSM approaches, specifically in larger tiary and Quaternary deposited coarser materials, as a result study areas with a diversity of landscape features (Lamich- of the erosion of the Andes mountains, along most of the hane et al. 2019), capturing most of the variance of the soil river edges (Goseen 1964), forming banks that do not store characteristics (Sindayihebura et al. 2017). Other advantages much SOC. We found low SOC contents in riparian areas, of RF include the modelling of high dimensional data, which in continental dunes, and in sandy deposits that are charac- would be difficult to afford with EK; the estimation of inter- terized by low water availability, coarser soil textures, and nal errors and the variable importance assessment; as well low-density vegetation (seasonal savanna). In riparian areas as a lower effort for data preparation and modelling time located specifically in silty-sandy banks, within the first (Breiman 2001; Camera et al. 2017; Grimm et al. 2008). 30 cm prevail dry conditions during a long part of the year, However, Ashtekar and Owens (2013) mentioned that a dis- which affects the accumulation of SOC. Aeolian erosion by advantage of this method could be that it only works with the Trade winds has acted as an essential factor in shaping soil-landscape relationships at specific locations and does the regional topography, particularly in the formation of con- not take into consideration extensive differentiation across tinental dunes and sandy deposits (FAO 1965). Intermedi- the landscape. Additionally, Hengl et al. (2015) stated that ate SOC contents were associated with Moriche palm and RF is effective only within the range of values used in the hyperseasonal savanna vegetation and loamy textured soils. training data, which provides irregular predictions for obser- vations outside of the range (Terra 2017). DSM Models Assessment EK, in turn, is based mainly on qualitative and quantita- tive criteria provided by the knowledge of an expert, that We compared two different DSM approaches (EK and RF) takes into account soil-landscape relationships for a specific for predicting SOC content at 0-10 cm and 10-30 cm. The region based on experience, analysis from sampling data, models performed comparably to map SOC variability in the and knowledge of interactions between covariates that influ- study area and gave a similar spatial prediction pattern. The ence the property distribution (Ashtekar and Owens 2013; models both predicted high SOC contents in dip features, Menezes et al. 2013; Ngunjiri et al. 2019). The use of EK where soil textures are finer, flooded conditions last longer, led to a better control and implementation of inference rules and vegetation has higher biomass production rates, and that explain soil-landscape relationships, even with limited low SOC values on the continental dunes with coarser soils, field data (Menezes et al. 2018). EK approaches have been low-density vegetation, and no seasonal flooding. However, satisfactorily applied by other authors (Akumu et al. 2015; with validation against field data, we found that the spatial Ashtekar et al. 2014; da Silva et al. 2016; Minai et al. 2021; predictions obtained with the EK represented field condi- Tsakiridis et al. 2019; Zhu et al. 2010). As evidence of its tions better. A key major difference between the modelling good performance, this approach has been widely used for approaches was associated with the riparian areas along Ari- the development of soil surveys by the US Department of poro and Meta rivers. The RF model predicted much higher Agriculture (USDA) (Zhu et al. 2001). On the other hand, values for these areas than the EK model, which showed EK is time consuming in terms of data selection, analysis, predictions that were more consistent with field observa- and processing. It requires knowledge about the mapping tions. Additionally, EK showed lower prediction errors than area and soil-landscape relationships, which may be limited RF. Both models were sensitive to the presence of outliers. for larger areas. RF underestimated the highest values and overestimated the lowest, while EK underestimated the highest and lowest val- Drivers of SOC Spatial Variability ues and gave overall predictions around the mean values. Lower performance in RF was attributed intrinsically Similar to other tropical landscapes studies: (Diwediga to the methodology. The algorithm estimates the values et al. 2017; Hamzehpour et al. 2019; Yang et al. 2007), in for each pixel, based on an average of the results of all the the Casanare flooded savannas, the vegetation coverage, 1 3 65 Page 14 of 20 Wetlands (2023) 43:65 topography, and soil texture were found to be the most influ- to generate better SOC-related models, especially for grassland ential factors of the SOC distribution. The combination of ecosystems with similar conditions around the world; improve these 3 features leads to different environmental conditions those already developed, generate effective policies for SOC that consequently influence the spatial distribution of SOC, sequestration that contribute to mitigate climate change; as as evidenced in this study. well as to improve soil health, biodiversity conservation and Vegetation coverage not only controls the inputs of organic ecosystem services (Gray et al. 2015). Similar landscapes material but also protects soils from erosion and conserves soil in South America lack carbon studies and their importance moisture (Lamichhane et al. 2019; Weil and Brady 2017). Sev- in terms of carbon accumulations is often understated. The eral studies have represented the vegetation coverage in spatial application of these study’s methodolgies in other areas would models across the landscape employing categorical land cover help improve our understanding of their importance and target maps or continuous vegetation indices derived from remote protection actions in priority areas. sensing (Gray et al. 2015; Minai et al. 2021; Mondal et al. 2017; Wan et al. 2019; Wang et al. 2020). Land cover maps and vegetation indices help to classify vegetation in terms of SOC Stocks Distribution types, biomass production, residues, and soil textures based on changes in leaf colors and leaf density (Motohka et al. 2010). The Casanare flooded savannas landscape holds important Specific combinations of these characteristics delimit areas SOC stocks (55.07 Mt of C), which may be affected by land from higher carbon fixation such as the semiseasonal savanna use changes and unsustainable management practices. While to lower carbon fixation such as the seasonal savanna. Veg- the pan-tropical wetlands map of Gumbricht et al. (2017) etation indices also show where to expect more plant water predicted the presence of peat deposits in this landscape, stress and less soil moisture content (Elhag and Bahrawi 2017), we found only mineral soils with high SOC contents. Nev- being relevant not only for soil moisture variability but also ertheless, the study area represents only around 4.43% of for SOC spatial modeling (Welikhe et al. 2017). Topography the flooded Llanos del Orinoco landscape, which accounts plays an imperative role in SOC storage across the landscape approximately 15 Mha (Hamilton et al. 2002). The numbers (Patton et al. 2019; Wiesmeier et al. 2019). It widely describes from our study would suggest that the total carbon stock water flow paths, water accumulation, soil erosion, and sedi- across this vast ecosystem is on the order of 109 tons of C, mentation, and consequently influences the soil moisture, with about 20 to 30% stored in flooded savanna ecosystems. microbial activity, and SOC accumulation (Sørensen et al. The expansion of agriculture and introduction of new 2006; Yun et al. 2019). Topographic covariates that indicate land use practices could result in significant wetland loss, potential moist/dry areas such as topographical wetness index contributing to the already high levels of loss in South improves SOC predictions (Adhikari et al. 2014; Lamichhane America and increasing CO2 emissions into the atmos- et al. 2019; Minasny et al. 2013; Siewert 2017; Taghizadeh- phere (Guo and Gifford 2002; Laganière et al. 2010; Pow- Mehrjardi et al. 2017; Wiesmeier et al. 2013). ers 2005; Zhou et al. 2019). Considered as the last agri- Additionally, as identified in several studies (Funes et al. cultural frontier in Colombia, recent developments in the 2019; Hounkpatin et al. 2018; Yang et al. 2007; Zinn et al. Colombian Llanos have seen the expansion of industrial 2005), the results also showed that texture was also key in crops (African oil palm, sugar cane, rice, and maize) as explaining the distribution of SOC. Coarser soil textures were well as increases in grazing area, which have replaced gal- correlated negatively with SOC, while finer textures showed lery forests, wetlands, and native vegetation (Usma and a positive relationship. Clay soils have a greater capacity to Trujillo 2011). This alters both the hydrological regime preserve or protect organic matter from microbial attack as and the composition and structure of plant communities well as to stabilize the organic matter through mineral–organic across the landscape (Lasso et al. 2010). In particular, matter binding. By promoting soil aggregations, silt and clay the wetlands in this landscape have been transformed into contents improve SOC protection (Six et al. 2002; Xie et al. minimally productive rice fields, due to the expansion of 2021), enhance soil moisture (Augustin and Cihacek 2016; agriculture in the region. Recurrent burning, a controver- Yang et al. 2016), water holding capacity, and plant available sial land preparation technique to enhance pasture produc- water. These conditions reduce or avoid SOC losses by soil tivity, has devastated native palm ecosystems and forests, erosion, enhance plant biomass production, and, consequently, which have been colonized by grasses. In addition, pol- conserve or accumulate SOC. On the other hand, sandy soils lution by fertilizers and the demand for water from crops do not physically or chemically protect SOC, they favor nutri- threaten the natural state of flooded ecosystems, Moriches, ent losses by leaching affecting plant density, both of which and other aquatic communities in the savannas, putting resulting in low SOC levels (Hamzehpour et al. 2019). their continued existence at risk (Lasso et al. 2014). As seen in this section, the identification and understanding Land use change in the tropics and subtropics is glob- of factors that control SOC, as gained from this study, is crucial ally the second largest C O2 source to the atmosphere, after 1 3 Wetlands (2023) 43:65 Page 15 of 20 65 combustion of fossil fuels (Barrezueta-Unda et al. 2019; Lal 2011), which will contribute to avoiding unwanted impacts 2015; Van der Werf et al. 2009). The use of wetland soils for on society and maintaining the environmental integrity agriculture, native vegetation conversion, deforestation, and of South American landscapes. Much attention has been intensification of fires have largely contributed to the GHG paid recently to swamp forests, but other wetlands in South emissions (Margono et al. 2014). In Peru, Moriche palm eco- America are in similar situations. Therefore, it is important systems are under increasing human pressure (Hergoualc’h to increase the understanding of SOC sequestration and et al. 2017; Horn et al. 2012) and we see evidence of this in interactions, provide the appropriate assistance to farm- our study area. ers about land management practices and the subsequent According to Mora-Fernández et al. (2015), the depart- creation and application of conservation strategies such as ment of Casanare is one of the most transformed and less silvopastoral or agroforest systems (Amézquita 1999; Jadán protected in the country (Romero Duque et al. 2018; Usma et al. 2015), which have demonstrated a positive balance and Trujillo 2011). Land use changes, mainly caused by the of SOC storage (Silva-parra 2018), and the definition of accelerated expansion of agricultural activities, the increase global GHC emissions mitigation strategies. of wildfires and the mechanical modification of the water regimes, have been the most harmful and stressful human activities for the ecosystems of the floodplain savannas (Mora-Fernández and Peñuela-Recio 2013). These changes Conclusions cause the loss of its natural vegetation, the variation of the water dynamics, and the alteration of carbon cycle (Ostle In this study, the status of SOC content and stocks at 0-30 cm et al. 2009). Although the study area represents only 15% soil depth in the Casanare flooded savannas was assessed. of the total area of the department, the intensive pressure Results demonstrated the presence of carbon levels lower of human development could result in the reduction of its than expected, but high relative to other South American SOC stocks and the release of important amounts of GHG flooded savannas. The spatial variability and driving fac- into the atmosphere (Ostle et al. 2009; Schreier et al. 1994; tors of SOC content were evaluated through the implemen- Sharma and Sharma 2022). At regional level, the impact of tation of two DSM approaches: Expert Knowledge (EK) a large-scale land use conversions of the flooded Llanos del and Random Forest (RF). Although both DSM approaches Orinoco ecosystem area (15 Mha) could transform this area performed very well, EK was considered slightly superior in a future source of important global emissions if correct to predict SOC across the study area. Vegetation coverage, decisions are not taken regarding the land management of topography and soil texture were the most relevant factors in the region. explaining the spatial variability of SOC content, attributing Acting on the commitments of the Paris Agreement, higher SOC contents to dip features, finer soil textures, and several countries have updated their GHG emissions miti- high-density vegetation. gation strategies, aiming to reduce GHG emissions and Our study showed very high carbon stocks in the study keeping a rise in global average temperature below 2 °C area and a high potential of GHG emissions associated (Carvajal-Agudelo and Andrade 2020). Given the impor- with intensive agricultural development in the region, in tance of wetlands for local climate resilience and the sig- addition to loss of large areas of wetlands. SOC stock nificant carbon impacts, the protection and restoration of quantification indicated the importance of maintaining wetlands should be an integral part of national actions or even improving the carbon sequestration in the area to address the climate change problem and support low to avoid carbon losses and, consequently, increasing CO2 emissions, climate resilient development. South America is emissions to the atmosphere. Inappropriate land and soil considered as the wettest continent in the world, counting management of these ecosystems could hamper the efforts with 22% of its area covered by wetlands (Junk et al. 2013). of Colombia and the global community on carbon seques- Nevertheless, conservation of these important ecosystems tration and reduction of C O2 emissions. Additionally, is lagging in comparison to other regions due to the slow these findings reinforce the view that wetlands, as the development of definitions, delineations, and classifica- Casanare flooded savannas and similar areas elsewhere tions of their wetlands (Prahalad and Kriwoken 2010). The in South America, are relevant carbon storage ecosys- national inventories are not updated, and the ecological tems that must be considered of worldwide interest and importance of these ecosystems is widely undocumented. protection. Thus, the development of accurate SOC and wetland maps plays an essential role, since an adequate inventory and Supplementary Information The online version contains supplemen- tary material available at https://d oi.o rg/1 0.1 007/s 13157-0 23-0 1705-3. monitoring of these carbon-rich ecosystems will lead to better ecological characterization, quantification of SOC Acknowledgments We acknowledge funding from the CGIAR stocks and assessment of GHG emissions (Page et  al. Research Program on Forests, Trees, and Agroforestry (CRP-FTA) 1 3 65 Page 16 of 20 Wetlands (2023) 43:65 with financial support from the CGIAR Fund Donors. We are grate- Alin A (2010) Multicollinearity. Wiley Interdisciplinary Reviews: ful for the financial support from the Bezos Earth Fund project Computational Statistics. 2:370–374. https:// doi. org/ 10. 1002/ Using genetic diversity to capture carbon through deep root sys- wics. 84 tems in tropical soils. We also wish to extend our thanks to Cesar Allam AS, Bassioni HA, Kamel W, Ayoub M (2020) Estimating the Botero who was fundamental during the soil sampling and the standardized regression coefficients of design variables in day- farmers who granted us access to their farms. lighting and energy performance of buildings in the face of mul- ticollinearity. Solar Energy 211:1184–1193. https:// doi. org/ 10. Author Contributions Javier M. Martín-López: Conceptualiza- 1016/j. solen er. 2020. 10. 043 tion, Field work, Methodology, Software, Validation, Formal analy- Amézquita E (1999) Propiedades físicas de los suelos de los Llanos sis, Investigation, Resources, Data Curation, Writing - Original Draft, Orientales y sus requerimientos de labranza. Palmas 20:28–30 Visualization. Ashtekar JM, Owens PR (2013) Remembering knowledge: an expert Louis V. Verchot: Conceptualization, Methodology, Investigation, knowledge based approach to digital soil mapping. Soil Horizons Formal analysis, Writing - Review & Editing, Supervision, Funding 54:0. https:// doi. org/ 10. 2136/ sh13- 01- 0007 acquisition. Ashtekar JM, Owens PR, Brown RA, Winzeler HE, Dorantes M, Libo- Christopher Martius: Writing - Review & Editing, Supervision, hova Z, Dasilva M, Castro A (2014) Digital mapping of soil Funding acquisition. properties and associated uncertainties in the llanos orientales, Mayesse da Silva: Conceptualization, Methodology, Validation, south america. GlobalSoilMap: Basis of the Global Spatial Soil Investigation, Writing - Review & Editing, Supervision, Project coor- Information System - Proceedings of the 1st GlobalSoilMap dination, Funding acquisition. Conference, 367–372. https:// doi. org/ 10. 1201/ b16500- 67 Augustin C, Cihacek LJ (2016) Relationships between soil carbon and Funding This work was supported by the CGIAR Research Program soil texture in the Northern Great Plains. Soil Science 181:386– on Forests, Trees, and Agroforestry (CRP-FTA) with financial sup- 392. https:// doi. org/ 10. 1097/ SS. 00000 00000 000173 port from the CGIAR Fund Donors and the Bezos Earth Fund project Barbosa RI, Silva dos Santos JR, Souza da Cunha M, Pimentel TP, Using genetic diversity to capture carbon through deep root systems Fearnside PM (2012) Root biomass, root:shoot ratio and below- in tropical soils. ground carbon stocks in the open savannahs of Roraima, Brazil- ian Amazonia. Australian Journal of Botany 60:405. https:// doi. Data Availability The datasets generated during and/or analyzed dur- org/ 10. 1071/ BT113 12 ing the current study are available from the corresponding authors on Barreto JS, Armenteras D (2020) Open data and machine learning reasonable request. to model the occurrence of fire in the ecoregion of “Llanos Colombo–Venezolanos.” Remote Sensing 12:3921. https:// doi. Declarations org/ 10. 3390/ rs122 33921 Barrezueta-Unda S, Velepucha-Cuenca K, Hurtado-Flores L, Jara- Competing Interests The authors have no relevant financial or non- millo-Aguilar E (2019) Soil properties and storage of organic financial interests to disclose. carbon in the land use pasture and forest. Revista de Ciencias Agrícolas 36:32–45 Barthelmes A, Joosten H (2018) Lineamientos para inventarios de Open Access This article is licensed under a Creative Commons Attri- turberas tropicales a fin de facilitar su designación como sitios bution 4.0 International License, which permits use, sharing, adapta- Ramsar. Nota informativa no9. Secretaría de La Convención de tion, distribution and reproduction in any medium or format, as long Ramsar Gland, Sui as you give appropriate credit to the original author(s) and the source, Behrens T, Schmidt K, Ramirez-Lopez L, Gallant J, Zhu A-X, provide a link to the Creative Commons licence, and indicate if changes Scholten T (2014) Hyper-scale digital soil mapping and soil were made. The images or other third party material in this article are formation analysis. Geoderma 213:578–588. https:// doi. org/ included in the article's Creative Commons licence, unless indicated 10. 1016/j. geode rma. 2013. 07. 031 otherwise in a credit line to the material. If material is not included in Blake GR (1965) Bulk Density. Methods of Soil Analysis C.A. - the article's Creative Commons licence and your intended use is not University 0/ Minnesota St. Paul, Minnesota 17. https:// doi. permitted by statutory regulation or exceeds the permitted use, you will org/ 10. 2134/ agron monog r9.1. c30 need to obtain permission directly from the copyright holder. To view a Boixadera J (2003) Hydromorphic and clay-related processes in copy of this licence, visit http://c reati vecom mons.o rg/l icens es/b y/4.0 /. soils from the Llanos de Moxos (Northern Bolivia). CATENA 54:403–424. https://d oi.o rg/1 0.1 016/S 0341-8 162(03)0 0134-6 Bouyoucos GJ (1936) Directions for making mechanical analyses of soils by the hydrometer method. Soil Science 42:225–230. References https:// doi. org/ 10. 1097/ 00010 694- 19360 9000- 00007 Breiman L (2001) Random forests. Machine Learning 45:5–32. Adhikari Kabindra, Alfred E. Hartemink, Minasny Budiman, Bou https:// doi. org/ 10. 1023/A: 10109 33404 324 Kheir Rania, Greve Mette B., Greve Mogens H. (2014) Digital Buriticá Mejia N (2016) Sabanas inundables de la orinoquía colom- Mapping of Soil Organic Carbon Contents andStocks in Den- biana -documento resumen-. Humboldt: 1–19 mark. PLoSONE 9(8). https:// doi. org/ 10. 1371/ journ al. pone. Cabrera-Amaya DM, Giraldo-Kalil LJ, Rivera-Díaz O, Castro-Lima 01055 19 F (2020) Riqueza, composición y distribución de las plantas Aitken R (1977) Wilderness areas in Scotland. Ph. D. thesis, University vasculares en sabanas y bosques ribereños de la cuenca baja of Aberdeen del río Pauto (Casanare-Colombia). Revista de La Academia Akumu CE, Johnson JA, Etheridge D, Uhlig P, Woods M, Pitt DG, Colombiana de Ciencias Exactas, Físicas y Naturales 44:1018– McMurray S (2015) GIS-fuzzy logic based approach in modeling 1032. https:// doi. org/ 10. 18257/ racce fyn. 1188 soil texture: using parts of the Clay Belt and Hornepayne region Chmura GL, Anisfeld SC, Cahoon DR, Lynch JC (2003) Global in Ontario Canada as a case study. Geoderma 239–240:13–24. carbon sequestration in tidal, saline wetland soils. Global https:// doi. org/ 10. 1016/j. geode rma. 2014. 09. 021 Biogeochemical. Cycles 17(4, 1111). https:// doi. org/ 10. 1029/ 2002G B0019 17 1 3 Wetlands (2023) 43:65 Page 17 of 20 65 Callens A, Morichon D, Abadie S, Delpey M, Liquet B (2020) Using Goseen D (1964) Geomorfología de los Llanos Orientales. Revista de La Random forest and Gradient boosting trees to improve wave Academia de Las Ciencias Exactas, Físicas y Naturales 12:129–140 forecast at a specific location. Applied Ocean Research 104. Gray JM, Bishop TFAA, Wilson BR (2015) Factors controlling soil https:// doi. org/ 10. 1016/j. apor. 2020. 102339 organic carbon stocks with depth in Eastern Australia. Soil Sci- Camera C, Zomeni Z, Noller JS, Zissimos AM, Christoforou IC, ence Society of America Journal 79:1741–1751. https://d oi.o rg/ Bruggeman A (2017) A high resolution map of soil types and 10. 2136/ sssaj 2015. 06. 0224 physical properties for Cyprus: a digital soil mapping optimi- Grimm R, Behrens T, Märker M, Elsenbeer H (2008) Soil organic zation. Geoderma 285:35–49. https:// doi. org/ 10. 1016/j. geode carbon concentrations and stocks on Barro Colorado Island - rma. 2016. 09. 019 Digital soil mapping using Random Forests analysis. Geoderma Carvajal-Agudelo BN, Andrade HJ (2020) Carbon capture regarding 146:102–113. https:// doi. org/ 10. 1016/j. geode rma. 2008. 05. 008 biomass from rural land use systems near the municipality of Gumbricht T, Roman-Cuesta RM, Verchot L, Herold M, Wittmann F, Yopal, Casanare, Colombia. Captura de Carbono En Biomasa Householder E, Herold N, Murdiyarso D (2017) An expert sys- de Sistemas de Uso Del Suelo, Municipio de Yopal, Casanare, tem model for mapping tropical wetlands and peatlands reveals Colombia. 24:13–22 South America as the largest contributor. Global Change Biology Castillo-Figueroa D, Martínez-Medina D, Rodríguez-Posada ME, 23:3581–3599. https:// doi. org/ 10. 1111/ gcb. 13689 Bernal-Vergara S (2019) Structural differences in mammal Guo LB, Gifford RM (2002) Soil carbon stocks and land use change: assemblages between savanna ecosystems of the Colombian a meta analysis. Global Change Biology 8:345–360. https:// doi. Llanos. Papéis Avulsos de Zoologia 59:11. https:// doi. org/ 10. org/ 10. 1046/j. 1354- 1013. 2002. 00486.x 11606/ 1807- 0205/ 2019. 59. 14 Hamilton SK, Sippel SJ, Melack JM (2004) Seasonal inundation pat- Costanza R, de Groot R, Sutton P, van der Ploeg S, Anderson SJ, terns in two large savanna floodplains of South America: the Kubiszewski I, Farber S, Turner RK (2014) Changes in the global Llanos de Moxos(Bolivia) and the Llanos del Orinoco(Venezuela value of ecosystem services. Global Environmental Change and Colombia). Hydrological Processes 18:2103–2116. https:// 26:152–158. https:// doi. org/ 10. 1016/j. gloen vcha. 2014. 04. 002 doi. org/ 10. 1002/ hyp. 5559 Cotrufo MF, Wallenstein MD, Boot CM, Denef K, Paul E (2013) The Hamilton SK, Sippel SJ, Melack JM (2002) Comparison of inunda- Microbial Efficiency-Matrix Stabilization (MEMS) framework tion patterns among major South American floodplains. Journal integrates plant litter decomposition with soil organic matter sta- of Geophysical Research Atmospheres 107:LBA 5-1-LBA 5-14. bilization: do labile plant inputs form stable soil organic matter? https:// doi. org/ 10. 1029/ 2000J D0003 06 Global Change Biology 19:988–995. https:// doi. org/ 10. 1111/ Hamzehpour N, Shafizadeh-Moghadam H, Valavi R (2019) Exploring gcb. 12113 the driving forces and digital mapping of soil organic carbon Craine JM, Brookshire ENJ, Cramer MD, Hasselquist NJ, Koba K, using remote sensing and soil texture. Catena 182. https:// doi. Marin-Spiotta E, Wang L (2015) Ecological interpretations of org/ 10. 1016/j. catena. 2019. 104141 nitrogen isotope ratios of terrestrial plants and soils. Plant and Haskins DM, Correll CS, Foster RA, Chatoian JM, Fincher JM, Soil 396:1–26. https:// doi. org/ 10. 1007/ s11104- 015- 2542-1 Strenger JM, Keys Jr JE, Maxwell JR, King T (1999) A geomor- da Silva MA, Silva MLN, Owens PR, Curi N, Oliveira AH, Candido phic classification system. Abstracts with Programs - Geological BM (2016) Predicting runoff risks by digital soil mapping. Society of America, 31:254 Revista Brasileira de Ciencia Do Solo 40:1–13. https://d oi.o rg/ Heitkamp F, Ahrends B, Evers J, Steinicke C, Meesenburg H (2020) 10. 1590/ 18069 657rb cs201 50353 Inference of forest soil nutrient regimes by integrating soil chem- de Menezes MD, Silva SHG, de Mello CR, Owens PR, Curi N (2018) istry with fuzzy-logic: regionwide application for stakeholders Knowledge-based digital soil mapping for predicting soil prop- of Hesse, Germany. Geoderma Regional 23:e00340. https:// doi. erties in two representative watersheds. Scientia Agricola org/ 10. 1016/j. geodrs. 2020. e00340 75:144–153. https:// doi. org/ 10. 1590/ 1678- 992x- 2016- 0097 Hengl T, Heuvelink GBM, Kempen B, Leenaars JGB, Walsh MG, Shep- de Menezes MD, Silva SHG, Owens PR, Curi N (2013) Digital soil herd KD, Sila A, MacMillan RA, De Jesus JM, Tamene L, Tondoh mapping approach based on fuzzy logic and field expert knowl- JE (2015) Mapping soil properties of Africa at 250 m resolution: edge. Ciência e Agrotecnologia 37:287–298. https://d oi.o rg/1 0. random forests significantly improve current predictions. PLoS 1590/ S1413- 70542 01300 04000 01 ONE 10:1–26. https:// doi. org/ 10. 1371/ journ al. pone. 01258 14 Diwediga B, Bao Q, Agodzo S, Wala K (2017) Potential storages and Hergoualc’h K, Gutiérrez-Vélez VH, Menton M, Verchot LV (2017) drivers of soil organic carbon and total nitrogen across river Characterizing degradation of palm swamp peatlands from space basin landscape : the case of Mo river basin ( Togo ) in West and on the ground: an exploratory study in the Peruvian Amazon. Africa. Ecological Engineering 99:298–309. https:// doi. org/ Forest Ecology and Management 393:63–73. https:// doi. org/ 10. 10. 1016/j. ecole ng. 2016. 11. 055 1016/j. foreco. 2017. 03. 016 Douglas DH (1994) Least-cost path in GIS using an accumulated Horn CM, Gilmore MP, Endress BA (2012) Ecological and socio- cost surface and slopelines. Cartographica: The International economic factors influencing aguaje (Mauritia flexuosa) resource Journal for Geographic Information and Geovisualization management in two indigenous communities in the Peruvian 31:37–51. https:// doi. org/ 10. 3138/ D327- 0323- 2JUT- 016M Amazon. Forest Ecology and Management 267:93–103. https:// Elhag M, Bahrawi JA (2017) Soil salinity mapping and hydrologi- doi. org/ 10. 1016/j. foreco. 2011. 11. 040 cal drought indices assessment in arid environments based Hounkpatin OKLL, Op de Hipt F, Bossa AY, Welp G, Amelung W on remote sensing techniques. Geoscientific Instrumentation, (2018) Soil organic carbon stocks and their determining fac- Methods and Data Systems 6:149–158. https://d oi.o rg/1 0.5 194/ tors in the Dano catchment (Southwest Burkina Faso). Catena gi-6- 149- 2017 166:298–309. https:// doi. org/ 10. 1016/j. catena. 2018. 04. 013 FAO (1965) Soil survey of the Llanos orientales, New York IGAC (2014) Estudio general de suelos zonificación de tierras: Depar- Funes I, Savé R, Rovira P, Molowny-Horas R, Alcañiz JM, Ascaso E, tamento de Casanare Herms I, Herrero C, Boixadera J, Vayreda J (2019) Agricultural IPBES (2019) Summary for policymakers of the global assessment soil organic carbon stocks in the north-eastern Iberian Peninsula: report on biodiversity and ecosystem services of the Intergovern- drivers and spatial variability. Science of the Total Environment mental Science-Policy Platform on Biodiversity and Ecosystem 668:283–294. https:// doi. org/ 10. 1016/j. scito tenv. 2019. 02. 317 Services. IPBES secretariat, Bonn, Germany.https:// doi. org/ 10. 5281/ zenodo. 26164 58 1 3 65 Page 18 of 20 Wetlands (2023) 43:65 Jadán O, Cifuentes M, Torres B, Selesi D, Veintimilla D, Günter S Minasny B, McBratney AB, Malone BP, Wheeler I (2013) Digital map- (2015) Influence of tree cover on diversity, carbon sequestration ping of soil carbon. Advances in Agronomy 118:1–47. https:// and productivity of cocoa systems in the Ecuadorian Amazon. doi. org/ 10. 1016/ B978-0- 12- 405942- 9. 00001-3 Bois et Forêts Des Tropiques 325:35–47 Minasny B, Rudiyanto, Sulaeman Y, Setiawan BI (2020) Open digital Jasiewicz J, Stepinski TF (2013) Geomorphons-a pattern recognition mapping for accurate assessment of tropical peatlands. Tropical approach to classification and mapping of landforms. Geomorphol- Wetlands Innovation in Mapping and Management- Proceedings ogy 182:147–156. https://d oi.o rg/1 0.1 016/j.g eomor ph.2 012.1 1.0 05 of the International Workshop on Tropical Wetlands: Innovation Jenny H (1941) Factors of Soil Formation: A System of Quantitative in Mapping and Management 2018:3–8. https://d oi.o rg/1 0.1 201/ Pedology. McGraw-Hill Book Company Inc 97804 29264 467-1 Junk WJ, An S, Finlayson CM, Gopal B, Květ J, Mitchell SA, Mitsch Mitra S, Wassmann R, Vlek PLG (2005) An appraisal of globa wetland WJ, Robarts RD (2013) Current state of knowledge regarding the are and its organic carbon stock. Current Science 88:25–35 world’s wetlands and their future under global climate change: a Mitsch WJ, Bernal B, Hernandez ME (2015) Ecosystem services of synthesis. Aquatic Sciences 75:151–167. https://d oi.o rg/1 0.1 007/ wetlands. International Journal of Biodiversity Science, Ecosys- s00027- 012- 0278-z tem Services & Management 11:1–4. https:// doi. org/ 10. 1080/ Köchy M, Hiederer R, Freibauer A (2015) Global distribution of soil 21513 732. 2015. 10062 50 organic carbon – Part 1: masses and frequency distributions of Mitsch WJ, Gosselink JG (2007) Wetlands, 4th edn. Wiley, Hoboken SOC stocks for the tropics, permafrost regions, wetlands, and the Mondal A, Khare D, Kundu S, Mondal S, Mukherjee S, Mukhopad- world. SOIL 1:351–365. https://d oi.o rg/1 0.5 194/s oil-1-3 51-2 015 hyay A (2017) Spatial soil organic carbon (SOC) prediction by Laganière J, Angers DA, Paré D (2010) Carbon accumulation in agricul- regression kriging using remote sensing data. Egyptian Journal tural soils after afforestation : a meta-analysis. Global Change Biol- of Remote Sensing and Space Science 20:61–70. https://d oi.o rg/ ogy 16:439–453. https://d oi.o rg/1 0.1 111/j.1 365-2 486.2 009.0 1930.x 10. 1016/j. ejrs. 2016. 06. 004 Lal R (2015) Restoring soil quality to mitigate soil degradation. 5875– Moore ID, Gessler PE, Nielsen GA, Peterson GA (1993) Soil attribute 5895. https:// doi. org/ 10. 3390/ su705 5875 prediction using terrain analysis. Soil Science Society of Amer- Lamichhane S, Kumar L, Wilson B (2019) Digital soil mapping algo- ica Journal 57:443–452. https://d oi.o rg/1 0.2 136/s ssaj1 993.0 3615 rithms and covariates for soil organic carbon mapping and their 99500 57000 20026x implications: a review. Geoderma 352:395–413. https:// doi. org/ Mora-Fernández C, Peñuela-Recio L (2013) Salud Ecosistémica de las 10. 1016/j. geode rma. 2019. 05. 031 sabanas inundables asociadas a la cuenca del río Pauto, Casanare, Lasso C, Rial A, Colonnello G, Machado-Allison A, Trujillo F (2014) Colombia. Yoluka ONG, fundación de investigación en biodiver- XI. Humedales de la Orinoquia (Colombia- Venezuela). Serie sidad y conservación, Fundación Horizonte Verde y Ecopetrol Editorial Recursos Hidrobiológicos y Pesqueros Continentales S.A, Bogotá, Colombia de Colombia. Instituto de Investigación de Recursos Biológicos Mora-Fernández C, Peñuela-Recio L, Castro-Lima F (2015) Estado del Alexander von Humboldt (IAvH). Bogotá conocimiento de los ecosistemas de las sabanas inundables en la Lasso C, Usma J, Trujillo F, Rial A (2010) Biodiversidad de la cuenca Orinoquia Colombiana. Orinoquia 19:253–271 del Orinoco - Bases científicas para la identificación de áreas pri- Motohka T, Nasahara KN, Oguma H, Tsuchida S (2010) Applicability oritarias para la conservación y uso sostenible de la biodiversi- of Green-Red Vegetation Index for remote sensing of vegeta- dad, Instituto de Investigación de Recursos Biológicos Alexander tion phenology. Remote Sensing 2:2369–2387. https:// doi. org/ von Humboldt, WWF Colombia, Fundación Omacha, Fundación 10. 3390/ rs210 2369 La Salle e Instituto de Estudios de la Orinoquia (Universidad Nahlik AM, Fennessy MS (2016) Carbon storage in US wetlands. Nature Nacional de Colombia). Bogotá, D. C., Colombia. https:// doi. Communications 7:1–9. https:// doi. org/ 10. 1038/ ncomm s13835 org/ 10. 1017/ CBO97 81107 415324. 004 Ngunjiri MW, Libohova Z, Minai JO, Serrem C, Owens PR, Schulze Liaw A, Wiener M (2003) Classification and regression by randomFor- DG (2019) Predicting soil types and soil properties with limited est. R News 3:18–22 data in the Uasin Gishu Plateau, Kenya. Geoderma Regional Margono BA, Potapov PV, Turubanova S, Stolle F, Hansen MC (2014) 16:e00210. https:// doi. org/ 10. 1016/j. geodrs. 2019. e00210 Primary forest cover loss in indonesia over 2000–2012. Nature Olson DM, Dinerstein E, Wikramanayake ED, Burgess ND, Powell Climate Change 4:730–735. https://d oi.o rg/1 0.1 038/n clima te227 7 GVNN, Underwood EC, D’Amico JA, Itoua I, Strand HE, Mor- McBratney AB, Mendonça Santos ML, Minasny B (2003) On digi- rison JC, Loucks CJ, Allnutt TF, Ricketts TH, Kura Y, Lamor- tal soil mapping. Geoderma. https:// doi. org/ 10. 1016/ S0016- eux JF, Wettengel WW, Hedao P, Kassem KR (2001) Terrestrial 7061(03) 00223-4 ecoregions of the world: a new map of life on Earth. BioScience Melton JR, Wania R, Hodson EL, Poulter B, Ringeval B, Spahni R, 51:933–938. https://d oi.o rg/1 0.1 641/0 006-3 568(2001)0 51[0933: Bohn T, Avis CA, Beerling DJ, Chen G, Eliseev AV, Denisov TEOTWA] 2.0. CO;2 SN, Hopcroft PO, Lettenmaier DP, Riley WJ, Singarayer JS, Ostle NJ, Levy PE, Evans CD, Smith P (2009) UK land use and soil Subin ZM, Tian H, Zürcher S, Brovkin V, van Bodegom PM, carbon sequestration. Land Use Policy 26:274–283. https:// doi. Kleinen T, Yu ZC, Kaplan JO (2013) Present state of global wet- org/ 10. 1016/j. landu sepol. 2009. 08. 006 land extent and wetland methane modelling: conclusions from a Page SE, Rieley JO, Banks CJ (2011) Global and regional importance model inter-comparison project (WETCHIMP). Biogeosciences of the tropical peatland carbon pool. Global Change Biology 10:753–788. https:// doi. org/ 10. 5194/ bg- 10- 753- 2013 17:798–818. https:// doi. org/ 10. 1111/j. 1365- 2486. 2010. 02279.x Minai JO, Libohova Z, Schulze DG (2021) Spatial prediction of soil prop- Patton NR, Lohse KA, Seyfried MS, Godsey SE, Parsons SB (2019) erties for the Busia area, Kenya using legacy soil data. Geoderma Topographic controls of soil organic carbon on soil-mantled Regional 25:e00366. https:// doi. org/ 10. 1016/j. geodrs. 2021. e00366 landscapes. Scientific Reports 9:1–15. https:// doi. org/ 10. 1038/ Minasny B, McBratney AB (2006) A conditioned Latin hypercube s41598- 019- 42556-5 method for sampling in the presence of ancillary information. Paul EA (2016) The nature and dynamics of soil organic matter: Plant Computers and Geosciences 32:1378–1388. https:// doi. org/ 10. inputs, microbial transformations, and organic matter stabiliza- 1016/j. cageo. 2005. 12. 009 tion. Soil Biology and Biochemistry 98:109–126. https://d oi.o rg/ 10. 1016/j. soilb io. 2016. 04. 001 1 3 Wetlands (2023) 43:65 Page 19 of 20 65 Peñuela L, Solano C, Ardila V, Galán S (2014) Sabana inundable y Sharma G, Sharma LK (2022) Assessment of soil carbon stock and ganadería, opción productiva de Conservación en la Orinoquia. important physicochemical properties in relation to land use Proyecto: “Fortalecimiento institucional y de política para incre- patterns in semi-arid region of Rajasthan, India. Journal of the mentar la conservación de la biodiversidad en predios privados Indian Society of Soil Science 70:191–203. https:// doi. org/ 10. en Colombia,” Sabana inundable y ganadería, opción produc- 5958/ 0974- 0228. 2022. 00019.6 tiva de conservación en la Orinoquia. Proyecto: “Fortalecimiento Shi X, Zhu A-X, Burt JE, Qi F, Simonson D (2004) A case-based institucional y de política para incrementar la conservación de reasoning approach to fuzzy soil mapping. Soil Science Society la biodiversidad en predios privados en Colombia.” Bogotá, of America Journal 68:885–894. https:// doi. org/ 10. 2136/ sssaj Colombia 2004. 8850 Powers JS (2005) Regional variation in soil carbon and 13C in forests Shiel D, Ladd B, Silva LCR, Laffan SW, van Heist M (2016) How are soil and pastures of northeastern Costa Rica Regional variation in soil carbon and tropical biodiversity related? Environmental Conser- carbon and d 13 C in forests and pastures of northeastern Costa vation 43:231–241. https:// doi. org/ 10. 1017/ S0376 89291 60000 11 Rica. https:// doi. org/ 10. 1007/ s10533- 004- 0368-7 Shimada M, Itoh T, Motooka T, Watanabe M, Shiraishi T, Thapa R, Prahalad VN, Kriwoken LK (2010) Implementation of the Ramsar Lucas R (2014) New global forest/non-forest maps from ALOS convention on wetlands in Tasmania, Australia. Journal of Inter- PALSAR data (2007–2010). Remote Sensing of Environment national Wildlife Law and Policy 13:205–239. https://d oi.o rg/1 0. 155:13–31. https:// doi. org/ 10. 1016/j. rse. 2014. 04. 014 1080/ 13880 292. 2010. 486697 Siewert M (2017) High-resolution digital mapping of soil organic car- Rice AL, Butenhoff CL, Teama DG, Röger FH, Khalil MAK, Rasmus- bon in permafrost terrain using machine-learning: a case study in sen RA (2016) Atmospheric methane isotopic record favors fossil a sub-artic peatland environment. Biogeosciences 15:1663–1682. sources flat in 1980s and 1990s with recent increase. Proceedings https:// doi. org/ 10. 5194/ bg- 15- 1663- 2018 of the National Academy of Sciences 113:10791–10796. https:// Silva-parra A (2018) Modeling soil carbon stocks and carbon dioxide doi. org/ 10. 1073/ pnas. 15229 23113 emissions (GHG) in production systems of Plain Altillanura. Romero Duque LP, Castro Lima F, Rentería Mosquera Á (2018) Con- Orinoquía 22:1–25. https:// doi. org/ 10. 22579/ 20112 629. 525 tribución al conocimiento de la vegetación de las sabanas de Sindayihebura A, Ottoy S, Dondeyne S, Van Meirvenne M, Van Casanare (Colombia). Revista U.D.C.A. Actualidad & Divul- Orshoven J (2017) Comparing digital soil mapping techniques gación Científica 21:197–205. https:// doi. org/ 10. 31910/ rudca. for organic carbon and clay content: Case study in Burundi’s v21. n1. 2018. 678 central plateaus. Catena 156:161–175. https://d oi.o rg/1 0.1 016/j. Roudier P, Hewitt AE, Beaudette DE (2012) A conditioned Latin catena. 2017. 04. 003 hypercube sampling algorithm incorporating operational con- Six J, Conant RT, Paul EA, Paustian K (2002) Stabilization mecha- straints. Digital Soil Assessments and Beyond - Proceedings of nisms of soil organic matter: implications for C-saturation of the Fifth Global Workshop on Digital Soil Mapping, 227–231. soils. Plant and Soil 241:155–176 https:// doi. org/ 10. 1201/ b12728- 46 Sjögersten S, Black CR, Evers S, Hoyos-Santillan J, Wright EL, Turner Sarmiento G (1984) The Ecology of Neotropical Savannas, Harvard BL (2014) Tropical wetlands: a missing link in the global carbon University Press, Cambridge, MA 02138. ISBN 0-674- 22460-4. cycle? Global Biogeochemical Cycles 28:1371–1386. https://d oi. Harvard University Press. https://d oi.o rg/1 0.4 159/h arvar d.9 7806 org/ 10. 1002/ 2014G B0048 44 74418 554 Soil Science Division Staff (2017) Soil Survey Manual. United States Sarmiento G (1983) The savannas of tropical America. In: F. Bour- Department of Agriculture Handbook No. 18 18, 1913–5. https:// liere, F. (Ed): Ecosystems of the World XIII. Tropical Savannas, doi. org/ 10. 1121/1. 19035 34 Elsevier, Amsterdam, pp. 245-288., in: Ecosystems of the World Sørensen R, Zinko U, Seibert J (2006) On the calculation of the topo- XIII. Tropical Savannas. pp. 245–288 graphic wetness index: evaluation of different methods based Saunois M, Stavert AR, Poulter B, Bousquet P, Canadell JG, Jackson on field observations. Hydrology and Earth System Sciences RB, Raymond PA, Dlugokencky EJ, Houweling S, Patra PK, 10:101–112. https:// doi. org/ 10. 5194/ hess- 10- 101- 2006 Ciais P, Arora VK, Bastviken D, Bergamaschi P, Blake DR, Szabó B, Szatmári G, Takács K, Laborczi A, Makó A, Rajkai K, Pász- Brailsford G, Bruhwiler L, Carlson KM, Carrol M, Castaldi S, tor L (2019) Mapping soil hydraulic properties using random- Chandra N, Crevoisier C, Crill PM, Covey K, Curry CL, Etiope forest-based pedotransfer functions and geostatistics. Hydrology G, Frankenberg C, Gedney N, Hegglin MI, Höglund-Isaksson L, and Earth System Sciences 23:2615–2635. https:// doi. org/ 10. Hugelius G, Ishizawa M, Ito A, Janssens-Maenhout G, Jensen 5194/ hess- 23- 2615- 2019 KM, Joos F, Kleinen T, Krummel PB, Langenfelds RL, Laruelle Taghizadeh-Mehrjardi R, Neupane R, Sood K, Kumar S (2017) Arti- GG, Liu L, Machida T, Maksyutov S, McDonald KC, McNorton ficial bee colony feature selection algorithm combined with J, Miller PA, Melton JR, Morino I, Müller J, Murguia-Flores machine learning algorithms to predict vertical and lateral dis- F, Naik V, Niwa Y, Noce S, O’Doherty S, Parker RJ, Peng C, tribution of soil organic matter in South Dakota, USA. Carbon Peng S, Peters GP, Prigent C, Prinn R, Ramonet M, Regnier P, Management 8:277–291. https:// doi. org/ 10. 1080/ 17583 004. Riley WJ, Rosentreter JA, Segers A, Simpson IJ, Shi H, Smith SJ, 2017. 13305 93 Steele LP, Thornton BF, Tian H, Tohjima Y, Tubiello FN, Tsu- Terra DF (2017) Accuracy and uncertainty of digital soil mapping ruta A, Viovy N, Voulgarakis A, Weber TS, van Weele M, van approaches to extract and transfer soil information from reference der Werf GR, Weiss RF, Worthy D, Wunch D, Yin Y, Yoshida area. Dissertação (Mestrado Em Ciência Do Solo) Universidade Y, Zhang W, Zhang Z, Zhao Y, Zheng B, Zhu Qing, Zhu Qiuan, Federal de Lavras 114 Zhuang Q (2020) The global methane budget 2000–2017. Earth Tsakiridis NL, Theocharis JB, Panagos P, Zalidis GC (2019) An evo- System Science Data 12:1561–1623. https:// doi. org/ 10. 5194/ lutionary fuzzy rule-based system applied to the prediction of essd- 12- 1561- 2020 soil organic carbon from soil spectral libraries. Applied Soft Schreier H, Shah PB, Lavkulich LM, Brown S (1994) Maintaining soil Computing Journal 81:105504. https:// doi. org/ 10. 1016/j. asoc. fertility under increasing land use pressure in the Middle Moun- 2019. 105504 tains of Nepal. Soil Use and Management 10:137–142. https:// Usma JS, Trujillo F (2011) Biodiversidad del departamento Casanare. doi. org/ 10. 1111/j. 1475- 2743. 1994. tb004 74.x Identificación de ecosistemas estratégicos. WWF - Colombia, Bógota, Colombia 1 3 65 Page 20 of 20 Wetlands (2023) 43:65 Valbuena R, Hernando A, Manzanera JA, Görgens EB, Almeida indicators at various scales. Geoderma 333:149–162. https://d oi. DRAA, Mauro F, García-Abril A, Coomes DA (2017) Enhanc- org/ 10. 1016/j. geode rma. 2018. 07. 026 ing of accuracy assessment for forest above-ground biomass esti- Xie E, Zhang Y, Huang B, Zhao Y, Shi X, Hu W, Qu M (2021) Spati- mates obtained from remote sensing via hypothesis testing and otemporal variations in soil organic carbon and their drivers in overfitting evaluation. Ecological Modelling 366:15–26. https:// southeastern China during 1981–2011. Soil and Tillage Research doi. org/ 10. 1016/j. ecolm odel. 2017. 10. 009 205:104763. https:// doi. org/ 10. 1016/j. still. 2020. 104763 Van der Werf GR, Morton DC, Defries RS, Olivier JGJ, Kasibhatla PS, Xu J, Morris PJ, Liu J, Holden J (2018) PEATMAP: Refining estimates Jackson RB, Collatz GJ, Randerson JT (2009) CO2 emissions of global peatland distribution based on a meta-analysis. Catena from forest loss. Nature Geoscience 2:737–738. https:// doi. org/ 160:134–140. https:// doi. org/ 10. 1016/j. catena. 2017. 09. 010 10. 1038/ ngeo6 71 Yang R, Liu F, Zhang G, Zhao Y, Li D, Yang J, Yang Fei, Yang Fan (2016) Vega LF, Nunes da Cunha C, Rothaupt K-O, Moreira MZ, Wantzen Mapping soil texture based on field soil moisture observations at a KM (2014) Does flood pulsing act as a switch to store or release high temporal resolution in an oasis agricultural area. Pedosphere sediment-bound carbon in seasonal floodplain lakes? Case 26:699–708. https:// doi. org/ 10. 1016/ S1002- 0160(15) 60078-9 study from the Colombian Orinoco-Llanos and the Brazil- Yang Y, Mohammat A, Feng J, Zhou R, Fang J (2007) Storage, patterns ian Pantanal. Wetlands 34:177–187. https:// doi. org/ 10. 1007/ and environmental controls of soil organic carbon in China. Biogeo- s13157- 013- 0495-9 chemistry 84:131–141. https:// doi. org/ 10. 1007/ s10533- 007- 9109-z Vepraskas MJ, Craft CB (2016) Wetland Soils: Genesis, Hydrology, Yun J, Chen X, Liu S, Zhang W (2019) Effects of temperature and Landscapes, and Classification, 2nd Editio. ed, CRC Press. CRC moisture on soil organic carbon mineralization.  IOP Conference Press, Boca Raton. doi:https:// doi. org/ 10. 1201/ b18996 Series: Materials Science and Engineering 562. https:// doi. org/ Vu DH, Muttaqi KM, Agalgaonkar AP (2015) A variance inflation 10. 1088/ 1757- 899X/ 562/1/ 012085 factor and backward elimination based robust regression model Zeraatpisheh M, Ayoubi S, Jafari A, Tajik S, Finke P (2019) Digital for forecasting monthly electricity demand using climatic vari- mapping of soil properties using multiple machine learning in a ables. Applied Energy 140:385–394. https:// doi. org/ 10. 1016/j. semi-arid region, central Iran. Geoderma 338:445–452. https:// apene rgy. 2014. 12. 011 doi. org/ 10. 1016/j. geode rma. 2018. 09. 006 Wan Q, Zhu G, Guo H, Zhang Y, Pan H, Yong L, Ma H (2019) Influ- Zhou Y, Hartemink AE, Shi Z, Liang Z, Lu Y (2019) Land use and ence of vegetation coverage and climate environment on soil climate change effects on soil organic carbon in North and North- organic carbon in the Qilian Mountains. Scientific Reports 9:1–9. east China. Science of the Total Environment 647:1230–1238. https:// doi. org/ 10. 1038/ s41598- 019- 53837-4 https:// doi. org/ 10. 1016/j. scito tenv. 2018. 08. 016 Wang X, Zhang Y, Atkinson P, Yao H (2020) Predicting soil organic Zhu A-X (1997) A similarity model for representing soil spatial infor- carbon content in Spain by combining Landsat TM and ALOS mation. Geoderma 77:217–242 PALSAR images. International Journal of Applied Earth Obser- Zhu A-X, Hudson B, Burt J, Lubich K, Simonson D (2001) Soil map- vation and Geoinformation 92:102182. https://d oi.o rg/1 0.1 016/j. ping using GIS, expert knowledge, and fuzzy logic. Soil Science jag. 2020. 102182 Society of America Journal 65:1463–1472. https:// doi. org/ 10. Wantzen KM, Couto EG, Mund EE, Amorim RSS, Siqueira A, Tiel- 2136/ sssaj 2001. 65514 63x börger K, Seifan M (2012) Soil carbon stocks in stream-valley- Zhu A-X, Qi F, Moore A, Burt JE (2010) Prediction of soil properties ecosystems in the Brazilian Cerrado agroscape. Agriculture, Eco- using fuzzy membership values. Geoderma 158:199–206. https:// systems and Environment 151:70–79. https:// doi. org/ 10. 1016/j. doi. org/ 10. 1016/j. geode rma. 2010. 05. 001 agee. 2012. 01. 030 Zhu A-X (1999) A personal construct-based knowledge acquisition Weil R, Brady N (2017) The Nature and Properties of Soils, 15th edn. process for natural resource mapping. International Journal of Pearson, United States of America, Columbus Geographical Information Science 13(2):119–141. https:// doi. Welikhe P, Quansah JE, Fall S, Elhenney WMc (2017) Estimation of org/ 10. 1080/ 13658 81992 41382 soil moisture percentage using LANDSAT-based moisture stress Zhu Q, Peng C, Ciais P, Jiang H, Liu J, Bousquet P, Li S, Chang J, index. Journal of Remote Sensing & GIS 6(2). https://d oi.o rg/1 0. Fang X, Zhou X, Chen H, Liu S, Lin G, Gong P, Wang M, Wang 4172/ 2469- 4134. 10002 00 H, Xiang W, Chen J (2017) Interannual variation in methane Whiting GJ, Chanton JP (2001) Greenhouse carbon balance of wet- emissions from tropical wetlands triggered by repeated El Niño lands: methane emission versus carbon sequestration. Tellus B: Southern Oscillation. Global Change Biology 23:4706–4716. Chemical and Physical Meteorology 53:521–528. https://d oi.o rg/ https:// doi. org/ 10. 1111/ gcb. 13726 10. 3402/ tellu sb. v53i5. 16628 Zinn YL, Lal R, Resck DVSS (2005) Texture and organic carbon rela- Wiesmeier M, Prietzel J, Barthold F, Spörlein P, Geuß U, Hangen E, Reischl tions described by a profile pedotransfer function for Brazilian A, Schilling B, von Lützow M, Kögel-Knabner I (2013) Storage and Cerrado soils. Geoderma 127:168–173. https://d oi.o rg/1 0.1 016/j. drivers of organic carbon in forest soils of southeast Germany (Bavaria) geode rma. 2005. 02. 010 - implications for carbon sequestration. Forest Ecology and Manage- ment 295:162–172. https:// doi. org/ 10. 1016/j. foreco. 2013. 01. 025 Publisher's Note Springer Nature remains neutral with regard to Wiesmeier M, Urbanski L, Hobley E, Lang B, von Lützow M, Marin- jurisdictional claims in published maps and institutional affiliations. Spiotta E, van Wesemael B, Rabot E, Ließ M, Garcia-Franco N, Wollschläger U, Vogel HJ, Kögel-Knabner I (2019) Soil organic carbon storage as a key function of soils - a review of drivers and 1 3