Abstract
Introduction
In the last couple of decades, many studies have aimed to reconstruct the Quaternary paleo-landscape architecture and stratigraphy of the North Sea Basin (Andresen et al., 2022; Bailey et al., 2020; Cohen et al., 2014; Cotterill et al., 2017a, 2017b; Coughlan et al., 2018; Fitch et al., 2005; Fleischer et al., 2022; Gaffney and Fitch, 2022; Gaffney et al., 2020; Hepp et al., 2012, 2017, 2019; Huuse and Lykke-Andersen, 2000; Lohrberg et al., 2020; Phillips et al., 2017; Prins and Andresen, 2021; Spinney, 2008; van Heteren et al., 2014; Winsemann et al., 2020). The repeated glacial-interglacial and post-glacial reshaping of the Quaternary landscape of the North Sea resulted in complex and heterogenous morphological structures. These structures which are now paleo-archives of past environmental changes can be related to various geological processes combined with major climatically-driven relative sea-level changes. In general, the Quaternary North Sea sub-surface is dominated by glacio-fluvial sands which are difficult to date and distinguish from each other due to strong similarities in material and grain characteristics (Fleischer et al., 2022). Thus, the identification of stratigraphic elements such as buried tunnel valleys, river channels, moraines, dunes, peats, and paleosols are important to unravel the Quaternary stratigraphic frameworks of the North Sea (Coughlan et al., 2018; Emery et al., 2020; Fleischer et al., 2022).
While this study focused on the southern part of the German North Sea sector and within a specific stratigraphic interval (uppermost 50 m), tunnel valleys and other comparable geologic structures were formed during all the Quaternary glaciation periods in the North Sea and are dominant in the upper strata of the sub-sea floor. The abundance, distribution, orientation, size, morphology, and inferred age of buried tunnel valleys across the North Sea are well documented in literature (Coughlan et al., 2018; Hepp et al., 2012; Huuse and Kristensen, 2016; Kristensen et al., 2007; Lohrberg et al., 2020; Lutz et al., 2009; Ottesen et al., 2020; Prins and Andresen, 2019; Prins et al., 2020; Stewart et al., 2013; van der Vegt et al., 2012). These buried tunnel valleys, formed at ice sheet margins, were active under direct glacial conditions and their infills played a key role in reconstructing the offshore Quaternary stratigraphy of the North Sea (Coughlan et al., 2018; Fleischer et al., 2022).
The Elbe Palaeovalley (EPV; Figure 1) is a drowned river valley located in the southeastern North Sea with a width of about 40 km, a length of approximately 210 km and a valley base of about 65 m below the present-day sea level. This valley extends from the southwest of Helgoland in a northwest-southeast direction through the whole German North Sea sector with an expression in the form of a shallow depression on the modern sea floor. Özmaral et al. (2022) and Papenmeier and Hass (2020) gave a comprehensive seismo-stratigraphic classification of the EPV infill which were ground-truthed with core data. The Özmaral et al. (2022) classification sub-divided the EPV infill into five major intervals as shown in Figure 2. These intervals from old to young are named Incised Channel Infills (FU), Transparent Lenses (LU), Oblique Unit (OU), Sigmoid Unit (SU) and Transparent Unit (TU). The dominant OU, SU, and TU units are of marine origin and the EPV itself incised Older Pleistocene Sediments denoted as BU. The LU units occurred as an in-channel aggradation in which two forms were identified: an inclined or dipping LU1 unit and an upper or shallower LU2 unit (Figure 2). The EPV thus played a major role in the late Pleistocene and early Holocene paleo-drainage network of northwest Europe as the meltwater and sediment flux from these areas were collected by it through a network of paleo-rivers (Ehlers et al., 2004; Pisarska-Jamroży, 2015).

Map of the German North Sea sector showing the study area (black square) in relation to the Elbe Palaeovalley as well as the location of the sub-bottom profiler section crossing the Elbe Palaeovalley (brown solid line).

Sub-bottom profiler section from the Elbe Palaeovalley with seismo-stratigraphic classifications (modified after Özmaral et al., 2022). The inset section details the western part of this sub-bottom profiler section alongside the Transparent Lenses 2 (LU2) and its relationship with other units. The Transparent Lenses 2 (LU2) is considered to be linked to the Palaeo-Ems deposits studied herein based on lithological correlation and radiocarbon ages from basal peat. See Figure 1 for location of the sub-bottom profiler section.
Recent studies of these paleo-valleys and other comparable natural archives have improved the understanding of the impact of relative sea-level changes on coastal landscapes and re-organization of the drainage system of the North Sea (Andresen et al., 2022; Blum and Törnqvist, 2000; Blum et al., 2013; Bourillet et al., 2003; Hepp et al., 2017; Meijer, 2002; Roberts et al., 2018; van Heteren et al., 2014). Extensive high quality geophysical, drill core, paleo-environmental, and sedimentological datasets from scientific and commercial cruises have continued to improve the understanding of paleo-landscape structures in the North Sea. The result of one of such research efforts is a discovery of a tributary to the EPV, the now submerged paleo-extension of the modern Ems River, which revealed the setting and evolution of the early Holocene landscape (Hepp et al., 2012). A small section of the river is seismo-acoustically and sedimentologically already well documented (Figure 3; Hepp et al., 2017, 2019) and a detailed deduction of its overall course is key in understanding the re-organization of the drainage system of the area since the Last Glacial Maximum (LGM). Based on this study, preliminary deductions from geophysical profiles obtained during the HE499 cruise in 2017 revealed that the linear course earlier defined for the Palaeo-Ems system (PES) in the southern sector (yellow dash line in Figure 3, Hepp et al., 2017) is an artifact owing to an initial less dense seismic coverage. Furthermore, toward the downstream part of the river course, the nature of the river is not fully understood as studies by Hepp et al. (2017) revealed that the channel fans out into about three sub-channels forming a delta at the western flank of the EPV (Figure 3).

