Chronological Assessment of the Balta Alba Kurgan Loess-Paleosol Section (Romania) – A Comparative Study on Different Dating Methods for a Robust and Precise Age Model

Loess-paleosol sequences (LPSs) are important terrestrial archives of paleoenvironmental and paleoclimatic information. One of the main obstacles for the investigation and interpretation of these archives is the uncertainty of their age-depth relationship. In this study, four different dating techniques were applied to the Late Pleistocene to Holocene LPS Balta Alba Kurgan (Romania) in order to achieve a robust chronology. Luminescence dating includes analysis of different grain-size fractions of both quartz and potassium feldspar and the best results are obtained using fine-grained quartz blue‐stimulated and polymineral post-infrared infrared-stimulated luminescence measurements. Radiocarbon (14C) dating is based on the analysis of bulk organic carbon (OC) and compound-specific radiocarbon analysis (CSRA). Bulk OC and leaf wax-derived n-alkane 14C ages provide reliable age constraints for the past c. 25–27 kyr. CSRA reveals post-depositional incorporation of roots and microbial OC into the LPS limiting the applicability of 14C dating in older parts of the sequence. Magnetic stratigraphy data reveal good correlation of magnetic susceptibility and the relative paleointensity of the Earth’s magnetic field with one another as well as reference records and regional data. In contrast, the application of paleomagnetic secular variation stratigraphy is limited by a lack of regional reference data. The identification of the Campanian Ignimbrite/Y-5 tephra layer in the outcrop provides an independent time marker against which results from the other dating methods have been tested. The most accurate age constraints from each method are used for two Bayesian age-depth modeling approaches. The systematic comparison of the individual results exemplifies the advantages and disadvantages of the respective methods. Taken as a whole, the two age-depth models agree very well, our study also demonstrates that the multi-method approach can improve the accuracy and precision of dating loess sequences.

Loess-paleosol sequences (LPSs) are important terrestrial archives of paleoenvironmental and paleoclimatic information. One of the main obstacles for the investigation and interpretation of these archives is the uncertainty of their age-depth relationship. In this study, four different dating techniques were applied to the Late Pleistocene to Holocene LPS Balta Alba Kurgan (Romania) in order to achieve a robust chronology. Luminescence dating includes analysis of different grain-size fractions of both quartz and potassium feldspar and the best results are obtained using fine-grained quartz blue-stimulated and polymineral postinfrared infrared-stimulated luminescence measurements. Radiocarbon ( 14 C) dating is based on the analysis of bulk organic carbon (OC) and compound-specific radiocarbon analysis (CSRA). Bulk OC and leaf wax-derived n-alkane 14 C ages provide reliable age constraints for the past c. 25-27 kyr. CSRA reveals post-depositional incorporation of roots and microbial OC into the LPS limiting the applicability of 14 C dating in older parts of the sequence. Magnetic stratigraphy data reveal good correlation of magnetic susceptibility and the relative paleointensity of the Earth's magnetic field with one another as well as reference records and regional data. In contrast, the application of paleomagnetic secular variation stratigraphy is limited by a lack of regional reference data. The identification of the Campanian Ignimbrite/Y-5 tephra layer in the outcrop provides an independent time marker against which results from the other dating methods have been tested. The most accurate age constraints from each method are used for two Bayesian age-depth modeling approaches. The systematic comparison of the individual results exemplifies the advantages and disadvantages of the respective methods. Taken as a whole, the two age-depth models

INTRODUCTION
Loess-paleosol sequences (LPSs) are widely spread across Central Eastern Europe ( Figure 1) and provide archives of paleoclimate information crucial for the reconstruction of the Pleistocene environmental evolution (e.g., Marković et al., 2015). In addition to semi-continuous dust input, the formation of LPSs is closely linked to changes in climate and environmental conditions, which are controlled by interglacial-glacial cycles (Gallet et al., 1996;Kohfeld and Harrison, 2003;Svensson et al., 2008;Újvári et al., 2010;Obreht et al., 2017;Zeeden et al., 2018b). Beyond reconstructing Milankovitch-scale climatic variations, LPSs also allow sub-millennial climatic fluctuations to be identified (Moine et al., 2017;Újvári et al., 2017;Zeeden et al., 2018c;Veres et al., 2018). In LPSs of the last glacial cycle, the imprint of past climatic fluctuations is preserved as a succession of weakly developed paleosols or pedogenetically overprinted loess layers, which often gradually turn into adjacent loess beds. These imprints of orbital to millennial climatic variations provide the backbone for chronostratigraphic correlations (at least) on hemispheric scales and interpretation of deposits across environmental boundaries (e.g., Rousseau et al., 2020). In paleoclimatology, this approach was used even before absolute dating techniques became available (e.g., Emiliani, 1955). Nevertheless, to unambiguously link the characteristic paleoclimatic imprints in sedimentary archives to a precisely defined chronostratigraphic interval, time anchors (tie points) are essential (Govin et al., 2015). In their absence, the chronologies of precisely dated archives (e.g., ice core records, speleothems) and even chronologies derived using astrochronological approaches remain speculative (e.g., Lowe and Walker, 2014). Hence, one of the main obstacles for the investigation and interpretation of climatic and environmental information preserved in sedimentary archives is the uncertainty of the age-depth relationship. Most analytical dating techniques applicable to LPSs have age uncertainties that exceed millennial-to-centennial-scale resolution. Thus, the exact timing and duration of individual short-term climatic fluctuations can be difficult to establish (e.g., Zeeden et al., 2018c). However, good age control of these events is of particular importance in understanding broad-scale atmospheric circulation patterns, climatic feedback mechanisms and possibly offsets in timing between different regions and deposit types. Therefore, in this study we combine absolute dating techniques (optically stimulated luminescence, radiocarbon) with correlative dating (tephrochronology, relative magnetic paleointensity, magnetic susceptibility correlation) and integrate the results in a robust absolute age model.
Optically stimulated luminescence (OSL) dating is the most commonly used absolute dating technique in loess chronological research (Roberts, 2008). Sediment accumulation is dated directly by measuring a signal that builds up over time in buried quartz and feldspar grains and resets during particle transport. The use of different grain-sizes and mineral extracts and measurement protocols has been well established (Murray and Wintle, 2000;Murray and Wintle, 2003;Timar-Gabor et al., 2011;Buylaert et al., 2012). The ubiquity of suitable sample material in loess, favorable conditions for signal resetting during eolian transport and a dating range covering at least the last interglacial-glacial cycle are the main advantages of this methods for dating LPSs (e.g., Lauer et al., 2016;Bösken et al., 2017;Timar-Gabor et al., 2017;Bösken et al., 2019;Constantin et al., 2019;Lomax et al., 2019;Fenn et al., 2020a). However, the age uncertainties of luminescence methods are on the order of thousands of years, and significant age offsets have also been reported between different mineral fractions dated as well as different protocols applied. In comparison, the age uncertainty of radiocarbon ( 14 C) analysis of carbonaceous and organic materials is substantially lower allowing centennial-scale resolution chronologies in LPS records (e.g., Song et al., 2012;Újvári et al., 2014c;Song et al., 2015). The method uses the radioactive decay of naturally occurring 14 C in organic matter or carbonates to calculate the time since the last equilibration with atmospheric CO 2 or other carbon pools (e.g., dissolved inorganic carbon in aquatic environments). Due to the half-life of 14 C and analytical limitations, maximum ages of 50-60 ka can typically be resolved with this approach for regular-sized (∼1 mg) samples. Since the atmospheric 14 C content has changed in the past, 14 C ages are calibrated against reference curves to obtain absolute ages (e.g., Reimer et al., 2020). Nonetheless, the application of 14 C analysis to loess is frequently limited by the absence/scarcity of suitable remains such as wood, charcoals, bones and carbonate fossils (e.g., Haesaerts et al., 2010;Pigati et al., 2013;Újvári et al., 2014c;Moine et al., 2017). Moreover, 14 C ages obtained on bulk organic carbon (OC) can be biased by relocation of organic matter in loess, which may result in an underestimation of sedimentation ages (e.g., Hatté et al., 2001). In contrast to OSL and 14 C analyses, magnetic stratigraphy is an indirect and correlative dating approach and includes age uncertainties that are not easily quantified. The method can be used for sediments accumulated beyond the dating limit of radiocarbon and luminescence techniques and is best applied to (several meter) long sequences (Heil et al., 2010;Song et al., 2018a;Zeeden et al., 2018b;Pfeifer et al., 2020;. In general, age control is achieved by "wiggle-matching", i.e., tying significant features of paleo-and rock magnetic parameters to independently dated master curves (Maher, 1999;Hambach et al., 2008). The signals used are typically the magnetic susceptibility and, to a lesser extent, the direction and/or relative paleointensity (RPI) of the Earth's magnetic field. If unequivocal correlations are not possible, potential correlations need to be verified using other dating techniques, which may result in an iterative process. While magnetic stratigraphy is a powerful chronostratigraphic tool for LPSs, it should be noted that ages can be younger than the deposition of the respective LPS (Liu et al., 2015), because soils develop on and into accumulated sediments post-depositionally.
The motivation for this study was to compare the results of the individual dating methods and obtain a robust absolute chronology for the LPS Balta Alba Kurgan (in Romanian: Balta Albȃ Kurgan; hereafter BAK; Figure 1). In this LPS, a well constrained tephra layer is included and forms an independent time marker. We apply luminescence dating, radiocarbon dating and magnetic stratigraphy. In addition to the standard (bulk OC) 14 C analysis, we also test the applicability of compound-specific radiocarbon analysis (CSRA) of leaf wax lipids for age determination of loess, a method that has been rarely applied to LPSs (Häggi et al., 2014). We also use an extended set of OSL dating protocols to validate the OSL age assignment, including one for coarse-grain feldspar that is typically not applied to loess in this region. Furthermore, we combine the results with statistical age-depth modeling to obtain a robust chronology of the BAK loess profile. Our approach illustrates the potentials and limitations of each individual method for dating loess sequences and may help to gain a better understanding of the advantages and disadvantages of the individual methods in dating of LPSs. Our observations and conclusions are highly relevant for future studies that focus on the interpretation of loess sequences to study the paleoclimatic and paleoenvironmental evolution in southeastern Europe and elsewhere.

REGIONAL SETTING AND STUDY SITE
The Lower Danube Basin (LDB) in southeast Romania (Figure 1) is characterized by vast Pleistocene sedimentary deposits of considerable thickness, which are mostly of eolian and fluvial origin. The deposits are incised by broad valleys and alluvial plains of the Danube and its major tributaries. These broad alluvial valleys acted as the major dust source for the coarse components of the loess formed during the Pleistocene (Buggle et al., 2008;Jipa, 2014). The loess deposits of the LDB reach thicknesses of more than 50 m. The thickness and the mean grain size of these loess deposits decrease with distance from potential source areas, i.e., especially the river beds (Smalley and Leach, 1978;Smalley et al., 2009;Jipa, 2014). The most intensively studied LPSs within the LDB are found on the Dobrogea plateau (Figure 1; e.g., Buggle et al., 2009;Buggle et al., 2013;Fitzsimmons and Hambach, 2014;Obreht et al., 2017), along the shores of the Black Sea (e.g., Necula et al., 2015) and along river Frontiers in Earth Science | www.frontiersin.org January 2021 | Volume 8 | Article 598448 cuts of the Danube and its tributaries across the Danube Plain (Constantin et al., 2014;Obreht et al., 2017;Antoine et al., 2019; Figure 1). These sections often exhibit several loess-soil couplets representing glacial-interglacial-cycles. The forelands of the Eastern Carpathians, however, are severely underrepresented in loess research, although the area is covered by vast loess deposits. The BAK LPS (45°16'40''N, 27°17'27''E, ca. 46 m above mean sea level, Figure 1) is located in this area, approx. 40 km east of the foothills of the Carpathian Bend. It is exposed at a road cut of road DC13 dissecting the northern slope of a loess plateau between the rivers Râmnicul Sȃrat and Buzȃu. The Carpathian foothills to the West comprise thick Neogene deposits. The loess plateau frames a basin that includes the alkaline lake Balta Alba, located approximately 2 km east of the section. Lake Balta Alba formed during the Holocene due to neotectonic subsidence in a paleo-channel of the Râmnicul Sȃrat (ILEC, 2019). The sampled BAK sequence is ∼15.4 m thick. Several large pitgrave burial mounds (kurgans) are preserved in the vicinity of the sampling site. The majority of these tumuli in the Lower Danube area most likely date to the Late Neolithic-Early Bronze Age transition (∼3 ka B.C.; Frînculeasa et al., 2015). Our sampled profiles are located close to the margin of a kurgan (Figures 2, 3). The uppermost 0.78 m of the profile formed after kurgan construction as indicated by the preservation of a nonploughed topmost horizon of a chernozem. Between 0.78 and 1.55 m depth, the sediments in the BAK profile consist of disturbed anthropogenic infill of the grave mound. Between 1.55 m and 2.10 m depth, the in-situ lower part of the primary Holocene chernozem (S0, nomenclature according to Marković et al., 2015) overlies bioturbated loess with a leopard skin pattern. From 2.10 m to 7.40 m depth, the sediment is a typical loess including pedogenetically overprinted intervals between 4.30-5.92 m and 6.40-6.95 m. A strongly bioturbated paleosol occurs between 7.40 m and 8.06 m depth. At 8.06-8.10 m depth a macroscopically visible tephra layer has been identified as laterally continuous loose pinkish-whitish patches of coarse silt, resembling other tephra occurrences within the Lower Danube loess Obreht et al., 2017;Zeeden et al., 2018c). Below the tephra, up to 8.70 m depth, the loess is carbonate-rich as well as bioturbated and strongly overprinted by pedogenesis. Underneath, a dark reddish paleosol with two to three distinct pedogenetic horizons is present down to 9.85 m depth. Other (weakly developed) paleosols can be found between 10.0 m and 10.2 m depth, as well as between 10.8 m and 11.4 m depth and are intercalated with carbonate-rich loess. The underlying loess reaches down to 13.2 m depth. A clay-rich, dark greyish to blackish paleosol resembling vertisol features with vertical cracks lies underneath, atop a carbonate-rich loess, which becomes increasingly reddish towards a depth of around 14 m. The paleosol between 14.0 m and 14.5 m depth exhibits chernozem-like characteristics. It overlies a carbonaterich loess horizon, which turns into loess with sub-mm sized pores filled with organic material at 14.7 m depth. At the base of the section, waterlogged, hydromorphic loess and sands are present.

