Abstract: Limited information is available about the use of intelligent system such as Artificial Neural Networks (ANN) to determine the most affecting factors on variability of soil organic carbon fractions (SOC) in the landscape scale. Therefore, this study was conducted to estimate SOC fractions by topographic attributes, selected soil properties and Normalized Vegetation Index (NDVI) data using ANN models. A total of 108 samples from surface soils (0-10 cm depth) were collected and various physical soil organic fractions were determined. The developed ANN models could explain 78-91% of the total variability in SOC fractions in the site studied. Sensitivity analysis using ANN models developed showed that NDVI as indication of vegetation cover was the most important factor for explaining variability of SOC fractions at the site. Furthermore, soil properties such as clay, silt and calcium carbonate and some topographic attributes which indirectly affect the total SOC content, also significantly influence the variability of SOC fractions. In overall, the results showed that the ANN models provide reliable prediction of SOC fractions by considering the NDVI, soil properties and terrain attributes.
INTRODUCTION
Soil Organic Carbon (SOC) plays a vital role in crop growth in natural ecosystems and at the same time is influenced by land use, soil type, climate and vegetation (Loveland and Webb, 2003; Onweremadu, 2008). It is also one of the important factors affecting soil quality, sustainability of agriculture, soil aggregate stability and plant production (Freixo et al., 2002; Loveland and Webb, 2003). Moreover, the intensity of global warming in the future is directly influenced by SOC cycle (Lal, 2004).
The SOC pool is strongly affected by land use changes and management strategies. SOC is derived from surface input of plants as well as roots and associated turnover of mycorrhizal hyphae (Lal, 2004; Lorenz et al., 2008). Change in land use such as clear or partial cutting of forest influences both quantity and quality of soil organic matter (SOM). It has frequently reported that the soils in natural ecosystems had higher SOM content, aggregate stability and saturated hydraulic conductivity as compared to their agricultural counterparts (Saviozzi et al., 2001; Puget and Lal, 2005).
Organo-mineral interactions protect SOC against biological degradation. The mineralogy and size distribution of the mineral fraction also affect SOC protection (Baldock and Skjemstad, 2000; Schmidt and Kogel-Knabner, 2002). Generally, clay and silt contents are positively correlated with SOM concentration and their amounts in the soil directly contribute to the accumulation of silt- and clay-protected SOC fractions (Six et al., 2002; Zinn et al., 2005; Bouajila and Gallali, 2008).
Knowledge about the formation and stabilization of soil aggregates in natural and disturbed ecosystems is necessary to address a variety of environmental concerns. These affairs are ranging from the fate and transport of hazardous wastes to the potential C sink strength of terrestrial ecosystem. The SOC and aggregates mutually protect each other. SOC is efficiently influenced the soil aggregation (Six et al., 1999). On the other hand, SOC is physically protected and stabilized between and within micro-aggregates (Six et al., 1999, 2000, 2002).
Evaluation of the total SOC pool provides little information concerning biochemical stability and duration of C stored in soils (Puget et al., 2005; Lorenz et al., 2008). Hence, the quantification of SOC in different fractions provides valuable information on functionality of the SOC pool (Bouajila and Gallali, 2008). Direct measurement of soil organic pool in different fractions is time consuming and laborious and costly affair (Yerokun et al., 2007). Therefore, the use of indirect techniques such as Digital Terrain Modeling (DTM) and Remote Sensing (RS) is highly desirable at the landscape scale.
Several soil scientists identified topography as one of the pedogenic factors (Florinsky et al., 2002) which significantly influences the spatial distribution of soil moisture, temperature and organic matter (Florinsky et al., 2004). Thus, quantitative information on the topographic characteristics has been employed in the form of Digital Terrain Models (DTMs). Since, the prediction soil properties using DTM describes the relationships between soil and topographic attributes at a point in the landscape (Moore et al., 1993; Bell et al., 1994; Thompson et al., 1997). Quantitative topographic data is often used in soil studies including for the modeling and prediction of soil properties. Because of a high spatial variability in soil properties within the landscape (Huggett, 2003), the use of indirect prediction approaches as the alternative methods have been widely used, for example, in DTM modeling to predict soil properties at a point in the landscape in a cost effective way. Moreover, RS has been widely used to evaluate surface SOM at landscape scale (Chen et al., 2000).
Artificial Neural Network (ANN) is a mathematical tool which has been inspired by biological neural networks and is a popular tool in the classification, prediction and recognition-based problems. It has widely been employed for prediction and modeling in environmental and biological concerns (Abdalla and Deris, 2005; Chayjan and Moazez, 2008; Fallah-Ghalhary et al., 2009; Dastorani et al., 2010).
To our knowledge, little attempt has been made to determine important environmental factors at the landscape scale that control SOC pool variability in primary particles and aggregate sizes using ANN. Therefore, this study was conducted: (i) to predict SOC pools in primary particles and aggregate sizes using soil, topographic and Normalized Difference Vegetation Index (NDVI) at the landscape scale and (ii) to determine soil, topographic and NDVI attributes that most explain the variability of the different fractions of SOC in a semiarid region in western Iran.
MATERIALS AND METHODS
Site description: This study was conducted in hilly region of upland Lordegan watershed located in western Iran (Fig. 1). The study area is located within 50°12′to 50°37′ E longitude and 31°58′ to 32°03′ N latitude.
| Fig. 1: | Location of the study area in western Iran |
The mean elevation of the area is approximately 1860 m a.s.l. The mean annual temperature and precipitation at the site are 15°C and 600 mm, respectively.
The hill slopes of the study area have been developed by extensive dissection of sedimentary Quaternary deposits. The soils of the study area were predominantly classified as Fine loamy, mixed, thermic, Typic Haploxerolls and Fine mixed, thermic, Typic Calcixerepts (Soil Survey Staff, 2006).
In total, 108 samples were collected from different land uses and slope position to capture all environmental variability in the studied area. Soil samples were collected from 0-10 cm depth using an auger; three sub-samples per 1 m2 area were made in to one composite sample to reduce micro-variability.
Laboratory analysis: To determine Particulate Organic Matter (POM), sub-samples of whole soil (2 mm) were dispersed in distilled water with Particulate Organic Matter (POM) removal at high-energy sonication (12500 J) for 13 min to complete aggregate breakdown and dispersion (Edwards and Bremner, 1967). The dispersed sample was passed through a 0.053 mm sieve. The material left on the sieve (>0.053 mm) was dried at 50°C in a ventilated oven and gently ground to pass a 0.053 mm sieve for determining POM. The silt- and clay-sized organo-mineral particles in suspension were separated considering the Stocks' law and using sedimentation and siphoning method (Bronick and Lal, 2005) and the mineral-associated organic C was measured in the clay and silt fractions which are hereby termed as organo-clay and organo-silt in this paper.
To determine the SOC pool in aggregate size fractions at first aggregates were separated using wet sieving method. A 100 g soil sample having intact aggregates which passed through a 4.75 mm sieve was capillary-wetted to matric suction of 30 kPa. The soil water content at 30 kPa was determined on a separate batch of aggregates using a pressure plate apparatus (Klute, 1986). The aggregates were separated into three sizes (i.e. 4.75-2, 2-0.25 and 0.25-0.053 mm) for chemical analysis (Cambardella and Elliott, 1993). The aggregates were wet-sieved in water for 30 min with vertical stroke of 1.3 cm and speed of 30 strokes min-1. The measured SOC in the three mentioned sizes are described hereby as macroaggregate-SOC, mesoaggregate -SOC and macroaggregate-SOC, respectively.
The SOC content was measured by the wet-oxidation method (Nelson and Sommers, 1982) in different fractions. Percentages of clay and silt were measured using a hydrometer method (Gee and Bauder, 1986) and the sand fraction was measured by sieving. Calcium carbonate equivalent (CCE) was measured by the Bernard’s calcimetric method (Black et al., 1965). Soil bulk density was measured using the core method (Blake and Hartge, 1986).
Digital terrain analysis and NDVI calculation: The elevation data were used to create 3x3 m Digital Elevation Models (DEM) using ILWIS (ITC, 1997). Then, primary and secondary topographical indices were generated from the DEM using ILWIS software and DIGEM software (http://www.geogr.uni-goettingen.de/pg/saga/digem). Topographic indices included primary and secondary indices. Primary indices were calculated directly from the DEM included elevation, slope, aspect, specific catchment area, profile, plan and mean curvatures. Plan curvature (PLANC) is curvature of the corresponding normal section which is tangential to a contour.
Vertical or profile curvature (PROFC) is curvature of corresponding normal section which is tangential to a flow line. Mean curvature (MEANC) is the average of normal section curvature (Wilson and Gallant, 2000). Secondary indices calculated from the combinations of the primary indices included Wetness Index (WI) which is the ratio of specific catchment area to slope gradient and indicates the spatial distribution of zones of surface saturation and soil water content in landscape. Wetness index was calculated using Eq. 1:
| (1) |
To reflect the erosive power of the terrain, Stream Power Index (SPI) was calculated using the Eq. 2:
| (2) |
The remote sensing data were used to build the model in this study included the Landsat ETM bands with spatial resolution of 30x30 m. The acquisition date of the image was 22 June 2010. The subset image covering the study area was corrected geometrically using the landform map of Iran 1:25000 scale as the reference. All image processing was performed by using ILWIS software.
The NDVI was calculated as the reflectance ratio from near-infrared (NIR) and red channel (R) of satellite or airborne sensors of Land sat satellite (ETM+) of 2002 as follows (Eq. 3):
| (3) |
Artificial neural network modeling: The most popular network used in engineering problems for nonlinear mapping was probably multilayer perceptron (MLP) with Back-Propagation (BP) learning rule (Haykin, 1994). A feed-forward back-propagating ANN structure was used to develop SOC fractions. Supervised learning uses known outputs to train the ANN and is more commonly used than unsupervised learning. Back propagation is a form of supervised learning where the error rate was sent back through the network to alter the weights to improve prediction and decrease the error (Kaul et al., 2005). The standard algorithm was based on the delta learning rule (Rumelhart and McClelland, 1986). For designing the ANN, MATLAB, software package was used. The topographical attributes, NDVI and selected soil properties were used as the input data for the three categories with SOC fractions as the target data in six output node (Table 1).
For designing the artificial neural network, the measured field data were used. The number of available data collected for this study was 108. The data set was shuffled; 64 of them were used for the learning process, 22 sets were used for testing and the remaining 22 sets were used for verification, respectively. The data sets for learning, testing and verification processes were selected randomly at different points on the landscape to avoid bias in estimation. In the modeling process, standardized variables were used and calculated as follows:
| (4) |
where, Xs is standard value, Xi is actual value, Xmean is the arithmetic mean of total value, Xmax is the maximum value and Xmin is the minimum value.
| Table 1: | Inputs and outputs categories of variables to establish the ANN models |
| BD: Bulk density; CCE: Calcium carbonate equivalent; TOC: Total organic carbon | |
The training process was performed using the BP in two steps, forward pass and backward pass (Levenberg- Marquardt training rule). In the forward step, an output pattern was presented to the network and its effect propagated through the network layer by layer. Then, the final computed output of the network was compared with the target output. In this step, a performance function (i.e. mean square error, MSE) was calculated and then second step of the BP algorithm was started by back propagation of the network error to the previous layer using the gradient-descent technique, the weights were adjusted to reduce the network error. This process was continued until the allowable network error was obtained. In many problems, a second hidden layer does not produce a large improvement in performance and varying the number of hidden neurons in the hidden layer is sufficient (El-Din and Smith, 2002). The number of hidden layers, the number of neurons in hidden layer and the number of iteration (Epoch) were selected by calibration through several test runs and trial and error. The best function for network was tansigmoid. Based on the R2 value of regression between the measured and predicted outputs, the number of neuron in hidden layer, iteration and finally, the best model was selected.
The Root Mean Square Error (RMSE), Mean Estimation Error (MEE) and coefficient of determination (R2) between the measured and the estimated values were used to evaluate the performance of models. The RMSE and MEE (Degroot, 1986) are as denoted below:
| (5) |
| (6) |
where, S(xi) denotes the predicted value, M(xi) is measured value and n is the total number of observations. The R2 also shows the degree to which two variables are linearly related to.
In order to identify the most important factors explaining the variability of SOC fractions, sensitivity analysis was done using the StatSoft method (StatSoft, 2004). A sensitivity ratio was calculated by dividing the total network error when the variable was treated as non variable by the total network error when the actual values of the variables were used. A ratio >1.0 implied that the variable made an important contribution to the variability in the property; and the variable with higher ratio was more important (StatSoft, 2004).
RESULTS AND DISCUSSION
Descriptive statistics: The descriptive statistics of the SOC fractions in surface soil samples (0-10 cm depth) for the studied area are given in Table 2. All selected variables followed normal distribution according to the Kolmogrov-Smirnov test. This was also confirmed by the values of skewness (Table 2) which varied from -1 to +1. The Coefficient of Variation (CV), as an index of variation of heterogeneity, was used. Among the SOC fractions, the highest CV was ascribed to POM (117%) and the lowest to SOC-clay (38.88%) (Table 2). Overall, almost all SOC fractions showed high variation in the studied region. It is likely that high variability in the SOC fractions, attributed to diversity of land uses in the studied area with different organic matter input varying in both quantity and quality, as well as to landscape position. It seems that the variability associated with SOC fractions depends on landscape position, causing differential movement of water at different positions in the landscape which leads to soil redistribution in different parts of the landscape Afshar et al. (2010).
| Table 2: | Descriptive statistics of selected soil physical, chemical and magnetic properties of surface (0-30 cm) soil samples at the study site in western Iran (N = 108) |
| Min: Minimum; Max: Maximum; SD: Standard deviation; C.V: Coefficient of Variation; SOC: Soil organic carbon; POM: Particulate organic matter; BD: Bulk density; CCE: Calcium carbonate equivalent; TOC: Total organic carbon | |
| Table 3: | Summary of the best structure and optimum parameters for the ANN models used for predicting selected soil properties at the site studied |
| SOC: Soil organic carbon; POM: Particulate organic matter | |
ANN modeling: For predicting SOC fractions in the selected hilly region, best structure of the ANN for each parameter was ascertained (Table 3). Each of the trained structures had 20 input nodes in three categories of soil, RS and topographic properties and six output nodes including SOC fractions (Table 3). The hidden-layer nodes were optimized 56, 45, 46, 51, 58 and 48 and the optimum iteration learning rates based on trial and error at 8000, 9000, 10000, 11000, 9000 and 8000 for POM, SOC-clay, SOC-silt, SOC-macro, SOC-meso, SOC- micro, respectively (Table 3).
The ANN model for POM resulted MEE,-0.012 and RMSE, 0.01, respectively. Also, the ANN models for SOC-clay, SOC-silt resulted MEE, 0.005 and -0.07 and -0.11 and 0.09, respectively (Table 4). The ANN models for SOC-macro, SOC-meso, SOC- micro resulted MEE, 0.05, -0.003 and 0.005 and RMSE, -0.12, 0.03 and 0.09, respectively. The ANN models developed for simulating SOC fractions explained 88, 91, 84, 78, 79 and 81% of the variability in the POM, SOC-clay, SOC-silt, SOC-macro, SOC-meso, SOC- micro, respectively, at the site studied (Table 4). The normalized predicted data versus normalized observed data for testing data set for different SOC fractions are illustrated in Fig. 2. The positive significant (p<0.05) correlation coefficients (r) of 0.93, 0.95 and 0.92 between the observed and the predicted POM, SOC-clay and SOC-silt were established which are presented in Fig. 2(a-c), receptively.
| Table 4: | Results of the sensitivity analysis (ratios) of the final ANN model used for predicting SOC fractions in the area studied |
| SCA: Specific catchment area; MEANC: Mean curvature; PLANC: Plan curvature; PROFC: profile curvature; SPI: Stream power index; STI: Sediment transport index; WI: Wetness index; CCE: Calcium carbonate equivalent; SOM: Soil organic matter; RMSE: Root mean square; RSP: Relative stream power; Shdrelif: Shaded relief; NDVI: Normalized difference vegetation index; BD: Bulk density; TOC: Total organic carbon; MEE: Mean estimation error; R2: Coefficient of determination | |
Furthermore, the scatter plots of observed and predicted of SOC-macro, SOC-meso, SOC- micro are presented in Fig. 2d-f, respectively. As it is seen in Fig. 2d the relationship between normalized observed data and testing data for SOC-macro was significant as 0.01 probability level with a coefficient of determination 79%. Also, the significant relationship resulted for SOC-meso and SOC-micro with coefficients of determination 79 and 81%, respectively (Fig. 2e-f).
Overall, the ANN models developed for predicting the SOC fractions in the present study by incorporating of NDVI created by ETM-Landsat and terrain attributes, explained 78-91% of the total variability in SOC fractions within the landscape. A part of the unexplained variability is probably due to other factors that affect the variability of SOC and may also contributed to uncertainty of remote sensing data especially due to accordance to unreal-time data. Moreover, as reported by other researchers (Kaul et al., 2005; Salazar et al., 2010) it is important to compare the results by the ANN models with those obtained by other statistical approaches. Hence the learning rate, number of hidden layer Lal (2004), number of hidden nodes and the training tolerance need to be determined accurately for developing models to predict SOC fractions.
| Fig. 2: | Scatter plots displaying relationships between standardized measured and estimated value of the SOC fractions from 0-10 cm surface soil layer in the studied area using ANN modeling. (a) POM, (b) SOC-clay; (c) SOC-silt; (d) SOC-macro; (e) SOC-meso; (f) SOC-micro |
However, the performance of the ANN models as compared with other approaches has greater realistic chance in SOC fractions prediction, especially when complex non-linear relationships exist among various factors. In such cases, the correlation study may provide inaccurate and even misleading results about the relationships (Liu et al., 2001).
The use of ANN modeling with additional hill slopes with greater variability in terrain attributes should help broaden the usefulness of the ANN-based SOC fractions prediction. In this regards, combining of terrain attributes with remotely sensed data with higher spectral and ground resolutions and real time remote sensing data could provide precise predictions.
ANN application has functional characteristics and provides many advantages over the other modeling approaches such as linear regression models. The most important advantage of using the neural network approach is that the network trained to find the relationships and the lack of them is assumed beforehand. Also, the other powerful attributes of ANN models include their flexibility and adaptively which play important role in material modeling (Liu et al., 2001; Kaul et al., 2005). It appears that the ANN approach may be sufficiently valid in predicting the SOC fractions using soil, RS and topographic attributes in the area studied. A reason for these findings can be attributed to the nonlinear relationship between soil and topographic attributes and the SOC fractions and the ANN technique can estimate these relations using nonlinear functions.
Determination of important factors explaining variability in SOC fractions: The relative importance of terrain attributes, remote sensing data and soil properties using sensitivity analysis based upon coefficients of sensitivity of the selected ANN model for estimating the SOC fraction, is presented in Table 4. The variables with high values made important contribution to the variability in SOC fractions.
The NDVI was identified as the most and first factor among the 20 input variables for explaining the variability in all SOC fractions in the study area. The NDVI, a remote sensing index, indicates the green cover on the land surface and displays a well documented relationship with crop and vegetation productivity and land use effect on SOC pools (Li et al., 2001; Pettorelli et al., 2005). Podeh et al. (2009) in a study in the Mazandaran forest computed the NDVI as major index to explore the spatial and temporal dynamics of land use cover. Land use can affect the distribution of SOC in different fractions due to total soil organic carbon, microbial activity, animal and root activities and presence of fungal hyphae (Kay, 1988; Lal, 2004; Lorenz et al., 2008). The abovementioned variables were significantly different among natural Quercus forest, disturbed forest and cultivated soils in the study area. In forest soils, probably the presence of polysaccharides and monosaccharides led to the formation of macroaggregates with higher SOC pool (Larre-Larrouy et al., 2004).
Among the soil properties, TOC was identified as the most important factor influencing SOC pools in different fractions with relative coefficient of sensitivity varying from 3.00 to 3.33 for different fractions (Table 4). Obviously higher TOC leads to higher SOC in different fractions. Furthermore, clay and silt contents were identified important variables for SOC-clay and SOC-silt fractions (Table 4). Clay and silt contents are generally positively correlated with SOM concentration and their amount in soil directly contributed to the accumulation of silt- and clay-protected SOC (Six et al., 2002). Clay content and calcium carbonate as the binding agents were also identified as the controlling factors for the variability in SOC concentration in macro- and mesoaggregates with relatively high coefficient of sensitivity. Clay content has a vital role in aggregation and subsequently affecting physical protection of SOC in macroaggregates. A study of Franzluebbers and Arshad (1996) on soil organic carbon pools affected by tillage practices in Canada reported strong relationships between clay content and macroaggregate size under different land uses. They indicated that soils containing 20-69% clay had higher POM compared to the soils with lower clay content. Bulk density and gravel content had lower contribution in explaining the variability in SOC fractions. Moreover, clay fraction of soil protects soil organic carbon by lower porosity with lower oxidation rate of SOC and also by surface adsorption of SOC on the super-active surfaces of clays (Christensen, 1992; Balabane and Plante, 2004).
Among the topographic attributes, WI, ProfC, Slope, SPI, STI and Shaded were identified as the most important factors that can be used in modeling of SOC fractions at the site studied (Table 4). All these factors, indirectly affected the vegetation density, microbial activity, total soil organic carbon, clay content and CCE also affected SOC fractions. For example, shaded relief indirectly influences TOC and SOC fractions. Carter et al. (1998) showed that shading of forest was one of the main factor in the degradation of SOC because of lower temperature and therefore, accumulation of SOC. Soil properties are significantly influenced by topography in hilly regions. Reicosky et al. (2005) and Afshar et al. (2010) showed that silt and clay were eroded from the convex slopes and transported to concave positions.
The results of sensitivity analysis also showed that the hydrological properties of landscape such as profile curvature, stream power index, wetness index and plan curvature that are related to moisture distribution over the landscape, are the most important factors that influence SOC fractions. Also, some of terrain attributes such as slope and sediment transport index which are related to erosion processes, influence the SOC content in the study area.
CONCLUSIONS
The designed ANN models were able to establish the relationship between the terrain attributes, soil properties and remote sensing data with SOC fractions content. The developed models were able to explain a great deal of total variability of different SOC fractions in the studied site. Sensitivity analysis results showed that NDVI as indicator of vegetation coverage, TOC, CCE and clay content were the most important factors explaining the SOC fractions. Among the terrain attributes such as slope, STI and SPI that related to erosion were the most important factors that influence SOC.