Soil Use Legacy as Driving Factor for Soil Erosion under Conservation Agriculture

Introduction

Soil erosion remains a major challenge worldwide, being the greatest threat for sustainable soil management and subsequent food production (Rickson et al., 2015). Soil erosion, its control and remediation practices are related to nine out of 15 sustainable development goals defined by the United Nations (FAO, 2019). This highlights the worldwide preoccupation about soil loss and its immediate and indirect consequences fueling the sustainability discussion. Changes in land use and cover generated by the disturbance of natural grasses are the main factor leading to accelerated soil erosion (Borrelli et al., 2017). In Uruguay, land use changes occurred since the early 2000s in form of agricultural intensification shifting towards continuous, annual cropping systems and simplified rotations which led to fresh water pollution and soil quality deterioration (Carrasco-Letelier and Beretta-Blanco, 2017; Ernst et al., 2018). Furthermore, Uruguay is affected by climate change, mainly shown by increasing precipitation rates and extreme weather events that directly alter water erosion risks (Munka et al., 2007; PNUD (Programa de las Naciones Unidas para el Desarrollo) Uruguay, 2007). The combination of deteriorated soils with more intense rainfall events increases the risk of soil erosion compromising soil and water quality.

Soil organic carbon (SOC) losses have been reported as a result of agricultural intensification increase in Uruguay (Beretta-Blanco et al., 2019; Grahmann et al., 2020) which could be related to incessant soil losses by erosion although soil conservation practices were implemented by law. SOC is the main determinant of aggregate stability of Mollisols (Novelli et al., 2013; Rubio et al., 2019), the most important agricultural soil type in Uruguay. Added to this, water erosion preferentially removes the light organic fraction of low density and therefore contributes to SOC pool depletion (Lal, 2003) and atmospheric CO₂ emissions (Lal, 2019). Soil erosion is strongly coupled with other nutrients’ runoff like nitrogen (N) which has consequences for freshwater and marine ecosystems leading to contamination, algae growth and overall biodiversity decline (De Vries et al., 2013). The particulate organic matter (POM) and potentially mineralizable N (PMN) were found to be sensitive indicators to detect the vulnerability for water erosion and overall soil degradation (Wander et al., 1998; Fabrizzi et al., 2003). Other parameters to monitor soil quality decrease and nutrient losses via erosion are total, fixed and volatile suspended solids and the carbon (C) and N enrichment ratios of the transported sediments (Palis et al., 1990; Holz and Augustin, 2021).

Conservation agriculture with its three principles of minimum soil tillage, crop rotation and residue retention is an often reported management system to effectively reduce soil erosion (Ernst and Siri-Prieto, 2009; Willett et al., 2019). To minimize soil erosion, soil cover must be coupled with a stable soil structure. As the basic unit of soil structure, stable aggregates are key for soil quality, soil fertility, and resistance to degradation (Le Bissonnais, 1996). A stable soil structure is crucial for maintaining soil porosity, gas exchange, water infiltration, erosion resistance, and SOC sequestration (Six et al., 2004; Morris et al., 2019). A reduction in tillage intensity increases soil structural stability, improves aggregate distribution and enhances soil porosity and therefore infiltration, thus decreasing soil erodibility, improving soil fertility and improving overall agronomic productivity (Azooz and Arshad, 2011; Verhulst et al., 2011; Palm et al., 2014). Straw residues on the soil surface act like barriers that prevent the soil from receiving directly the high kinetic energy of raindrops during heavy rainfalls (Turtola et al., 2007), hence reducing surface soil dispersion and crust formation. In addition to this direct effect of residue retention on soil infiltration, conservation agricultural practices also contribute indirectly through improved structural stability, bulk density and pore structure (Kumar and Goh, 1999; Zhang et al., 2007; Ranaivoson et al., 2017).

Long-term experiments are important to understand how conservation agriculture modifies the dynamics in biochemical and geophysical processes over decades, which in consequence affect soil erodibility (Richter et al., 2007; Johnston and Poulton, 2018; Grahmann et al., 2020). As for surface runoff and erosion determination, surface runoff plots that count with a predefined, surrounded area and water collection tanks, so called “Wischmeier” plots are required (Kinnell, 2016). Experimental facilities of long-term surface runoff plots are used to document the effects of conservation agriculture on soil erosion providing a valuable, but seldom available experimental data set.

Since 2013, it is mandatory for Uruguayan farms larger than 50 ha to present a land use and management plan that ensures an estimated average soil erosion rate smaller than a previously defined soil-specific tolerance threshold (MGAP, 2013). The application of this law, which was updated lastly in 2018, aims to promote conservation agriculture practices reducing the risk of soil erosion in cropping systems based on the Revised Universal Soil Loss Equation (RUSLE) (Pérez-Bidegain et al., 2018; Zurbriggen et al., 2020). The RUSLE model was calibrated and validated for Uruguay and available in an adapted software tool (EROSION 6.0) (Renard et al., 1997; García Prechac et al., 2017; Pérez-Bidegain et al., 2017). Hereof, it is important to emphasize that the current K factors of the Uruguayan RUSLE, which indicates how susceptible a soil is to water detachment, was calculated by Puentes (1981) based on soil characteristics determined in a soil sampling campaign conducted in 1976 (MGAP-DSA, 1976).