Details of the study area north of the Juist Island showing the location of the geophysical profiles available for this study together with the interpreted and hypothetical course of the Palaeo-Ems system from Hepp et al. (2017). The key geophysical profiles discussed in this study are labeled with their respective code and shown as brown solid lines. The purple stars show the locations of the sediment cores (Core 09-3 and Core 09-6) used in this study.
Understanding the paleo-drainage networks of the North Sea is crucial for various reasons. In pre-historic times, coasts, and river banks were sought-after settlement areas, because hunting and fishing offered a good livelihood there. At the same time, the flowing waters represented a natural network of paths through marshland, moorland, swamps, and sometimes mark boundaries or territories (Warnke et al., 2014). Furthermore, the locations of the rivers are often preserved thereby holding relicts of the former landscapes and human habitation (Bailey et al., 2017, 2020; Coles, 1998; Flemming et al., 2017). A detailed knowledge of the paleo-drainage stratigraphic infill is also important for various offshore geotechnical and engineering installations including but not limited to cable routing and windfarm development (BSH, 2021; Weihrauch et al., 2009).
This study thus aimed at understanding the drainage systems of the paleo-landscape of the North Sea shelf by answering the following distinct but closely related research questions. Firstly, what are the exact drainage pathways of the Palaeo-Ems system (PES) at the western margin of the Elbe Palaeovalley as well as its southern connection to the Modern Ems River? Secondly, what is the Palaeo-Ems system/Elbe Palaeovalley spatial and stratigraphic relationship? What is the PES spatial relationship with other glacial and post-glacial structures within the study area? The answers to these research questions will contribute to the morpho-stratigraphic reconstruction of the paleo-landscape which developed after the LGM within the German sector of the North Sea.
Materials and methods
Data
The study area is located north of the East Frisian Islands Norderney and Juist, within the southern part of the German North Sea sector (Figure 1). It extends between northings 5,955,000 to 6,000,500 m and eastings 351,700 to 396,700 m (WGS84 UTM Zone 32N Projection) and covers an area of about 1800 km2 (Figure 3). The dataset available for this study consists of high-resolution sub-bottom profiler (SBP) sections (sediment echo sounder) acquired during the 2013 HE403 and 2017 HE499 cruises using the research vessel Heincke. Additional grids of closely spaced Boomer and core data from commercial surveys carried out between 2007 and 2010 were sourced in areas with now restricted access due to ongoing windfarm operations, as complementary dataset (Figure 3).
Geophysical data
A total of about 380 lines of both sediment echo-sounder (SES) and Boomer data constitutes a dense grid of high-resolution 2D acoustic profiles available for this study (Figure 3). The SES profiles were collected using a ship-mounted sediment echosounder Innomar SES2000 medium-100 model with an acoustic power greater than 247 dB (~5.5 kW). Differential Global Positioning System data were also obtained along each of the SES track lines. The SES data provide a very high resolution of the uppermost 15–20 m below the sea floor. The SES2000 is a parametric echo sounder system in which signals of new frequencies are created by the interaction of two transmitted signals of different frequencies owing to non-linearities in the sound propagation. On board of the research vessel Heincke used during the HE499 expedition, the adopted SES2000 medium-100 model was operated in multi-frequency mode and secondary frequencies of 4, 8, and 12 kHz were generated at alternate times. The acquired acoustic data were digitally converted to the standard SEGY data format using the SesConvert tool and merged with differential GPS navigation. A total of 116 east-west profiles in a north-south progression which are between 200 to 250 m apart were acquired and due to extreme wave and weather conditions at the time of the survey, north-south profiles and sediment cores required for ground truthing were not collected. About 80 SES profiles acquired during the HE403 expedition in 2013 were also made available for this study. These profiles were acquired in a similar manner as those from the HE499 expedition.
The Boomer data from the commercial survey were obtained using a single-channel seismic system operated using an energy of 200 J and a 3 m single-channel streamer with eight hydrophones. The Boomer data produced a high-resolution image of the upper 15–30 m with a very high vertical resolution of about 20–50 cm. The data which were acquired in both the north-south and east-west directions were recorded using the NWC software (Nautik Nord GmbH) with a shot rate of 3.5 shots per second and a 16-bit resolution. A grid comprising of about 184 profiles with a mesh width of approximately 200 m was acquired in total. Conventional seismic processing procedures were carried out specifically to meet the needs of high-resolution shallow-water SBP data with special emphasis on water-column static correction. Where applicable, water-column height variations caused by tides, weather and current were accounted for using procedures outlined by Abegunrin et al. (2020) prior to interpretation.
The processed datasets were then interpreted in a standard workstation. The SBP data were analyzed using the IHS Kingdom Advanced software for conventional seismic and stratigraphic interpretations. Key stratigraphic surfaces were identified based on detailed information observed from SBP sections which are marked by changes in impedance contrast. The identified surfaces were then mapped across the survey area based on amplitude strength, character, and continuity. Channels and other seismic-scale structures defined on the SBP sections were also identified and mapped. Due to the different acquisition parameters and degree of resolution, the SBP data was gridded using the cubic spline gridding algorithm while the Boomer data was gridded using the flex gridding algorithm. The gridded surfaces were then used in generating isopach maps for intervals of interest. The identified seismo-stratigraphical units and structures were then tied to core samples at core sites by linearly converting the two-way travel time to depth based on an average sound velocity of 1550 m/s, a value fitting well for the uppermost 50 m below the sea floor based on previous studies (Lutz et al., 2009; Winsemann et al., 2020; Hanno Keil, personal communication, February, 2021). This enables a direct integration of SBP and core data in ascertaining and correlating boundaries and lithology with seismo-stratigraphical units and structures. The sinuosity index (SI) of each of the meandering loop of the channel course was obtained by taking the ratio of the length between the inflection points of each meandering loop and the straight-line distance between these inflection points. The sinuosity index was then plotted against distance to obtain the sinuosity index chart which was sub-divided into highly sinuous (SI: >1.5), moderately sinuous (SI: 1.5–1.1) and straight (SI: <1.5) according to Charlton (2008) classification. The channel gradient was also calculated by dividing the vertical difference in the depth of the channel base between two points by its horizontal distance.
Core Data
Two sediment cores from commercial surveys carried out in 2009 within the vicinity of the SBP traverses, both with overall recovery >90%, were made available for this study. Due to commercial restrictions, the cores selected for this study were coded 09-6 and 09-3 (Figure 3). These cores, which are used for ground-truthing the geophysical data used in this study, were recovered from water depth between 29 and 31 m using the University of Bremen operated Geo-corer 6000 system. The coring system yielded sediment cores with a recovery rate between 90% and 97% for all encountered lithologies. This core was segmented into 1 m sections, immediately opened, photographed, visually described and documented on board and later further characterized in onshore laboratories. These two cores alongside five other cores from the PES course are the basis of Hepp et al. (2019) publication.
Results
Deductions from sub-bottom profiler sections: Mapped structures and seismo-stratigraphic units
Three distinct geological structures namely buried tunnel valleys, lacustrine-fill unit, and straight channel structure were identified alongside the PES over the SBP data coverage. These structures were mapped and described in order to avoid misinterpretation and false mapping which is typical in areas with complex generations of geologic structures. The mapped structures are separated by distinct seismic reflectors marked by major depositional or erosional boundaries. On representative SBP sections presented in Figure 4a to f, the base and infills of various generations of buried tunnel valleys (purple line), near horizontal lacustrine-fill unit (green line), channels of the PES (red line), other previously unidentified straight channels (blue line) and depositional units (black and yellow lines) were identified. Figure 4a to f also depicts the seismo-stratigraphic relationship of these structures. Due to the intensive erosional activity in the study area, glacial and post-glacial structures often appear at similar stratigraphical depths. All identified structures are completely buried with no expression on the modern sea floor.

