Reconstructing rates and patterns of colluvial soil redistribution in agrarian (hummocky) landscapes

Humans have triggered or accelerated erosion processes since prehistoric times through agricultural practices. Optically stimulated luminescence (OSL) is widely used to quantify phases and rates of the corresponding landscape change, by measuring the last moment of daylight exposure of sediments. However, natural and anthropogenic mixing processes, such as bioturbation and tillage, complicate the use of OSL as grains of different depositional ages become mixed, and grains become exposed to light even long after the depositional event of interest. Instead, OSL determines the stabilization age, indicating when sediments were buried below the active mixing zone. These stabilization ages can cause systematic underestimation when calculating deposition rates.


Introduction
Soils are constantly changing by continuous or periodical redistribution and mixing of soil material. Humans have triggered and accelerated these processes since prehistoric times by agricultural practices (van Andel et al., 1990;Favis-Mortlock et al., 1997;Dreibrodt et al., 2010;Dotterweich, 2008;Dotterweich et al., 2014). Over time, human-induced erosion and deposition rates in hilly areas increased and larger areas of land were subjected to these soil redistribution processes (Dreibrodt et al., 2010;Kołodyńska-Gawrysiak et al., 2018;Kappler et al., 2019). Currently, soil erosion is one of the major threats to global soil fertility and food security (Pimentel, 2006;Amundson et al., 2015). Although historical erosion may have contributed substantially to total soil erosion budgets (Enters et al., 2008), recent human-induced erosion rates far exceed those found in natural landscapes (Montgomery, 2007) and oftentimes largely exceed rates of soil production through weathering (Alewell et al., 2015). Improved soil conservation measures are urgently required to challenge land degradation and to sustainably conserve ecosystem services (Bouma, 2014). The design of efficient measures will benefit from insights into rates of natural and human-induced soil redistribution, and nature-based approaches may even seek to use those processes and the resulting changes in landforms. Therefore, quantitative data on past soil transport and deposition rates are required.
Human-induced erosion occurs either indirectly, by exposing soil to wind and water (e.g. after deforestation), or directly by working the soil with ploughs. Vertical soil mixing processes are also intensified by humans. In a natural setting, bioturbation is the main mixing process, with the most active mixing in the top few centimetres of the soil (e.g. burrowing animals, root growth and decay) (Gabet et al., 2003). Bioturbation often shows an exponential decline of mixing rates with depth (e.g. Reimann et al., 2017). With the introduction of ploughs, the thickness of the most active mixing layer increased to a few decimetres and mixing intensified. The resulting combination of lateral and vertical transport processes constitutes a complex system in which the rates of lateral and vertical transport of soil material are not easily distinguished (Johnson et al., 2014). Especially in colluvial settings, where tillage-related transport plays a dominant role, the constant vertical mixing complicates the determination of depositional ages through the application of conventional geochronological methods (e.g. by disturbing the stratigraphic order). This disturbance demands additional steps in the geochronological analysis.
Optically stimulated luminescence (OSL) dating is often used to establish depositional ages and rates of colluvial sediments by measuring the OSL signal of grains (Fuchs and Lang, 2009). When buried, soil minerals such as quartz or feldspars build up an OSL signal from incoming ionizing radiation from the surrounding soil matrix and cosmic rays. When the grain is exposed to daylight, this OSL signal is reset, referred to as bleaching. OSL thus measures the last exposure of a grain to daylight: the moment of bleaching.
The application of OSL on colluvial sediments comes with several difficulties. First, burial age can be wrongly estimated due to insufficient bleaching of some grains during sediment transport, or due to entrainment of younger or older grains in the sediment by vertical mixing processes. However, Fuchs and Lang (2009) argue that there is only a small chance of insufficient bleaching in colluvial sediments, because the sediments are transported at the surface during many small events, providing sufficient opportunity for bleaching. Moreover, the constant vertical mixing during ploughing and through bioturbation also increases the chances of bleaching by overturning the soil and might be an even more dominant process for the bleaching of colluvial sediments than lateral transport (Berger and Mahaney, 1990). When considering an eroding hillslope, there are thus three mechanisms of mixing and bleaching affecting the OSL chronology: (1) pre-bleaching on the source location by vertical mixing (Reimann et al., 2017); (2) bleaching during lateral transport by water or tillage (Fuchs and Lang, 2009); and (3) post-bleaching on the target location by vertical mixing .
Second, while pre-bleaching and bleaching during transport reduce the chance of insufficient bleaching of the colluvial deposits, post-depositional bleaching prevents OSL signal buildup after the colluvium is deposited. The reworking of the topsoil through natural (bioturbation) or anthropogenic (ploughing) processes continuously exposes grains to light after their initial deposition and burial. In natural, non-eroding systems ( Figure 1A), the bleaching by bioturbation leads to an OSL age that decreases exponentially towards the surface (Johnson et al., 2014;Román-Sánchez et al., 2019). When such landscapes are exposed to ploughing and erosion, the agedepth patterns change substantially. Ploughing either draws furrows that fill up with sediments by water flow or successive plough activity (ard ploughs) (Nyssen et al., 2011;Pavelka et al., 2017), or tilts entire sods of soil material (mouldboard ploughs) (Andersen et al., 2016). Current plough activity with mouldboards homogenizes the plough horizons in just a few years (e.g. Schimmack et al., 1994;Schuller et al., 2004). Plough layers created by ards are also well mixed, because the furrows fill up with material from the entire range that is reworked (Lewis, 2012). Soil reworking by ploughing thus causes continuous resurfacing of material in the plough layer, and induces bleaching of a part of the grains. As these grains get mixed over the entire plough layer, the original age-depth profile gets overprinted ( Figure 1B) (Huisman et al., 2018). As a further complication, not all surfaced sediments are always bleached, because part of the sediments will be enveloped in large clumps or smaller aggregates (Šarauskis et al., 2008). This Figure 1. Conceptual diagram comparing expected OSL age-depth relations at different landscape positions in an idealized system for: (A) natural systems; (B) agrarian systems. The natural system is expected to be stable (negligible erosion and deposition), with low OSL ages near the surface due to bioturbation. When such natural systems are converted to agrarian systems, the OSL age-depth relations will change accordingly. The non-eroded (hill-top plateau) position in B indicates how ploughing affects the OSL age in the top layer. The erosion (upper slope) and deposition (lower slope) positions indicate how the OSL age-depth relations change relative to the time when the fields were ploughed, but no erosion had yet occurred.
[Colour figure can be viewed at wileyonlinelibrary.com] 2409 RATES AND PATTERNS OF SOIL REDISTRIBUTION IN AGRARIAN LANDSCAPES insufficient bleaching implies that the idealized age-depth profiles in Figure 1B will not be that clear in the field and additional actions are required to extract the age of the plough layer and buried colluvium.
The effects of mixing and erosion by ploughing on OSL ages differ per landscape position ( Figure 1B). On non-eroding positions, ploughing will affect only the top layer of the soil, where the soil is actively reworked. This will cause a discontinuity in the age-depth profiles. On eroding positions, the original soil profile will be truncated and the OSL age of deeper layers gets reset when they are incorporated in the plough horizon. On positions with deposition, the supply of sediments leads to burial of colluvium below the active plough layer. In this stabilized colluvium, the soil is no longer actively reworked and the OSL age can build up. This implies that the burial age measured with OSL in colluvium indicates the moment of stabilization of a colluvial layer instead of the moment of deposition. Mixing depths by tillage change through time following agricultural development, and the depths where sediments get stabilized change accordingly. Therefore, there is an increasing bias in deposition rates when these are calculated with the stabilization ages and corresponding depths. We suggest that more realistic deposition rates can be calculated by correcting for this post-depositional bleaching effect, using the mixing depths at the moment of stabilization.
Our study site is the colluvial infilling of a kettle hole in the Uckermark region, northeastern Germany. In this lowland setting, we expect complex spatiotemporal patterns of colluviation due to very short transport distances and a history of increasingly intensive land use. Kettle holes are the centres of closed depressional catchments, which act as traps for the eroded sediments from the slopes. This means that kettle holes contain the complete geo-archive of the catchment, storing a very local and direct signal of soil redistribution. Kettle hole catchments in northeastern Germany nowadays are heavily influenced by erosion, mainly due to recent tillage (Li et al., 2002;Sommer et al., 2008;van Oost et al., 2009). The determination of accurate deposition rates requires dense vertical sampling of the colluvial soil to obtain a detailed chronology (e.g. Kołodyńska-Gawrysiak et al., 2018). However, current studies often ignore lateral spatial variation in a kettle hole by sampling only one location in the centre of the depressions, where the colluvial layer is often thickest. Due to the complex topographical evolution of kettle holes under erosion (van der Meij et al., 2017), we expect that there is considerable spatial variation in deposition rates, implying that a spatial sampling scheme is required to quantify the colluvial infilling history and hence erosion history.
The objectives of this study are (1) to tailor an OSL dating approach to quantitatively study soil redistribution in agrarian landscapes, solving the difficulties OSL dating faces in colluvial settings and (2) to reconstruct spatial and temporal patterns of colluviation in an agrarian kettle hole. We first describe how we analysed OSL measurements from colluvial settings where partially bleached and rejuvenated samples occur. Using this method, we derived stabilization ages for five densely sampled locations in the kettle hole catchment. Second, we determined the deposition depths corresponding to the different stabilization ages, using a reconstruction of plough depths in the study region. We used these deposition depths to calculate accurate deposition rates. Third, we propose a conceptual model describing the complex spatial and temporal colluvial infilling of agrarian kettle holes in hummocky landscapes. The presented methodology addresses the previously mentioned difficulties when dating colluvial settings and shows new possibilities for using OSL in intensively managed landscapes.