While the use of simulation models to predict soil erosion has been widely applied in the country, the actual state of soil quality was not considered yet. As the supposedly stable K factor was defined in the 1980s and was calculated considering soil texture, total soil carbon content and two classification factors associated to soil structure and permeability, this might lead to an underestimation of erosion in today´s degraded soils and does not consider the dynamics of soil aggregate stability, the history of land use and the subsequent SOC losses at the regional scale (Alewell et al., 2019). The deterioration of soil properties caused by their past use and management, and particularly tillage, can have a great impact on erosion (Beniston et al., 2015) and has rarely been studied in a remediation context. Soil management practices, like deep ploughing or subsoil loosening lead to irreversible soil structure changes, which in turn affect soil quality over many decades (Schneider et al., 2017) and cause persistent soil legacy effects of previous land use. This may also account for the cover management factor (C factor), which is directly influenced by the vegetation type, growth stage of the vegetation and root mass, not only from a yearly perspective, but over a longer period of time (Gyssels et al., 2005; Panagos et al., 2015a). In the current study, the history of soil use and crop management intensity for each experimental unit was transformed to a newly developed index of agricultural intensification (IAI) combining the number of years under tillage activities and crop production with the length of agricultural production periods (intensification sequence index, Caviglia and Andrade, 2010; Martinez et al., 2020).

The objectives of this study were 1) to quantify surface runoff, soil and nutrient losses in a continuous cropping system under conservation agriculture implemented in long-term surface runoff plots and 2) to determine the variability in physical, chemical and hydraulic soil properties due to soil use and management history. We hypothesized that current conservation practices reduce soil losses to a minimum and assumed that differences in physical and chemical soil parameters between plots were existent and generated by differences in historical land use management which in turn explains the variability of soil losses between the plots nowadays. Furthermore, this paper provides starting points to discuss an update of the erodability factor K that might be affected by previous tillage and crop management through effects on SOC and soil structure.

Materials and Methods

Experimental Site

The study was conducted at INIA La Estanzuela, Uruguay (lat. 34°20′S, long 57°41′W, 66 m a.s.l.), in a soil classified as Vertic Argiudoll (USDA Soil Taxonomy; Soil Survey Staff, 2014). The site has a warm temperate climate with annual rainfall of 1,125 mm and an annual reference evapotranspiration of about 1,180 (± 143) mm over 55 years (Grahmann et al., 2020). Maximum monthly mean temperature ranges from 15°C in the winter to 28°C in the summer, and the minimum monthly mean temperature varies from 6°C in the winter to 17°C in the summer (Baethgen et al., 2021). Available water holding capacity was 92.7 mm in 0–56 cm depth for undisturbed soil (Hill et al., 2008). In 0–15 cm soil depth, the particle size fraction contains 178 g/kg sand, 450 g/kg silt, and 372 g/kg clay and soil organic carbon (SOC) averaged 2.93% (± 0.09) in 0–7.5 cm soil depth for the experimental site in 2019.

Long-Term Wischmeier Plots

The long-term experiment was installed in 1984 when six Wischmeier runoff plots were built with a size of 22.1 m length and 3.5 m width along an agricultural field´s slope side. The difference in height between top and bottom end of each plot was measured by means of an automatic level and level staff in 2020, and the resulting slope ranged between 4.3 and 4.7%. Each plot was completely surrounded by a metal border of 10 cm height of which 5 cm were inserted into the soil. The lower side ended up in a funnel system which transported the runoff water to a first cement tank of 750 L capacity embedded in the soil. During extreme rainfall events, the overflow system led through a fractionator, with one seventh (1/7) of water passing to a second tank of 1000 L capacity. The total runoff water storage capacity captured up to 100 mm runoff events in each plot.

A schematic overview of the historical treatments and land use changes for each experimental period and plot can be found in Figure 1. Since the beginning of the experiment, four experimental periods were distinguished. In the first experimental period (1984–1997), half of the plots were managed with conventional tillage using a chisel plough to 10–15 cm soil depth. Additionally, two different crop rotations were investigated: continuous cropping and crop-pasture rotations. The treatments were randomly distributed. The data of this experimental period were used to validate the USLE model for contrasting soil use systems (García-Préchac, 1992). Between 1998 and 2008, the experiment was discontinued due to financial constraints, leaving the plots without any use under perennial grass cover. Since 2008, the second experimental period, four of the plots were recovered to test a soybean (Glycine max)-fallow rotation under no-till management focusing on the assessment of the environmental impact of agrochemicals. In the third experimental period, cover crops rotating with soybean were incorporated in half of the plots from 2012 to 2014 and the research focus was set on phosphorus runoff (Lizarralde et al., 2015). In the fourth experimental period which started in 2015, two more plots were taken into operation with continuous crop rotation of soybean – cover crop (Black oats (Avena strigosa))—soybean – cover crop (Ryegrass (Lolium multiflorum))—soybean – fallow – maize (Zea mays)—wheat (Triticum aestivum L.)—soybean – cover crop (Rye (Secale cereale) and vetch (Vicia villosa))—soybean. Since 2015, all six plots were used and managed uniformly (“uniformity trial”) with best soil conservation practices using no-till, residue retention and cover crops.

FIGURE 1. Experimental design of the Wischmeier plots and their history of soil use and management since the experimental onset in 1984.

Sowing and harvest operations were conducted with experimental agricultural machinery, pesticide applications were conducted manually according to best management practice. The cover crop mix of rye and hairy vetch was sown in May 2019 and roller-crimped in November 2019, while the previous cover crops were terminated with glyphosate. Fertilizer application in form of solid urea and diammonium phosphate during the last six study years was identical in each of the plots and was applied according to previous soil tests (Table 1).

TABLE 1. Planting and harvest dates and average fertilizer rates per year and crop since uniform management started in 2014 (CC: cover crop).

Index of Agricultural Intensification