Uninterpreted (left) and interpreted (right) representative sub-bottom profiler sections showing the identified structures and their inter-relationships within the study area. See the brown solid lines in Figure 3 for locations of all the sub-bottom profiler sections. The Palaeo-Ems deposits, which are the target of this study, are younger than the buried tunnel valley and lacustrine-fill unit deposits and older than the straight river channel deposits.
Two prominent generations of buried tunnel valleys trending approximately in a NE-SW and NNW-SSE direction were identified and mapped on SBP sections (Figure 4f). These structures are characterized by steep dipping flanks and are more V-shaped. The internal structure of these buried tunnel valleys exhibits homogenous infills and sharp boundaries with the underlying and adjacent strata. The buried tunnel valleys are depicted by seismic reflector R1 (purple color in Figure 4c, d, and f) and the valley floor generally exhibits an undulating and overdepeened course. Seismic facies analyses of the NE-SW buried tunnel valley infills as seen in Figure 4f revealed steeply inclined low- to moderate-amplitude, semi-continuous to discontinuous reflectors. The NNW-SSE buried tunnel valley have a steeply inclined moderate- to high-amplitude, semi-continuous to continuous internal reflectors (Figure 4c, d and f) and it is sometimes incised by the PES (Figure 4c). Both buried tunnel valley generations do not show any evidence of internal structuring or erosional horizons. The lithology and age of this structure cannot be ascertained due to the lack of sediment cores and geological ground-truth data.
Furthermore, west of the PES, a shallow depression structure was identified as an adjoining structure to the PES (green line in Figure 4a–c) and also as a discrete structure (green line in Figure 4d). Both forms of this structure have the same approximate erosive depth (about 15 m), architecture and infill patterns (Figure 4a–d). This structure is referred to as lacustrine-fill unit due to the characteristic near-horizontal seismic reflectors of its infill. On SBP sections (Figure 4a–d), the base is marked by seismic reflector R2 (green color) characterized by an uneven and continuous surface which have a sharp and distinct contact with the substrate. The seismic facies of the infills are defined by high-amplitude, continuous reflectors at the base overlain by a characteristic low- to medium-amplitude, semi-continuous to continuous reflectors (Figure 4d). The stratigraphic infill of this structure was also not deduced due to no cores being drilled through it.
On SBP sections, the PES, which is the focus of this study, was mapped as a filled incision with an asymmetric geometry characterized by both steep and gentle flanks in areas where the SBP traverses cut its course perpendicularly (Figure 4a and c). Reflection terminations marked these channel flanks which cut into underlying sediments and sedimentary structures. The channel outline is depicted by seismic reflector R3_1 (red line in Figure 4a–d) and represents an uneven erosional stratigraphic surface truncating underlying basal seismic unit. Seismic facies analyses of the PES infills are characterized by high-amplitude and continuous reflectors compared to the substrate in which the channel incised. The infill reflectors are undisturbed and parallel to each other (Figure 4a–d). In some areas, the uppermost infill of the channel shows a characteristic secondary incision by younger structures (Figure 4a and b). Based on the resolution of the SBP data, four seismo-stratigraphical units (Unit II, III, IV, V) were deduced from the base of the channel to the sea floor, while Unit I underlies the base of the channel (Figure 5a). In chronological order, Unit I which is characterized by low-amplitude and discontinuous reflectors (Figure 5a) corresponds to fine- to medium-grained sand (Figure 5b). Unit II, which is defined by medium- to high-amplitude, semi-continuous to continuous reflectors (Figure 5a), is composed of peat with some intercalated sand layers (Figure 5b). This unit marks the base of the Palaeo-Ems system. Unit III is defined by medium-amplitude, continuous reflectors (Figure 5a) with a combined lithology of silt, fine sand and clay (Figure 5b). This unit constitutes the main PES stratigraphic infill and the top of PES is marked by seismic reflector R4. Unit IV is predominantly fine sand seen on SBP sections as medium amplitude, semi continuous to continuous reflectors (Figure 5a and b). Unit VI is characterized by high-amplitude, continuous reflectors (Figure 5a) and the lithology of this unit is fine sand with some pockets of medium-grained sand (Figure 5b). These seismo-stratigraphic units were also identified and correlated across the two core samples used in this study (Figure 5b) with the major differences being in their thickness and depth of occurrence.