Study area
The study area is located in the Uckermark district in the hummocky ground moraine landscape of Brandenburg, northeastern Germany (Figure 2). After retreat of the Weichselian ice cap (~20 ka ago) (Lüthgens et al., 2011), melting of dead ice formed an irregular pattern of kettle holes, or depressional catchments (Andersson, 1998). Most of these small catchments are not connected to the regional surface water drainage network and act as closed systems. The varying sizes and hydrological properties of the kettle holes are important factors for geo-and biodiversity in the hummocky landscape (Kalettka and Rudat, 2006). The current climate of the Uckermark region is characterized as sub-continental, with a negative water balance (potential evapotranspiration > rainfall). Average annual rainfall in the period 1996 to 2017 was 479 mm (Grünow meteorological station) (DWD Climate Data Center, 2018b). The highest rainfall amount and intensity occurred in the summer months. In winter, precipitation mainly fell as snow, which can trigger surface runoff when melting. Average annual temperature in the same period was 9.1°C, with 0.6 and 17.7°C as average winter and summer temperature (DWD Climate Data Center, 2018a). The study site is part of the experimental CarboZALF-D field , which represents a landscape laboratory with a research focus on the feedbacks between erosion processes and carbon dynamics (Rieckh et al., 2015;Aldana Jague et al., 2016;Gerke et al., 2016;Hoffmann et al., 2018). Before the establishment of the CarboZALF-D experimental site in 2009, an extensive coringbased on a GIS analysis of soil forming factorstook place in the surrounding Quillow catchment (165 km 2 ) to assure representativeness of the soil pattern at the regional scale. Closed depressions make up~12% of the arable land in the Quillow catchment (Deumlich et al., 2010), of which 50% are covered by colluvial sediments.
With the introduction of agriculture in the area, the natural vegetation made way for agricultural use. This agricultural use triggered soil erosion and is thus essential for understanding erosion dynamics and landform evolution (de Alba et al., 2004). The habitation development and land use history in the immediate vicinity of the kettle hole were difficult to reconstruct in great detail, because local archaeological data are limited. However, the main trends in land use, with a focus on changes in plough use, during the Holocene in northeastern Germany could be relatively well reconstructed using a literature survey, including published work on comparable lowland settings in neighbouring parts of Germany, Poland and southern Scandinavia. We provide this reconstruction in Document S1. The outcomes of the reconstruction are used as an approximation of the most likely land use development at the study site. Key elements for our study are the timing of the adoption and expected soil impact of the ard plough and mouldboard plough, and the moment that the kettle hole was artificially drained (Table I).

