Spatial variability of soil chemical properties of Moso bamboo forests of China

Regassa Terefe · Kun-yong Yu · Yangbo Deng ·Xiong Yao · Fan Wang · Jian Liu

Abstract This study investigates the spatial variability of soil organic matter (SOM), soil organic carbon (SOC) and pH in the upper 20-cm layer and 20–40 cm layer in Moso bamboo ( Phyllostachys pubescens Pradelle) forests using a geostatistics model. Interpolation maps of SOM, SOC, and pH were developed using ordinary kriging (OK) and inverse distance weighted (IDW) methods. The pH, SOC, and SOM of the two soil layers ranged from 4.6 to 4.7, from 1.5 to 2.7 g kg −1 and from 20.3 to 22.4 g kg −1 , respectively. The coef ficient of variation for SOM and SOC was 29.9–43.3%while a weak variability was found for pH. Gaussian and exponential models performed well in describing the spatial variability of SOC contents with R 2 varying from 0.95 to 0.90. The nugget/sill values of pH are less than 25%, which indicates a strong spatial correlation, while the nugget/sill values of SOC and SOM fall under moderate spatial correlation. Interpolation using ordinary kriging and inverse distance weighted methods revealed that the spatial distribution of SOM, SOC, and pH was inconsistent due to external and internal factors across the plots. Regarding the crossvalidation results, the ordinary kriging method performed better than inverse distance weighted method for selected soil properties. This study suggests that the spatial variability of soil chemical properties revealed by geostatistics modeling will help decision-makers improve the management of soil properties.

Keywords Cross-validation · Geostatistics · Inverse distance weighted · Ordinary kriging · Semi-variance

Introduction

Soils are critical resources for the growth and development of plants. Improvement of soil properties plays a signif icant role on the restoration of vegetation as key factors of structure and function for a degraded environment. Factors contributing to declining soil fertility are poor soil management,land use changes, exploitation of marginal lands, erosion,and nutrient leaching (Lemenih and Kassa 2014). For plant growth and the transformation of nutrients in soil ecosystems, soil properties provide a suitable nutrient environment and chemicals (Tripler et al. 2006; Bhunia et al. 2018a).Several studies have been applied to dif ferent soil types for the evaluation of spatial variability of soil properties under various management practices (Jin and Jiang 2002; Rodríguez et al. 2009; Najaf ian et al. 2012; Shit et al. 2015;Bhunia et al. 2018a). Identify the spatial distribution of soils with dif ferent properties for a specif ic area will create better conditions for farming and environmental management practices through conserving soil fertility (Cambardella et al. 1994; Marklein and Houlton 2012).

The measurement of soil properties over a wide area under dif ferent agro-ecological practices is time- consuming and expensive. Geostatistical techniques are one of the most ef fective methods for analyzing the spatial distribution of soil properties (Liu et al. 2014). Using a geostatistical method, specif ically for spatial interpolation, the values of soil properties can be estimated for locations where ground observations are not available (Goovaerts 1999; Pebesma 2006; Reza et al. 2010). Dif ferent interpolation techniques,such as ordinary kriging (OK), co-kriging (CoK), inverse distance weighting (IDW) and artif icial neural networks,have been developed (Zhao et al. 2010; Li and Shao 2014;Gao et al. 2016; Ghorbanzadeh et al. 2019).

John et al. ( 2019) used the OK method to study the spatial variability of dif ferent soil properties in Ovia North-East,Nigeria. Bhunia et al. ( 2018b) reported the assessment of spatial variability of nitrogen, pH, phosphorus, potassium,electrical conductivity and organic carbon in lateritic soils of West Bengal, India using the OK method. Yang et al.( 2016) reported on the spatial distribution of SOC and total nitrogen for farmlands in northern China using this method and Yao et al. ( 2019) studied the spatial variability of soil organic carbon in Changting city, Fujian province,China using both the OK and IDW methods. In recent years,researchers have undertaken numerous experiments on spatial variability of soil properties using various techniques.However, due to soil variability, ecological diversity, forest species distribution, and other related factors, not all areas have been covered. The nature of soil is dynamic, which requires further spatial variability evaluation under dif ferent agro ecology systems, soil types, and forest species.