Information regarding crop rotation and tillage intensity for each plot since 1984 were used to estimate an index of agricultural intensification (IAI), a way of classified scoring for historical soil use intensity. The IAI is estimated as the product of two indices, one accounting for cropping sequence intensification (ISI) and the second one for tillage intensification (ITI). Land use intensity is estimated using the index of sequence intensification (ISI) and is calculated on a yearly basis to account for the number of crops with active plant growth during the study period (Caviglia and Andrade, 2010; Novelli et al., 2013). The index enlarges with an increasing number of crops planted per year. To account for the reported beneficial effects on soil quality and structure through pasture incorporation into crop rotation (Studdert et al., 1997; Ernst et al., 2018), a coefficient of 0.5 was applied during the first experimental phase to the cropping sequence length when pastures were included (Wischmeier and Smith, 1978).

Tillage intensity was evaluated using a specifically developed index of tillage intensification (ITI). This index evaluates the total number of tillage events and weights them based on the C tillage factors provided by Panagos et al. (2015b). A weighting factor of 1 was applied to conventional tillage (conform to the C factor for tilled bare soil in the RUSLE, Renard et al., 1997); 0.25 for no-till and 0.01 for broadcast sowing of cover crops by hand causing minimal soil disturbance.

The IAI was calculated for each of the four experimental periods and accounted for fallow time when certain plots were discontinued (Table 2). The discontinuation period of all plots between 1998 and 2008 was not included in the IAI determination. The second and fourth experimental phases have identical IAI for all plots as they were managed uniformly over the respective periods. Table 2 demonstrates that the ISI had a higher weight than the ITI due to intensive cropping with up to two crops per year. To strengthen index interactions and its hazard effects on soil quality, ISI and ITI were multiplied (IAI, Table 2) as undertaken in the Soil Quality Rating (SQR) procedures presented by Mueller et al. (2007). Otherwise using sums, NT_Mix and TIL_CROP would have led to a similar IAI, although their soil use history was disparate. Especially during the first experimental period between 1984 and 1997, IAI differed between plots with lowest IAI in plot 4 and 6 and highest in plot 3 and 5.

TABLE 2. Land use history expressed as index of agricultural intensification per plot and experimental period (ISI: index of sequence intensification, ITI: index of tillage intensification, IAI: index of agricultural intensification).

Due to the limited number of long-term Wischmeier plots, two plots with similar land use history and IAI were grouped into three treatment groups (Table 3). The no-till treatment group NT_Mix had differences in historical cropping intensity with crop-pasture rotation in plot 4 and continuous cropping in plot 6.

TABLE 3. Treatment groups according to their multiplier IAI criteria based on the summed experimental period averages of the index of sequence intensification (ISI) and the index of tillage intensification (ITI).

Continuous Measurements

