Application of Wischmeier-Smith, EPIC, M-USLE, and WEPP Methods for Determination of Erodibility Factor (K) of Soil

ABSTRACT


INTRODUCTION
Agricultural development in Indonesia has increased the economy and people's welfare, but has had an impact on the decreasing land quality (degradation of land resources) (Septiyanti, 2022). Degradation of land resources is increasing as a result of population growth (Wahyunto & Dariah, 2014). Land degradation that occurs in Indonesia is caused by erosion accelerated by human activities which results in a decrease in the physical, chemical and biological quality of the soil, compaction of the soil surface, and a decrease in soil infiltration capacity so that the volume of surface runoff increases and has an impact on increasing river discharge and flooding (Kurnia et al., 2010).
Factors that cause and affect the amount of erosion are climate, soil, topography, vegetation cover crops, and human activities. Soil is an important factor in determining the amount of soil erosion. Each soil type has a different soil sensitivity to erosion. The sensitivity of soil to erosion, also known as the soil erodibility factor, is the property of whether or not a soil is easily destroyed by the force of falling raindrops or the force of surface runoff which is the interaction between the physical and chemical properties of the soil. Soil properties that affect soil erodibility are structure, texture, organic matter, and soil permeability (Satriawan & Fuady, 2015).
Calculation of soil erodibility values can be determined by laboratory analysis and using simulations. Soil erodibility values in an area can show varying values due to the use of different methods to calculate soil erodibility. Based on an international scale, the methods used to determine the erodibility of tropical soils vary widely, such as using the EPIC equation for the cases of Thailand (Liu et al., 2021) and India (Masroor et al., 2022) and using the Wischmeier-Smith equation for the case of India (Majhi et al., 2021) and Malaysia (Yusof et al., 2021). Regions with subtropical climates such as in China use the modified USLE (USLE-M) and WEPP methods (Wang et al., 2015) as well as the EPIC, Shirazi and Torri equations (Guo et al., 2022). The majority of research in Indonesia to determine the soil erodibility factor (K) uses the Wischmeier-Smith equation such as research by Ashari (2013), Kusumandari (2014), Agustina & Dewi (2020), Joniardi et al. (2020), and Djufri et al. (2021). In Indonesia, not many have used the EPIC, USLE-M, and WEPP methods to determine soil erodibility factors, while other tropical regions have used these methods. Therefore, it is important to conduct this research to apply and evaluate four methods of determining soil erodibility factors, namely Wischmeier and Smith, EPIC, USLE-M and WEPP, for tropical soils in Indonesia.

Time and Location
This research was conducted in September 2021 -June 2022. Soil sampling was carried out on moors, gardens, and shrubs located in Nglanggeran Village, Patuk District, Gunung Kidul Regency, D.I. Yogyakarta (Figure 1). The location of this research was geographically located at 110°32'20.101" East Longitude and 7°51'14.818" South Latitude. The climate in Nglanggeran Village was included in the tropical monsoon climate, where due to changes in wind direction, this village has a rainy season and a dry season. This climate normally got an excessive amount of rainfall during the rainy season and usually in the form of frequent thunderstorms.  The altitude was between 200 -700 masl, the average temperature was 23 -27 °C, the average rainfall was 3024 mm/y, the slopes were steep (> 45%), and the soil type was Latosol. Geologically, Nglanggeran is composed of andesitic lava igneous rock, volcanic sandstone, pumice breccia, shale and andesitic volcanic breccia.
Nglanggeran Village consists of 5 hamlets, namely Karangsari, Doga, Nglanggeran Kulon, Nglanggeran Wetan, and Gunung Butak. Soil sampling was carried out on disturbed and undisturbed top soil layer on three different land uses, namely shrubs, gardens, and moors (dry land field). Each type of land use was performed at 3 different locations as repetition. Map of soil sampling locations was presented in Figure 2. The condition of the land for sampling was described in Table 1. Table 1. Condition of the land surface for soil sampling based on different land use Soil analysis was carried out at the Soil Laboratory of the Faculty of Agricultural Technology (UGM) and the LaPitaya Laboratory, Yogyakarta. The erosion simulation was carried out in a green house located on Jalan Matraman, Maguwoharjo Village, Depok District, Sleman Regency, D.I Yogyakarta.

Soil Sampling
Soil sampling was carried out in a disturbed and undisturbed manner. Disturbed soil samples were taken using a shovel in the topsoil layer and undisturbed soil samples

