Tracing erosion rates in loess landscape of the Trzebnica Hills (Poland) over time using fallout and cosmogenic nuclides

Loess landscapes are highly susceptible to soil erosion, which affects soil stability and productivity. Erosion is non-linear in time and space and determines whether soils form or degrade. While the spatial variability of erosion is often assessed by either modelling or on-site measurements, temporal trends over decades to millennia are very often lacking. In this study, we determined long- and short-term erosion rates to trace the dynamics of loess deposits in south-western Poland. We quantified long-term (millennial) erosion rates using cosmogenic (in situ 10Be) and short-term (decadal) rates with fallout radionuclides (239+240Pu). Erosion processes were studied in two slope-soil transects (12 soil pits) with variable erosion features. As a reference site, an undisturbed soil profile under natural forest was sampled. The long-term erosion rates ranged between 0.44 and 0.85 t ha−1 year−1, whereas the short-term erosion rates varied from 1.2 to 10.9 t ha−1 year−1 and seem to be reliable. The short-term erosion rates are up to 10 times higher than the long-term rates. The soil erosion rates are quite consistent with the terrain relief, with erosion increasing in the steeper slope sections and decreasing in the lower parts of the slope, while still maintaining high values. Soil erosion rates have increased during the last few decades owing to agriculture intensification and probably climate change. The measured values lie far above tolerable erosion rates, and the soils were found to be strongly imbalanced and exhibit a drastic shallowing of the productive soils horizons.