Moso bamboo (Phyllostachys pubescensPradelle) is one of the most common, useful and dominant forest species in the Tianbaoyan Nature Reserve, Fujian Province. This area mountainous, with soils highly exposed to degradation,leaching, and serious erosion from intensive and extensive human activities (Du et al. 2015). Consequently, intensive soil management with fertilization, site reclamation, and regular understory removal is increasingly common in southern China, especially in the main bamboo- producing provinces of Zhejiang and Fujian (Zhou et al. 2006). However, there is insuf ficient information on spatial variability of important soil chemical properties in response to the environment and to cover species. There are def iciencies in the interpolation map for the management of soil chemical properties in the Reserve using the dominant species, Moso bamboo. A better understanding of the spatial distribution of soil chemical properties with specif ic species is required for sustainable land management, for the determination of inputs and for improving productivity of the Moso bamboo forests. Consequently, improving knowledge of the spatial distribution of soil organic matter (SOM), soil organic carbon (SOC) and pH is essential for sustainable land management and maintaining or improving Moso bamboo forest productivity. Spatial variability information of soil properties depends on the specif ic soil studies, but the soil of Moso bamboo forests in the Tianbaoyan Nature Reserve lacks such information and needs to be determined. Therefore, this research evaluates the spatial variability of SOM, SOC, and pH in the upper 20-cm layer and the 20–40 cm layer in Moso bamboo forests and examines the spatial distribution of SOM, SOC, and pH using a geostatistical technique.

Materials and methods

The study was carried out in Moso bamboo forests at the Tianbaoyan Nature Reserve of Yong’an, Fujian Province,China, at 25° 55′–25° 58′ N latitude and 117° 31′–117° 33′E longitude. The area has remnants of Daiyunshan Mountain ecosystems, with low to middle altitudes ranged from 580 to 1605 m a.s.l. (Liu et al. 2010). The climate is subtropical monsoon with average temperatures of 23 °C, a minimum temperature of − 11 °C, and a maximum temperature of 40 °C. Average annual precipitation is approximately 2000 mm. The dominant soil types are Oxisols. The major species areP. pubescens, among which the mixed tree speciesChoerospondias axillaris(Roxb.) andLiriodendron chinense(Roxb.) belongs to high dense plantations in the Moso bamboo forest. To improve stand productivity and economics value, application of fertilizer, removal of the understory, and timely harvesting has been widely practiced in these forests. Fertilization was carried out over the whole study area using the protocol established by the local Forest Bureau (Fig. 1).

Soil sampling and analyses

One hundred and thirty-four 25.8 m × 25.8 m sample plots were established and soil samples were collected from each plot in August 2018. In the center of each plot, a soil prof ile was established to a depth of 40 cm in two layers, the upper 20-cm and the 20–40 cm. Bamboo roots are mostly distributed in the upper 40-cm layer for the purpose of this study to estimate the spatial distribution of soil chemical properties. From each plot, soil samples taken from the same layers were bulked together. Geographical coordinates for all collected samples were recorded using a portable GPS. The samples were registered, numbered, and pretreat for analysis.All samples were air-dried, impurities removed, and sieved through 2 mm and 0.2 mm mesh sieves before analysis.

Soil carbon was determined with a CN elemental analyzer (Vario Max, Elementar Co., Ltd., Langenselbold,Hesse, Germany) using a 0.1 g sample passed through a0.2 mm sieve and ignited for 15 min at 1000 °C. Soil organic matter was determined using oxidation by heating with K2 Cr2 O7solution in a concentration of H2 SO4(0.1 g f ine soil K2 Cr2 O7, 5 ml concentrated H2 SO4, 0.1 g f ine soil passed through a 0.2 mm sieve, digested for 5 min at 250 °C and titrated with 0.2 mol L−1FeSO4mixed with concentrated H2 SO4) (Tang et al. 2017; Yu et al. 2019). Soil pH was measured in a water and potassium chloride (1 M KCl) suspension in a 1:2.5 (soil: liquid ratio) potentiometrically using a glass-calomel combination electrode (Bielders et al. 2002).

Fig. 1 Location of the study area and sample plots in Fujian Province, China

Statistical and geostatistical analyses