Daily and 10 min record precipitation data were provided by the Research Unit of Climate and Geographic Systems, GRAS, INIA, Uruguay to calculate rainfall intensity parameters during 2017–2019 (available at http://www.inia.uy/gras). Meteorological data were obtained from an automated weather station (Campbell Scientific, Inc, Logan, UT, United States) equipped with pluviograph and pluviometer located approximately 700 m from the experimental site. According to Sasal et al. (2010), daily rainfall events were classified into three groups: small (< 40 mm), intermediate (41–69 mm) and large (>70 mm). After each rainfall event, the amount of precipitation and the volume of runoff water per plot were recorded. The water height in the tank was recorded and the volume of runoff water was calculated using previously calibrated conversion factors for each tank and plot. The calibration was carried out yearly, filling the tanks with a rising, precise amount of water and simultaneously recording the water table height in each of the six tanks. The runoff coefficient was determined for each rainfall event as the ratio of runoff per rainfall.

Historical data on water runoff, sediment losses, and SOC were continuously measured since the installation in 1984 and available for the present study. Between 2017 and 2019, water samples were taken after each runoff event to evaluate the amount of fixed (mineral material), volatile (organic material) and total suspended solids with the gravimetric method Nr. 2540 of the American Public Health Association (APHA et al., 2012). For this, the runoff water was homogenized in the tank and a maximum 1,000 ml water sample was taken; for small rainfall events, a smaller sample proportion was available. In the laboratory, 100 ml sample were added to previously washed and dried glass microfiber filters <0,2 µm (934-AH Whatman™; Little Chalfont, United Kingdom) using a Büchner funnel and a vacuum pump with an adjusted pressure between 10 and 20 mbar. Total solids were determined by drying the filter at 105°C for 4 h and weighed afterwards with an analytical precision balance. In a subsequent step, the same filter was incinerated at 550°C for 2 h, and after cooling to room temperature, weighed a second time to determine fixed suspended solids. Volatile suspended solids were calculated as the difference between total and fixed solids. Method quality of each batch was controlled by using a total of three control samples of a 25 and 100 mg/L kaolin solution and distilled water blanks.

For the soil loss (ERO) calculations on a hectare basis per runoff event, the following equation was used.

where TSS is the amount of total suspended solids in mg per liter, RUN is the total amount of water in the tank in liter and A is the Wischmeier plot area in m².

In four occasions of high runoff water collection during spring 2019 (9^th of September, 2^nd of October, 11^th of October, 13^th of October), sediments were collected through 10 L runoff samples and water evaporated gently over 48 h at 50°C. Sediment samples were subsequently analyzed for total carbon (C) and N_tot by LECO. Data were used to calculate the sediment nutrient enrichment for C and N as the ratio of nutrient concentration in eroded sediment to that of soil samples taken in October 2019 in 0–7.5 cm soil depth (Palis et al., 1990).

Crop yield was determined manually for two subsamples per plot, cutting two crop rows per 4 linear meters.

Soil Measurements

Field Sampling

Historical sampling campaigns for SOC were carried out with a soil auger taking one composite sample per plot in 0–20 cm soil depth in the years 1984, 1986, 1987, 1989–1996, and in 0–15 cm soil depth in the years 2011, 2014–2016.

On 27^th of December 2018, all plots were sampled for one composite soil sample at three depths (0–7.5, 7.5–15, 15–30 cm) and on 25^th of October 2019, all plots were sampled for two composite samples at five depths (0–7.5, 7.5–15, 15–30, 30–45, 45–60 cm).

Soil Chemical Analysis

SOC and total N (N_tot) were analyzed by dry combustion at 900°C followed by infrared detection (LECO Truespec; Wright and Bailey (2001). Before 2011, SOC samples were analyzed with the Tinsley method (heated dichromate/titration), therefore SOC values obtained since 2011 were converted using a previous determined factor of 0.81 (Grahmann et al., Forthcoming 2022; in press).

Soil pH was determined potentiometrically (1:2.5 soil/distilled water suspension; Beretta-Blanco et al., 2014). Phosphorus was measured colorimetrically by the Bray-1 method using a 1:10 (w/v) soil/solution ratio and an extraction time of 5 min (Bray and Kurzt, 1945). Cation exchange capacity was analyzed by extracting the exchangeable cations Ca, Na, Mg and K with 1 M ammonium acetate (NH₄OAc) at pH 7, and reading the extract by atomic emission (K and Na) according to Jackson (1964) and atomic absorption (Ca and Mg) or by atomic emission with ICP-OES equipment.

Soil Physical Properties

Undisturbed soil samples were taken for the determination of bulk density in December 2018 and December 2019 for two depth increments (0–7.5, 7.5–15 cm). Three sites located at the upper, middle and lower part of each plot were selected. Samples were taken using a soil probe with a cylinder volume of 98.2 cm³ and samples dried for 48 h at 105°C. Penetration resistance was measured from 0 to 80 cm soil depth with a hand held Penetrologger (Eijelkamp, Giesbeek, Netherlands, cone base area 1 cm²) at 4 days (13th of September, 19th of September, 30th of October and 23rd of December) in 2019. Resistance measurements were executed according to the manufacturer’s instructions by applying the electronic penetrometer together with a datalogger, allowing for immediate storage and processing of the data in the datalogger. In parallel, soil moisture was measured volumetrically in 6 cm soil depth.

Three soil blocks were sampled in April 2019 in the upper, middle and lower part of each plot (15 cm × 15 cm × 20 cm). From those samples, a visual evaluation of soil structure (VESS) was conducted in the field as described by Guimarães et al. (2011). The samples were then air dried for the determination of four indicators of aggregate stability in the laboratory. The stability of aggregates to quantify the main processes of aggregate breakdown associated with water stress were assessed as proposed by Le Bissonnais (1996). This methodology evaluates three mechanisms of aggregate breakdown: 1) Aggregate slaking due to fast wetting (Treatment 1); 2) Differential swelling (Treatment 2) and 3) Mechanic breakdown (Treatment 3). Three 5 g samples of 3–5 mm aggregates were taken for the determination of aggregate stability for each of the treatments. For treatment 1, the aggregates were immersed in 250 ml of distilled water for 10 min. For treatment 2, the aggregates were capillary rewetted on a tension table at 3 cm tension for 30 min before their immersion in water. For treatment 3, the aggregates were rewetterd in ethanol, then immersed in water (to avoid fast wettering) and mechanically agitated 10 times using a Feodoroff agitator. Each one of these treatment samples was gently transferred to a 50 μm sieve, previously immersed in ethanol. Fragments greater than 50 μm were oven-dried at 40°C and dry-sieved through seven different sieve sizes (3,000, 2000, 1,000, 500, 250, 100, and 50 μm). The aggregate stability for each treatment was expressed as the normalized mean weight diameter for each individual treatment, and the mean aggregate diameter of all treatments (MWD_{LeBissoinnais} ) is evaluated.

The fourth evaluated indicator corresponds to the USDA wet aggregate stability (Aggregates_USDA; method 1B1b2a1 in USDA-NRCS, 2004), that follows a disruption of an initially 3 g air dried sample of 1–2 mm soil aggregates by submerging and wet sieving in distilled water through a 0.5 mm sieve. A following dispersion in sodium hexametaphosphate solution (Calgon) is performed for determining >0.5 mm sand mass. Sand mass is subtracted from both the initial mass (i.e., 3 g) and the mass retained after sieving in water. Aggregate stability is computed as the ratio between the latter and former dry sand-free soil masses.

Soil Hydrological Measurements

Infiltration rates were determined by double-ring infiltrometers (Eijelkamp, Giesbeek, Netherlands) in May 2019 in five selected sites per plot according to the Eijkelkamp user manual (Eijelkamp, 2018). Two rings were installed in non-disturbed rows and three rings in the planting rows. Philip’s infiltration equation was applied for steady state infiltration rates (Philip, 1957). Surface saturated hydraulic conductivity (K_s) was estimated through the measurement of steady-state infiltration, solving Wooding equation (Wooding, 1968), using the White and Sully method (White and Sully, 1987; Logsdon and Jaynes, 1993; Angulo-Jaramillo et al., 2000).

Soil Biological Properties

Particulate organic matter (POM) fractions of 212 μm, 53 µm and <53 µm (mineral-associated organic matter, MAOM) were separated by an adapted method of Cambardella and Elliott (1992) using 6.66 g of dry, 2 mm sieved soil and 20 ml of calgon (5% sodium hexametaphosphate). After drying, each fraction was analyzed for SOC.