(a) Boomer profile M273 with stratigraphy of Core 09-6 and seismo-stratigraphic interpretation. (b) Correlation of seismo-stratigraphic units across the two selected drilled cores used in this study. See Figure 3 for locations of the Boomer profile and cores.
In addition, another structure, straight channel structure, marked by seismic reflector R5 were also identified and mapped within the study area (blue line in Figure 4a, b, and e). The base of this structure is an uneven erosional stratigraphic surface that cuts underlying structure and older sediments. In most cases, the course of this structure occurred as discrete channels (Figure 4e) while in some areas, they incised the upper parts of older structures such as the PES (Figure 4a and b). The typical internal pattern of the infill is a transparent unit defined by low- to medium amplitude, continuous and sometimes dipping internal reflectors (Figure 4e). The lithology and age of these straight channel structure cannot be deduced accurately because of lack of sediment cores.
Downstream changes along the reconstructed Palaeo-Ems system
An overview of the mapped Palaeo-Ems system along its paleo-flow direction is presented in the interpreted isopach map shown in Figure 6b. The northern course of the documented PES channel, west of the EPV (Figure 3), as well as the southern course toward the coast are hereafter referred to as the northern and southern extensions respectively. Within the southern extension, the southernmost area (about 20 km north of the Juist Island) was characterized by intense erosional processes and all the SBP sections from this area are mostly featureless. From about 20 km north of the Juist Island, the PES was identified as incisions on SBP sections and the overall mapped course spans over 100 km in length (Figure 6b). This channel course runs from south to north and the mapped course presented in Figure 6b depicts the last open phase of the river which was filled during subsequent marine transgression. Overall, the channel shows evidence of lateral migration and abandoned meanders that is, oxbow lakes (Figure 6b).

Isopach maps of different parts of the study area based on available profiles for this study. Isopach map (a) with interpretation (b) of the reconstructed Palaeo-Ems system course. Note: the black arrows in (a) depict the limits between the Boomer data and the sub-bottom profiler sections and both data are gridded using different gridding algorithms. P1 and P2 in (b) represents major sub-paths of the Palaeo-Ems course. Isopach map (c) with interpretation (d) of the buried tunnel valley courses, which are older and stratigraphically below compared to the Palaeo-Ems. Isopach map (e) with interpretation (f) of the straight river channel courses, which are younger and stratigraphically above compared to the Palaeo-Ems.
The course of this reconstructed ancient PES can be divided into three main tracts from upstream to downstream: (1) the meandering, single channel tract defined by a high sinuosity index (2) the moderately sinuous, single channel part marked by an abrupt decrease in sinuosity index (3) the straight part within the delta area with low sinuosity index (Figure 7a–c). The sinuosity index of each of the meandering loops within the upstream tract is relatively high and variable ranging from about 1.6 to about 3.0 with an average value of 2.3. The moderately sinuous part is less meandering with sinuosity index between 1.3 and 1.5 and an average value of 1.4 (Figure 7a). The straight part has a relatively low and less variable sinuosity index with values being constantly less than 1.1. Figure 7b also depicts a gradient plot obtained from the depth of the deepest point of the seismically resolved base of the channel across profiles. From this plot, it was possible to deduce the trend and gradient of the channel in a south-north progression. The channel has a shallow gradient and dips gently at about 25 km (gradient of about 0.04%) before flattening out at about 37 km (Figure 7b). In the northern extension, the PES can be sub-divided into two major pathways (P1 and P2), each subsequently sub-divided into at least three (P1) and two (P2) sub-paths (Figure 6b). The major pathway P1 (black stars) appeared to be stratigraphically deeper than the major pathway P2 (green stars) as seen in Figure 7b. The channel width increases from about 1100 m to about 1700 m in a south-north progression. However, as seen in Figure 4d, the width especially within the southern extension could be less than 300 m.