To maintain a normal distribution of measurement data, outliers are removed before descriptive statistics. Outlier data refers to the extreme values caused by measurement error or by any other way those results arise in the data distribution.It is def ined as extreme values that were caused by errors in soil sampling and analysis. Extreme values were considered as any value is not within three standard deviations of the mean (± 3 s, whereand s are the average value and the standard deviation for each soil variable) (Rosemary et al.2017; Yao et al. 2020). Statistical data such as minimum,maximum, average, standard deviation, skewness, kurtosis and coef ficient of variation are determined. Therefore,spssv.21 software was used for descriptive statistical analysis. The spatial distribution of soil properties in un-sampled areas was evaluated by geostatistics (Goovaerts 1999).Semi-variogram describes the spatial correlation between data values. Ordinary kriging (OK) was performed in GS+v10.0, and interpolation maps of selected chemical properties using OK and inverse distance related (IDW) methods were generated in ArcGIS v10.2.2.

Ordinary kriging (OK)

A best f it theoretical model of semi-variogram is used in the OK and IDW methods to forecast the distribution of soil properties (Saito and Goovaerts 2000; Pang et al. 2011;Yao et al. 2020). The experimental semi-variogram is calculated using the known equation as indicated by Schöning et al. ( 2006). There were three models used for this study:Gaussian (G), exponential (E), and spherical (S) to explore the spatial continuity of soil distribution and various soil properties (Goovaerts 1999). Based on the highest values of the determinant coef ficient (R2) and the minimum values of the residual sum of squares (RSS), our best-f itted model was selected according to Guan et al. ( 2017).

Nugget (C 0), sill (C 0+C), and range (A) are the three parameters acquired from the selected models used for the determination of the spatial composition of soil chemical properties. The ratio of the nugget to sill is used to describe the degree of spatial correlation of soil properties and the range is used to ref lect the spatial autocorrelation distance of soil properties (Ye et al. 2018; Yao et al. 2020).C 0represents the size of the error, but is generally not zero because of unavoidable errors in the process of sampling and measurement.C 0+Cis the semi-variance value of the semi-variogram

when it reaches the maximum range withh, representing the overall characteristics of regionalized variables, and is the spatial structure value and the sum of nuggets.

Cross-validation model

The ef ficiency and error of the soil properties prediction map were determined by the cross-validation method (Zirlewagen and von Wilpert 2004; Bhunia et al. 2018b). The semi-variance function of geostatistics is used to f it the spatial variability of soil properties and to estimate unknown points within the f itting range (John et al. 2019; Yao et al.2019). For the process of cross-validation, the root means square error (RMSE) and the mean error (ME) of the model were also calculated. As the value of the RMSE approaches zero, it indicates the accuracy of the prediction of the model.

where,nrepresents the number of points,PiandMiare the predicted and measured values of the observations.

Results

Statistical description of pH, SOC and SOM

The statistical summaries of soil chemical properties in the upper 40-cm layers are shown in Table 1. The pH in the Moso bamboo forest increased from 4.1 to 5.2. Soil organic carbon (SOC) went from 0.5 g kg−1to 4.5 g kg−1.The minimum value of soil organic matter (SOM) was 1.8 g kg−1and the maximum 46.3 g kg−1. The coef ficient of variation (CV) was the lowest for pH (4.2%) and was highest for SOM (43.3%). There are minor dif ferences in soil chemical properties because of erosion, and forest degradation aggravates the loss of organic matter and organic carbon. The distribution of SOM, SOC, and pH varied from slightly negatively skewed (skewness ≤ –0.95)to moderately positively skewed (skewness > 0.74). As the soil depth increased, the variation of SOC and SOM decreased (Table 1).

Geostatistical analyses of SOC, SOM and pH

With soil organic carbon, the determination of coef ficient(R2) is the highest at 0.95 and 0.90 at both soil depths,respectively. The R2value of soil organic matter at both depths was the lowest at 0.15 to 0.54, while the R2value of pH was moderate. The range of variation in soil indices was moderate with dif ferences from 1.55 to 37.26 m. A nugget value is commonly regarded as the degree of spatial variability resulting from random factors, e.g., sampling errors (Goovaerts 1999). In this study, nugget values in soil pH, ranged from 0.002 to 0.009, SOM values from 18.3 to 47.7 in dif ferent soil layers.