Potentially mineralizable nitrogen (PMN) was measured in three occasions in 2019 by anaerobic incubation of fresh soil samples over 7 days and subsequent determination of the produced ammonium (NH₄) by colorimetry (Waring and Bremner, 1964; Bundy and Meisinger, 1994). Analyzed samples were from identical origin as samples used to determine aggregate stability (see above). Details for POM and PMN analysis are reported in (Fabrizzi et al., 2003).

Soil Erosion Modeling and Statistical Analysis

The Revised Universal Soil Loss Equation (RUSLE) is the mostly used, widespread model to predict average annual soil losses resulting from rainfall erosion of cropland (Renard et al., 1997; Tang et al., 2015). The estimation of soil erosion with this model is defined by six factors: R is a runoff-rainfall erosivity factor; K is a soil erodibility factor; LS is a topographic factor combining slope length (L) and slope steepness (S); C is a cover-management factor and P is a supporting practices factor. RUSLE was used to model the corresponding annual long-term erosion losses for the three treatment groups. The RUSLE factors R (rainfall erosivity), L (slope length factor) and S (slope steepness factor) did not change within the experimental site. The same applied for the P factor (support practice factor) as all plots were managed identically. The 30-years R factor for this region was 4248 MJ mm/ha/h/yr (Pérez-Bidegain et al., 2017) and the K factor was 0.023 t ha h/MJ/ha/mm (Puentes, 1981). The three treatment groups differed in the C factor (cover-management factor), which ranged between 0.06 and 0.3 (Clerici and García Préchac, 2001; Hill et al., 2008).

A widespread weakness of long-term runoff plots is their unreplicated nature (Packer et al., 1992; Williams et al., 2009; Ramos-Scharrón and Figueroa-Sánchez, 2017) which can result in an unbalanced experimental design during their long-term history which is why in this study, plots were grouped as pairs per treatment. However, repeated measures over time for continuous variables and pseudo-replicated point measurements within each plot did not allow an analysis of variance. Therefore, with acknowledgement of pseudoreplication (Davies and Gray, 2015), data were analyzed using descriptive statistics, paired “t” tests (Microsoft Office Excel 2010), a general linear model (SAS Version 9.4; PROC GLM) (SAS Institute, Cary, United States, NC) and correlation analysis to evaluate this comprehensive data set.

Results

Evolution of Soil Organic Carbon

SOC decreased continuously in the first experimental period until 1996 and tended to slightly increase in the study years. In June 2015, SOC reached a high level in all treatment groups after soybean harvest. In April 2015, Black oats were sown broadcast and high root biomass was present during sampling. A sharp decrease in all treatments was recorded in December 2018 when maize and wheat were cropped within the same calendar year. TIL_CROP had lowest SOC content in most of the measured years and showed a sharp SOC depletion during the first experimental period. TIL_CP was 6 years longer under fallow and had the lowest IAI, showing higher SOC in most sampling events after conservation agriculture was applied in all plots (Figure 2).

FIGURE 2. Topsoil carbon content (SOC in %) for each treatment group (fallow period between 1997 and 2010 was excluded) since the onset of the experiment in 1984 until 2019.

Soybean yields ranged from 1,322 kg/ha in 2014–5,472 kg/ha in 2017. NT_Mix and TIL_CROP were most intensively cropped, having three times more grain extracted compared with TIL_CP over the last 10 years (Supplementary Table S1). In the last experimental period between 2015 and 2019, lowest crop yields were obtained in plot 3 (TIL_CROP), and for soybean also in plot 1 (TIL_CP) and highest yields were observed in plot 2 (TIL_CP) and plot 6 (NT_Mix).

Soil Chemical and Biological Properties

SOC was the lowest in both years and most depths in TIL_CROP and increased for all treatments in the second year (Supplementary Table S2). The opposite happened for N_tot which was lower in 2019, having lowest content for both years in TIL_CROP in the first 30 cm soil depth. Soil pH was slightly acidic in the top soil for all treatments and increased with decreasing soil depth, in many occasions pH was lowest in NT_Mix. Particulate organic matter (POM) averaged 0.08 g/100 g soil for coarse POM (>212 µm) and 0.09 g/100 g soil for fine POM (53–212 µm), while mineral-associated organic matter (MAOM, <53 µm) averaged 2.49 g/100 g soil. Significantly highest POM was measured in the TIL_CP with lowest IAI and 6 years longer fallow period. Potentially mineralizable nitrogen (PMN) averaged 47, 37 and 42 mg/kg in 0–7.5 cm and 3, 10 and 16 mg/kg in 7.5–15 cm soil depth for the sampling dates 25^th October, 11^th November and 29^th November 2019, respectively. Lowest PMN was measured on two occasions in TIL_CROP (data not shown).

Soil Physical and Hydrological Properties

The average infiltration rate for all treatments was 27.5 mm/h. Infiltration rate was highest with NT_Mix and lowest in TIL_CROP, the same was true for saturated hydraulic conductivity (K_s) which averaged 22.4 mm/h. Both parameters were significantly lowest in plot 3 and highest in plot 6 (p = 0.0011 and 0.0018, respectively; Table 4). Bulk density averaged 1.35 g/cm³ in both depths (0–7.5, 7.5–15 cm) in 2018 and 1.22 g/cm³ in 7.5 and 1.26 g/cm³ in 15 cm soil depth in 2019 and was lowest in both years in plots with lowest IAI in TIL_CP and increased with higher IAI due to historical tillage intervention and longer continuous cropping. Penetration resistance in 0–15 cm soil depth was 1,482, 1,382 and 1,646 kPa for TIL-CP, NT_Mix and TIL_CROP, respectively averaged over four measurement days. It was significantly highest in three out of four measurement days in TIL_CROP (p = 0.0059).