Shrubs 1
The land has grass vegetation and shade due to large trees and bamboo Shrubs 2 The land has vegetation litter and shade due to large trees Shrubs 3 The land has vegetation litter and shade due to large trees Garden 1 The land with cocoa trees with vegetation litter was cleaned Garden 2 The land with cocoa trees with vegetation litter was not cleaned Garden 3 The land with cocoa trees with vegetation litter was cleaned Moors 1 The land was planted with sweet potatoes and cassava Moors 2 The land was planted with corn and cassava Moors 3 The land was planted with cassava and porang were taken using a ring sampler and hammer at a depth of 5-10 cm. The soil samples taken were put in plastic and labeled. Disturbed soil samples are used to test soil properties, including texture, structure, specific gravity, and organic matter content. Meanwhile, undisturbed soil samples were used for testing soil volumetric weight properties.

Soil Analysis
The soil properties included soil structure, texture, permeability, and soil organic matter. Soil structure was determined quantitatively by calculating the average diameter of the soil using a digital microscope. From the results obtained, the soil structure code was determined according to Arsyad (2010) using Table 2. Soil permeability was analyzed using the constant head permeameter method using Equation 1 (Hubbert, 1957).
where k is soil permeability (cm/h), Q is the amount of water flowing in each measurement (ml), t is time (h), L is thickness of soil sample (cm), and A is surface area of soil sample (cm 2 ). Table 3 provided soil permeability code according to Arsyad (2010). Table 2. Soil structure code according to soil structure type Table 3. Soil permeability code according to soil permeability class Soil texture was analyzed using the pipette method. The principle of the pipette method was to take soil samples from suspension with a depth of 20 cm for the first pipette to determine the content of the silt and clay fraction in the soil and 5 cm for the second pipette to determine the content of the clay fraction. Then the soil texture criteria was determined using the soil texture triangle from the percentage of sand, silt and clay fractions obtained. The organic matter content was tested at the LaPitaya Laboratory. (1)

Erosion Simulation
The steps for measuring sediment mass data and surface runoff using a rain simulator are as follows: 1. Conducted rain simulator testing. Rain intensity and rain distribution values were obtained based on trials using 3 erosion plots which were arranged systematically under a rain simulator within 15 minutes. The discharge and water that comes out were regulated to get the desired value of rain intensity and distribution of rain. A sketch of the rain simulator system was presented in Figure 3. 2. Plots measuring 32 x 45 x 20 cm were made conical with clear plastic so that no rainwater flows into the shelter. 3. Erosion plots were set with a slope according to the soil sampling location. 4. Rain simulator was operated for 30 minutes. 5. Soil sediment and running water were collected in a plastic cup. Sampling was carried out every 3 minutes in 30 minutes so that there were 10 repetitions on each plot. 6. The mass of soil and water were measured by separating soil and water sediments using filter paper. 7. The amount of soil loss and surface runoff rate were determined from the soil mass and water acquired in step 6.  Wischmeier & Smith (1978) proposed an equation to calculate soil erodibility factor as presented in the following: where K is soil erodibility factor (t.h.MJ −1 .mm −1 ), M is dimension of soil particle size, P s is the percentage of very fine sand and silt contents, P c is the percentage of clay content, a is the content of organic matter (%C × 100/58), b is the number of soil structure code (Table 2), and c is soil permeability class (Table 3). (2) (3)

EPIC Equation
The erodibility factor in the EPIC (Erosion Productivity Impact Calculator) method is calculated using the parameters of soil particle distribution and soil organic matter content using the equation (Sharpley & Williams, 1990): where K is soil erodibility factor (t.h.MJ −1 .mm −1 ), S is percentage of sand content, M is percentage clay content, F is the percentage of silt content, E is (1 -S/100), and C is percentage of carbon organic content.

M-USLE Equation
Erodibility factor from M-USLE (Modified Universal Soil Loss Equation) consider the effect explicit runoff for an event (Ke) and was calculated using Equation (5) (Kinnell et al., 2018): where K e is soil erodibility factor (t.h.MJ -1 .mm -1 ), A is the rate of soil loss (ton/ha), Q R is runoff ratio, E is kinetic energy of raindrops (MJ/ha), I 30 is maximum rain intensity for 30 min (mm/h), LS is factor of land slope and slope length, C is crop cultivation factor, and P is factor of soil conservation practiced.

WEPP Equation WEPP (Water Erosion Prediction Project) equation calculate erodibility factor (Ki)
according to the following relation (Kinnell, 1993): where K i is soil erodibility factor (t.h.MJ -1 .mm -1 ), D i is erosion rate, I is rainfall intensity, Q is runoff ratio, and S f is land slope factor.