Data collection and preparation
Data collection For OSL dating, two quantities need to be determined. The first quantity is the equivalent dose (D e ); the estimation of radiation dose received by quartz grains since the last bleaching event. The second quantity is the dose rate (DR); the yearly absorbed radiation dose, based on radionuclide contents of the deposit, sample depth and water content.
We took OSL samples from five locations, with varying colluvium thickness (78-110 cm), topographic position and underlying material (Figures 2 and 3). The soils were described following FAO guidelines (FAO, 2006;IUSS Working Group WRB, 2015) and the German soil classification system (KA5) (Ad-hoc-AG Boden, 2005). At two sites on the fringes of the depression (P2, P3), soil pits were dug and OSL samples were collected horizontally from the pit walls using PVC tubes of  5 cm diameter. In the centre of the depression, samples were obtained at three sites (BP5, BP6, BP8) by coring using vertical PVC cores of 5 cm diameter. For each of these three locations, three cores were sampled. The first core was used for profile descriptions. The second and sometimes third core was used to extract material for D e and DR determination from 5 cm increments. The vertical distance between samples was mostly 10 cm, with some larger distances up to 18 cm.

Lab analysis
The samples were opened and prepared in subdued amber light conditions at the Netherlands Centre for Luminescence dating (NCL, Wageningen University & Research). The material was carefully selected from the centre of the cores to prevent possible mixing and contamination effects due to smearing along the sample tube edges. The 180-250 μm sand fraction was obtained by wet sieving, and chemically treated with HCl (10%) and H 2 O 2 (10%) to remove carbonates and organic matter, respectively. The quartz fraction of the samples was obtained using density separation and then etched using 40% HF for 45 min to remove remaining feldspar contamination and the outer alpha-exposed rim of the quartz grains. The remaining quartz grains were sieved again over a 180 μm sieve to avoid grains that were heavily affected by the etching.
The identification of unbleached or rejuvenated grains optimally requires D e measurements on single grains, to avoid averaging of the OSL signal (Duller, 2008). As previous research in this region showed that less than 5% of the quartz grains carry an OSL signal that fulfilled all quality requirements (Lüthgens et al., 2011), we decided to use very small aliquots (1 mm) containing only 15-20 grains as a proxy for single-grain quartz D e determinations. We applied the SAR protocol of Murray and Wintle (2003). On average, 31% of our measured aliquots passed standard rejection criteria mainly associated with the small aliquot size (i.e. few OSL emitting quartz grains).
The DRs were determined from the remaining soil material in the soil cores. The samples were dried overnight at 105°C to determine gravimetric moisture content. LOI at 450°C was used to estimate organic matter content. The remaining material was ground and mixed with wax in a 70/30 sediment/wax ratio and moulded into 1 or 2 cm thick pucks, depending on the available material.
The pucks were used to measure radionuclide activity concentrations for K-40 and several nuclides of the U and Th decay chains using high-resolution broad-range gamma-ray spectrometry. The measured activity concentrations were converted to infinite matrix beta and gamma dose rates using the conversion factors of Guérin et al. (2011). DRs were then calculated taking into account the attenuation effects of organic matter and moisture (Madsen et al., 2005), as well as the grain-size-dependent attenuation of beta dose (Mejdahl, 1979), the site and depthdependent contribution of cosmic dose (Prescott and Hutton, 1994), internal alpha dose rate (0.010 ± 0.005 Gy ka À1 ) (Vandenberghe et al., 2008) and the contribution of gamma radiation from different soil layers (Aitken, 1985).