(a) Measured sinuosity index of each of the meandering loop of the Palaeo-Ems system showing a palaeo-seaward decrease in sinuosity and sub-division in different tracts. (b) Gradient of the Palaeo-Ems system reconstructed from the base of the channel deposits along the paleo-river course (i.e. the black line in Figure 6b), as well as the base of the Transparent lenses (LU) of the Elbe Palaeovalley. Note: the purple vertical lines represent ranges while the circles denote mean values. The base of the channel was better imaged from the Boomer profiles due to lower frequency and deeper penetration depth. (c) Sub-division of the Palaeo-Ems system along its paleo-flow direction based on sinuosity, gradient and channel sub-divisions.
Buried tunnel valley and straight channel structure morphology
The morphology of the two prominent buried tunnel valley generations within the study area is depicted in the interpreted isopach map presented in Figure 6d. The NE-SW valley is larger in size than the NNW-SSE valley. The NE-SW valley has a width of up to 800 m and a length of about 2.5 km while the NNW-SSE valley has an approximate mapped length of about 7.8 km and a width of about 300 m. The width and thickness of the NNW-SSE valley infill progressively become smaller from north to south (Figures 4c, d, f, and 6d). While the NE-SW valley was captured within the southernmost part of the study area at a depth below sea level of about 20 m, the NNW-SSE valley transcends the study area in an approximate north-south direction with an approximate depth below sea level of about 12 m in the south to more than 25 m further north. These depth values are based on the seismically resolved base of each of the valley generations and both valley generations extend beyond the area mapped in this study. Furthermore, the NNW-SSE valley shows a cross-cutting relationship with the uppermost part of the NE-SW valley as seen in Figure 4f. These buried tunnel valleys have overdeepened courses which is evident from the thickness of its infills and the thalweg lack a general dip.
The straight channel structure is a channel network that runs parallel to each other (Figure 6f). From a stratigraphic point of view, this structure represents a relatively young system. This channel network has a general NW-SE trend and extends beyond the study area. Within the study area, the mapped length of individual channels within the channel network varies from 2 to 4 km with width ranging between 180 and 220 m. The sinuosity index of each of the channel is about 1.01. The individual paleo-channel within this straight channel structures do not show any form of cross cutting relationship with each other (Figure 6f).
The Palaeo-Ems system stratigraphic relationship with the Elbe Palaeovalley and other geologic structures
For the first time, the PES was also analyzed in detail with respect to its spatial and chronological interaction with the older Elbe Palaeovalley (EPV). The PES constitutes one of the major tributaries toward the southwestern part of the EPV (Figure 1) where the river mouth fans out into two major and at least five minor pathways (Figure 6b). Detailed seismo-stratigraphic interpretation of profiles from the PES in the study area identified units that most likely belong to the LU2 unit of Özmaral et al.’s (2022) classification for the infills at the western flank of the EPV (Figure 2). The LU2 of Özmaral et al. (2022) and that identified in this study show consistency in terms of geometry, depth of occurrence, and seismic facies characteristics. On low-noise SBP sections, the base of both the PES and LU2 unit of the EPV within the study area are inseparable and their infills form a continuous unified depositional system (Figure 8). The acoustic unit of this unified depositional system is defined by medium-amplitude, semi-continuous reflectors at the base overlain by a thin acoustically transparent interval. The upper boundary is represented by low-amplitude, semi-continuous reflectors (Figure 8). On some of the interpreted profiles, the unified PES and LU2 seismic unit is directly overlain by a transparent interval (Figure 8). Generally, the gradient of the base of the PES was observed to be shallower than the gradient of the base of the LU2 (Figure 7b).

Sub-bottom profiler section showing the unified Palaeo-Ems system/Elbe Palaeovalley morphological relationship. See Figure 3 for location.
The spatial and stratigraphic relationship of the PES with other identified structures presented in Figure 4 is shown in the interpretative model in Figure 9. As seen in Figure 4a to c, younger structures show evidence of re-use of older structure’s pathways independent on time slices: The PES transcends parts of NNW-SSE valley (Figure 4c) while the straight channel structure incised both the lacustrine-fill unit and PES (Figure 4a and b). The oldest identified structures are the buried tunnel valleys bounded by seismic reflectors R1 (purple) and R3_1 (red) as seen in Figure 9. This was followed by the infill unit of the lacustrine-fill unit defined by reflectors R2_1/R2_2 (green) and R3_1/R4 (red/black). Next is the PES channel (base reflector R3_1; red) in which the base incised the earlier formed structures. As seen in Figure 9, the PES incised the eastern part of the lacustrine-fill unit. The PES also incised more than 15 m of the NNW-SSE buried tunnel valley over an approximate length of about 4 km. Reflector R3_2 (blue) marks the base of LU2 unit of the EPV and it is also a time equivalent with reflector R3_1 (red). Evidence for this include similar age, infill lithology and depositional settings (Table 1). Reflector R4 (black) is a key stratigraphic surface representing a period of non-deposition/erosion that separate older structures from younger ones (Figure 9). Reflectors R4 (black) and R6 (yellow) defined the Stratified unit whose areal extent is limited to the western area. The top of LU2 unit is bounded by reflector R5 (light blue) which also marks the base of the straight channel structure (Figure 9). All identified structures are overlain by laterally deposited units bounded by reflectors R4/R6 (Stratified unit), R5/R6 (Transparent unit) and R6/SF (Mobile North Sea sand) respectively. The sea floor (SF) consistently dips northward with depth to sea floor ranging from less than 15 m in the south to more than 35 m in the north.