Data Analysis
Statistical data analysis used to compare the four methods of determining soil erodibility was to use one way ANOVA. One way ANOVA analysis was used to determine whether there was a statistically significant difference between the means of two or more independent groups. The hypothesis used in this analysis was as follows: H0: There is no significant difference between the four soil erodibility methods H1: There is significant difference among the four soil erodibility methods At the significance level α = 5%, H0 was rejected if F Count ≤ F Tabel or if Sig. < α. In other words, H1 was accepted if F Count > F table .

Soil Properties
Data collection and measurement of soil properties in this study included soil structure, soil texture, soil volumetric weight or bulk density (ρ b ), soil particle density (ρ s ), soil permeability (k) and soil organic matter (OM). The results of measuring soil properties are shown in Table 4.

Table 4. Soil properties in Nglanggeran Village
Based on Table 4 it is known that soil structure in the three land uses have a very fine granular. This soil structure is able to absorb water into the soil properly so as to reduce surface water runoff and prevent soil erosion by surface water runoff. Soil texture in the land use of moors and gardens is clay while in the shrub is clay loam. Clay -textured soils have a higher percentage of clay as compared to silt and sand. Soils with a high clay content have more micro-pores than macro-pores, making it difficult for water to pass into the soil. This causes a slow rate of soil permeability. Whereas clay loam textured soils have almost the same percentage of sand, silt and clay fractions so they have a higher permeability rate than clay textured soils due to the balanced comparison of micro and macro pores in the soil.
The value of bulk density ρ b in this study ranged from 1.10 -1.16 g/cm 3 , this figure corresponds to the critical value for healthy agricultural soil, which is less than 1.20 g/ cm 3 for clayey soils (Brouwer & Jenkins, 2015). The greater the bulk density value causes the permeability rate to decrease and vice versa. The value of bulk density or soil density is also influenced by the amount of organic matter. The higher the organic matter, the lower the bulk density value because adding organic matter to the soil can increase the amount of soil pore space and form a crumbly soil structure so that it will reduce the bulk density of the soil and vice versa (Saputra et al., 2018). The more organic matter content contained in the soil, the smaller the particle density value (Hanafiah, 2005). The particle density ρ s values in this study ranged from 2.27-2.46 g/ cm 3 . Soil particle density is influenced by mineral content.
Soil permeability values in this study ranged from 0.08 -1.16 cm/hour (slow class). Many factors affect soil permeability, especially texture, structure, aggregate stability, porosity, pore size distribution and organic matter content (Mulyono et al., 2019). The impact of tillage can also affect the permeability capacity of the soil, compaction by rain, animals and heavy equipment can drastically reduce the ability of the soil to absorb water and close the pores of the soil (Mulyono et al., 2019). The value of organic matter in this study ranged from 2.08 -4.48%. Organic matter can improve soil aggregate, increase water holding capacity, aeration pores and infiltration rates and facilitate root penetration thereby increasing land productivity and crop yields.

Land Slope
Slope measurements in this study used an abney level and obtained the degree of slope shown in Table 5 the land use of moors, gardens and shrubs are 16.50 (30%), 14.89 (27%), and 25.17 (47%), respectively. The steeper the land, the soil will easily erode so that the soil erodibility factor becomes high. This is because the slope of the land will determine the magnitude of the velocity and volume of runoff water which affects the amount of the erosion value.