Sampling
For our multidisciplinary dating approach, four subsections (BAK1 to BAK4) were defined along the BAK sequence (Figures 2, 3). The composite profile comprises BAK1 (0-8.50 m), BAK2 (8-12 m) and , which are directly adjacent to one another (thus forming a vertical sampling line) and located slightly to the South of the center of the kurgan. BAK4 is situated approximately 10 m to the South of the main profile (BAK1-3) and spans 3.66 m. Prior to sampling, the outcrops were thoroughly cleaned to obtain fresh and unweathered sediments. Samples for tephrochronological analyses were collected from the visible volcanic ash layer.
Samples for OSL dating were taken using 20 cm long metal tubes. Additionally, a sample of the surrounding 30 cm of sediment was taken for radionuclide analysis (i.e., approximate penetration depth of γ-radiation). In total, six samples were analyzed in this study. For 14 C analysis, twelve samples of 6 cm thickness (∼20 kg of sediment each) adjacent to OSL samples were taken from the BAK1 profile. Oriented paleomagnetic samples were obtained from subsections BAK2 and BAK4 ( Figure 2). These subsections allow for cross-correlation with an overlap of approximately 85 cm, which includes the tephra layer. Sampling was performed by inserting square brass tubes of 2 cm edge length through an orientation frame into the wall. Full spatial orientation was provided by means of magnetic compass measurements combined with an inclinometer to define the dip of the frame. The resulting sediment prisms were subsequently carefully pushed from the tubes into cubic plastic boxes of 8 cm 3 volume to retain orientation. The vertical space between the center of the samples is ∼2.1 cm. Overall 304 specimens were obtained: 120 from the BAK2 profile and 184 from the BAK4 profile. Twelve additional samples from the BAK1 and BAK2 subsections were selected for mineral magnetic analyses.

Tephra Geochemical Analysis
Bulk samples collected from the volcanic ash bed identified at BAK were subjected to glass-shard geochemical analyses at the Bayerisches GeoInstitut, University of Bayreuth, Germany. Measurements were performed on thin sections of volcanic glass shards containing sediment embedded in epoxy using a JXA8200 microprobe (JEOL, JP) and single-grain, wavelengthdispersive electron microprobe analysis (WDS-EPMA). Measurements were performed using an accelerating voltage of 15 keV and a 6 nA beam current. Peak counting times were 10 s for Na, 30 s for Si, Al, K, Ca, Fe and Mg, 40 s for Ti and Mn and 60 s for P. Precision is estimated at <1-6% (2σ) and 10-25% (2σ) for major and minor element concentrations, respectively. Geochemical results are summarized in Supplementary Table  6 and an image of the analyzed glass shards is provided in Supplementary Figure 15 of Supplementary Section 4.