Classification of the Palaeo-Ems system/Elbe Palaeovalley stratigraphic infill relationship.
LGM: Last Glacial Maximum; R: major seismic reflectors; SF: sea floor.
Age date was obtained from the basal peat.
Discussion
Integrated sub-bottom profiler (SBP) and core data interpretation was used to describe the overall seismo-acoustic, sedimentologic characteristics, and morphological relationships of the Palaeo-Ems system (PES) in relation to other identified geologic structures within the study area. These interpretations revealed a comprehensive picture of the PES evolution which is herein discussed.
The PES is present largely as a buried meandering channel with no morphological expression on the modern sea floor. The meandering part of the channel is highly to moderately sinuous and the meandering nature was attributed to low velocity of the flowing river and lithological composition coupled with the shallow south-north gradient. While channel gradient and sinuosity are independent morphologic elements (Miller, 1988), we observed that the meandering course of the PES is characterized by a low gradient and high sinuosity index which became less meandering (moderately sinuous) as it approaches the downstream area (Figure 7). As the PES approaches the downstream area toward the western flank of EPV, the river course sub-divided into several sub-paths (Figure 6b), forming a delta which is characterized by straight, smaller channels with sinuosity index that is slightly greater than 1 and near zero gradients.
The base of the PES channel as imaged from high-resolution Boomer profiles corresponds to Unit II in Figure 5b and represents a peat horizon. All core samples taken along the channel course described by Hepp et al. (2017) penetrated the basal peat with thickness varying between 40 and 160 cm and average of about 80 cm (Table 1). Radiocarbon ages obtained from the base and top of the peat section provides evidence for the onset of peat growth at about 11.5 cal ka BP and its termination around 9.7 cal ka BP (Hepp et al., 2019; Table 1). The basal peat consists of well-preserved plant fibers, seeds, and pollens and thus represents the pre-marine phase of the PES evolution. Typically, the basal peat is overlain by clay, silty clay, and sand of varying thickness denoted by units III, IV, and V (Figure 5). A 10–20 cm shell agglomerate was also identified between units IV and V at some locations by Hepp et al. (2019). Hepp et al. (2019) reported an upward replacement of freshwater diatoms by brackish and coastal planktonic diatoms within Unit III. This unit represents the brackish phase in the PES evolution. Units IV and V are deposited under full marine conditions characterized by almost complete absence of freshwater and brackish diatoms. The channel infills are thus deposited in an environment that changed rapidly from freshwater to brackish and finally to marine conditions over a period of about 200 years as a result of relative sea-level rise (~2.5 m based on Vink et al., 2007 sea level curve) during the early Holocene (Hepp et al., 2019). In general, the PES channel incised sub-strata Pleistocene glacial and glacio-fluvial sediments/structures (Unit I).
An interesting observation along the PES course is the change in its fluvial style from upstream to downstream. The channel course changed from a meandering style in the upstream part to a more linear and then distributive (i.e. deltaic) style in the downstream part (Figure 6b). We link this downstream trend to natural processes of a fluvial system approaching a standing body of water and thus, this might be controlled by backwater hydrodynamic effects (Colombera et al., 2016; Ganti et al., 2016).
The conceptual model presented in Figure 9 depicts the stratigraphic relationships between the various structures within the study area based on deductions from the seismic expressions in Figure 4a to f. The buried tunnel valleys (Figure 4c, d, and f) are widespread geologic features in northern Europe (Lutz et al., 2009) and basically two prominent generations with NE-SW and NNW-SSE orientations were identified within the study area (Figures 4f and 6d). The widths and thicknesses of the buried tunnel valleys which became progressively smaller in the southern part of the study area were attributed to subsequent erosional processes. These buried tunnel valleys are comparable in size, length, orientation and morphology to other Pleistocene buried tunnel valleys mapped and described by Lutz et al. (2009), Stewart et al. (2013) and Lohrberg et al. (2020) in the southern North Sea and adjacent areas. These buried tunnel valleys are the oldest identified structures and composed of mainly glacio-fluvial sands (Coughlan et al., 2018; Hepp et al., 2012) incising Neogene and late Palaeogene sediments (Lohrberg et al., 2020; Winsemann et al., 2020). While the exact processes responsible for the origin and formation of these tunnel valleys, which are widespread in the North Sea and northern Europe, remains contentious in the literature (Cofaigh, 1996; Huuse and Lykke-Andersen, 2000; van der Vegt et al., 2012; Winsemann et al., 2020), they are generally believed to occur roughly parallel to the direction of ice flow and serve as key conduits for meltwater discharge underneath the continental ice sheets (Özmaral et al., 2022; Stewart et al., 2013). The base of the buried tunnel valleys is undulating and lack a general dip across profiles, these characteristic features clearly distinguished tunnel valleys from fluvial channels which have a hydrodynamic flow and gradient (Lonergan et al., 2006). The PES incised about 15 m into the buried tunnel valleys (Figure 4c) emphasizing a relatively younger age (early Holocene) for the PES and possibly re-use of older buried tunnel valley pathways by the PES in some parts. Such re-use of older structures by younger glacial valleys or rivers which are usually expressed on seismic sections as complex valley infill architecture of multiple cut and fill structures are well documented in the North Sea (Lutz et al., 2009; Prins and Andresen, 2019; Stewart et al., 2013 amongst others). The low gradient of the PES channel suggests a low velocity for the river with a moderate discharge resulting in channel cut-offs preserved as oxbow lakes. Further northwest of the PES area, similar fluvial structures are known from the western part of the Dogger Bank as well as the north-eastern part of the EPV (Emery et al., 2020; Fitch et al., 2005; Hepp et al., 2017). These structures exhibit close morphological similarities and time of development as the PES. However, the advancing sea level of the North Sea resulted in an earlier drowning of these more northern river systems.
West of the PES channel, the lacustrine-fill unit was also identified and mapped as a characteristic laterally stratified sediment succession. Careful examination of the margin of this structure revealed that the PES cuts the eastern side of this structure (Figure 4a–c). This again shows that the PES channel is relatively younger than the lacustrine-fill unit. Seismic facies analyses of the PES infill revealed a minimum of a two-phase infill made up of peat-clay/silt sequences (Hepp et al., 2019). Reflector R4, a key stratigraphic surface represents a period of erosion and/or non-deposition (hiatus). The straight channel structure is the last identified structure within the study area interpreted as younger fluvial river channel networks. The characteristic progradational infill represents incisions infill of older sediments and sedimentary structures. In fact, the uppermost infill of the PES and lacustrine-fill unit are incised by these fluvial river channels (Figure 4a and b). This indicates that the fluvial river channels are younger than the buried tunnel valleys, lacustrine-fill unit and the PES which are all present within the study area. The fluvial river channels also occur as discrete structures with a general NW-SE trend (Figures 4e and 6f). All identified and mapped geological structures within the study area are overlain by laterally deposited fine sand and silt separated from the overlying fine to medium mobile sand by reflector R6. This mobile sand is widely distributed throughout the North Sea (Zeiler et al., 2000) and are referred to as the Upper Marine Formation by Coughlan et al. (2018).
Although studies by Papenmeier and Hass (2020) and Özmaral et al. (2022) detailed the evolution and stratigraphic infill of the EPV, this study for the first-time shed light on the PES/EPV morphologic and stratigraphic relationship. The EPV, a now mostly buried geomorphological shallow trough, was the main drainage pathway for meltwater from the southern margin of the Scandinavian ice sheet as well as rivers draining the North European plain (Figge, 1983; Özmaral et al., 2022; Papenmeier and Hass, 2020). The PES on the other hand has been shown to be one of the major tributaries of this approximately 210 km long EPV situated in the south-eastern North Sea (Hepp et al., 2017). While the PES played a major role in the drainage network of the hinterland areas, its drainage pathways and morpho-stratigraphic relationship at the western margin of the EPV remain unknown. Based on this study, a detailed interpretation of high-resolution SBP sections west of the EPV revealed that the PES spatially sub-divided into two major pathways (P1 and P2) with each pathway branching out into smaller paths as it approaches the western rim of the EPV (Figure 6b). Here, the branching of PES formed a delta with linear courses as opposed to the meandering course in its upstream area. This linear nature was attributed to a reduced velocity of the flowing river, erodibility nature of the substrate and/or a sudden change in slope or topography.
While the EPV itself is older than the PES, detailed analysis of SBP sections from the western margin of the EPV revealed that the PES formed a unified depositional system with the EPV. In fact, the early infill sequence of the PES was observed to coalesce with the LU unit of the Özmaral et al. (2022) EPV classification. Sediment cores revealed the presence of peat at the base of the PES (Figure 5) and EPV LU2 unit (Figure 2; Özmaral et al., 2022) which are deposited in a wetland setting. Radiocarbon dates by Hepp et al. (2019) and Özmaral et al. (2022) on these basal peats revealed similar ages in the range of 9.6–11.5 cal ka BP (Table 1). These similar ages corroborate the assumption that both the PES and LU unit are contemporaneous. In addition, the similar infill lithology of silty clay and/or clayey silt and sand also points toward a simultaneous deposition (Table 1). This implies that the PES constitutes not only part of the paleo-drainage system that emptied into the EPV but also contributes to the shaping of the second phase (LU unit) of the EPV evolution at the earliest.
The PES /EPV relationship can thus be explained by understanding their infill history and morphological patterns. The first phase of the EPV formation which is represented by the FU unit (Figure 2; Özmaral et al., 2022) is a classic ice marginal valley infill formed as a result of water collected from the retreating ice. This FU unit formed a braided river system at the base of the EPV with characteristic morphological terraces during periods of low sea level (Özmaral et al., 2022). Post-glacial relative sea-level rise led to a retreat in coastline and the EPV transformed from a valley to an in-channel aggradation deposition observed as LU unit (Özmaral et al., 2022). In-channel aggradation as a result of sea-level rise (Martin et al., 2009; Törnqvist, 1994) can result in super-elevation of flowing river relative to the surrounding floodplain (Mohrig et al., 2000; Slingerland and Smith, 2004). The EPV drainage networks and adjacent areas thus adapted to the new hydrodynamic regime during the post-glacial relative sea-level rise (Emery et al., 2020; Özmaral et al., 2022) and maintained elevation relative to the sea level. This also possibly aided the changes in Palaeo-Ems River course resulting in different sub-paths as observed at the western flank of the EPV (Figure 6b). The different PES pathways formed a deltaic situation west of the EPV which is evident in the near-zero gradient seen on the gradient plot of the sub pathways (green and black symbols in Figure 7b). Morphological terraces initiated during the incision stage of the EPV are accompanied by valley widening during the subsequent valley infill stages (Özmaral et al., 2022). Such back-cutting and widening of the EPV are believed to be responsible for the subsequent capture of the Palaeo-Ems River by the EPV with the infill of the PES corresponding with the last phase of LU unit (LU2) of Özmaral et al. (2022) classification. The steeper gradient of LU (yellow and purple symbols in Figure 7b) compared to that of the PES is also believed to be as a result of the back-cutting and widening of the EPV. The infill of both the PES and LU unit represents a last stage deposition prior to a full marine condition thus indicating a relatively young age for the joining event. Both the EPV and Palaeo-Ems systems were later drowned as a result of the fast-rising relative sea level overwhelming the adaptation capabilities of the joint drainage system. Similar modern-day analog system displaying a delta filling phenomenon within a larger valley is also visible around the Rio de la Plata near Buenos Aires, Argentina. Here, the Paraná River delta was captured by the valley of the Rio de la Plata. While the PES/EPV and Paraná River/Rio de la Plata systems differ in various aspects, including their origin and climatic conditions, both systems show similarities in terms of their dimensions and morphology. We thereby believe that this modern analog can at least partially explain the evolution of the Palaeo-Ems River and its stratigraphic relationship with the EPV during the Quaternary.
Furthermore, given the differing topographic and major river inputs of both the PES and paleo-drainage networks archived in the Dogger landscape, comparisons were made between both systems as they are both transgressed by the early Holocene relative sea-level rise in the North Sea. Fitch et al. (2005), Bailey et al. (2017), Gearey et al. (2017), Hepp et al. (2017), Emery et al. (2020), Gaffney and Fitch (2022) amongst others gave a comprehensive seismo-stratigraphic description of river channel networks from the Dogger Bank area. These river channels have similar evolution, morphological and sedimentological characteristics as the PES. Core analyses of one of the paleo-channel described by Gearey et al. (2017) revealed a phase of basal peat formation during the early Holocene (age: 9.2–10.2 cal ka BP), upper core silt unit dated 7.0 cal ka BP and top core marine shell-rich sands. Similar to the PES, the Dogger Bank River channel networks are interpreted to be fluvial in origin and their stratigraphic infills are deposited in sedimentary environments that also changes from freshwater to brackish and then to full marine (Emery et al., 2020; Gearey et al., 2017). Both river systems thus experienced similar paleo-climatic conditions during the early Holocene prior to marine transgression at about 8 cal ka BP.
Conclusion
This study provides a new insight into the overall course of the PES as well as its stratigraphic relationship with other geologic structures, particularly the EPV, using high-resolution SBP sections. The overall morphology of the PES was dictated by the topography of the area, underlying structures, intensity of fluvial activities as well as the relative timing of these events. In a south-north progression, the PES occurred in most parts as a meandering river and became more linear toward the western margin of the southern part of the EPV. Over the course of its length, the river showed various abandoned pathways which supported intense fluvial activities at the time of formation. In addition, the PES incised and re-used the pathway of the older structure in some areas. At the western margin of the southern EPV, the river branched into two major and at least five minor sub-paths forming a large delta. It was deduced, based on available data for this study, that the younger PES connected to the older pre-existing EPV structure and formed a unified depositional system with the last stage of the EPV LU interval (LU2). This indicates that the PES played a key role in shaping the early phases of the EPV evolution. Other mappable geologic structures include buried tunnel valleys, lacustrine-fill unit and fluvial river channels. In terms of its stratigraphic infill, the PES is made up of a polyphased peat and clay/silt sequence. Above this PES infill are laterally deposited silt and fine sand which are subsequently overlain by the fine- to medium-grained mobile sands of the Upper Marine Formation.
The relatively low gradient of the PES was considered as the key factor responsible for its drowning as the river was unable to adjust to the prevailing relative rise in sea level which completely flooded the paleo-landscape before 8.2 cal ka BP. Nevertheless, the PES, a drowned predecessor of the modern Ems River, played a major role in the paleo-drainage network of the hinterland areas and constituted one of the major tributaries in the southern part of the EPV. Together with other tributaries of the EPV, the PES formed part of the drainage system that developed in the area since the Last Glacial Maximum and it is also significant to the early Holocene coastal landscape development of the German North Sea sector.