Introduction
The increase in soil erosion is a direct consequence of agricultural exploitation and threatens soil stability, quality and its productivity (Rickson 2014;Guzmán et al. 2015;Alewell et al. 2017;Golosov et al. 2021). Long-term and intense erosion removes topsoil from the upper part of a slope and deposits the eroded material at the toe of the slope, thus leading to irreversible changes in the natural structure of those soils and their corresponding horizons (Świtoniak et al. 2016;Zádorová and Penížek 2018;Golosov et al. 2021). One of the materials most susceptible to erosion are loess deposits (Licznar et al. 1981;Yang et al. 2006;Zhang et al. 2018;Poręba et al. 2019). Loess materials are widely distributed across the world (Muhs 2013;Schaetzl and Attig 2013;Pasquini et al. 2017). In Europe, the most extensive deposits on Earth cover an area from France to Russia, having developed during the Last Glacial period (Haase et al. 2007;Lehmkuhl et al. 2020). Productive soils have developed on the loess, including Chernozems, Pheozems and Luvisols (Altermann et al. 2005;Gerlach et al. 2012;Labaz et al. 2018;Kabała et al. 2019;Loba et al. 2020). As a consequence, since the Neolithic period, loess areas have been deforested and transformed into arable lands, which has given rise to intense soil erosion processes (Altermann et al. 2005;Gerlach et al. 2012;Poręba et al. 2019).
So far, 137 Cs has predominantly been used as the isotope tracer for soil erosion in loess landscapes (Yang et al. 2006;Poręba et al. 2011Poręba et al. , 2015Poręba et al. , 2019. 137 Cs is a fallout radionuclide (FRN) that is distributed globally as a result of, for example, nuclear weapons testing in the 1960s and nuclear power plant accidents, and so it enables the investigation of erosion rates over the last few decades (Alewell et al. 2014(Alewell et al. , 2017Zollinger et al. 2015;Meusburger et al. 2016). However, it is characterised by a short half-life of 30.17 year, and recent estimations have indicated that more than 70% of the global 137 Cs has disappeared due to its radioactive decay (Xu et al. 2015). At sites exhibiting erosion and a subsequent additional loss of 137 Cs, its detectability becomes increasingly difficult. Moreover, the Chernobyl accident in 1986 caused a spatially inhomogeneous distribution of 137 Cs in Central and Western Europe, which has caused some limitations in its application (Arata et al. 2016a;Alewell et al. 2017). As a result, 239+240 Pu isotopes are currently more often used due to their longer half-lives (24,110 and 6561 years, respectively) and because Pu is absent from the volatile fraction released by the reactor accident of Chernobyl (Arata et al. 2016a).
Long-term erosion rates (over several millennia) can be determined using cosmogenic nuclides, such as 10 Be (Hidy et al. 2010;Zollinger et al. 2015;Calitri et al. 2019). Based on their origin, two types of 10 Be are distinguished: (1) meteoric 10 Be, which is constantly produced in the upper atmosphere by the spallation of nitrogen and oxygen by cosmic rays and is deposited at the Earth's surface by rainfall (Graly et al. 2010;Willenbring and von Blanckenburg 2010;Wyshnytzky et al. 2015), and (2) in situ 10 Be, which is directly produced in the crystal lattice of quartz by the interaction of cosmogenic rays (Hidy et al. 2010). For loess landscapes, there is a lack of studies using meteoric or in situ 10 Be to determine denudation rates. Meteoric 10 Be has been applied to studying palaeoclimatic variations (Gu et al. 1997;Zhou et al. 2015), estimating palaeoprecipitation (Sartori et al. 2005), establishing time scales for loess deposition (Chengde et al. 1992) and in tracking the translocation of 10 Be in Luvisols (Jagercikova et al. 2015). In most cases, 10 Be and 239+240 Pu or 137 Cs have been used separately for calculating erosion rates (Arata et al. 2016a;Waroszewski et al. 2018;Musso et al. 2020). Recent studies, however, have combined these two types of isotopes to compare medium-and long-term erosion rates.
For instance, Zollinger et al. (2015), Calitri et al. (2019) and Jelinski et al. (2019) applied both types of isotopes to compare erosion processes for the last 50-60 years with long-term rates in order to crystallise the effects of anthropogenic pressures and climate change on soil processes.
In this study, we made a first attempt to apply in situ 10 Be and 239+240 Pu to quantify soil erosion in a loess landscape, with special focus on (1) documenting long-and short-term erosion rates, (2) cross-checking recent and past erosion rates to determine how much erosive processes have intensified in recent decades and (3) verifying whether isotopic methods are suitable tools for studying erosion processes in the south-western Polish loess belt.

Study area
The investigation area was located in the Trzebnica Hills, south-western Poland (Fig. 1). The subglacial glaciotectonic disturbances in this region are estimated to have come from the Odra Glaciation (Saalian-Drenthe, marine isotope stage [MIS] 8/6) (Krzyszkowski 1993;Jary 1996). In general, the local geology is dominated by Quaternary deposits (loess, glacial till, fluvioglacial sediments), but Neogene deposits (clays, sands and gravels) occur locally in small areas (Pachucki 1952;Dyjor 1970;Dyjor and Kościówko 1982;Jary 1996). Loess, which is the youngest deposit, is distributed across the Trzebnica Hills as a layer with varying thickness (Jary 1996). The outcrop in Zaprężyn provides a record of loess dynamics and pedogenesis from the Last Interglacial (Eemian, MIS 5e) to the Upper Pleniweichselian (MIS2) (Jary and Ciszek 2013).
The soils in this region, mostly developed from the loess deposits, are characterised by subsoil clay illuviation (Luvisols) or the presence of a dark humic topsoil (Chernozems and Pheozems) (Licznar et al. 1988;Zmuda et al. 2009;Kabała and Marzec 2010;Glina et al. 2014) and have been affected by denudation processes due to agricultural use since the Neolithic period (Poręba et al. 2011).
The native vegetation is represented by oak-hornbeam forest, although most of the land is used for agriculture due to the productivity of the soils (Anioł-Kwiatkowska 1998). The main crops cultivated are wheat, beetroot and corn. According to the Köppen-Geiger climate classification, the Trzebnica Hills experience warm summers and a humid, continental climate. The mean annual temperature is 8 °C, with average temperatures in the coldest and warmest months being − 3 °C (January) and 17 °C (July), respectively (Bac and Rojek 2012). The annual average precipitation is ca. 600 mm (Bac and Rojek 2012).

Field sampling
The sites sampled were in the southern part of the Trzebnica Hills, close to the village of Wysoki Kościół. Two transects along two slopes with pronounced erosion features were samples in order to track the multidirectional character of soil erosion (Table 1, Fig. 1). The two slope catenas started at the top of the local hill, sharing sample WK1 as the starting profile for both transects, and with the first transect (WK1, WK2-WK7) trending to the south and the second (WK1, WK8-WK12) to the north-east. The range of slope inclination was similar for both transects and reached the highest values with 12-13° at the shoulder position (Table 1). The altitudinal differences along the first catena (WK1-WK7) are about 36 m and along the second catena (WK1-WK12) 20 m. For each transect, the soil profiles were described according to the guidelines for soil description (Food and Agriculture Organization [FAO] 2006) and classified according to the FAO-World Reference Base (WRB) system (International Union of Soil Sciences Working Group WRB 2015). Undisturbed soil samples were taken, using stainless steel rings (100 cm 3 ), from all the soil horizons for bulk density analysis. About 1 kg of bulk soil material was sampled for physicochemical and geochemical analyses. To obtain enough quartz (0.25-0.60 mm fraction) for in situ 10 Be analysis, 8-9 kg of soil material were taken at depths of 20 cm from the surface to 100 cm. For the 239+240 Pu analyses, samples were taken every 5 cm from the surface to a depth of 40 cm. However, the Ap horizon (0-20/0-25 cm) was mostly homogenous due to ploughing, so only one sample was taken from the first 20/25 cm. In addition, a reference soil profile was sampled in a nearby forested area. To overcome the large sampling number, a large amount of sample material (about 1 kg) was taken per depth increment and homogenised. Prior to the measurements, the soil samples were dried, crushed and sieved (2 mm mesh).

Soil properties
The particle-size distribution was measured using a sieving (sand fraction) and hydrometer method to determine the silt and clay fractions (van Reeuwijk 2002). The pH of the samples was measured potentiometrically (in deionised water) in a 1:2.5 suspension (Kabała et al. 2016). The hydrolytic acidity was extracted using a 1-M sodium acetate (CH 3 COONa) solution and potentiometric titration, while the exchangeable ions (Ca 2+ , Mg 2+ , K + , Na + ) were extracted using a 1-M ammonium acetate (CH 3 COONH 4 ) solution at pH 7(van Reeuwijk 2002). Calcium carbonate (CaCO 3 ) was measured conductometrically using a Scheibler apparatus. The soil organic carbon (SOC) content was determined by dry combustion at 550 °C and the non-dispersive infrared absorption of CO 2 , using a Ströhlein CS-mat 5500 analyser (prior to the analysis, the calcium carbonates were removed, if present). The bulk density was measured using the dry weight method (van Reeuwijk 2002).
The geochemical composition was determined using X-ray fluorescence (XRF). Approximately 5 g of soil was milled to 50 μm in a tungsten carbide disc swing mill (Retsch® RS1, Germany). The milled samples were then weighed into plastic cups using a Prolen® foil and measured using an energy-dispersive, helium-flushed XRF To estimate the degree of weathering of the soil layers, the Chemical Index of Alternation (CIA) (Nesbitt and Young 1982) was used (molar ratios):

Determination of cosmogenic, in situ 10 Be
10 Be was extracted from pre-cleaned quartz grains from the 0.25-0.60 mm fraction using an isotope dilution method that followed the modified protocol of von Blanckenburg (Kohl and Nishiizumi 1992;von Blanckenburg et al. 1996) at the University of Zurich. The 10 Be/ 9 Be ratios were measured using a Tandy accelerator mass spectrometer (AMS) at ETH Zurich (Christl et al. 2013), normalised to the ETH AMS standard S2007 N ( 10 Be/ 9 Be = 28.1 × 10 −12 nominal) and calibrated to ICN 01-5-1 ( 10 Be/ 9 Be = 2.709 × 10 −11 nominal) (Nishiizumi et al. 2007), both associated with a 10 Be half-life of 1.387 ± 0.012 Myr.

Determination of the 239+240 Pu activity in soils
The soil samples were prepared and measured for plutonium isotope analysis according to the method of Ketterer et al. (2004). Concentrations of 239 Pu and 240 Pu were measured relative to the 242 Pu spike using an Agilent 8800 Triple quad ICP-MS spectrometer equipped with an ICA Apex-IR nebuliser, which allows the detection of ultra-trace-elemental concentrations down to the parts-per-quadrillion level. The resulting concentrations were converted into the combined activity of 239+240 Pu, corrected to the preparation blanks and normalised to the standard reference material IAEA-447 (IAEA-CU-2009-03 2012. The reproducibility of the laboratory preparation was checked using randomly-picked duplicates of the samples.

Calculation of soil redistribution rates
The long-term erosion rate was determined based on the assumption of the secular equilibrium of 10 Be production and decay in the upper 80 cm of the soil. We used the CRONUS online calculator to convert 10 Be concentrations into an erosion rate (Balco et al. 2008; https:// hess. ess. washi ngton. edu). The calculator includes all production channels of 10 Be by secondary cosmic rays through the specified soil depth. For the short-term denudation rates, two conversion models were applied: 1. The profile distribution model (PDM) of Walling and Quine (1990) and Zhang et al. (1990), which was developed to convert FRN inventories into soil redistribution rates for undisturbed sites. A recent study by Calitri et al. (2019), however, showed that it can also be used for agricultural soils: where E is the erosion rate (t ha −1 year −1 ), t is the year of sampling, t 0 is 1963 (year of thermonuclear weapons testing), X is the % reduction in total inventory with regard to the local reference value and h 0 is a profile shape factor; and. 2. Modelling Deposition and Erosion rates with Radio-Nuclides (MODERN) (Arata et al. 2016a, b), which models the FRN depth profile using a stepwise function. For each increment, inc returns a value Inv inc , which is the total inventory of the sampling site, measured for the whole depth profile, d (cm). MODERN provides results in centimetres (cm) of soil loss or gain and can be calculated in t ha −1 year −1 using the equation: where Y is the soil erosion or deposition rate (t ha −1 year −1 ), x * is the soil loss or gain returned in centimetres from MODERN, xm is the mass depth (kg m −2 ), d is the total depth increment considered at the sampling site, t 1 is the sampling year and t 0 is the reference year.

Statistical analysis
To test the statistical association between the chosen soils and relief characteristics, the Pearson's correlation coefficients were calculated using Statistica 13 software. (

Morphology and soil properties
The soils along the two studied transects predominantly developed from loess deposits that had a thickness of up to 1.5 m ( Table 2). In some cases, the thin loess mantle (shallowed by erosion) was found to be underlain by a glaciofluvial substrate (WK2) and/or a calcareous glacial till (WK7, WK8, WK9). These lithic discontinuities were clearly recognisable in the grain size distribution and in part in the bulk density ( Table 2). The loess deposits had a silt loam texture with a dominance of coarse silt, while the glaciofluvial and glacial sediments had a sandy loam, loam or clay loam texture. In general, the loess mantle was completely decalcified along the first transect (except for the glacial till in WK7), while the profiles WK8 and WK9 of the second transect contained some small amounts of calcium carbonate (0.02-5.80%). As a result, the presence of calcium carbonate is linked to a high base saturation and relatively high pH ( Table 2). All soils were characterised by a low SOC content, oscillating between 0.09 and 0.88%, and by a very low nitrogen content (0.02-0.08%).f Almost all the soils revealed a clear morphological differentiation related to clay illuviation. Complete Luvisols, with E and Bt horizons, were only found at the top of the studied hill (WK1) and in the mid-slope section (WK4-WK6), however. The other soil profiles had strong features related to erosion, such as the absence of an E horizon and the incorporation of the Bt into the Ap horizon (WK2, WK3, WK10, WK11), or the deposition of eroded fine-grained material covering the remnants of an older A horizon (WK6; Ab horizon at a depth of 45 cm) and the formation of a thick colluvium (nearly 60 cm in depth) in the toe slope (WK12). In profile WK7, a clear stratification of sand and silt was detected over the glacial till. The sandy material revealed a diagonal lineation. The morphology of the soils suggests a strong differentiation caused by significant geomophodynamic processes acting on the slopes. This led to the shallowing of the loess deposits and exhumation of glaciofluvial and glacial sediments.

Geochemical composition of the soils
The content of major oxides in the soils was relatively similar in both transects (Table S1), oscillating in the range of 61.1-87.0% for SiO 2 , 5.6-11.8% for Al 2 O 3 , 1.0-3.7% for Fe 2 O 3 and 1.7-2.8% for K 2 O. The CaO content mainly ranged from 0.4 to 0.7%, whereas the Ca content was significantly higher in horizons bearing carbonate nodules. The Zr and Hf content often indicates the presence of aeolian material (De Vos and Tarvainen 2006), and they were used to discriminate these from glacial sediments. The average Zr and Hf contents in the loess layer were 302 and 8.9 mg kg −1 , respectively (McLennan 2001). The lowermost horizons of WK2, WK7, WK8 and WK9 had lower Zr and Hf contents, with a maximum of 256 and 6.9 mg kg −1 , respectively. This reflects the lithological discontinuity observed in the field. The relatively high variability of Hf and Zr in soil profiles WK5 and WK6 highlights the dynamics of slope processes and/or loess sorting along the slope (Table S1). The reference soil profile (WK0) situated in the forested area was characterised by similar values for the major oxides, Zr and Hf as in the studied soils, which confirm their origin as loess deposits. The majority of the soils revealed CIA values above 50 (Table S1), corresponding to a low degree of weathering. Only the samples from profiles WK8 and WK9 and the lowermost horizon of WK7 had CIA value of below 50, representing less-weathered material (Nesbitt and Young 1982).

In situ 10 Be and erosion rates
The average content of in situ 10 Be in the soils ranged from 0.9 to 1.5 ( × 10 5 ) atoms g −1 (Table 3). In both transects, the soils situated on the upper part of the slope had slightly lower concentrations than the soils in lower/toe slope positions. The calculated long-term erosion rates were similar for both transects. In profile WK1 (the starting point for both transects), the long-term erosion rate was 0.46 t ha −1 year −1 . In the mid-slope position, the erosion rate increased, oscillating between 0.56 and 0.85 t ha −1 year −1 , while in the toe slope, it decreased again down to 0.44-0.50 t ha −1 year −1 . The site with the highest erosion rate was WK9, with the lowest erosion rate measured in WK6 (Fig. 2).

Activity of 239+240 Pu and erosion rates
The plutonium activity in the soils was generally low, ranging from 0.002 to 0.520 Bq kg −1 (Table 4, Fig. 3). The isotopic ratio of 240 Pu/ 239 Pu (Table 4) can be used to determine the origin of the plutonium in soils (Alewell et al. 2014(Alewell et al. , 2017. The ratio refers to Northern Hemisphere mid-latitude weapons testing fallout (Alewell et al. 2017) when its value ranges between 0.14 and 24 (typically around 0.18), whereas ratios of 0.37 to 0.41 indicate influence from the Chernobyl accident (Alewell et al. 2017). The mean ratio of 240 Pu/ 239 Pu in the investigated soils was 0.18, indicating no influence from the Chernobyl fallout. The erosion rates varied considerably depending on the model applied (Fig. 2, Table S2). The PDM gave higher erosion rates, ranging from 1.4 to 16.9 t ha −1 year −1 , whereas MODERN gave lower rates of between 1.2 and 10.9 t ha −1 year −1 . In general, the soils of the second transect experienced higher erosion rates

Correlation analysis
A positive correlation was detectable between the soil erosion rates (obtained from 10 Be and 239+240 Pu) and the slope gradient (Table 5), indicating, not surprisingly, higher erosion with a steeper slope. The negative correlation between the CIA values and erosion rates suggested that slope processes led to the removal of the topsoil and exposed substrates in deeper horizons that were less weathered or even had a different lithology.

Long-and short-term erosion rates in loess deposits
The calculated long-term erosion rates were consistent with the terrain relief. In both transects, erosion increased in the mid-slope position, which was steeper (WK2, WK3, WK8, WK9), and mostly decreased in the lower, flatter parts of the slope. Erosion processes still occurred in the toe slope positions. However, the erosion rate was lower than in the upper parts of the slope (Fig. 2). This cross-check and agreement with independent data (Table 6) indicated that the assumption of secular equilibrium is plausible.
To determine the short-term erosion rates based on 239+240 Pu, two conversion models were applied. The PDM showed considerably higher values (from 1.4 to 16.9 t ha −1 year −1 ) than the MODERN (1.2 to 10.9 t ha −1 year −1 ). The PDM assumed that the activity of 239+240 Pu had an exponential decay function within the soil profile (Arata et al. 2016b;Calitri et al. 2020). The MODERN, however, did not make any assumptions about the 239+240 Pu distribution along the profile (Arata et al. 2016b) and thus may more precisely simulate the behaviour of FRN (e.g., under ploughing activities) (Alewell et al. 2017). Therefore, the results obtained using this model were deemed to be more reliable. Similar to the long-term rates determined by 10 Be, the mid-slope positions (WK3,4,5,8,9 and 10) exhibited, in general, the highest erosion rates. Although site WK4 exhibits a relatively thick soil profile, the present-day erosion rates are high. This indicates that such high erosion rates must be a recent process, as otherwise such a thick profile could not persist. This mid-slope position has a relatively steep slope and is, therefore, particularly susceptible to erosion. In general, the measured erosion rates are high but are frequently measured or modelled for arable land in Europe (Cerdan et al. 2010).
While in situ 10 Be does not form any inorganic or organomineral complexes nor colloids (because this nuclide is directly produced in the crystal lattice of quartz by the interaction of cosmogenic rays), it also does not move along the soil profile. Pu, however, accumulates in soils and sediments through atmospheric deposition. Consequently, part of Pu might be lost with time through pedogenetic processes or leaching. Pu has, however, an extremely low solubility and a high affinity to organic matter (Alewell et al. 2017) or is strongly adsorbed onto clay particles. In undisturbed soils, essentially, the entire inventory is concentrated in the top 30 cm (Ketterer et al. 2011;Iurian et al. 2015). Because the isotopes of 239+240 Pu are strongly adsorbed, the redistribution of these isotopes has occurred as a result of physical particle movements such as erosion (Alewell et al. 2017). Luvisols are encountered in the investigation. These soils are characterised by clay illuviation. It might be that a small part of Pu migrated along the soil profile due to the translocation of clay particles. This seems rather unlikely, given the fact that the content of Pu strongly decreases with soil depth, except in profile WK12. In this soil, however, no clay translocation was observed. Therefore, the Pu content is primarily affected by soil redistribution.

Soil erosion rates in loess landscapes
The long-term erosion rates of the studied pedons (0.44-0.85 t ha −1 year −1 ), calculated using in situ 10 Be, were in a good agreement with values from other studies of loess landscapes using different methods (Table 6). In Germany, Dreibrodt et al. (2010Dreibrodt et al. ( , 2013 determined the erosion in the Early Bronze Age (0.3-0.6 t ha −1 year −1 ) and of Fig. 2 Short-and long-term soil erosion rates along the studied toposequences. The thickness of the loess mantle is not to scale. Drillings were made down to 2 m the Late Neolithic to Bronze Age (0.4-0.5 t ha −1 year −1 ) by relating the mass of the deposited sediments to the that of the catchment area and/or by reconstructing slope cross sections. Gillijns et al. (2005) and Kołodyńska-Gawrysiak et al. (2018) analysed closed depression catchments in Belgium and Poland, respectively. The obtained erosion rate from 430 AD to today was estimated to be 2.1 t ha −1 year −1 for sites in Belgium (Gillijns et al. 2005), whereas those of the selected sites in Poland varied between 0.24 and 0.27 t ha −1 year −1 from the Late Vistulian (Weichselian) up to today (Kołodyńska-Gawrysiak et al. 2018).
The short-term erosion rates of the studied soils (1.2-10.9 t ha −1 year −1 ), using 239+240 Pu as a tracer, differed from the data presented so far (Table 6). Erosion rates were also calculated for the Trzebnica Hills using USLE . These calculated rates ranged between 2.5 and 4.3 t ha −1 year −1 -although these values lie within the range determined by 239+240 Pu, they were in general a factor of 2 lower. Other studies on erosion rates in Polish loess landscapes have shown a wide variability. Poręba et al. (2019Poręba et al. ( , 2015 determined soil erosion rates of 4.9 to 39.9 t ha −1 year −1 using 137 Cs and 2.2 to 30.7 t ha −1 year −1 using 210 Pb ex . Rejman et al. (2008) used sediment yields over a 10-year period, obtaining erosion rates from 0.4 to 95.0 t ha −1 year −1 , whereas Święchowicz (2016) and Rejman and Brodowski (2010) determined with sediment yields soil losses of 47.3 t ha −1 year −1 and 98.2 t ha −1 year −1 , respectively. Rafalska-Przysucha and Rejman (2015) assessed erosion rates of 24.3 t ha −1 year −1 using  sedimentary archives from closed depressions. The discrepancies between these results might be explained by the different methods used and/or by the landscape relief. USLE involves mathematical modelling only and includes parameters that depend on precipitation and crop rotation, which change over time (Kaszubkiewicz. et al. 2011). Therefore, a direct comparison with FRN-derived erosion rates may not be practicable.
The determination of erosion rates using 137 Cs in Central and Western Europe may be difficult because two different sources exist (nuclear weapons testing in the 1950s and 1960s and the Chernobyl accident). When using this approach, the proportions of these two sources must be determined (Alewell et al. 2017;Poręba et al. 2019). These conditions were met by the investigation of Poręba et al. (2019Poręba et al. ( , 2015. The erosion rates determined with 210 Pb ex were much higher than in our study (Poręba et al. 2019). This apparent discrepancy to our erosion rates is related to a different terrain relief and land-use intensity. Rejman et al. (2008), Rejman and Brodowski (2010) and Święchowicz (2016) measured very high erosion rates, which was likely predominantly due to an almost complete absence of vegetation on their plots. Closed depressions as an archive of soil erosion refer to Modern Times; however, in research presented by Rafalska-Przysucha and Rejman (2015), this covers a longer time span (188 year) than 239+240 Pu analysis. Our data lie within the range of those from other loess areas in Europe. In Germany, the erosion rates range mostly between 2 and 10 t ha −1 year −1 when using the CORINE database (Cerdan et al. 2010). When using other methods such as the reconstruction of the slope cross section or by relating the mass of the deposited sediments to the catchment area, then the erosion rates are normally in the range of 3.2 and 13.3 t ha −1 year −1 (Dreibrodt et al. 2010(Dreibrodt et al. , 2013. In Belgium, analysis of the sediments in a closed depression catchment showed erosion rates from 5.5 to 9.8 t ha −1 year −1 (Gillijns et al. 2005), the erosion rates from sediment yields ranged from 0.5 to 7.9 t ha −1 year −1 (Verstraeten and Poesen 2001;Evrard et al. 2008), whereas with the SEDEM model (sediment delivery model) and 137 Cs, soil losses ranging from 2.0 to 5.0 t ha −1 year −1 and 3.0 t ha −1 year −1 were determined, respectively (Van Rompaey et al. 2001;Van Oost et al. 2003). The relatively narrow range of erosion rates for soils on loess (Europe) may indicate the starting point of their agricultural use.

Cross-check of long-and short-term erosion rates
Applying isotopes that cover different time ranges along two transects has enabled a cross-check of the long-and medium-term erosion rates (Zollinger et al. 2015;Jelinski et al. 2019) but has also provided an insight into the multidirectional course and intensity of erosional processes. The short-term erosion rates are distinctly higher (up to 10 times) than the long-term rates. This finding fits very well with the results presented in Kołodyńska-Gawrysiak et al. (2018), proving that, since prehistoric times, soil erosion rates (0.39-0.57 t ha −1 year −1 ) have increased by a factor of almost 10 compared to the time span between the Middle Ages and Modern Times (3.7-5.9 t ha −1 year −1 ).
Moreover, the results from in situ 10 Be and 239+240 Pu analyses show the intensity of erosion processes along hillslopes. The pedons located in mid-slope positions experienced generally higher erosion rates. This is especially evident in the second transect where the soil morphology exhibits features that are strongly linked to erosion, such as the absence of an E horizon and the incorporation of the Bt into the Ap horizon. Although the soil WK1 apparently looks undisturbed, the CIA values (Table S1) indicate that the profile is disturbed (the CIA value should decrease with increasing soil depth which is, however, not the case here).
Usually, in the toe slope position, eroded material is deposited and, consequently, accumulation is expected (Henkner et al. 2017). At our sites, however, erosion was still measurable here. This highlights that intense erosion occurred along the studied slopes. During heavy rain, the material was eroded from the entire length of the slope, even from sites where the slope gradient was decreased (Fig. 4). The erosion rates have increased during the last few decades  (2015) even at sites with a low slope gradient, such as WK7 and WK12 (Table 1). The morphology of site WK12 indicated the accumulation of colluvial material, but the recent erosion rates were 2.5 times higher than the long-term rates. At site WK7, this increase was almost by a factor of 20. The fact that the soil at WK7 experienced higher erosion rates than at WK12 may be due to the protective effect of redcurrants (Ribes spicatum Robson) for a few years at the latter site, which may have improved soil resilience. A further explanation is that the position of profile WK7 is in the marginal zone of the hill. There, periodic surface runoff may have increased the erosion rates. In general, the high rates of shortterm erosion are related to the intensification and mechanisation of agriculture (Foucher et al. 2014;Kopittke et al. 2019;Poręba et al. 2019), although agricultural activities alone may be responsible for the increase erosion rates. Climate change might be an additional factor causing higher erosion rates because it is giving rise to drier soils, fewer rainfall events, but increasing event intensity (Routschek et al. 2014;Zollinger et al. 2015;Zádorová and Penížek 2018). Kundzewicz and Matczak (2012) also noted an increase in rainfall intensity in Poland. The effect of climate change on soil erosion can, however, not be further quantified. Alewell et al. (2015) posited that tolerable soil erosion rates must be less than or equal to the soil production rates, otherwise the soil would start to degrade. Because soil production rates strongly decrease with the age of a soil, also the tolerable erosion rates (as a soil destructive process) decrease with time. Most of the European soils in the lowlands have an age of > 10 kyr. Tolerable erosion rates for soils in alpine climates and having a surface age of > 10 kyr are between 0.5 and 1 t ha −1 year −1 . Also, Verheijen et al. (2009) showed that tolerable erosion rates for European soils should be less than 1 t ha −1 year −1 . Mediterranean to alpine soils show all after 10 kyr strongly reduced formation rates (< 1 t ha −1 year −1 ; Egli et al. 2014). Consequently, tolerable erosion rates in the range of 0.5 to 1 t ha −1 year −1 are also applicable to our investigation area. Only the long-term erosion rates in our investigated area were in this range or below. However, the short-term erosion rates were 10-22 times higher and thus far beyond any tolerable rates. Therefore, the present-day soil loss considerably exceeds soil production and will lead to significant soil thinning and reduce its productivity .

Conclusions
The long-term erosion rates of loess landscapes calculated using in situ 10 Be were in agreement with the results from sedimentary archives of closed depressions and slope cross-section reconstructions. The short-term erosion rates determined using 239+240 Pu differed from previously published results in Poland that used other approaches ( 137 Cs, USLE, sediment yields, closed depressions). The main reasons for these discrepancies are the differences in conceptual approach and local terrain variability, such as topography, land-use modifications over time, etc. The recent erosion rates determined for the investigated sites using 239+240 Pu are comparable to those from loess sites in Germany, Belgium and France. The short-term erosion rates are up to more than one order of magnitude higher than the long-term rates. This indicates the strong recent impact of agriculture and probably also climate change. The current soil erosion values are far above any tolerable erosion rate. Therefore, the soil is strongly imbalanced, resulting in drastic soil thinning. From a methodological point of view, the use of in situ 10 Be and 239+240 Pu can provide insights into the temporal evolution of soil erosion rates, although loess soils have problematic characteristics, such as clay migration, that may affect their determination.
for the photographic documentation. In addition, we are grateful to two unknown reviewers for the helpful comments on an earlier version of the manuscript. We are also grateful to the owners of the fields: Lech Bienias, Kazimierz Kwaśniewski, Jerzy Placha, Czesław Waszuk, for enabling conducting these research.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.