Luminescence Dating
Polymineral fine-grained samples (4-11 µm) were prepared for luminescence dating using conventional preparation techniques (Frechen et al., 1996). A separate batch of the fine-grained fraction was etched in 34% H 2 SiF 6 for one week to extract fine-grained quartz. Laboratory treatment of sand-sized samples included sieving to isolate the 63-90 µm fraction, chemical treatment with HCl (10%), H 2 O 2 (10%) and Na 2 C 2 O 4 (0.01N) to remove carbonates, organic material and clay, respectively, as well as density separation (ρ 2.58 g cm −3 , ρ 2.62 g cm −3 and ρ 2.68 g cm −3 ) to isolate potassium feldspar and quartz. We used an automated Risø TL/OSL DA 20 reader equipped with a calibrated 90 Sr beta source for equivalent dose determination. All samples were prepared and measured in the Cologne Luminescence Laboratory (CLL, University of Cologne, DE). Frontiers in Earth Science | www.frontiersin.org January 2021 | Volume 8 | Article 598448 5 Blue light emitting diodes at 470 nm (FWHM 20 nm) were used for the stimulation of the quartz samples and the luminescence signal was detected through a Hoya U-340 filter (blue stimulated luminescence (BSL); Murray and Wintle, 2000;Murray and Wintle, 2003). The initial 0.8 s of the signal minus a background of the last 4 s were used. Laboratory experiments included preheat plateau tests using a series of preheat temperatures between 180°C and 280°C and dose recovery tests. For the dose recovery test, the samples were illuminated with blue light emitting diodes for 100 s at room temperature and a laboratory dose in the range of the natural dose was given to the samples. Equivalent doses were calculated using the arithmetic mean.
Post-infrared infrared stimulated luminescence measured at 290°C (p-IR IRSL 290 or pIRIR 290 ; Thiel et al., 2011) and infrared stimulated luminescence measured at 50°C (IR 50 ; Wallinga et al., 2000) were used to measure the polymineral fine-grained samples and the sand-sized potassium feldspar samples. Stimulation was carried out with infrared diodes (870 nm, FWHM 40 nm) and the signals were detected through an interference filter (410 nm). The initial 4 s of the signal minus a background of the last 20 s were used. Laboratory experiments included prior-IR stimulation temperature tests  using a range of temperatures from 50°to 220°C and dose recovery tests. For the dose recovery tests, the samples were illuminated for 24 h in a Hönle SOL2 solar simulator and a laboratory dose in the range of the natural dose was applied to the samples. Equivalent doses were calculated using the arithmetic mean. Fading tests (Auclair et al., 2003) were carried out for all samples for the IR 50 measurement protocol. IR 50 ages were corrected for fading following the approach of Huntley (2006).
For the determination of the environmental dose rate, radionuclide concentrations of 238 U, 232 Th and 40 K were measured using high-resolution gamma-ray spectrometry. The dose rate was calculated with the open access web-based program DRAC v.1.2 (Durcan et al., 2015), using conversion factors reported by Guérin et al., 2011, alpha-efficiency values of 0.10 ± 0.014 (pIRIR 290 ; Schmidt et al., 2018), 0.07 ± 0.02 (IR 50 ) and 0.04 ± 0.01 (BSL) and a water content of 25 ± 10% for samples with clay and fine silt contents of about 25-30% and 15 ± 10% for samples with clay and fine silt contents of about 20%. Additionally, the gravimetric water content was determined by weighing before and after drying the dosimetry samples in the oven at 50°C. The internal beta dose rate contribution of the polymineral and the feldspar samples was calculated assuming a potassium content of 12.5 ± 0.5% (Huntley and Baril, 1997). The cosmic dose rate was calculated following Prescott and Hutton (1994). We refer to the Supplementary Section 1 for further details.

Radiocarbon Analysis
Sediment samples for 14 C analysis were oven-dried at 40°C. For bulk OC 14 C analysis, root remains were removed by handpicking prior to subsequent chemical treatment with a standard acid-alkali-acid (AAA) protocol. The purified bulk OC was converted to graphite using an automated graphitization system and analyzed at the CologneAMS facility at the University of Cologne, DE (Rethemeyer et al., 2019). Due to the low OC content of the samples, which produced 0.30 to 0.63 mg C (less than 1 mg C needed for the routine measurement), we determined the background of this specific suite of samples with size matched standards at 42,700 ± 590 14 C years (F 14 C 0.0049 ± 0.0004, fraction modern, Reimer et al., 2004).
For CSRA, ∼1 kg of the dried sediment was ground in an agate mortar. The total lipid extract (TLE) was sequentially extracted by ultrasonication using methanol (MeOH), MeOH: dichloromethane (DCM, 1:1, v/v) and DCM: hexane (Hex, 1:1, v/v). The TLE was saponified with 0.5 M KOH in MeOH and water (9:1, v/v) at 80°C for 2 h. After addition of water, neutral lipids (including n-alkanes) were removed by liquid-liquid phase separation using DCM. For extraction of n-alkanoic acids, the remaining KOH mixture was acidified to pH 1 and the acidfraction was extracted with DCM. n-Alkanes were eluted over SiO 2 columns with Hex. n-Alkanoic acids were methylated with MeOH and 12 M HCl (95:5, v/v) using MeOH of known 14 C isotopic composition (F 14 C MeOH 0.0002) to form fatty acid methyl esters (FAMEs).
Individual n-alkane and n-alkanoic acid homologs were isolated using a 7680 Agilent gas chromatograph (GC, Agilent Technologies, United States) equipped with a CIS 4 injection system (Gerstel, DE) and coupled to a preparative fraction collector (PFC; Gerstel, DE). Trapping of individual compounds was achieved at room temperature. For subsequent processing, individual compounds were flushed from the glass traps with 1 ml DCM. The purity of each isolated compound was determined by GC-FID (Agilent 7890B, Agilent Technologies, United States) and compounds with purities >98% were processed further. Two molecular standards C 18 n-alkane (Fluka, Prod. No. 74691-5g, Lot. 0001448903 with an F 14 C value of <0.0008) and squalane (C 30 isoprenoid; Fluka, Prod. No. 85629-50 ml, Lot. 0001418796 with an F 14 C value of 1.0187 ± 0.0033) were also processed to detect extraneous carbon (blank C) introduced during sample handling and PCGC purification.
Compounds were transferred into tin capsules for combustion and 14 C analysis with the EA-GIS (elemental analyzer coupled to a gas injection system) periphery at CologneAMS . Due to the small and variable sample sizes of the lipid biomarker compounds used for 14 C analysis (9-42 µg C, see Supplementary 2, Supplementary Table 2), the extraneous C introduced during combustion in tin capsules and analysis with the EA-GIS set-up was also quantified. Following the procedure described in Hanke et al. (2017) a constant contamination of 1.0 ± 0.1 µg C with an F 14 C value of 0.44 ± 0.05 was determined (see Supplementary Section 2.2 for further details). To assess if the 14 C age of a sample is distinguishable from background, the massdependent detection limit was calculated following Hanke et al. (2017). For example, the mass-dependent detection limit of a sample with a mass of 16 µg C (Col6363) was F 14 C 0.015, corresponding to ∼38,000 14 C yr BP (Supplementary Figure 7). F 14 C values of all compound samples were corrected for the PCGC processing blank (F 14 C 0.25 ± 0.017 and 0.6 ± 0.3 µg C) Frontiers in Earth Science | www.frontiersin.org January 2021 | Volume 8 | Article 598448 6 and FAMEs for the addition of one methyl group added during transesterification using mass balance.

Magnetic Stratigraphy
For age determination using magnetic stratigraphy, the variations of the relative paleointensity (RPI) of the Earth's magnetic field recorded in the loess sequence and the magnetic susceptibility of the sediment were used. Magnetic measurements were conducted in different laboratories. A list of the devices used in the individual labs is given in Supplement 2 (Supplementary Table 4). The measurements performed include determination of the mass and the magnetic susceptibility (given as mass specific magnetic susceptibility χ) and frequency dependence of magnetic susceptibility (κ FD ) of the oriented samples. The measurement of the natural remanent magnetization (NRM) was followed by progressive alternating field (AF) demagnetization in 12 steps up to 80 mT. The data was plotted in orthogonal vector diagrams (Zijderveld, 1967) and analyzed using the program Remasoft 3.0 (Chadima and Hrouda, 2006). The characteristic remanent magnetization (ChRM) was isolated using principal component analysis (Kirschvink, 1980) embedded in the Remasoft 3.0 software (Chadima and Hrouda, 2006). For most samples the ChRM was determined between the 15 mT and the 40 mT AF step. In every second sample, an anhysteretic remanent magnetization (ARM) was imparted by superposition of a 100 µT DC field and 100 mT alternating field. Subsequently the ARM was demagnetized in a 25 mT alternating field. The RPI was determined as suggested by Levi and Banerjee (1976), i.e., dividing the remanence of the sample measured after 25 mT AF demagnetization by the remanence of the ARM after the 25 mT AF demagnetization. The value of 25 mT was chosen to ensure magnetic cleaning of viscose remanence, whilst a reasonable amount of magnetization is still preserved. Furthermore, mineral magnetic investigations were performed on 12 selected samples to ascertain whether the sediment had preserved paleomagnetic information over geological time scales. Besides determination of the hysteresis parameters, the temperature dependence of the magnetic susceptibility (κ TD ) was measured.
To provide chronostratigraphic constraints by magnetic stratigraphy we correlate the magnetic properties of the BAK site to those of independently dated reference records. We visually align our RPI data with the upgraded version of the high-resolution global paleointensity stack GLOPIS-75 (Laj et al., 2004) and GLOPIS-GICC05 (Laj and Kissel, 2015). Additionally, we compare the data with the RPI values gained from the archaeological site of Poiana Ciresului in northeastern Romania ( Figure 1; Zeeden et al., 2009;Zeeden et al., 2011). To reduce noise, the RPI of the BAK profile is smoothed with a three-point running average. Please note that the original age model of Poiana Ciresului (Zeeden et al., 2009;Zeeden et al., 2011;Zeeden and Hambach, 2020) is based on an RPI correlation for which the Laschamp geomagnetic excursion is set at 38.5 ka BP instead of 41 ka BP (Laj and Kissel 2015). We therefore show an adjusted age model in this study that correlates the four tie points used in Zeeden et al. (2009) with GLOPIS-GICC05. We also compare the magnetic susceptibility data of the BAK pmag composite to the NGRIP δ 18 O record (Andersen et al., 2004) and to different loess records from the region ( Figure 1). As for the RPI correlation, the individual tie points result from wiggle matching of prominent fluctuations. For the time interval of interest, we use the data of the loess sites Rasova (Figure 1; Zeeden et al., 2016;Zeeden et al., 2018c;Zeeden et al., 2019) and Vlasca in Romania (Figure 1; Obreht et al., 2017). Please note that the magnetic susceptibility of Vlasca below the Campanian Ignimbrite/Y-5 tephra layer includes new data for which the age model is preliminary only. Additionally, we correlate BAK to the magnetic susceptibility record of the Batajnica site in Serbia (Figure 1; Marković et al., 2009;Basarin et al., 2014) for which OSL ages have recently been published (Avram et al., 2020). An overview of the data used for the age models at the respective loess sites used for correlation with the BAK site is shown in Supplementary  Table 3.
The age models resulting from the paleomagnetic dating approaches arise from linear interpolation between age tie points. The age-depth relationship was linearly extended to include the top and the bottom of the sequence.