TABLE 4. Treatment averages for physical and chemical soil properties measured in 2019 (GLM-SAS: capital letters indicate evidence of significant different treatment groups (p < 0.05); MWD-mean weight diameter, SOC-soil organic carbon, PMN-potentially mineralizable nitrogen, POM-C-particulate organic matter carbon).

Treatments had no statistically significant effect on Aggregates_USDA and the VESS index. Nevertheless, TIL_CP had higher Aggregates_USDA and an “intact” soil structure score (VESS index 2 with high aggregate porosity) whereas NT_Mix and TIL_CROP obtained VESS scores of three and more (firm structure with low aggregate porosity). For measured aggregate stability according to USDA and Le Bissoinais, mean weight diameter of soil aggregates was lowest in TIL_CROP, having the highest IAI (Table 4).

Most pronounced differences in soil quality parameters caused by different soil use and management history were found for infiltration rate (Infil_mm/h), saturated hydraulic conductivity (Ks_mm/h) and PMN in 7.5 cm soil depth (Figure 3). Although not statistically substantiated, relative differences in VESS and SOC in 7.5 cm soil depth between treatments groups were recognized.

FIGURE 3. Spider diagram of relative deviations (in %) from the population mean for soil parameters assessed in 2019.

Soil Erosion, Runoff and Nutrient Losses

Yearly average runoff was 2.7 times greater in 2018 than 2017 while the runoff coefficient was 2.5 times greater in 2018 (Table 5).

TABLE 5. Rainfall classification and average annual runoff information for the study period.

Yearly average soil losses amounted 233, 805 and 139 kg/ha for the three respective study years with highest erosion rates in 2018 in all treatment groups. Yearly differences can mostly be explained by climatic variability. The higher erosion was due to more intense rainfall events and about 100 mm more rainfall in 2018 compared with the other two studied years. Also, the number of high, intensive rainfall events above 40 mm was the highest in 2018 (6), but only 2 and 3 events were recorded in 2017 and 2019, respectively (Table 5). Overall, rainfall occurred in about a third of each calendar year. A power regression and moderate correlation was found between rainfall and surface runoff for 98 events during 2017 and 2019 (R² = 0.53; y = 0.0013x^2.1045). Few rainfall events above the 1:1 line caused higher runoff than the actual amount of rainfall of the corresponding event and was influenced by the actual state of soil moisture saturation and hence, previous rainfall events (data not shown). We found overall higher erosion rates in the NT_Mix and lowest in the TIL_CP to sustain our hypothesis that a smaller IAI score is condensed in soil erosion rates due to historical intensive cropping and tillage (Figure 4). As soil and crop management was similar in the last five study years, measured differences in runoff and erosion between plots were elicited by soil management legacy or due to natural spatial variability between plots. The coefficient of variation (CV) among Wischmeier plots was high and between 52% (8^th October) and 120% (3^rd September) in 2017 for soil erosion events of more than 5 kg/ha. In 2018, CV increased between 23% (16^th December) to 145% (29^th May) and in 2019, CV between plots was between 45% (26^th July) and 74% (30^th January) for soil erosion events with more than 5 kg/ha average soil loss.

FIGURE 4. Accumulated soil loss per plot and year (kg/ha).

The estimated long-term soil erosion rate with RUSLE was 1,000 kg/ha, 1,600 kg/ha, and 2,100 kg/ha for NT_Mix, TIL_CP, and TIL_CROP, respectively.

In 2019, volatile solids corresponded to one third of total soil losses, representing mainly the removal of straw residues and decomposing material on the soil surface (Table 6). The share of fixed material increased with increasing total soil loss. Moderate to high correlations were found between the amount of lost volatile solids and maximum rainfall intensity (R² = 0.64 for NT-Mix; R² = 0.70 for TIL_CP; R² = 0.45 for TIL_CROP) and between the fixed and volatile solid ratio and the respective runoff for each treatment group (R² = 0.70 for NT_Mix; R² = 0.60 for TIL_CP; R² = 0.45 for TIL_CROP).

TABLE 6. Fixed and volatile solids in kg/ha/day for selected soil erosion events in 2019, sorted from shortest to longest rainfall duration per rainfall event.

Peaks of intense runoff events mostly coincided with high rainfall, however not all recorded rainfall events over 70 mm resulted in erosion (eg., May 2019, Figure 5). The plots were permanently covered and cropped during the monitoring period. Not for all reordered erosion events, rainfall intensity was the relevant factor. There was only one exceptional erosion event during the study period in NT_Mix in July 2018, leading to more than 350 kg/ha soil loss within 1 month. On 30^th of June, wheat was planted and hence soil was disturbed by the planter disks. Due to heavy rainfall the days before and after planting, wheat emergence was reduced and crop establishment was poor and soil was not covered properly. Although it was an exceptional erosion event in NT_Mix, the overall magnitude of erosion was much lower than the tolerable threshold.

FIGURE 5. Monthly pattern of soil erosion for the respective treatment groups in relation to climatic variables of monthly rainfall and highest monthly 10-min rainfall intensity.

The enrichment ratio showed that SOC was 6.8 times higher in the sediment eroded in TIL_CROP than contained in the topsoil, this carbon enrichment ratio was even higher in NT_Mix and lowest with the lowest IAI in TIL_CP (Table 7). The nutrient loss through runoff was slowed down with lowest IAI which in turn is related to the higher SOC in the topsoil observed in the latest study years (Figure 2). Highest sediment enrichment for C and N in NT_Mix did not lead to lowest SOC content in the top soil and is explained by the overall management effect of straw retention and no-tillage over several decades in this treatment group.