To compare the performance of each interpolation approach, AME, ME, RMSE and R2were calculated(Table 2). ME values are below 0.05 and in the range of 0.001 to 0.034. In addition, SOC, SOM, and pH values were close to zero, indicating that the OK and IDW models were relatively unbiased. With regards to the RMSE prediction, the values are close to 1, which is good for the interpolation of each soil property. However, all indices express that the results are consistent for the predictable maps of soil properties from OK and IDW.

Table 1 Summary of descriptive analysis of SOC, SOM, and pH

Table 2 Accuracy assessment of OK and IDW methods for SOC,SOM, and pH

Discussion

SOC, SOM and pH in Moso bamboo forests

Soil chemical properties are essential for the growth and development of dif ferent plants and crops. This study evaluated the dispersion of variables by using the value of the coef ficient of variation (CV). In soil science, the CV can be divided into three categories. When it is ˂ 10%, the variability is weak; a CV between 10 and 100%, variability is moderate; and when the CV is100%, the variability is strong (Fang et al. 2016). Our f indings indicated that the CVs varied from 4.2% to 43.3%. The lowest was with pH,which means the pH of the study area soils had weak variability. For the other soil chemical properties, the CVs values indicated moderate variability.

From the upper 20 cm layer, pH ranged from 4.1 to 5.1, while in the 20–40 cm layer, it ranged from 4.1 to 5.2. According to these results, we noted that as the depth increases, pH also increases (Table 1). The pH of the study area categorized as highly acidic, which might af fect the growth and productivity of Moso bamboo forests. Most of the farmers in the study area they prefer to apply fertilizer for two main reasons: (1) to neutralize the soil; and, (2) to enhance soil fertility to increase productivity. Soils with high clay and organic matter contents are more able to adapt to a drop or rise in pH compared to sandy soils. Soil organic matter content also changes under dif ferent management activities. Soil organic carbon is important for agricultural production because it improves soil structure and cation exchange capacity, and water-retention capacity, thus having positive impacts on soil fertility (Liu et al. 2019). The results here show that organic matter and carbon levels were lower in the study soil because of the higher pH. As the soil becomes more acidic, the availability of organic matter and carbon decreases. Soil pH also signif icantly af fects SOC contents because it regulates nutrient availability, organic matter turnover, and a series of soil processes (Kemmitt et al. 2006).

Spatial variation of SOC, SOM and soil pH

Our study used the semi-variograms model for the analyses of selected soil chemical properties (Fig. 2). The highest determination of coef ficients (R2) was recorded for SOC for both soil depths 0.95 and 0.90, respectively (Table 3).Having high R2values indicates that the model has a good semi-variance structure and can better ref lect the spatial distribution. Soil pH and SOM had average R2values of 0.63 in 0–20 cm, 0.74 in 20–40 cm, 0.54 in 0–20 cm, and 0.15 in 20–40 cm, respectively. The determination of coef ficient values of all soil properties revealed a good semi-variance structure along with models at both depths except for SOM in the second layer. In comparing models, the highest R2(0.95 and 0.9) values were obtained under the Gaussian and exponential models for both soil depths in SOC (Table 3).

The results indicate that the nugget value was lower in pH and higher in SOM in both layers, mainly due to the larger magnitude of soil organic matter. These positive variations highlight the positive ef fects of nugget values on the inherent variability of nutrient sampling errors and any other environmental factors. In addition, the nugget values decreased as soil depth increased. The spatial variability of soil properties is expressed by the sill values (Liu et al. 2014), and these ranged from 0.04 to 95.5, and decreased with depth similar to nugget values.