Age Modeling
To combine the age information of the individual dating methods, we used two different Bayesian age-depth modeling approaches. The Bchron software (Parnell et al., 2008;Parnell, 2018) runs in the R environment (R Core Team, 2020). The BChron age-depth modeling approach relies on an increasing sediment column and a gamma distribution as prior information for a sediment accumulation model. The ages and their uncertainties are assumed to be correct. We also apply the ADMin approach (Zeeden et al., 2018a), which is originally tailored to luminescence ages. Here, we apply this model to the combination of ages from different dating methods. Initially, it reconstructs a probability density for the random and systematic parts of uncertainty, which represent unsystematic random effects and a possibly systematic offset, respectively. It then uses the random part for modeling and recombines the model results with systematic uncertainty. This model defines no priors about sedimentation rates, but resamples ages and uncertainty 'as is' and uses combinations, which represent chronological order.
After careful inspection of the results, individual age data for each method were selected as described below (Luminescence dating, Radiocarbon Analysis and Magnetic Stratigraphy). 14 C ages are calibrated using the IntCal20 dataset (Reimer et al., 2020) and reported using the 2σ uncertainty ranges ( Table 2). When multiple 14 C ages were included from individual layers, they were combined by averaging calibrated age ranges. An uncertainty of 10% is assigned to the ages resulting from magnetic stratigraphy in order to provide a conservative estimate. Reported OSL ages are given with 1σ uncertainty, which is considered as such in the age model. The tephra layer identified as the Campanian Ignimbrite/Y-5 (see Tephrochronology) is embedded in the age model using its reported 40 Ar/ 39 Ar age of 39.85 ± 0.14 ka (Giaccio et al., 2017). Please note, calibrated radiocarbon ages refer to AD 1950, while there is no reference date for luminescence dating and magnetic stratigraphy. For modeling, we did not adjust ages to a common reference point in time, since the offset of c. 70 years is well within the uncertainty ranges of the respective methods and also falls within the temporal resolution of our sampling intervals (which is ∼100 yr for magnetic stratigraphy and ∼200 yr for OSL and 14 C dating). Therefore, adjusting the age scales in this study would imply a degree of precision that is not appropriate. The R script for age-depth modeling and data interpolation is available in the Supplements 5 and 6.
Typical IR50 and pIRIR290 signals are depicted in Supplementary Figure 1. First IR stimulation temperature tests showed no dependency of the equivalent dose on the stimulation temperature, which indicates signal stability Supplementary Figure 2B). Dose recovery tests indicate that laboratory given doses were recovered within <5% of unity for the IR 50 and pIRIR 290 signals of potassium feldspars and polymineral fine-grained samples BAK1-3 (7.47 m), BAK1-2 (7.93 m) and BAK2-7  Figure 3). We measured fading rates between 1.30 ± 0.55% g 2days and 3.86 ± 0.52 g 2days for the IR 50 signals and between −2.01 ± 0.37% g 2days and −0.06 ± 0.06% g 2days for the pIRIR 290 signals. Mean equivalent doses are lower for IR 50 measurements ranging from 64.7 ± 3.2 Gy to 93.0 ± 4.7 Gy for IR 50 and from 108.9 ± 5.4 Gy to 149.5 ± 7.5 Gy for pIRIR 290 ( Table 1). Dose rates are between 3.40 ± 0.16 Gy/ka and 3.98 ± 0.23 Gy/ka for IR 50 and between 3.49 ± 0.14 Gy/ka and 4.25 ± 0.2 Gy/ka for pIRIR 290, respectively. This results in fading corrected IR 50 ages between 22.6 ± 2.1 ka and 43.9 ± 3.9 ka (Figure 4). Post-IRIR 290 ages for the corresponding samples are between 28.2 ± 2.5 ka and 44.4 ± 3.4 ka (Figure 4). The pIRIR ages were not fading corrected due to the low values measured. The calculated ages of both sample sets are in good agreement for 50% of the samples and results differ for those samples with low fading rates (∼1.3% g 2days ) as well as for sample BAK1-3 (7.47 m depth; Table 1).

Radiocarbon Dating
Twelve bulk OC 14 C ages were obtained for subsection BAK1 (Figure 3). In the following, conventional 14 C ages are reported. The uppermost sample collected from the Holocene soil at 1.90 m depth provided an age of 6,440 ± 60 14 C yr BP (Table 2; Figure 5). Up to 7.05 m depth, 14 C ages increase with depth and appear in stratigraphic order (28,400 ± 280 14 C yr BP at 7.05 m). Below 7.05 m depth, 14 C ages range from 27,300 ± 220 to 29,100 ± 220 14 C yr BP and are either reversed or agree in consecutive depth intervals ( Figure 5).
Compound-specific 14 C ages of n-alkanes and n-alkanoic acids were analyzed for four depths (Table 2; Figure 5). In general, the 14 C ages of n-alkanoic acids increase with increasing molecular weight (chain length). At 3.35 m depth, n-C 16 alkanoic acid has an age of 2,430 ± 110 14 C yr BP, while n-C 26 and n-C 28+30 alkanoic acids have 14 C ages of 14,300 ± 230 and 15,100 ± 230 14 C yr BP, respectively. The n-C 29+31 alkanes are slightly older than the n-C 28+30 alkanoic acids (15,950 ± 170 14 C yr BP). At 4.77 m depth, one 14 C age could be obtained for n-C 24+26 alkanoic acids (18,700 ± 190 14 C yr BP). Concurrent n-C 29 and n-C 31 alkanes are older (21,600 ± 290 14 C yr BP and 22,000 ± 310 14 C yr BP, respectively) than the n-C 24+26 alkanoic acids. At 6.35 m depth, a total of 5 compounds could be analyzed. Analogous to the sample from 3.35 m depth, 14 C ages increase with increasing chain length: n-C 16 alkanoic acid has an age of 3,190 ± 115 14 C yr BP, n-C 24 alkanoic acid (17,850 ± 230 14 C yr BP) and n-C 26 alkanoic acid (22,500 ± 370 14 C yr BP) are younger than n-C 28+30 alkanoic acids (25,000 ± 350 14 C yr BP). The value obtained for the n-C 29+31 alkanes shows the highest 14 C age (29,200 ± 570 14 C yr BP). At 7.47 cm depth, n-C 24+26 alkanoic acids (23,100 ± 290 14 C yr BP) are younger than n-C 29+31+33 alkanes (28,400 ± 420 14 C yr BP) in agreement with the age patterns observed in the overlying samples.

Magnetic Stratigraphy
The paleomagnetic samples of subsections BAK2 and BAK4 were combined into a composite profile by aligning their magnetic susceptibility curves (Supplement 3, Supplementary Figure 8). The depth scale of this BAK pmag composite corresponds to the depth scale of the main profile (BAK1-3). Overall, 236 oriented samples are included, of which the upper 109 and lower 127 originate from the BAK4 and the BAK2 profile, respectively. The magnetic properties are reported and discussed in Supplementary Section 3 (Mineral Magnetic analyses and paleomagnetic evidence).
In contrast to the absolute dating techniques used in this study, the ages obtained via correlation of the magnetic evidence of the BAK site with reference records (Supplementary Table 3) are a result of correlative interpretation (see Magnetic Stratigraphy). Because of the scarcity of high resolution paleomagnetic records from the region for the relevant time range, the directional data do not contribute to the magnetic stratigraphy in this study. However, directional data of the ChRM of the BAK pmag composite are available in the Supplementary Section 3.3, including a comparison with the inclination and the declination determined at the Poiana Ciresului site.