Rainfall
Precipitation measurements (in this case rain) were carried out using AWS (Automatic Weather Stations) and the results are shown in Table 6. Based on Table 6, it is known that the total rainfall value for January 2022 to April 2022 is 1088.40 mm and the average rain intensity is 29.28 -35.04 mm/day. According to the BMKG (2021) rain intensity of 20 -50 mm/day is categorized as moderate rain. Rain intensity values obtained from AWS is used as rain intensity in the erosion simulation process using a rain simulator tool. while the EPIC method only uses soil texture and organic matter parameters. Based on Table 7, the erodibility value using the Wischmeier-Smith method is lower than the EPIC, due to differences in soil structure and permeability parameters used in the calculations. The pattern of soil aggregate arrangement (structure) cannot be fully reflected in soil texture and organic matter content. Natural factors such as compaction by rainwater impacts and human or anthropogenic factors such as tillage practices and compaction by human traffic and transportation will determine soil ability in flowing water through either infiltration or percolation. When soil structure and permeability are included as parameters in the erodibility calculation, it will affect the value, whether it is relatively higher or lower depending on the nature of both. The Wischmeier-Smith method produces relatively smaller erodibility values than the EPIC method. If viewed from the soil structure, it should produce a K factor higher than EPIC. The very fine granular structure has a poor movement of water masses in the soil and this is reflected in its permeability, which is relatively slow. The impact is that the soil is easily saturated, the bonds between soil particles/aggregates are easily disrupted so that it is easily eroded or has low erodibility. Further research needs to be conducted to examine these findings. Soil erodibility calculations using the M-USLE and WEPP methods are based on measurements of erosion rates and surface runoff using a rain simulator. The resulting erodibility factor (K) ranges from 2×10 -4 to 10 -2 . This values are smaller than the K factor resulting from the Wischmeier-Smith and EPIC methods. There are factors of rainfall intensity and land slope that affect the calculation of the K factor by measuring erosion rates and surface runoff with rain simulation. In this study, the smallest rain intensity obtained from the rain simulator during the calibration process was 227.49 mm/h. On the other hand, the rain intensity obtained from 44 field data is 1.36 mm/h. Thus, the calibration of the rain simulator has not been successful and this condition will affect the measurement of erosion rate and surface runoff for each land use sample. With a much larger value of rainfall intensity, a much larger K factor should be obtained. However, the M-USLE and WEPP calculations yield a smaller K factor. Because the erosion plot which contains disturbed soil samples is not able to produce full volume so there is a room for inundation above the soil surface. This puddle is like a layer that protects the soil from rainwater blows, so the soil is more resistant to erosion energy (low K factor).
To find out whether there are differences in soil erodibility values with the four methods, a statistical test was carried out, namely homogeneity test and one way ANOVA with the results presented in Tables 8 and 9. Based on Table 8, the significance value is 0.105 and greater than 0.05 so that the data is homogeneous and can be followed by a one-way ANOVA test.

Land use
Value of soil erodibility using different methods Wischmeier-Smith EPIC M-USLE WEPP Moors 1.6 x 10 -1 2.9 x 10 -1 3.0 x 10 -2 1.0 x 10 -2 Gardens 1.2 x 10 -1 3.3 x 10 -1 1.0 x 10 -3 2.0 x 10 -4 Shrubs 2.2 x 10 -1 3.2 x 10 -1 1.0 x 10 -3 1.3 x 10 -3 Table 8. Results of homogeneity test of variance Table 9. Results of one way ANOVA test for soil erodibility value Based on Table 9, it is known that the significance value is 0.000 and less than 0.05 so it can be concluded that the average of the four soil erodibility methods is significantly different. This study focuses on the implementation and evaluation of the calculation of the erodibility factor, so that the difference in values is sufficiently evaluated with the one way ANOVA test. The evaluation of the soil erodibility model is reviewed by comparing the erodibility factors from the literature having similar conditions to the current study sites. Ashari (2013) stated that the erodibility factor of the soil in Nglanggeran Village for land use of moors (dry land), gardens and shrubs using the Wischmeier-Smith method is 0.17 -0.34. Soil erodibility factors in our study using the Wischmeier-Smith equation have values in that range so that the results obtained are in accordance with the literature. For the EPIC, M-USLE, and WEPP methods, no one has ever used these methods to determine the soil erodibility factor in Nglanggeran Village, so there is no reference to soil erodibility values as a comparison. The soil type approach is used to compare the soil erodibility factors of the four methods used. The type of soil in Nglanggeran Village is Latosol soil. According to Puslitbang Pengairan Bandung, Latosol soil has an erodibility value of 0.075 (Fahliza et al., 2013). In this study it was identified that there were no soil erodibility calculation results that were close to the reference values of Latosol soils. Therefore, the application of the EPIC, M-USLE, and WEPP methods for Nglanggeran tropical soil needs to be reevaluated, especially the technical aspects of measuring erosion rates and surface runoff.

CONCLUSIONS
The results of applying different methods to determine soil erodibility factors in tropical soil of Nglanggeran Village are 1.2×10 -2 to 2.2×10 -1 using the Wischmeier-Smith, 2.9×10 -1 to 3.3×10 -1 with EPIC, 10 -3 to 3×10 -4 using M-USLE, and 2×10 -4 to 10 -1 using WEPP method. Based on the literature (reference values from similar studies at the same site, or from soil types), the K factor using Wischmeir-Smith method is within the reference range. Of the four methods used in this study, all of them can be used to determine the K factor of tropical soils. However, the M-USLE and WEPP methods still have deficiencies in the process of simulating erosion rates and surface runoff which will have an impact on the calculation of the K factor.