The spatial correlation is represented by the value of the nugget/sill ratio. If the ratio is < 25%, it is strong, between 25% –75% it is moderate and > 75% it shows a weak spatial correlation (Fang et al. 2012; Wang et al. 2016). The values in this study of the nugget/sill ratios showed a strong correlation with pH. However, moderate spatial correlations were found for both SOM and SOC. A strong spatial correlation is attributed to inherent soil properties such as parent materials, texture, topography, and vegetation as described by (Fu et al. 2014). External activities such as fertilizer applications and management practices primarily regulate the spatial variability of soil and indicate a weak spatial correlation (Bhunia et al. 2018b). Moderate spatial correlation is controlled by both internal and external factors. There was a strong spatial correlation for pH at both depths with 17.4% in the upper 20 cm and 6.0% in the 20–40 cm layer, while for soil organic matter, the upper 20 cm layer showed a moderate spatial correlation at 50.0% and a strong spatial correlation for the 20–40 cm layer of 22.1%. For soil organic carbon at both depths, there was a moderate spatial correlation of 35.8% in the upper 20 cm layer and 27.6% in the 20–40 cm layer (Table 3). The strong spatial correlation for pH might be due to the inf luence of inherent factors. In the case of soil organic carbon, the inf luence was from both inherent and external factors. This might be related to the application of fertilizers and to management practices used for the enhancement of the productivity of Moso bamboo forests.The semi-variogram model is used for the estimation of spatial distribution values of distance and inverse function of spatial covariance (Bofana and Costa 2017). Thus, the variability of values increases with increasing distance between point locations, and it is also the same for spatial autocorrelation and semi-variogram (Nayanaka et al. 2010; Bofana and Costa 2017).

Cross-validation of OK and IDW

Fig. 2 Semi-variograms with f itted models at 0–40 cm depths for SOC, SOM, and pH. The parameters are given in Table 3

Table 3 Semi-variogram indicators of best f itted models for SOM, SOC, and pH contents

The cross-validation process used to test the semi-variogram model is validated by using the OK and IDW methods using the sample plots referring to the nearest samples(Tang et al. 2016; Bhunia et al. 2018b). The prediction values for OK and IDW were plotted against the measured values to evaluate the interpolation performances displayed (Fig. 3). Except for SOC, the OK method performed better for SOM and pH showing higher coef ficient of determination (R2). In this study, the cross-validation values generally were found in the middle range, are similar to other studies conducted in dif ferent agro-ecologies.

Spatial interpolation of SOC, SOM, and pH

Fig. 3 Cross-validation of OK and IDW methods interpolation of SOC, SOM and pH at 0–20 cm and 20–40 cm depths

The OK and IDW interpolation methods were used to predict the spatial variability of SOC, SOM, and pH (Fig. 4 a–l).Based on the interpolation of geostatistics, the position points are estimated by the measured sample values. The essence of this method is to estimate the weighted average of the local estimation, even if the unknown points distributed on the maps are estimated with known measured data points(Sigua and Hudnall 2008; Bhunia et al. 2018b). This study is in agreement with other studies that the number and area coverage of sampling points have an impact on OK and IDW interpolation accuracy (Liu et al. 2010; Bhunia et al. 2018b;Yao et al. 2019). The more sampling points the more precise the interpolation results. Therefore, it is better to use more sampling points than fewer.

The interpolation maps of SOM, SOC, and pH revealed modest inconsistency distributions. The distribution of SOM was relatively higher in the north and lower in the south, whereas an average was found in the central part of the area (Fig. 4 a–d). Comparing the ordinary kriging and inverse distance weighted methods for SOM showed some variation in interpolation maps. For OK, the prediction maps show the distribution of the point samples concentrated in the north and northwest, whereas, with the IDW method, the concentration was found in the central part of the area. In the case of SOC, the predictions of samples for both methods were almost the same (Fig. 4 e–h). The distribution of pH indicated strongly acidic according to the USDA soil pH classif ication, and is similar to the f ield data.

Fig. 4 Spatial interpolation of SOM, SOC and pH, using IDW and OK

Overall, the distribution of SOM, SOC, and pH was inconsistent due to external and internal factors. Among the main factors, the ef fect of leaching, forest degradation, severe erosion and strongly acidic soils are some of the main factors.Although high pH mainly related to intensive farming activities and undulating topography, soils are exposed to leaching and erosion which accelerates the removal of the topsoil (Liu et al. 2006).

Conclusions

This study provides spatial distribution maps for soil organic matter, soil carbon and pH in a Moso bamboo forest. There were inherent factors af fecting the pH with both inherent and external factors af fecting soil organic carbon. These are due to human interventions through the application of fertilizers and management practices. There was a better cross-validation performance in producing the spatial distribution maps for soil organic matter and pH in comparing IDW to OK.

This f inding has resulted in some changes in soil properties which showed modest inconsistent distribution throughout the sample plots. The inconsistency of soil properties distribution might be the result of internal and external factors. However, to determine the ef fects of these variables on spatial variability of soil properties, additional research is needed. This study provides valuable information for the sustainable development of Moso bamboo forest in the study area.