Data analysis
We took several steps to determine the stabilization ages and deposition depths of our samples. These steps are summarized in Figure 4 and explained in detail below. The complex analysis is required as the measured D e distributions result from several soil mixing and transport processes ( Figure 1).
Age modelling for determining stabilization ages The OSL age is calculated by dividing the palaeodose, which is the best estimate of the dose received during burial, by the DR. An age model is required to compute the palaeodose from the measured D e distribution of a sample. The most commonly used age models are the central age model (CAM) and the minimum age model (MAM) (Galbraith et al., 1999;. The CAM assumes the OSL sample to be unaffected by partial bleaching and/or post-depositional mixing, and is similar to a weighted mean of the measured D e s. If some part of the grains has an age older than the event of interest (e.g. due to partial bleaching or mixing), then part of the measured D e distribution should be ignored for palaeodose estimation. The MAM provides a means to do this (Galbraith et al., 1999;Galbraith and Roberts, 2012) and the bootstrapped version of the model (Cunningham and Wallinga, 2012) has been shown to provide robust results for small aliquots of both well-bleached and heterogeneously bleached quartz (Chamberlain et al., 2018). A drawback of the MAM is that it is very sensitive to D e outliers at the low end of the distribution.
We expect D e distributions for our samples to be highly overdispersed, mainly due to mixing of grains with different depositional ages, and inclusion of grains that were surface exposed after stabilization. We are interested in the doses representing the stabilization of the sediment below the active mixing zone. However, some grains exhibit higher D e as they were not surfaced by mixing, and some exhibit lower D e as a consequence of sporadic deep bioturbation beyond the active mixing zone (e.g. beetle burrows or roots). This requires additional processing steps, which are outlined in the following sections. From the D e distributions, we iteratively calculated the two standard deviations range of the mean and removed the D e s which fell outside that range window. This process was repeated until no more outliers were identified. The process removed outliers both at the low and high end of the distribution; removal of low outliers is required, while removal of high outliers has no effect on the results obtained in the next step.
Bootstrapped MAM The bootstrapped MAM was used to derive likelihood functions of palaeodoses from the D e distribution (Cunningham and Wallinga, 2012). We used the logged version of the MAM, except when there were D e s with negative values. In those cases, we used the unlogged version. The bootstrapped MAM requires information on the expected over-dispersion of the D e distribution in the absence of mixing or heterogeneous bleaching (i.e. the over-dispersion that would be obtained for a well-bleached unmixed sample of these sediments). Assuming that some of our samples are well-bleached and unmixed, we can obtain the minimum over-dispersion as the best estimate of σ b for our set of samples (n = 32), by applying the bootstrapped MAM to the corresponding dataset of over-dispersion values of our samples calculated through the CAM (Chamberlain et al., 2018). We expect no over-dispersion in our dataset of over-dispersions. Therefore, we set σ b = 0 for determining the σ b of our dataset. This approach robustly selects the σ b from the lowest overdispersions present in our measurements. For our samples this approach provided a σ b of 15 ± 1%, which is reasonable for this aliquot size and this type of sediment. This minimum scatter observed in well-bleached and unmixed D e distributions is most likely caused by unaccounted experimental uncertainties and beta dose rate heterogeneity in the surroundings of the sample. The distribution of palaeodoses was divided by the DR to obtain the distribution of bootstrapped MAM ages.
The age-likelihood distributions should contain all uncertainties in DR estimates, and all unshared errors, both random and systematic. The unshared errors consist of random errors and the unshared systematic error (USS). The USS is a systematic error inside one sample, but a random error between samples (Rhodes et al., 2003). We calculated the unshared errors by summing the errors resulting from determination of moisture and organic matter, and the errors associated with the selected grain size, the cosmic dose determination and the counting statistics of the gamma spectrometer. These errors were inserted in the uncertainty during bootstrapping. Additional errors (shared errors) are associated with the calibration of the beta source, the estimation of internal alpha radiation (Vandenberghe et al., 2008) and the conversion of radionuclide activity to DR (Guérin et al., 2011). These shared errors (2.4%) were added to all presented ages with their uncertainty after the OxCal analysis (see next section).
Bayesian age modelling For samples with known stratigraphic relationships, information on the order can be combined with the age likelihoods to further refine ages and identify inconsistencies using Bayesian approaches. We used the online program OxCal 4.3 (Bronk Ramsey, 2009) to constrain the prior age-likelihood functions (bootstrapped MAM ages) with the stratigraphy into posterior age-likelihood functions, which represent the stabilization ages (Rhodes et al., 2003;Bronk Ramsey, 2008;Cunningham and Wallinga, 2012). We selected the sequence-deposition model, which only uses the relative stratigraphic position of the samples. Other, more complex models can also process information on the absolute depth of the samples, the shape of the depositional curve and the number of depositional events in a chronology (Bronk Ramsey, 2008). However, these extra data and parameters are difficult to estimate for our dataset and therefore, we only used relative stratigraphic information.

Deposition depths and rates
Deposition depths. The ages resulting from OxCal represent the moment that the sediments were buried below the active soil mixing zone in which they were more or less continuously mixed. The actual deposition age is thus overprinted by postbleaching of the sediments due to ploughing (Figure 1). To derive correct deposition rates, we need to correct for this postbleaching in the active soil mixing zone. We focus here on soil mixing by ploughing rather than natural bioturbation, because ploughing is currently the dominant soil mixing process. We cannot accurately estimate the deposition age at the sample depth, but we can estimate the level of the soil surface at the moment of stabilization, using the reconstructed plough depths (Table I). We refer to this level as the deposition depth. The deposition depth equals the sample depth, minus the mixing depth at the time of stabilization ( Figure 5). The deposition age corresponding to this deposition depth is identical to the stabilization age at the stabilization (i.e. the sample depth). Over a period with constant plough depth, deposition rates calculated using either stabilization or deposition depths and ages will be identical. However, when an increase in plough depth occurs between two ages, the rates calculated from deposition depths and ages will show that the actual deposition rates were actually higher than those calculated with conventional stabilization depths and ages. Often, the current soil surface is included in the calculation of deposition rates from stabilization ages and depths. This will also lead to incorrect results if the depths of the dated samples are not corrected for postbleaching. Sediments on the surface are not yet stabilized, so the surface represents a deposition age and depth of zero. Therefore, it should not be compared with stabilization depths of dated subsurface sediment layers. The steep red curve in Figure 5 shows the resulting overestimation in deposition rates.
In some cases, uncertainty in stabilization ages and introduction dates of the different plough types may prevent attribution to one historical plough regime, causing two possible sets of deposition depths for one stabilization depth. Likewise, the un- Figure 5. Conceptual overview of the relation between ages and depths from the moments of stabilization and deposition, see Table I for plough codes. The stabilization ages can be derived through OSL and OxCal analysis. By subtracting the plough depth from the stabilization depths, the corresponding deposition depths are calculated. The uncertainty in plough depth reconstruction propagates into the deposition depths, visualized by the grey uncertainty band. The red line shows the incorrect calculation of deposition rates when comparing stabilization ages and depths of the samples with the deposition age and depth of zero at the soil surface. This conceptual figure does not show the uncertainties associated with the ages, merely the uncertainty associated with the ploughing depths. [Colour figure can be viewed at wileyonlinelibrary.com] 2413 RATES AND PATTERNS OF SOIL REDISTRIBUTION IN AGRARIAN LANDSCAPES certainty in plough depth propagates into deposition depths ( Figure 5). To account for these uncertainties in the final deposition rates, we used a Markov chain Monte Carlo (MCMC) approach where we iteratively resampled stabilization ages (n = 5000) from their likelihood curves and plough regimes to calculate the corresponding deposition depths. We assumed uniform distributions for the introduction ages and mixing depths of the ploughs. Ultimately, the resulting collection of 5000 deposition depths per sample was used to calculate a final probability function of deposition depths. The deposition ages, which are identical to the stabilization ages, and deposition depths were then used to calculate the deposition rates.
Deposition rates. During our analysis, we derived likelihood functions of the deposition ages and depths. These likelihood functions are not normally distributed and can therefore not be presented as an average age with standard deviation. Instead, we used the notation used in radiocarbon dating (e.g. Van der Plicht, 1993) showing the most probable age or depth range (one sigma, 68.3%) for each sample. In the case of bi-modal likelihood functions, the one-sigma interval may comprise two distinct age ranges.
From these age [a (a)] and depth [d (cm)] ranges, we calculated deposition rates [dep.rate (cm a À1 )] between two samples. We assumed the one-sigma age and depth ranges to be normally distributed and used the average value of these ranges as the age or depth of the samples used in Equation (1). In case of two distinct ranges caused by bi-modal distributions, we assumed that the age or depth ranged from the overall minimum to the overall maximum of both ranges. The average of this new range was used as input to Equation (1). This approach recognizes the larger uncertainty associated with the bi-modal distributions, but uses an average value which may not fall in the one-sigma probability range.
To calculate the uncertainty of the deposition rates [ε dep. rateAtoB (cm a À1 ), using Equation (2)], we used half of the onesigma probability range to represent the one-standard-deviation uncertainty of the ages and depths [ε a (a) and ε d (cm)]. We assigned a standard deviation of 1 cm to the sample depths, to account for the 5 cm depth ranges in sampling.

Field and lab results
The five sampled locations varied in their horizonation and properties. A detailed overview of the soils is given in Figure S1. In brief, locations P2 and P3, on the slopes, contain a layer of colluvium 75 and 78 cm thick, respectively. The colluvium is underlain with soil developed from glacial till. Both profiles display macro-and micromorphological signs of clay migration. P2 shows an eluvial and illuvial horizon in the colluvium above the fossil Ah. P3 shows only an eluvial horizon in the colluvium, indicating that the clay migrated to the soil below, which is documented by clay cutans along pores in thin sections (data not shown). There is no former surface horizon recognizable above the eluvial horizons. Locations BP5, BP6 and BP8 were sampled from the centre of the depression and are underlain by peat. These profiles contain a layer of colluvium of 93, 91 and 68 cm and show redoximorphic (gleyic) properties below 71, 61 and 30 cm, respectively, without evidence of clay migration. Table II shows the details of the OSL sample preparation and analysis. Soil organic matter content of the colluvium ranges between 2.1 and 3.6% (average 2.1%). Soil organic matter content in the peaty buried surfaces was much higher (7, 13 and 40%). Total DR for the colluvial samples ranged between 2.04 ± 0.09 and 2.65 ± 0.10 Gy ka À1 . For the two peaty samples with significantly higher sample moisture and organic matter content (NCL-7317060 and -061) the total DR is much lower (1.25 ± 0.11 and 1.85 ± 0.13 Gy ka À1 , respectively). The DR and radionuclide measurements show no correlation with age, depth or location. Combined with the absence of sedimentologically different layers in the colluvium, the fine texture (sandy loam) and the constant reworking of the sediments after deposition, we don't expect any substantial extra heterogeneity in beta DR that is not accounted for by incorporation of σ b in the age modelling.

Age results
We summarized the age results of every step in our analysis (Figure 4) in age-likelihood plots ( Figure S2 and Figure 6). The different curves demonstrate how the age likelihood of a sample changes through the steps in the analysis. The obvious young and old outliers were removed from the measurements during the outlier removal phase, resulting in narrower age distributions. The maximum age of the old outliers is much higher for the fringe positions (P2 and P3, 10-50 ka) than the depression positions (BP5-8, 1-5 ka, with one outlier BP8, 10 ka). This difference in maximum age of the outliers either indicates that during the transport of material from the fringes to the depression, previously partially or unbleached grains were bleached, or that there was no new uptake of old grains from the soil below by deep bioturbation. Some samples have very young outliers (e.g. NCL-7317039, -065, -150, Figure S2). These outliers are sporadically present in the samples and do not follow a stratigraphic order for the different locations when present.
The bootstrapped MAM selected the younger D e s from the measured D e distributions. This results in a narrow age distribution, reflecting the younger part of the age-likelihood distribution after removal of outliers (Figures 6B, C and E). Adding stratigraphic constraints using OxCal results in some alterations in the likelihood functions, compared to the MAM ages. In some cases the OxCal-based stabilization ages are similar to the bootstrapped MAM ages ( Figures 6B, C and E), while for other samples there is a considerable difference between MAM and stabilization ages. The stabilization ages can either be younger (Figures 6D and G) or older ( Figure 6F) than the MAM ages.
Following the increasing plough depth over time (Table I), the differences between stabilization depths and deposition depths increase with time ( Figure 7). For samples with stabilization ages ranging over multiple periods of reconstructed land use, the deposition depths show multiple depth ranges (e.g. the top of BP6, Figure 7).

Deposition rates
The rates calculated with the stabilization and deposition depths (stabilization rates and deposition rates, respectively) are similar, but represent different processes: stabilization rates represent the speed at which sediments get stabilized below the active mixing layer and deposition rates represent the actual process of deposition on top of the soil. The two rates differ between samples where a change in plough regime has occurred: the deposition rates are higher. Stabilization rates between the uppermost samples and the soil surface are much higher than the deposition rates, with the most extreme rate being the one near the surface of location P3. The stabilization rates show an increase from 0.005 to 0.787 cm a À1 (157-fold increase), with one outlier: a highly uncertain rate of 3.098 cm a À1 (sample NCL-7317070). For the deposition rates, which are corrected for deposition depths, the range is from 0.008 to 0.788 cm a À1 (99-fold increase), with the same very uncertain sample with a high rate of 2.168 cm a À1 (Figure 8). The relative uncertainties of the deposition rates increased by an average of 7% when the deposition depths were used, which reflects the uncertainty associated with the reconstructed plough regimes. At the fringes of the depression (P2 and P3), colluviation started at 2900-2200 BCE and 700 BCE-300 CE, respectively. In the central depression, this process started much later. BP5 and BP8, on the edges of the peat body, started receiving sediments a few hundred years ago; 1480-1600 CE and 1710-1800 CE. BP6, close to the centre of the depression, has the youngest onset of colluviation, 1820-1870 CE. This shows a gradual infilling of the depression from the sides. Deposition rates have increased substantially over time (Figure 8). The rates at the onset of colluviation, from P2 and P3, were on the order of 0.01 cm a À1 . These rates have increased by two orders of magnitude up to~1 cm a À1 in recent times.

Discussion
The benefit of advanced age and rate modelling in colluvial settings OSL dating has been widely applied in colluvial settings (Eriksson et al., 2000;Fuchs and Lang, 2009;Dreibrodt et al., 2010;Kappler et al., 2018;Kołodyńska-Gawrysiak et al., 2018). Colluvial deposits often undergo active soil formation (e.g. bioturbation) and land use (e.g. tillage), and this is likely to induce post-bleaching of sediments after deposition ). In contrast, transport and mixing processes can also promote post-depositional uptake of older grains. Consequently, colluvial sediments can contain non-bleached grains and post-bleached grains in their corresponding D e distributions, which is clearly demonstrated in our samples ( Figure 6 and Figure S2). The young outliers are only sporadically present, and their ages do not follow a stratigraphic order. These grains thus represent a more recent mixing process than post-depositional tillage (e.g. bioturbation) (Wilkinson et al., 2009). To avoid any influence of these grains on the MAM, these outliers were removed before applying the MAM. One Table II. Age and rate results of the OSL analysis. Stabilization ages and deposition depths are presented in one-sigma ranges. Deposition rates are presented as average and one-sigma uncertainty. The rates are calculated between the sample and the next higher sample or the soil surface. Details of the OSL sample preparation and intermediate age results during the analysis are provided in Table S1   Location Depth ( Figure 6. Selection of age-likelihood plots for a few representative samples. Data for all samples are provided in Figure S2. The plots on the left (A-D) are from location P2 and the plots on the right (E-G) are from location BP5. Note that for plot A, a linear X axis is used to allow visualization of the differences around the zero-year mark. The curve representing the measurements corresponds to the D e distribution. The other curves correspond to ages at various steps in our analysis (see Figure 4) could argue that young and old outliers even each other out, allowing the use of a simple CAM to calculate the age. However, that would be true only coincidentally, and in fact we found different amounts of rejuvenated and unbleached outliers in most of our D e distributions, invalidating such an approach. Therefore, in colluvial settings with a high probability of post-depositional overprinting of the OSL ages (e.g. by intensive agrarian land use), a more advanced OSL data analysis is inevitable. Our analytical framework that deals with these young and old outliers combined an outlier removal approach with a bootstrapped MAM, and constrained age distributions with stratigraphic information (OxCal). Finally, we corrected the stabilization depths for the post-bleaching effect using archaeological information, allowing us to calculate deposition rates. We described in the Introduction how ploughing mixes the entire plough horizon, so we expect a homogeneous distribution of bleached grains throughout the entire plough horizon. We illustrate the validity of our framework to extract the ages of these bleached grains with sample NCL-7317068 ( Figure 6A), which was taken from the plough layer of location P2. We expected this sample to have an age close to zero years, because of its constant reworking. The age will not be completely zero, because it takes some years for the recently bleached grains to be reworked through the plough layer (cf. Schimmack et al., 1994). The measured D e distribution of this sample corresponds to a one-sigma age range of À225 to 389 a (mean: 82 a). This wide distribution is partly caused by the occurrence of non-bleached grains and shows that grains in the active mixing zone are not all completely bleached. The different analysis steps (Figure 4) gradually resulted in decreasing uncertainty and diminishing OSL age towards zero years. After outlier removal the one-sigma age range was À25 to 109 a (mean: 42 a), analysis through the MAM resulted in ages of 14-23 a (mean: 18.5 a) and after Bayesian modelling the final age estimate was 13-22 a (mean: 17.5 a). This finding shows that using a MAM and OxCal to calculate stabilization ages from heterogeneous colluvial samples substantially decreases age uncertainty and provides a reasonable age for a plough layer. The small deviation from the zero age may represent the time it takes for bleached grains to get reworked into the plough layer, or it can represent the lower dating limits of our framework due to methodological and statistical limitations (Madsen and Murray, 2009;Wallinga and Cunningham, 2013). The two patterns we see for sample NCL-7317068 are also visible in most other samples ( Figure S2 and Table S1): (1) the uncertainty of the OSL ages reduces substantially during all steps in the analysis and (2) the stabilization age often corresponds well to the mode of the density function of the measured age before outlier removal. The second point indicates that the samples contain sufficient grains that have been bleached during ploughing. There were also differences between stabilization ages and measured ageswith bi-modal distributions (sample NCL-7317040, -151) or stratigraphic discontinuities (sample NCL-7317041). Our framework favoured the younger, lower peak of the distributions, because we used a MAM, or because the younger ages differed the most from the ages of the samples above and below the sample. Such complex age distributions can thus result in a loss of age information following our framework. Nonetheless, the age we are interested in, the younger stabilization age, gets selected. Another limitation of OxCal is the changing of some ages for no apparent reason, such as the lowest sample at BP5 (NCL-7317060, Figure 6G). This was the lowest sample taken at this location, in the former surface horizon. There was no chronologic information present below this sample, so there was no clear reason to reduce the age of this sample. This effect did not occur in the samples from the former surface horizons of the other locations.
Only the near-surface sample NCL-7317068 allowed a formal validation, because we know the age of this sample should be close to zero. Additionally, the stratigraphic consistency of the measured age distributions, similar age patterns in most samples and the correlation of our ages to archaeological indices of land use change (see next section) give us confidence that our methodology successfully extracts the stabilization ages from the heterogeneous D e distributions. Our approach could be applied to other colluvial settings as well, after adjustments to local archaeological information.
By bootstrapping the MAM and calculating the deposition depths using an MCMC approach, we created a robust model to calculate stabilization ages and associated deposition depths. The ages and depths comprised the uncertainties resulting from the lab analysis as well as from the archaeological reconstruction, while avoiding disturbance by unbleached and rejuvenated grains. The results of our study clearly show the impact of the post-bleaching effect on the calculation of deposition rates (Figure 7). Deposition depths, corrected for this post-bleaching, are located substantially shallower in the soil than the conventional stabilization depths. Rates derived from these deposition depths indicate more intensive deposition before agricultural intensification than rates derived from stabilization depths. Also, the high deposition rates near the surface that would be obtained when using stabilization depths were shown to be significantly lower in reality. The correction for post-depositional mixing increased the uncertainty of the deposition rates. If the post-depositional mixing is ignored and stabilization depths and ages are used to calculate rates, the uncertainty would be lower. However, this approach would calculate rates of stabilization rather than the rates of deposition in which we are interested. The two types of rates are correlated, but represent a different process. Because we know that post-depositional mixing occurs in agrarian settings, and we know that the uncertain reconstructed ploughing depths affect the deposition rates, this uncertainty must be considered in the calculation of deposition rates. On top of that, the relative uncertainty of deposition rates is on average only 7% higher than the uncertainty of stabilization rates, indicating that the stabilization ages contribute the majority of uncertainty to the deposition rates.
OSL is not the only geochronological method that must take the effect of post-depositional mixing into account. Our reasoning also applies to radiocarbon dating, which dates organic soil particles such as pollen, organic matter and charcoal, and radionuclide dating (e.g. Cs-137) (Alewell et al., 2014). These particles are also mixed through the soil by bioturbation and tillage. The depths at which these particles are found therefore do not necessarily correspond to the layer at which they were deposited.
Spatial and temporal patterns of deposition Spatially varying onset of deposition Our sampling design with high spatial density revealed unexpected patterns in colluviation. The sediment ages can be grouped in two age groups: older sediments (>300 years) and recent sediments (<300 years). The onset of colluviation varied by several thousands of years over very small distances, in two distinct phases of mass movement (Figures 3, 7 and 9), as discussed below. This demonstrates that sampling locations for geochronological studies in colluvial settings should be selected with great care and, specifically, not solely placed at positions with the thickest colluvium.
An explanation for the two phases of deposition can be found in changes in the hydrological situation of these kettle holes. Originally, the depressions were filled with water and peat (kettle lakes, Figure 9A). With deforestation and the start of agriculture (Table I), slope sediments were exposed and started to move towards the depression, where they were deposited around the kettle lake (old sediments in Figure 9). Deforestation also led to wetter conditions in the kettle hole by reduced evaporation (Twine et al., 2004;Lamentowicz et al., 2008;Mao and Cherkauer, 2009), increasing the water level and peat growth simultaneously with sediment deposition. There are two plausible blocking mechanisms preventing the sediments from entering the central depression: (1) the rising groundwater level and resulting peat formation following deforestation ( Figure 9B1) and (2) a terrace-like landform, formed by (pre)historic ploughing around the depression ( Figure 9B2). These terrace-like shapes surrounding kettle holes can still be observed in currently forested kettle holes (van der Meij et al., 2017) as well as in a large number of agrarian kettle holes in the Uckermark region which are currently filled with water. Both blocking mechanisms may have co-existed and both would be related to the wet conditions inside the central depression.
In the 18th century, population density increased substantially after the Thirty Years' War (Hinze, 1981). The higher demand for agricultural land forced people to drain the wet depressions in the landscape. The construction of drainage ditches and tile drainage enabled the reclamation of these wet spots in the landscape. In the last centuries, agricultural drainage greatly reduced the amount and extent of water bodies in kettle holes (Kochanowska et al., 1998). In our land use reconstruction (Table I and Document S1), we used historical maps to show the construction of a ditch draining our kettle hole between 1787 and 1826 CE. Larger areas of peatlands in the Uckermark were drained even earlier (Behre, 2008;Kaiser et al., 2012Kaiser et al., , 2014. The lower water level made the central depression accessible for agricultural use. The start of artificial drainage coincides with the oldest colluvium ages at locations BP5 and BP6 (Figures 7  and 8). The sample from the former surface horizon in BP5 and the lowest samples from BP8, located at the edges of the colluvium, show ages which are~100 years older than the introduction of artificial drainage, indicating earlier deposition at the edges of the peaty depression. The drainage substantially increased deposition rates in the depression, forming a layer of young sediments, which originated from the fringe positions and higher upslope areas ( Figure 9C). Ultimately, the entire kettle hole became covered with young sediments ( Figure 9D). This conceptual model shows that spatial and temporal variation in deposition has to be considered in geochronological studies in colluvial settings. Geochronological approaches such as ours can prevent incorrect correlation of colluvial layers, by separating colluvial layers with different ages.
The interpretation above indicates that the fringes of the depression have two depositional phases, which are visible in their age-depth profiles (Figure 7). These phases correspond to a hiatus in soil horizonation in the colluvium. The old colluvium in P2 and P3 shows signs of clay migration (lessivation). In the recent sediments, these signs are not optically visible, indicating a different pedogenesis; closer to the centre of the depression, upward movement of water and gleyzation become more important. The absence of a former surface horizon between the two layers in P2 and P3 indicates either erosion of a part of the earlier colluvium, or incorporation of the former surface horizon into the newly added colluvium by ploughing. The presence of clay migration features can thus be used as a proxy for colluvium age, as these features are more frequently observed in old colluvial infillings of depressions, but not in younger ones (e.g. Kołodyńska-Gawrysiak et al., 2018). A consequence of possible erosion on the fringes is underestimation of deposition rates through the hiatus, because part of the old sediments would have been removed.
The two distinct phases of colluvial infilling we observed in our study area should be expected in other kettle holes underlain with peat as well, which are the dominant type of kettle holes in the Uckermark region (see 'Study area'). Another geochronological study covering four kettle holes in the Uckermark region showed substantially different stratigraphies and chronologies, but none of these locations had a substantial peat body in the subsurface (Kappler et al., 2018). Rather, the peat that was observed was interlaid with colluvial deposits. This regional diversity in kettle hole infilling again stresses the importance of understanding the stratigraphy in a study area first before taking geochronological samples. Each location might have experienced a different geomorphic and land use history, affecting the chronology of the deposits and the representability of the results for a larger area. Phases and rates of deposition As indicated before, two phases of colluvial soil redistribution could be distinguished in our kettle hole. The start of the oldest phase corresponds to the timing of (pre)historic settlement (Kappler et al., 2018). Colluviation at location P2, 2900-2200 BCE, started a few hundred years earlier than the first ages of Kappler et al. (2018), but the acceleration of colluvial dynamics between 1500 CE and now is detected in both study sites. The oldest ages at P2 coincide with the habitation increase and agricultural intensification in the Neolithic in the area (Jahns, 2000(Jahns, , 2001, and a general increase of luminescence ages derived from colluvium in northeastern Germany (Kappler et al., 2019). The start of colluviation at P3, 700 BCE-300 CE, coincides with settlement during the Iron Age and Roman Age (Jahns, 2001). The different upslope topography of positions P2 and P3 could also have caused the differences in onset of colluviation. P2 had a steeper upslope area, with a much shorter distance to the crest. The shorter transport distances and higher erosion potential could have caused this position to accumulate colluvium at an earlier stage than P3, where the transport distance was much larger.
The deposition rates at CarboZALF-D increased by two orders of magnitude between the start of colluviation and the present (Figure 8). Recent deposition rates are similar to those reported for other kettle holes (e.g. Frielinghaus and Vahrson, 1998). The increase in deposition rates coincides with an increase in population density and agricultural intensification from the 18th century onwards (Hinze, 1981;Dotterweich, 2008). The more intensive agricultural practices, together with larger, more connected fields (Dotterweich, 2008), increased the sediment flux towards the depression. The recent high deposition rates found in the central depression should, however, not only be attributed to this agricultural intensification, but also to the reclamation of the central depression. The reclamation created a large elevation gradient, which promoted transport of previously deposited sediments from the fringes into the depression (Figure 9). The deposition rate is thus dependent on landscape position and does not correspond one-to-one with increased erosion rates in the catchment.
The deposition rates we found are comparable to the erosion rates found in various Central European loess catchments, as compiled by Kołodyńska-Gawrysiak et al. (2018). It is, however, difficult to make a quantitative comparison between erosion and deposition rates from different studies, due to their dependence on local topography and land use history. This also complicates the evaluation of our improved methodology to calculate deposition rates from colluvial sediments; deposition rates calculated with deposition ages and depths show the same trends as the rates calculated with stabilization ages and depths, except for the youngest rates. Nonetheless, we believe that our approach yields more realistic results than conventional methods, because it considers the effects of incomplete bleaching and post-bleaching, which are processes that affect sediments in every setting where tillage plays a role in sediment mixing and transport.

Conclusions
Determining the timing and rates of erosion and deposition is difficult due to various mixing and transport processes. This hampers the linking of deposition histories to land use history, and the assessment of landscape responses to deforestation and ploughing, for example. OSL methods are often used to determine deposition histories, but the effect of postbleaching by tillage, for example, at the location of deposition is often not considered. We present a novel approach 2419 RATES AND PATTERNS OF SOIL REDISTRIBUTION IN AGRARIAN LANDSCAPES to include tillage effects in the calculation of stabilization ages of colluvial sediments and corresponding deposition rates. We applied this framework to 32 OSL samples to study the complex spatial and temporal infilling of a kettle hole in an agrarian hummocky landscape in northeastern Germany. Our main findings were: • Colluvial deposits show a large spread in the level of bleaching of the grains. However, OSL ages can be extracted using advanced age modelling techniques and corresponding deposition rates can be derived using an archaeological reconstruction of plough depths. • The combination of outlier removal, a bootstrapped MAM and stratigraphic correction substantially reduces the age uncertainty of colluvial deposits. If colluvial deposits are ploughed after deposition, a correction for post-bleaching by tillage is required to derive deposition ages and depths from the stabilization ages determined by OSL dating. Although this procedure results in slightly greater uncertainty, the deposition age estimates are more accurate. • The calculated stabilization ages generally coincide with the mode of the measured D e distributions, indicating that a substantial amount of grains are well bleached during postdepositional mixing.
• The investigated kettle hole shows a complex infilling, with the oldest colluvium at the fringes of the depression and the youngest colluvium in the centre of the depression and overlying the older colluvium at the fringes. This was probably as a consequence of changing hydrological situation. The start of deposition varied by thousands of years over a distance of tens of metres. • Deposition rates in the kettle hole increased by two orders of magnitude between the start of colluviation (~3000 BCE) and the present. The spatial and temporal deposition patterns depend heavily on landscape position and land use, such as artificial drainage. Therefore, a local increase in deposition rates does not necessarily imply increased erosion rates throughout the catchment. This study demonstrates that determination and interpretation of colluvial sediment ages and deposition rates is possible using state-of-the-art OSL dating methods. We recommend the use of numerical dating techniques in combination with thorough soil-geomorphological field surveys and historical analysis of land use to elucidate complex spatiotemporal patterns of landscape change.