Luminescence Dating
Dating of quartz (BSL) is often deemed to be the most reliable luminescence technique, because quartz signals bleach quickly during transport allowing the signal to build up from zero after deposition. In contrast, feldspar and polymineral IR 50 signals bleach at a lower rate, but might suffer from signal loss (anomalous fading) after deposition (Wintle, 1973). Feldspar and polymineral pIRIR signals show the slowest bleaching rates, but are supposedly less prone to anomalous fading . The signal loss can be corrected for, but fading correction requires various assumptions (Thomsen et al., 2008;Guérin et al., 2015). Due to the different bleaching rates, similar ages for three different signals of the same sample are a reliable indicator of complete resetting of the luminescence signal prior to deposition Klasen et al., 2018). Ideally, we can distinguish a partially bleached signal from a completely bleached signal by comparison of the fading corrected IR 50 with the pIRIR 290 signal and the BSL signal. We applied this approach to samples from three depth intervals. For sample BAK1-2 (4-11 µm; 7.93 m depth), the obtained ages are similar for all three protocols and we conclude that fading correction of the IR 50 signal was successful and the luminescence signal was completely reset prior to deposition. A complete reset would be expected for eolian sediments. For this sample, all three protocols appear to be equally suitable for dating. For sample BAK1-5 (4-11 µm, 6.35 m depth), BSL and pIRIR 290 ages are in good agreement, while the IR 50 age is younger. This is likely due to the low fading rate of 1.23 ± 0.8%. For sample BAK1-3 (4-11 µm, 7.47 m depth), BSL and pIRIR 290 ages agree well ruling out incomplete bleaching, while the IR 50 age is significantly older ( Table 1). For two of the three samples, the IR 50 luminescence signal in combination with fading correction results in different age estimates, which is most likely due to the IR 50 fading correction and the underlying assumptions included in the modeling approach (Thomsen et al., 2008;King et al., 2018). Therefore, for dating the BAK profile we regard the fading corrected IR 50 ages to be less robust than the pIRIR 290 ages obtained on the polymineral fine grain fraction. As outlined above, the BSL and pIRIR 290 protocols provide equally robust age estimates for the BAK profile. However, we focused on pIRIR 290 dating for the lowermost samples, because we did not want to rely solely on quartz dating for ongoing measurements of older samples (not presented here) as it was shown that quartz signals are problematic for doses above 150 Gy (fine-grain) and 250 Gy (coarse-grain), respectively (Avram et al., 2020). Nevertheless, the corresponding quartz and pIRIR 290 ages of samples BAK1-2, BAK1-3 and BAK1-5 are a good indicator that pIRIR 290 is a suitable dating protocol for the section.
Part of the uncertainty related to the luminescence age calculation is the estimation of the paleo-water content. Most studies use an average water content value for the entire period between sediment accumulation and sampling. For most samples, this period covers different paleoclimatic conditions, suggesting that sediment water contents are not constant. For the BAK samples we measured water contents between 1 and 30%, but it is possible that samples dried out during transport if sample containers were not completely air tight. Therefore, we use the grain size composition of each sample to estimate a paleo-water content (see Supplement 1, Supplementary Figure 4) of 25 ± 10% for samples with a fine silt and clay content above 25% and 15 ± 10% for samples below this value.
The fine-grain (4-11 µm) pIRIR 290 and BSL ages of BAK1-1 (8.23 m depth) and BAK1-2 (7.93 m depth), which bracket the CI/Y-5 tephra agree well with its established 40 Ar/ 39 Ar age of 39.85 ± 0.14 ka (Giaccio et al., 2017), while coarse-grain (63-90 µm) pIRIR 290 ages at 8.23 m and 7.93 m depth seem to overestimate the tephra age slightly (∼1 ka when 1σ error is considered; Figure 4). The older coarse-grain feldspar ages may be explained by contributions of coarse-grained tephra particles that were not removed completely during sample 2 | Radiocarbon isotope data and ages obtained on bulk organic carbon (bulk OC) and n-alkanes and n-alkanoic acids (alk. acids). Calibrated ages are given as 2σ uncertainty ranges. preparation. Dating tephra deposits using luminescence techniques is not always straightforward, even when the surrounding sediments are dated (cf. Bösken and Schmidt, 2020). Independent of the dated grain-size range, samples taken in close proximity to the tephra layer are particularly critical, since fresh geologic deposits are prone to relocation of radioelements leading to an incorrect ascertainment of dose rates (Krbetschek et al., 1994;Biswas et al., 2013;Bösken and Schmidt, 2020). In the Middle and Lower Danube Basins, loess layers embedding the CI/ Y-5 tephra have been often dated using luminescence dating methods (Constantin et al., 2012;Fitzsimmons et al., 2013;Veres et al., 2013;Anechitei-Deacu et al., 2014;Bösken et al., 2017;Obreht et al., 2017;Zeeden et al., 2018c): depending on the specific dose rate of the profiles, equivalent doses of loess layers above the CI/Y-5 tephra range from 114 ± 7 Gy to 200 ± 4 Gy for fine-grained quartz and 104 ± 6 Gy to 221 ± 12 Gy for fine-grained polymineral pIRIR 290 measurements, respectively, while values below the tephra layer range from 136 ± 6 Gy to 286 ± 9 Gy (fine-grained quartz) and from 172 ± 12 Gy to 193 ± 4 Gy (finegrained polymineral pIRIR 290 ), respectively. These values result in ages between 26-44 ka above and 39-57 ka below the tephra. In the BAK sequence, equivalent doses and ages are in the same age range as reported in those studies. This highlights the consistency of luminescence ages among the BAK site and other loess sites in the region. Given the fact that fine-grain pIRIR 290 results seem to provide the most accurate age around the tephra, we only use the fine-grain pIRIR 290 ages for our integrated age-depth model.

Radiocarbon Analysis
Sedimentary bulk OC contains a mix of organic matter (OM) from different sources with different 14 C isotope signatures that may either be deposited syn-depositionally or introduced postdepositionally. For example, OM reworked from older deposits (e.g., eroded from bedrock or soils in the source area of the dust) may have a 14 C isotope signal, which is older than the actual age of loess deposition or potentially radiocarbon-dead. In contrast, OM introduced after sediment deposition, e.g., leaching by percolating pore water, may be younger than the sediment (e.g., Hatté et al., 2001). Typically, aboveground biomass is preferably used for dating terrestrial sediments, since roots can extend up to several meters below the surface introducing modern C into deeper layers of the loess profile. In addition to macroscopic plant remains, epicuticular leaf wax lipids such as long-chain n-alkanes (n-C 27 to n-C 35 ) and n-alkanoic acids (n-C 24 to n-C 32 ) are preserved in sedimentary records. As done by us in the BAK profile, these lipids have previously been tested as target OC for loess 14 C chronologies (e.g., Häggi et al., 2014;Haas et al., 2017;Zech et al., 2017;Bliedtner et al., 2020). In the following we will discuss to what extent the 14 C ages of bulk OC and of leaf wax lipids provide a representative estimate of the actual age of sediment deposition in the BAK profile. Above 7.05 m depth, bulk OC 14 C ages successively increase with depth, suggesting that bulk OC may provide meaningful age information ( Figure 5). Below this depth, 14 C ages do not increase continually, which could be due to very high sedimentation rates in this part of the profile. However, an age reversal at 7.47 m depth indicates that mixing with younger OM (e.g., by bioturbation) likely influences 14 C ages in this part of the profile (Figures 2, 5). This interpretation is strongly supported by the presence of the independent age marker (CI/Y-5 tephra layer, 39.85 ± 0.14 ka) at ∼8 m depth, since the calibrated bulk OC 14 C age from below this depth yields 32.71-33.87 cal ka BP ( Table 2). Analogous to soils, loess deposits are not closed systems and the 14 C isotope signature of bulk OC or individual compounds is determined by their intrinsic turnover times determined by input vs. output flux (e.g., Trumbore, 2009). Soils are also prone to physical disturbance such as bioturbation by roots, which penetrate into deeper horizons and admix younger OM into older sediments (e.g., Zech et al., 2017). Similar to leaves, roots contain some proportion of n-alkanoic acids, with dominance of the n-C 24 homolog (e.g., Wiesenberg et al., 2012). In the BAK profile, n-C 24 alkanoic acid is less 14 C-depleted than bulk OC, suggesting the presence of "younger" root-derived material in the older BAK sediments ( Figure 5). Evidence for admixture of younger OC throughout the loess sequence also comes from the two 14 C ages obtained on n-C 16 alkanoic acids from 3.35 and 6.35 m depth. Both alkanoic acids are significantly younger than the other analyzed OC fractions from the same sediment samples with offsets of >25 kyr. Short-chain n-C 16 and n-C 18 alkanoic acids are produced by plants (leaves and roots, e.g., Wiesenberg et al., 2012) but also derive from microbial sources (e.g., phospholipid n-alkanoic acids in soils; Zelles, 1999). Admixture of OC from roots and microbes with "younger" 14 C isotope signatures has particularly large effects on relatively old samples and likely limits the applicability of bulk OC 14 C dating in loess deposits to sediments <25 14 C ka, which per mass balance are less biased by such contributions. Similar results from LPS sequences in Central Asia (Song et al., 2018b) indicate that this is likely characteristic for LPS deposits in general and not a feature inherent to the BAK site. For dating loess deposits, OM sub-pools that are less prone to admixture of younger material, likely provide a more accurate age assignment. The 14 C ages of the other epicuticular leaf wax lipids in BAK sediments are older than the age of the n-C 24 alkanoic acid ( Figure 5), with n-C 29 and n-C 31 alkanes consistently providing the oldest ages and being in good agreement with bulk OC 14 C ages. However, the systematic depletion in 14 C of n-C 29 and n-C 31 alkanes relative to the corresponding n-C 28+30 alkanoic acids in all investigated samples indicates that different biogeochemical processes affect 14 C signatures of individual compounds and compound classes. The n-alkanoic acids show a systematic increase in age with increasing chain length ( Figure 5), an observation that is well known from marine sediments (e.g., Kusch et al., 2010), but data from soils are limited (van der Voort et al., 2017). In addition to compound source variability, this observation is also considered an effect of selective degradation leading to relative 14 C enrichment of the short-chain n-alkanoic acid pool relative to the long-chain homologs (e.g., Matsumoto et al., 2007). The older ages of long-chain n-alkanoic acids relative to short-chain n-alkanoic acids as well as n-alkanes relative to long-chain n-alkanoic acids likely also reflect a fast turnover of short-chain n-alkanoic acids and slow turnover of n-alkanes in the loess (e.g., Wiesenberg et al., 2008;van der Voort et al., 2017). In addition, n-alkanes in sediments may also contain some proportion derived from geological ( 14 C-free) sources resulting in bias towards 14 C ages that are older than the actual sedimentation age (e.g., Bliedtner et al., 2020). In subsection BAK1, however, the chain lengths distribution of n-alkanes is similar to that found in grasses and herbs implying negligible contribution of n-alkanes from petrogenic sources (cf. Supplement 2). In summary, n-alkane ages closely correspond to bulk OC 14 C ages and we conclude that both n-alkanes and bulk OC reflect the best sediment age estimates in this study. For the age-depth model we therefore use n-alkanes and bulk OC 14 C ages. However, 14 C ages from 7.05 m depth and below are not considered in age-depth modeling, since these ages are biased by admixture of "younger" OC resulting in significant underestimation of actual sedimentation ages.