TABLE 7. Sediment enrichment ratio for carbon (SOC) and nitrogen (N_tot) in surface runoff samples (n = 4, SD = standard deviation).

Discussion

Historical Land Use and Current Soil Quality

During the first experimental period (1984–1997), which represented 55% of the total running time since the installation of the evaluated experiment, a generalized depletion of SOC was observed in all evaluated plots (Figure 2). This depletion was greater in TIL_CROP systems where SOC in 1996 represented 81% of the initial SOC content whereas in TIL_CP and NT_Mix, the average losses were 10 and 3%. When compared to no-till systems, intensive tillage can promote C mineralization by breaking soil aggregates and therefore increasing the access of soil microorganisms to previously protected C pools (Six et al., 1999, 2000). Additionally, tillage can promote soil erosion losses (Verhulst et al., 2010). Despite that, after 1997 all plots received the same tillage management until 2019, but SOC levels in TIL_CROP treatments remained below the C levels observed for the other two treatments during the evaluated period. The observed differences in near-surface SOC levels generated in early stages of the experiment can be related to the positive effects of historical no-tillage in NT_Mix on soil physical indicators, like penetration resistance, aggregate stability, bulk density, and available water capacity (Blanco-Canqui and Ruis, 2018), that were not eliminated with later changes in tillage (Table 4; Figure 3). Similar reverse effects on the decline in soil fertility and soil aggregation were found in a 20 year old grass pasture site with a prior history of cropping activity (Jones et al., 2016). There, the loss of SOC under pasture increased with greater years under cropping and soil aggregation and mineralizable N did not improve with perennial pastures. The same was true for TIL_CROP with highest cropping activity where most of the measured soil quality parameters did not recover at the same level as for TIL_CP or NT_Mix. The TIL_CP treatment which incorporated pastures in the rotation and had a longer discontinuation period, counted with only three soybean harvests and extracted 10 Mg of grain/ha. However, in the NT-Mix and TIL_CROP treatments 33 Mg of grain/ha were harvested over ten harvest years (Supplementary Table S1). As expected, the combination of lower C and nutrient extraction with the lower tillage intensities in TIL_CP maintain an overall higher soil quality when compared to more extractive and tillage intensive systems (Amsili et al., 2021). The observed soil degradation in TIL_CROP is consistent with that reported for continuous annual cropping systems in the Pampas region where a continuous nutrient extraction combined with relatively low mineral fertilizer inputs (Table 1) might explain the historical decline of SOC (Figure 2). The observed depletions in SOC can significantly deplete obtainable yields (Ernst et al., 2020; Rubio et al., 2021a). When short-term management strategies such as increasing crop fertilization and mechanical soil decompaction have been proven insufficient to remediate yield depletion (Ernst et al., 2020; Rubio et al., 2021b), our results indicate that 5 years of conservation agriculture do not compensate for the degradation of soils generated by previous land use history. Significant differences in several soil properties confirm our hypothesis that the soil management legacy still lingers on the current state of soil (Table 4). Further studies are needed to determine the amount of time needed to recover soil health and yield in degraded areas. The intensification of agricultural production, in the form of historical cropping sequence and tillage activity, was represented by the study´s treatment groups. However, it should be avoided to generalize that higher cropping intensity automatically leads to worse soil conditions and hence to higher erosion rates. This depends on the management practice at the cropping system level and on multiple options and combinations cropping intensification is realized (Mouratiadou et al., 2021). Permanent cover, growing roots and continuous N and C input may also lead to beneficial effects on soil quality which highly depends on the management (crop rotation, tillage, cover crops, etc) (Xiong et al., 2019; Cassman and Grassini, 2020). However, sustainable intensification as a condition to avoid soil deterioration was not given in the current study with simplified and soybean focused crop rotation in the second and third experimental period.

Historical Land use and Current Soil Erosion

The distinction between volatile and fixed suspended solids in water samples is done to determine total mineral soil loss and additional removal of organic material, most present as SOC, by water erosion. According to Lal (2019), there is a big lack of research in this area showing that particulate organic carbon or POM are related to losses of volatile solids. In 2019, higher POM was found in treatments with higher volatile solid losses. The annual losses of volatile solids in 2019 increased with decreasing IAI averaging 84 kg/ha in TIL_CP, 58 kg/ ha in NT-Mix and only 40 kg/ha in TIL_CROP showing the opposite pattern of POM (TIL_CROP < NT_Mix < TIL_CP). Highest losses of volatile solids in TIL_CP were also measured in 2017 (148 kg/ha in TIL_CP, 106 kg/ ha in TIL_CROP and 74 kg/ ha in NT_Mix), however no POM was assessed for the other 2 years.

Several studies reported that the soil removed by erosion is 1.3–5 times richer in organic matter than the remaining soil (Bagarello and Ferro, 2017) which was even higher in NT-Mix and TIL_CROP (Table 7). In line with our results, Bertol et al. (2007) analyzed organic carbon in runoff sediments under soybean cropping and found higher losses in no-till compared with conventional tillage, but overall C enrichment ratio was much smaller (between 1.00 and 1.17) than in the current study. Sediment ratios did not follow a certain treatment tendency for the selected rainfall events (Table 6). But higher SOC combined with improved soil structure indices in TIL_CP may lead to the assumption that SOC was tighter bound to MAOM and hence less prone to get lost as fixed solids (attached to sediments, Holz and Augustin, 2021).