Magnetic Stratigraphy
The paleomagnetic chronostratigraphic framework for the BAK pmag composite profile is essentially based on the correlation of the smoothed RPI values with the GLOPIS-GICC05 data set. The RPI of loess sequences is not widely used for age determination even though a number of studies gained reasonable results (e.g., Liu et al., 2005;Zeeden et al., 2009;Rolf et al., 2014;Li et al., 2020). At the BAK site, the basic requirement of homogeneity of magnetic properties is fulfilled as shown in the mineral magnetic assessment (Supplementary Section 3.3). The BAK pmag composite and GLOPIS-GICC05 show greater similarities in the upper than in the lower part, but allow for correlation of significant features in both cases ( Figure 6). However, various confounding factors introduce uncertainties when aligning different records. The record of the past variability of the intensity and direction of the Earth's magnetic field in sediments depends to a large extent on the accumulation rate and the complex remanence acquisition processes. Additionally, hiatuses, bioturbation and chemical overprinting may influence the record (Evans and Heller, 2001;Roberts and Winklhofer, 2004;Jin and Liu, 2010). Therefore, not all highs and lows of the RPI reference record may be similarly expressed in the BAK profile. In fact, only when considering the stratigraphic position and age of the CI/Y-5 tephra layer, six tie points (1-6, Table 3) can be determined for visual correlation of the two RPI records ( Figure 6). Additional tie points would be ambiguous. Thus, we limit the number of tie points, but acknowledge that this approach may average out possible variations in accumulation rate. The age-depth data resulting from the correlation ( Table 3) imply that the BAK pmag composite profile comprises two magnetic excursions, Mono Lake and Laschamp. The Mono Lake excursion (Channell et al., 2020) may correspond to the RPI low in the BAK profile between tie points 2 and 3 ( Figure 6). We decided not to tie this RPI low to the GLOPIS-GICC05 curve, because the sections above and below would be substantially stretched and clinched, respectively. It is still under debate whether the Mono Lake geomagnetic excursion occurs globally synchronous, asynchronous, or may be a series of excursions across a given time interval (Korte et al., 2019). From the correlation applied in Figure 6, the putative Mono Lake excursion is dated at 33.2 ka BP. Despite the offset to the Mono Lake excursion in GLOPIS-GICC05, dated at 34.15 ka BP, it is within the time interval associated with the excursion in other parts of the world (Channell et al., 2020). The second geomagnetic excursion is the Laschamp with an age of ∼41 ka (Channell et al., 2020). It is recorded only as a minor decrease of 3 | Tie points of RPI correlation (1-6, cf. Figure 6) assigned to the ages resulting from correlation with GLOPIS-GICC05. Tie points resulting from correlation of magnetic susceptibility data (A-H, cf. Figure 6) are shown as mean ages and their age variation. These values do not represent confidence intervals. the RPI in the BAK composite, but is supported by the stratigraphic position just below the CI/Y-5 tephra layer ( Figure 6). Thus, this excursion defines tie point 4 (8.22 m). Using this tie point, the BAK record can be directly linked to the RPI record of Poiana Ciresului. However, even though Mono Lake and Laschamp can be distinguished at the Poiana Ciresului site (Zeeden et al., 2009), the magnetic record is partly corrupted due to the high number of archeological layers (Nițu et al., 2019) and a visual side-by-side comparison is difficult. Because the RPI record of the Poiana Ciresului site does not extend below the Laschamp geomagnetic excursion, a comparison with BAK is not possible for the complete BAK pmag composite profile. In the lowermost part of the BAK pmag composite, the RPI record shows a similar trend matching two different sections of the GLOPIS-GICC05 record, i.e., at 51.9 ka and 46.9 ka. Considering the thickness of the BAK LPS and the lack of macroscopic evidence suggesting an increase of the accumulation rate such as sandy layers or a general increase in grain size, we tentatively FIGURE 6 | Correlation of the BAK pmag composite. The position of the CI/Y-5 tephra is indicated as a grey horizontal line. RPI of BAK vs. depth is correlated (blue lines) with GLOPIS-GICC05 vs. age and compared to the RPI vs. age of the Poiana Ciresului site. The tie points are shown by blue numbers (1-6, cf. Table 3). The positions of geomagnetic excursions are noted. The blue arrow is described in the text. The χ vs. depth curve of BAK is correlated (red lines) with the χ vs. age curves of the sites Rasova and Vlasca, the χ vs. depth curve of Batajnica (note different scale) and the NGRIP δ 18 O (Greenland interstadials are numbered). Tie points of the magnetic susceptibility correlation are labeled A to H (cf. Table 3). Examples of stylized representations of the trend of the magnetic susceptibility in the lower part are shown with purple lines next to Vlasca and BAK. Tie points labeled with question marks are ambiguous and not considered for age modeling.
Frontiers in Earth Science | www.frontiersin.org January 2021 | Volume 8 | Article 598448 favor the correlation that leads to an age of 51.9 ka instead of 46.9 ka (blue arrow in Figure 6). This solution results in a moderate decrease of accumulation that could be explained by hiatuses, while the alternative would require a sudden and extreme increase of sediment deposition rates. The second magnetic dating approach is based on correlation of the magnetic susceptibility. This technique was first proven to be useful by Heller and Liu (1986) and is widely applied since then. We correlate eight tie points of the magnetic susceptibility curve of the BAK pmag composite profile (labeled A to H from top to bottom) to selected loess records of the region (Figure 6), as well as to the NGRIP δ 18 O isotope record. The ages obtained using this method depend either on the age models of the reference records or the relationship of the fluctuations of the magnetic susceptibility to the interplay of interstadial-stadial events expressed in the NGRIP record. In both cases, the results should only be used with caution. On the one hand, the age models of the selected loess records are partly compiled from multiple dating approaches, including radiocarbon and luminescence dating (Supplementary Table 3). Given that we also aim at testing the feasibility of different dating methods, the use of these data may result in a vicious circle. On the other hand, the magnetic susceptibility is biased by non-climatic factors such as grain-size distribution and local dust influx. The characteristic features of the magnetic susceptibility curve of a section may therefore not coincide with those of the NGRIP record, although they may look similar. Furthermore, a delayed response of the loess properties to climate change may have occurred (Újvári et al., 2014a). Nevertheless, we apply the correlation of the magnetic susceptibility in order to provide a comparison to the RPI derived age model.
Tie point A (5.83 m) is a characteristic peak that is recognized in all reference records (Figure 6). At the BAK site, it is a horizon rich in organic matter, which could be a small pedogenic horizon. Because of its thinness, it might be a short-term regional signal; thus, we did not correlate the peak to the NGRIP record. The small peaks defined as tie points B (6.17 m) and C (6.39 m) correspond to Greenland interstadials (GI) 3 and 4 in the NGRIP curve. They are also easily identified in the magnetic susceptibility curve of Vlasca, but not unequivocally recognizable at Batajnica and Rasova. In contrast, tie points D (7.11 m) and E (7.61 m) are based on a striking resemblance of the shape of the BAK curve with Batajnica. The further correlation of D is based on the correlation of Batajnica with Rasova, which show clear similarities in the range between C and D. By comparison, the Vlasca record shows higher similarity to the NGRIP curve than to the BAK curve. However, the position of E in Batajnica can be correlated very well with Vlasca. The identification of tie point E in the Rasova profile is more complicated since characteristic features are missing.
The CI/Y-5 tephra is found between tie points E and F (8.55 m) in the BAK, Vlasca and Rasova profiles. The correlation of the curves in the profile section below this marker horizon is less obvious than above, which might be due to regional differences in the intensity of pedogenetic overprinting. A general trend consisting of multiple changes between moderate decreases and fast increases of the magnetic susceptibility towards larger depth can be recognized more or less clearly in all magnetic susceptibility records. Stylized representations of this trend are visualized with purple lines along the Vlasca and BAK curves in Figure 6. Based on this observation, the tie points F, G (9.36 m) and H (9.75 m) are defined. F and G are located within an increase of magnetic susceptibility, H is positioned in a small minimum above a stronger decrease of magnetic susceptibility, which is interpreted to reflect the beginning of GI 14. It should be noted from the tie lines that the age model resulting from the correlation of the site Rasova with NGRIP (Zeeden et al., 2016;Zeeden et al., 2018c;Zeeden et al., 2019) is different to the correlation we suggest in this study. Accordingly, the positions of the tie point within the Rasova record provide minimum age estimates for the magnetic susceptibility age model of BAK. In contrast, the polymineral pIRIR 225 ages of Batajnica in this part of the section appear slightly too old (shown in Figure 6) and will not be further discussed in this study.
The BAK age models derived from correlation of RPI and magnetic susceptibility include all tie points (Figure 7) that are not denoted with a question mark in Figure 6. All used tie points are equally weighted and considered in the integrating age model (Integrating Age Models). According to the age models based on magnetic susceptibility and RPI correlation, the sediment at the depth of the CI/Y-5 tephra was deposited 39.3 and 40.6 kyr ago, respectively. Due to the good agreement with the known age of 39.85 ± 0.14 ka (Giaccio et al., 2017), the models appear robust. Overall the models resulting from magnetic dating agree well. In the upper part, the age models diverge slightly, because tie points C, D and E do not correspond to tie points 2 and 3. As a consequence, the age models shift between 5.8 and 8.2 m depth with a maximum age offset of 2.3 ka. In general, ages that are younger than expected might be explained (to a certain degree) by a lock-in effect, which causes the sediment to carry the magnetic signal of a younger time interval than that of sediment deposition. The lock-in depths of paleomagnetic signals in LPSs strongly depend on the syn-depositional water availability as well as at the degree of diagenesis and pedogenesis (e.g., Liu et al., 2015) and are therefore difficult to estimate. However, the BAK profile is characterized by good preservation of subtle sedimentological structures such as the CI/Y-5 tephra layer and macroscopically visible layers of fine sands within unaltered loess as well as in interstadial soils. Thus, an overprinting of the primary depositional signal by postdepositional processes (bioturbation, pedogenesis) seems minimal. A conservative estimate of the potential temporal offset can be calculated assuming a lock-in depth of 0.1-0.15 m in unaltered loess and a slightly higher lock-in depth in interstadial soils. This leads to a potential temporal offset of ∼1 kyr for the average accumulation rate of the BAK pmag composite profile. However, the ages resulting from RPI correlation, which are also based on the magnetic properties, exclude a significant bias from a lock-in effect. Thus, the age discrepancies may rather arise from a mismatch during correlation. In the lower part of the BAK profile, tie point G causes a different shape of the age-depth plot of the susceptibility/ NGRIP based model (Figure 7). This difference can be explained by the different frequency of tie points. Irrespective, the age constraints provided by the magnetic stratigraphy approaches overall suggest rather consistent accumulation rates over time.

Integrating Age Models
Age-depth models are used for constraining the chronological age assuming that samples in stratigraphic order must have chronological order (e.g., Buck et al., 1996;Bronk Ramsey, 2009). Age-depth modeling generally allows for reducing uncertainty of the age estimates when uncertainties of individual ages overlap. This approach has been previously employed to loess deposits in the Danube catchment (e.g., Fenn et al., 2020a;Fenn et al., 2020b;Sümegi et al., 2020). In this study, we select the most accurate ages obtained with each method for age-depth modeling (for details see discussion of individual methods). Hence, age information for the upper half of the BAK profile is provided by 14 C data, while the lower part is covered by results of luminescence dating and magnetic stratigraphy ( Figure 8A). Unfortunately, while we were aiming at high data coverage allowing cross-method validation, robust OSL and radiocarbon ages overlap only at 5.59 and 6.35 m depth.
We apply two age-depth models, one of which is conservative in reducing uncertainty (ADMin; Figure 8C) and the other interpolates ages between the dated horizons (BChron; Figure 8B). Both models use ages and their respective 2-sigma uncertainty. We compare the model outputs in a similar fashion as the data to discuss the chronology of the BAK profile using dating and modeling approaches; model output is available in the supplementary materials ( Supplementary Tables 8-10). The BChron model ( Figure 8B) shows an age span of ∼6 to ∼54 ka BP for the uppermost 9-10 m of the BAK profile. In the upper half of the profile, the sedimentation model produces increasing uncertainties between individual ages. This is due to relatively few age control points in this section and the gammadistribution based sedimentation model inherent to this method. At ∼6 m depth, the model is clearly constrained by the 14 C ages and uncertainty is confined by the high precision of the 14 C ages. Around the CI/Y-5, uncertainty is strongly reduced by this welldated marker bed. For the oldest three samples, the model provides the best fit using the younger ages included in the 95% uncertainty interval. This is probably due to the sedimentation model, which forces the results towards a constant accumulation rate. The ADMin model ( Figure 8C) does not interpolate ages, but only produces model ages for dated depths. Up to ∼5.5 m depth, both the BChron and the ADMin model results are similar; below this depth, the ADMin model results include higher uncertainty between several 14 Cdated layers at ∼6 m depth. In general, since the random portion (not shared portion) of uncertainty is used in the model, it reduces the uncertainty of the age model less. In this study, the random component makes up an average of about 60%. Analogous to BChron, the ADMin age models best fit is shifted towards younger ages for the lowermost three ages, but not as much as the BChron model.
Overall, the BChron and the ADMin models are in good agreement ( Figure 8D). The age-depth models based on magnetic and luminescence data are in excellent agreement and are supported by the age of the CI/Y-5 tephra. Both models have their strengths. The BChron model interpolates ages and directly constructs a chronology for the whole sedimentary sequence. The ADMin model deliberately does not do this to avoid potential incorrect interpolation (cf. Telford et al., 2004;Trachsel and Telford, 2017). Because of their similarity we do not consider one model better than the other. Further constraints of the BChron and ADMin Bayesian models are based on a suite of assumptions about the sedimentation process itself and on the correctness of both ages and uncertainties. The uncertainties of magnetic and luminescence ages may be over-or underestimated, which may influence the age-depth models. Especially for magnetic dating, uncertainties are not easy to calculate, since they include uncertainties that are associated with each stage of paleomagnetic analysis (Tauxe et al., 2010;Heslop and Roberts, 2020). In case of the magnetic data set of the BAK site, we only estimate uncertainty, i.e., we have a component of random uncertainty (e.g., through proxy data) and also systematic effects from correlation with a possibly biased reference chronology. We acknowledge that such an uncertainty should not be used for modeling per se (Zeeden et al., 2018c), but in the absence of alternative quantifiable methods this is a common approach. In this case, we argue that the mean modeled age should be correct within ∼±10% of the real age, because this would account for FIGURE 7 | Age models derived from RPI (blue) and magnetic susceptibility/NGRIP correlation (red). Because different ages arise from the correlation of most tie points by magnetic susceptibility correlation, the mean value is given by the red line framed by a shaded area confined by the minimum and maximum ages. This area does not represent the confidence interval. The tephra layer is indicated in grey at 39.85 ka (Giaccio et al., 2017).
Frontiers in Earth Science | www.frontiersin.org January 2021 | Volume 8 | Article 598448 possibly systematic parts of uncertainty and possibly uncertainty sources which were not accounted for by radiometric dating. Another problem is the confined uncertainty at ∼6 m depth with an age of ∼30 ka, which may not be a realistic model outcome as it relies on the assumption that the 14 C age is unambiguous. Finally, due to the good agreement of the BChron and ADMin age models, we do not favor one model over the other. In both models, the composite section of the BAK LPS analyzed in this study comprises an age of ∼56 ka (with a difference of ∼2 ka). The age-depth trend reflects semi-constant sediment accumulation similar to existing records from the region (e.g., Tecsa et al., 2020;.

Evaluation of the Multi-Method Approach and Implications for Future Studies
Most studies of LPSs combine OSL dating with either 14 C analysis or magnetic stratigraphy to gain age control (e.g., Novothny et al., 2011;Song et al., 2012;Újvári et al., 2014a;Újvári et al., 2014b;Song et al., 2015;Song et al., 2018b;Zhang et al., 2018), while a  combination of 14 C analysis and magnetic stratigraphy is less commonly used (Zeeden et al., 2011). Only in a very few cases all three methods have been applied in concert (e.g., Zeeden et al., 2009;Li et al., 2020). However, none of these studies compared the results of different dating approaches systematically against each other or investigated the effects of age-depth models. In this respect, we provide a unique case study for comparison and integration of the three dating methods most commonly used in loess research. Above, we illustrate that each of the methods provides reasonable age constraints, even though all individually are subject to method-specific biases and limitations. In the following, we compare the results of the BAK section between ∼5.6 and ∼8.7 m depth in a stratigraphic context to evaluate the multi-method approach. To allow us to directly compare magnetic stratigraphy at those depths constrained by luminescence and/or radiocarbon ages, we use linear interpolation between tie points (the respective age values are marked with an asterix). The ages obtained using the different methods show a considerable spread of up to ∼15 ka at 6.35 m depth and >25 ka at 7.47 m depth ( Figure 9). However, some 14 C and OSL ages were identified to be biased by post-depositional processes (n-alkanoic acids with chain lengths <28) and uncertainty in the correction procedure for anomalous fading (IR 50 4-11 µm and 63-90 µm) and will, thus, not be discussed further. Ages obtained by IR 50 deviate substantially from the age constraints provided by magnetic stratigraphy, radiocarbon chronology and the other OSL ages (Figure 9). This indicates that IR 50 measurements are less suitable for many geochronological questionsat least at the BAK LPSdue to reduced accuracy. The remainder of the data cluster in variable groups over depth. At 5.59 m depth, bulk OC 14 C and OSL (pIRIR 290 4-11 µm) ages are in very good agreement (29.9 ± 0.7 ka and 28.2 ± 2.5 ka, respectively). The pIRIR 290 (4-11 µm) protocol was identified as the most reliable OSL protocol for this profile (setting), which is supported by the overlap with bulk OC 14 C ages. In contrast, the results of magnetic stratigraphy seem to be too young in the uppermost section. Down to 6.35 m depth, the magnetic stratigraphy correlation shows an offset of ∼2,000 years relative to the bulk OC 14 C ages. At 6.35 m depth, OSL (pIRIR 290 4-11 µm) and 14 C (n-C 28+30 alkanoic acids) ages as well as magnetic stratigraphy correlation (28.5* ka) range between 28.5 ± 2.5 ka and 30.6 ± 0.7 ka, while OSL (BSL 4-11 µm), 14 C (bulk OC and n-alkanes) and RPI correlation (30.4* ka) age estimates range between 30.3 ± 3.2 ka and 32.3 ± 1.8 ka. Because the latter group includes those ages identified as the most accurate for 14 C and magnetic stratigraphy and we regard BSL and pIRIR as equally suitable we favor the older age assignment here. Although the bulk OC 14 C age at 7.05 m depth (32 ± 0.8 ka) is in very good agreement with the age obtained by magnetic stratigraphy correlation (31.7* ka), it probably underestimates the sedimentation age due to the contribution of younger OC from overlying strata or roots. A bias towards young 14 C ages is also evident for the section below implying that 14 C ages do not reliably reflect sediment ages older than ∼30-32 ka BP in the BAK section. As a consequence, the age of magnetic stratigraphy tie point D (7.11 m) is likely also too young. In contrast, the RPI correlation shows a better fit with the overall trend and likely provides appropriate results. At 7.47 m depth and below, the ages obtained by pIRIR 290 (4-11 µm), BSL (4-11 µm) and magnetic stratigraphy overlap within the analytical uncertainty of the OSL ages. In contrast, the coarse-grain pIRIR 290 ages seem to overestimate the time of deposition at the BAK site. Around 8 m, the presence of the CI/Y-5 tephra supports the ages obtained by OSL and the magnetic susceptibility based magnetic stratigraphy. Please note that the RPI based magnetic stratigraphy uses the tephra to adjust the complete correlation; thus, we cannot conclude that the tephra supports the RPI based ages.
Direct comparison of the methods allows us to discern details that otherwise may not be recognized. The consistent results of certain OSL ages with the 14 C results and magnetic stratigraphy reveal that OSL may provide accurate ages if a suitable measurement protocol is chosen. As shown above, we suggest to refrain from using the conventional IR 50 protocol as this might lead to less precise ages. In future studies, we suggest comparing at least two different luminescence protocols and/or grain-size fractions. It might also be interesting to further investigate coarse-grain feldspar and single-grain measurements. Our 14 C analyses of bulk OC and leaf wax lipids support previous findings that show biogeochemical processes inherent to loess deposits may bias bulk or molecular level 14 C ages. For this specific LPS profile we conclude that sediments can be reliably 14 C-dated up to a maximal age of 30 ka. This age limit may be higher or lower in other LPS, depending on the proportion of younger (e.g., root-derived) OM relative to contemporaneous OM. With respect to magnetic stratigraphy, both magnetic susceptibility and RPI correlations seem to result in valuable age constraints. Whether preference should be given to one of these methods in future studies, depends on the aim of the study. It has been shown that the correlation to δ 18 O isotope records such as NGRIP alone can be quite complicated and that comparison with local records helps determining tie points. Thus, if additional records from the region are available, the use of magnetic susceptibility is less time-consuming when establishing chronologies. However, a higher availability of RPI and directional data of LPSs is desirable. The RPI would enable deducing the regional field behavior of the Earth's magnetic field and the directional data could potentially provide an additional tool for correlation. Finally, we recommend the use of a carefully selected modeling approach that fits both the specific conditions for a given site and the specific age data generated by the chosen methods. However, we do not advise a specific modeling approach, since the model selection criteria are typically subjective.
The overall good agreement of the individual age data obtained with different methods at BAK is a rare observation in LPS studies and highlights the significance of the study (cf. Li et al., 2020). Especially the use of the tephra in the BAK profile offered the unique possibility of comparing the dating methods against a well-defined time marker. We propose to use this approach in other geochronological Frontiers in Earth Science | www.frontiersin.org studies that test dating methods against each other. The results highlight that ages from LPSs should be obtained using at least two independent yet complementary dating methods to allow verification of individual results. However, higher data density and age validation using additional techniques seem to substantially increase the reliability of LPS age models. The increased accuracy of the BAK age-depth models may enable correlation of BAK paleoclimatic proxy information to archives with better age control such as speleothems, marine sediments or ice cores (e.g., Svensson et al., 2008;Govin et al., 2015;Staubwasser et al., 2018;Zhang et al., 2019). At best, a resolution beyond the capability of the discussed methods may be achieved, which may in turn allow comparison of paleoclimate trends across regions and archives. This work may assist future studies to find a suitable combination of dating methods and special techniques to obtain accurate and precise chronologies for loess-paleosol sequences. Due to the variability of LPS (e.g., differences in lithology, grain sizes, carbon contents, mineral compositions, vegetation influences), it remains to be tested whether our conclusions regarding the applicability of each method in the BAK profile can be transferred to other LPSs. We regard our current study as a first step to a generalization to facilitate obtaining high quality age models for LPSs.

DATA AVAILABILITY STATEMENT
All 14 C datasets presented in this study are included in the article and the Supplementary Material and are also uploaded in the PANGAEA data repository (Berg, 2020). The data sets of OSL (i.e., DRAC files for age calculation) are available in the CRC806 database (Nett et al. 2020). The magnetic stratigraphy data sets generated for this study are available in the PANGAEA data repository (Scheidt et al., 2020) and in the data repository of the CRC 806. The raw data of the mineral magnetic analyses supporting the conclusions of this article and the data of the sites Vlasca, Batajnica, Poiana Ciresului will be made available by the authors, without undue reservation. FIGURE 9 | Ages determined for BAK using the different methods. Outliers are not shown. For better visibility, n-alkane 14 C ages are not distinguished based on chain length and only the mean age of the tie points of the magnetic stratigraphy-based age model are shown. OSL data error bars reflect 1σ uncertainty. Calibrated 14 C ages are presented with 95% (2σ) confidence intervals. The age-depth relationship of the tephra layer is depicted as a dashed line.

AUTHOR CONTRIBUTIONS
DB, HB, FL, NK, MM and JR initialized the project. NK suggested the study. UH, JN, DV suggested the site and the investigation approach. SS coordinated the work at the manuscript. JN provided the OSL part, with help from NK, DB, HB, LO: JN (sampling, data interpretation, writing), NK (measurements, data interpretation, writing), DB (review and editing), HB (review and editing), LO (sample preparation, measurements, review and editing). SB provided the 14 C part with help from AS, SK, JR: SB (data interpretation, writing), AS (measurements, data interpretation, writing), SK (data interpretation, review and editing), JR (data interpretation, review and editing). SS provided the magnetic stratigraphy part with help from UH, CL, FM: SS (sampling, measurements, data interpretation, writing), UH (sampling, data interpretation, review and editing), CL (measurements, review and editing), FM (Master's thesis: measurements, preliminary data interpretation; review and editing). DV provided the tephrochronology part with help of UH: DV (data interpretation, writing), UH (measurements, review and editing). CZ provided the age modeling (calculation, interpretation, writing). SP provided the regional information with help of FL, MM: SP (sampling, field survey, interpretation, writing), FL and MM (interpretation, review and editing). The final language check was provided by SK.

FUNDING
The study was funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) -Project number 57444011 -CRC 806 "Our Way to Europe." DV's work on tephrochronology of BAK was supported through a fellowship provided by the Alexander von Humboldt Foundation.

ACKNOWLEDGMENTS
The study is part of the CRC 806 which aims at understanding routes and timing of anatomically modern humans' migration from Africa to Europe. We thank all colleagues who make this CRC successful. We thank A. Zander of the Cologne Luminescence Laboratory for measuring the radionuclide concentrations of the dosimetry samples. We thank D. Haase for sharing the shapefiles of the loess distribution for Figure 1. This is IPGP contribution number 4178.