As discussed in the previous section, historical land use had a legacy effect on soil quality parameters, several of which reportedly may control soil erodibility (Victoria et al., 2001; Taleshian Jeloudar et al., 2018; Alaboz et al., 2021). However, the soil loss data of the present study do not support a corresponding effect of historical land use on soil erosion. An explanation to the relatively low erosion rates measured across the 3 years of study (Figure 4) is found in the no-tillage practices that were implemented in 2009, 12 years after the end of the initial experimental period (Figure 1). Cover crops were later consistently incorporated in all plots and possibly further improving soil conservation. Soil quality parameters might be responding to this management, although noticeable effects are largely limited to the surface soil layers (Liebig et al., 2004) and only gradually reach deeper soil. Another finding that supports this explanation is small differing ratio of fixed:volatile sediments across historic use treatments for several sampling days in 2019 (Table 6). Consequently, in spite of the latest conservation practices, soil quality parameters are still exhibiting the legacy effects, while erosion losses can be mitigated substantially through enhanced soil cover. As proposed recently by (Willett et al., 2019), no-tillage with residue cover is a functioning soil conservation practice to prevent nutrient losses as shown in the current study by reduced soil erosion and subsequent minimal C and N loss during the experimental period.

Infiltration rate, saturated hydraulic conductivity and PMN were the soil quality parameters more affected by the soil use history (Figure 3). Infiltration rate and saturated hydraulic conductivity are key soil hydrologic properties that affect soil water dynamics, hence soil water erosion (Toy et al., 2002). RUSLE allows for estimating long-term soil erosion average rates, hence its results cannot be compared with yearly measurements of sediment losses over the three study years. Nevertheless, RUSLE can be used as an indicator of sustainability. Erosion estimations were based on the nationally calibrated RUSLE/USLE model and assumed a scenario of continuing soil use and management in each of the treatment groups. This resulted in relatively high soil loss rates (i.e. between 1,000 and 2,100 kg/ha/yr), but still well below the tolerable soil loss threshold of 7,000 kg/ha/yr. The lowest modelled erosion was obtained for NT_Mix, as a consequence of using a lower C-factor to account for improved soil cover from remaining crop or pasture residues. Although it was not the purpose of this paper to validate the RUSLE model, the modelled erosion rate of 1,000 kg/ha/yr in NT_Mix performed averagely in the range between 71 and 1,478 kg/ha/yr for the measured yearly rates (Figure 4), but the model obtained much higher rates for the remaining treatments TIL_CP and TIL_CROP. However, Knapen et al. (2008) and Blanco-Canqui et al. (2009) showed that conservation tillage practices may also contribute to reduced erosion through a smaller erodibility factor K. Therefore, if the erodibility factor is modified in response to soil use (or soil quality parameters), modelled and observed erosion rates could perform better. We therefore propose that these measured soil legacy effects in soil properties should be reflected in changes of the RUSLE soil K factor, which cannot be static but dynamic due to ongoing, permanent processes that affect soil quality (Alewell et al., 2019). The runoff plots have existed since 1984, which implied considerably changes in treatments and different land use management over time. Overall, TIL_CROP showed a noticeable degradation pattern compared with the other two treatment groups which should be translated into a higher erodibility. Although Sasal et al. (2010) reported a linear regression between ISI and cumulative runoff, this was not confirmed with the applied IAI in the current study. Runoff was lowest in all 3 years in TIL_CP, but highest in NT_Mix whereas the treatment with the highest IAI, TIL_CROP, had medium runoff rates for all 3 years. Hence, increasing historical soil use and management intensity does not automatically lead to higher runoff rates, but rather depends on a complex combination of cropping and tillage management factors and subsequent soil process effects.

Conclusion

Over the study period, conservation agriculture controlled soil erosion and soil loss was minimal and far below the national threshold of 7,000 kg/ha. Historical soil use and management caused significant long-term effects on soil properties leading to adverse effects in soil quality in historically intensively tilled and cropped soils. We found significant differences, especially for soil physical parameters between the soil legacy treatments which were mainly due to the historical management conducted during the first experimental period from 1984 to 1997. Five years of conservation agriculture with year-around soil cover did not remediate soil degradation caused by continuous cropping and ploughing. The agricultural intensification index was sensible to detect cumulative treatment differences due to previous soil and land use. The legacy of soil use and management affected the resilience of current sustainable cropping systems. This study found evidence that soil erodibility is affected by cropping management changes over time, and requires an adjustment of the K factor when crop and tillage activity had historical modifications. However, it seems that surface cover during the study period offsets the importance of soil quality status regarding the proneness to soil loss. Conservation agriculture compensated for reduced soil quality in the short-term and protected TIL_CROP in order to avoid soil erosion.

Data Availability Statement

The raw data supporting the conclusion of this article will be made available by the authors on request, without undue reservation.

Author Contributions

KG and VR conducted the samplings and measurements, KG was responsible for evaluation and analysis of the data. KG did the writing of the original draft and the visualization, MP-B assisted with the RUSLE modelling and added to the discussion. JQ added to the discussion. KG and VR had the project idea, supervised and added to the discussion.

Funding

This study has been financially supported by the National Institute of Agricultural Research (INIA), Uruguay (Grant number N-23622).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank Jorge Sawchik for the re-installation of this experiment. We thank Wilfredo Mesa, Emiliano Barolín, Héctor Vergara, Gualberto Soulier, Leonardo Silva, Nicolás Leiva and Marcelo Schusselin for trial maintenance; Carolina Lizarralde for previous trial analysis and Henning Schwalm for the analysis of suspended solids in the laboratory.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fenvs.2022.822967/full#supplementary-material

References